ROMS
Loading...
Searching...
No Matches
mod_diags.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE mod_diags
3#ifdef DIAGNOSTICS
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! avgzeta Time-averaged free surface !
13! !
14! Diagnostics fields for output. !
15! !
16! DiaBio2d Diagnostics for 2D biological terms. !
17! DiaBio3d Diagnostics for 3D biological terms. !
18! DiaBio4d Diagnostics for 4D bio-optical terms. !
19! DiaTrc Diagnostics for tracer terms. !
20! DiaU2d Diagnostics for 2D U-momentum terms. !
21! DiaV2d Diagnostics for 2D V-momentum terms. !
22! DiaU3d Diagnostics for 3D U-momentum terms. !
23! DiaV3d Diagnostics for 3D V-momentum terms. !
24! !
25! Diagnostics fields work arrays. !
26! !
27! DiaTwrk Diagnostics work array for tracer terms. !
28! DiaU2wrk Diagnostics work array for 2D U-momentum terms. !
29! DiaV2wrk Diagnostics work array for 2D V-momentum terms. !
30! DiaRUbar Diagnostics RHS array for 2D U-momentum terms. !
31! DiaRVbar Diagnostics RHS array for 2D V-momentum terms. !
32! DiaU2int Diagnostics array for 2D U-momentum terms !
33! integrated over short barotropic timesteps. !
34! DiaV2int Diagnostics array for 2D U-momentum terms !
35! integrated over short barotropic timesteps. !
36! DiaRUfrc Diagnostics forcing array for 2D U-momentum terms. !
37! DiaRVfrc Diagnostics forcing array for 2D V-momentum terms. !
38! DiaU3wrk Diagnostics work array for 3D U-momentum terms. !
39! DiaV3wrk Diagnostics work array for 3D V-momentum terms. !
40! DiaRU Diagnostics RHS array for 3D U-momentum terms. !
41! DiaRV Diagnostics RHS array for 3D V-momentum terms. !
42! !
43!=======================================================================
44!
45 USE mod_kinds
46!
47 implicit none
48!
49 PUBLIC :: allocate_diags
50 PUBLIC :: deallocate_diags
51 PUBLIC :: initialize_diags
52!
53!-----------------------------------------------------------------------
54! Define T_DIAGS structure.
55!-----------------------------------------------------------------------
56!
58
59 real(r8), pointer :: avgzeta(:,:)
60
61# ifdef DIAGNOSTICS_BIO
62# if defined BIO_FENNEL
63 real(r8), pointer :: diabio2d(:,:,:)
64 real(r8), pointer :: diabio3d(:,:,:,:)
65# elif defined ECOSIM
66 real(r8), pointer :: diabio3d(:,:,:,:)
67 real(r8), pointer :: diabio4d(:,:,:,:,:)
68# endif
69# endif
70
71# ifdef DIAGNOSTICS_TS
72 real(r8), pointer :: diatrc(:,:,:,:,:)
73 real(r8), pointer :: diatwrk(:,:,:,:,:)
74# endif
75
76# ifdef DIAGNOSTICS_UV
77 real(r8), pointer :: diau2d(:,:,:)
78 real(r8), pointer :: diav2d(:,:,:)
79 real(r8), pointer :: diau2wrk(:,:,:)
80 real(r8), pointer :: diav2wrk(:,:,:)
81 real(r8), pointer :: diarubar(:,:,:,:)
82 real(r8), pointer :: diarvbar(:,:,:,:)
83
84# ifdef SOLVE3D
85 real(r8), pointer :: diau2int(:,:,:)
86 real(r8), pointer :: diav2int(:,:,:)
87 real(r8), pointer :: diarufrc(:,:,:,:)
88 real(r8), pointer :: diarvfrc(:,:,:,:)
89 real(r8), pointer :: diau3d(:,:,:,:)
90 real(r8), pointer :: diav3d(:,:,:,:)
91 real(r8), pointer :: diau3wrk(:,:,:,:)
92 real(r8), pointer :: diav3wrk(:,:,:,:)
93 real(r8), pointer :: diaru(:,:,:,:,:)
94 real(r8), pointer :: diarv(:,:,:,:,:)
95# endif
96# endif
97
98 END TYPE t_diags
99!
100 TYPE (t_diags), allocatable :: diags(:)
101!
102 CONTAINS
103!
104 SUBROUTINE allocate_diags (ng, LBi, UBi, LBj, UBj)
105!
106!=======================================================================
107! !
108! This routine allocates all variables in the module for all nested !
109! grids. !
110! !
111!=======================================================================
112!
113 USE mod_param
114# if defined DIAGNOSTICS_BIO && defined ECOSIM
115 USE mod_biology
116# endif
117!
118! Imported variable declarations.
119!
120 integer, intent(in) :: ng, lbi, ubi, lbj, ubj
121!
122! Local variable declarations.
123!
124 real(r8) :: size2d
125!
126!-----------------------------------------------------------------------
127! Allocate module variables.
128!-----------------------------------------------------------------------
129!
130 IF (ng.eq.1 ) allocate ( diags(ngrids) )
131!
132! Set horizontal array size.
133!
134 size2d=real((ubi-lbi+1)*(ubj-lbj+1),r8)
135!
136! Diagnostic arrays.
137!
138 allocate ( diags(ng) % avgzeta(lbi:ubi,lbj:ubj) )
139 dmem(ng)=dmem(ng)+size2d
140
141# ifdef DIAGNOSTICS_BIO
142# if defined BIO_FENNEL
143 IF (ndbio2d.gt.0) THEN
144 allocate ( diags(ng) % DiaBio2d(lbi:ubi,lbj:ubj,ndbio2d) )
145 dmem(ng)=dmem(ng)+real(ndbio2d,r8)*size2d
146 END IF
147 IF (ndbio3d.gt.0) THEN
148 allocate ( diags(ng) % DiaBio3d(lbi:ubi,lbj:ubj,n(ng),ndbio3d) )
149 dmem(ng)=dmem(ng)+real(n(ng)*ndbio3d,r8)*size2d
150 END IF
151# elif defined ECOSIM
152 IF (ndbio3d.gt.0) THEN
153 allocate ( diags(ng) % DiaBio3d(lbi:ubi,lbj:ubj, &
154 & ndbands,ndbio3d) )
155 dmem(ng)=dmem(ng)+real(ndbands*ndbio3d,r8)*size2d
156 END IF
157 IF (ndbio4d.gt.0) THEN
158 allocate ( diags(ng) % DiaBio4d(lbi:ubi,lbj:ubj,n(ng), &
159 & ndbands,ndbio4d) )
160 dmem(ng)=dmem(ng)+real(n(ng)*ndbands*ndbio4d,r8)*size2d
161 END IF
162# endif
163# endif
164
165# ifdef DIAGNOSTICS_TS
166 allocate ( diags(ng) % DiaTrc(lbi:ubi,lbj:ubj,n(ng),nt(ng),ndt) )
167 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng)*ndt,r8)*size2d
168
169 allocate ( diags(ng) % DiaTwrk(lbi:ubi,lbj:ubj,n(ng),nt(ng),ndt) )
170 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng)*ndt,r8)*size2d
171# endif
172
173# ifdef DIAGNOSTICS_UV
174 allocate ( diags(ng) % DiaU2d(lbi:ubi,lbj:ubj,ndm2d) )
175 dmem(ng)=dmem(ng)+real(ndm2d,r8)*size2d
176
177 allocate ( diags(ng) % DiaV2d(lbi:ubi,lbj:ubj,ndm2d) )
178 dmem(ng)=dmem(ng)+real(ndm2d,r8)*size2d
179
180 allocate ( diags(ng) % DiaU2wrk(lbi:ubi,lbj:ubj,ndm2d) )
181 dmem(ng)=dmem(ng)+real(ndm2d,r8)*size2d
182
183 allocate ( diags(ng) % DiaV2wrk(lbi:ubi,lbj:ubj,ndm2d) )
184 dmem(ng)=dmem(ng)+real(ndm2d,r8)*size2d
185
186 allocate ( diags(ng) % DiaRUbar(lbi:ubi,lbj:ubj,2,ndm2d-1) )
187 dmem(ng)=dmem(ng)+2.0_r8*real(ndm2d-1,r8)*size2d
188
189 allocate ( diags(ng) % DiaRVbar(lbi:ubi,lbj:ubj,2,ndm2d-1) )
190 dmem(ng)=dmem(ng)+2.0_r8*real(ndm2d-1,r8)*size2d
191
192# ifdef SOLVE3D
193 allocate ( diags(ng) % DiaU2int(lbi:ubi,lbj:ubj,ndm2d) )
194 dmem(ng)=dmem(ng)+real(ndm2d,r8)*size2d
195
196 allocate ( diags(ng) % DiaV2int(lbi:ubi,lbj:ubj,ndm2d) )
197 dmem(ng)=dmem(ng)+real(ndm2d,r8)*size2d
198
199 allocate ( diags(ng) % DiaRUfrc(lbi:ubi,lbj:ubj,3,ndm2d-1) )
200 dmem(ng)=dmem(ng)+3.0_r8*real(ndm2d-1,r8)*size2d
201
202 allocate ( diags(ng) % DiaRVfrc(lbi:ubi,lbj:ubj,3,ndm2d-1) )
203 dmem(ng)=dmem(ng)+3.0_r8*real(ndm2d-1,r8)*size2d
204
205 allocate ( diags(ng) % DiaU3d(lbi:ubi,lbj:ubj,n(ng),ndm3d) )
206 dmem(ng)=dmem(ng)+real(n(ng)*ndm3d,r8)*size2d
207
208 allocate ( diags(ng) % DiaV3d(lbi:ubi,lbj:ubj,n(ng),ndm3d) )
209 dmem(ng)=dmem(ng)+real(n(ng)*ndm3d,r8)*size2d
210
211 allocate ( diags(ng) % DiaU3wrk(lbi:ubi,lbj:ubj,n(ng),ndm3d) )
212 dmem(ng)=dmem(ng)+real(n(ng)*ndm3d,r8)*size2d
213
214 allocate ( diags(ng) % DiaV3wrk(lbi:ubi,lbj:ubj,n(ng),ndm3d) )
215 dmem(ng)=dmem(ng)+real(n(ng)*ndm3d,r8)*size2d
216
217 allocate ( diags(ng) % DiaRU(lbi:ubi,lbj:ubj,n(ng),2,ndrhs) )
218 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)*ndrhs,r8)*size2d
219
220 allocate ( diags(ng) % DiaRV(lbi:ubi,lbj:ubj,n(ng),2,ndrhs) )
221 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)*ndrhs,r8)*size2d
222
223# endif
224
225# endif
226!
227 RETURN
228 END SUBROUTINE allocate_diags
229!
230 SUBROUTINE deallocate_diags (ng)
231!
232!=======================================================================
233! !
234! This routine deallocates all variables in the module for all nested !
235! grids. !
236! !
237!=======================================================================
238!
239 USE mod_param
240
241# ifdef SUBOBJECT_DEALLOCATION
242!
243 USE destroy_mod, ONLY : destroy
244# endif
245!
246! Imported variable declarations.
247!
248 integer, intent(in) :: ng
249!
250! Local variable declarations.
251!
252 character (len=*), parameter :: myfile = &
253 & __FILE__//", deallocate_diags"
254
255# ifdef SUBOBJECT_DEALLOCATION
256!
257!-----------------------------------------------------------------------
258! Deallocate each variable in the derived-type T_DIAGS structure
259! separately.
260!-----------------------------------------------------------------------
261!
262! Diagnostic arrays.
263!
264 IF (.not.destroy(ng, diags(ng)%avgzeta, myfile, &
265 & __line__, 'DIAGS(ng)%DiaBio2d')) RETURN
266
267# ifdef DIAGNOSTICS_BIO
268# if defined BIO_FENNEL
269 IF (ndbio2d.gt.0) THEN
270 IF (.not.destroy(ng, diags(ng)%DiaBio2d, myfile, &
271 & __line__, 'DIAGS(ng)%DiaBio2d')) RETURN
272 END IF
273
274 IF (ndbio3d.gt.0) THEN
275 IF (.not.destroy(ng, diags(ng)%DiaBio3d, myfile, &
276 & __line__, 'DIAGS(ng)%DiaBio3d')) RETURN
277 END IF
278# elif defined ECOSIM
279 IF (ndbio3d.gt.0) THEN
280 IF (.not.destroy(ng, diags(ng)%DiaBio3d, myfile, &
281 & __line__, 'DIAGS(ng)%DiaBio3d')) RETURN
282 END IF
283
284 IF (ndbio4d.gt.0) THEN
285 IF (.not.destroy(ng, diags(ng)%DiaBio4d, myfile, &
286 & __line__, 'DIAGS(ng)%DiaBio4d')) RETURN
287 END IF
288# endif
289# endif
290
291# ifdef DIAGNOSTICS_TS
292 IF (.not.destroy(ng, diags(ng)%DiaTrc, myfile, &
293 & __line__, 'DIAGS(ng)%DiaTrc')) RETURN
294
295 IF (.not.destroy(ng, diags(ng)%DiaTwrk, myfile, &
296 & __line__, 'DIAGS(ng)%DiaTwrk')) RETURN
297# endif
298
299# ifdef DIAGNOSTICS_UV
300 IF (.not.destroy(ng, diags(ng)%DiaU2d, myfile, &
301 & __line__, 'DIAGS(ng)%DiaU2d')) RETURN
302
303 IF (.not.destroy(ng, diags(ng)%DiaV2d, myfile, &
304 & __line__, 'DIAGS(ng)%DiaV2d')) RETURN
305
306 IF (.not.destroy(ng, diags(ng)%DiaU2wrk, myfile, &
307 & __line__, 'DIAGS(ng)%DiaU2wrk')) RETURN
308
309 IF (.not.destroy(ng, diags(ng)%DiaV2wrk, myfile, &
310 & __line__, 'DIAGS(ng)%DiaV2wrk')) RETURN
311
312 IF (.not.destroy(ng, diags(ng)%DiaRUbar, myfile, &
313 & __line__, 'DIAGS(ng)%DiaRUbar')) RETURN
314
315 IF (.not.destroy(ng, diags(ng)%DiaRVbar, myfile, &
316 & __line__, 'DIAGS(ng)%DiaRVbar')) RETURN
317
318# ifdef SOLVE3D
319 IF (.not.destroy(ng, diags(ng)%DiaU2int, myfile, &
320 & __line__, 'DIAGS(ng)%DiaU2int')) RETURN
321
322 IF (.not.destroy(ng, diags(ng)%DiaV2int, myfile, &
323 & __line__, 'DIAGS(ng)%DiaV2int')) RETURN
324
325 IF (.not.destroy(ng, diags(ng)%DiaRUfrc, myfile, &
326 & __line__, 'DIAGS(ng)%DiaRUfrc')) RETURN
327
328 IF (.not.destroy(ng, diags(ng)%DiaRVfrc, myfile, &
329 & __line__, 'DIAGS(ng)%DiaRVfrc')) RETURN
330
331 IF (.not.destroy(ng, diags(ng)%DiaU3d, myfile, &
332 & __line__, 'DIAGS(ng)%DiaU3d')) RETURN
333
334 IF (.not.destroy(ng, diags(ng)%DiaV3d, myfile, &
335 & __line__, 'DIAGS(ng)%DiaV3d')) RETURN
336
337 IF (.not.destroy(ng, diags(ng)%DiaU3wrk, myfile, &
338 & __line__, 'DIAGS(ng)%DiaU3wrk')) RETURN
339
340 IF (.not.destroy(ng, diags(ng)%DiaV3wrk, myfile, &
341 & __line__, 'DIAGS(ng)%DiaV3wrk')) RETURN
342
343 IF (.not.destroy(ng, diags(ng)%DiaRU, myfile, &
344 & __line__, 'DIAGS(ng)%DiaRU')) RETURN
345
346 IF (.not.destroy(ng, diags(ng)%DiaRV, myfile, &
347 & __line__, 'DIAGS(ng)%DiaRV')) RETURN
348# endif
349# endif
350# endif
351!
352!-----------------------------------------------------------------------
353! Deallocate derived-type DIAGS structure.
354!-----------------------------------------------------------------------
355!
356 IF (ng.eq.ngrids) THEN
357 IF (allocated(diags)) deallocate ( diags )
358 END IF
359!
360 RETURN
361 END SUBROUTINE deallocate_diags
362!
363 SUBROUTINE initialize_diags (ng, tile)
364!
365!=======================================================================
366! !
367! This routine initialize all variables in the module using first !
368! touch distribution policy. In shared-memory configuration, this !
369! operation actually performs propagation of the "shared arrays" !
370! across the cluster, unless another policy is specified to !
371! override the default. !
372! !
373!=======================================================================
374!
375 USE mod_param
376# if defined DIAGNOSTICS_BIO && defined ECOSIM
377 USE mod_biology
378# endif
379!
380! Imported variable declarations.
381!
382 integer, intent(in) :: ng, tile
383!
384! Local variable declarations.
385!
386 integer :: imin, imax, jmin, jmax
387 integer :: i, idiag, j
388# ifdef SOLVE3D
389 integer :: itrc, iband, k
390# endif
391
392 real(r8), parameter :: inival = 0.0_r8
393
394# include "set_bounds.h"
395!
396! Set array initialization range.
397!
398# ifdef DISTRIBUTE
399 imin=bounds(ng)%LBi(tile)
400 imax=bounds(ng)%UBi(tile)
401 jmin=bounds(ng)%LBj(tile)
402 jmax=bounds(ng)%UBj(tile)
403# else
404 IF (domain(ng)%Western_Edge(tile)) THEN
405 imin=bounds(ng)%LBi(tile)
406 ELSE
407 imin=istr
408 END IF
409 IF (domain(ng)%Eastern_Edge(tile)) THEN
410 imax=bounds(ng)%UBi(tile)
411 ELSE
412 imax=iend
413 END IF
414 IF (domain(ng)%Southern_Edge(tile)) THEN
415 jmin=bounds(ng)%LBj(tile)
416 ELSE
417 jmin=jstr
418 END IF
419 IF (domain(ng)%Northern_Edge(tile)) THEN
420 jmax=bounds(ng)%UBj(tile)
421 ELSE
422 jmax=jend
423 END IF
424# endif
425!
426!-----------------------------------------------------------------------
427! Initialize module variables.
428!-----------------------------------------------------------------------
429!
430 DO j=jmin,jmax
431 DO i=imin,imax
432 diags(ng) % avgzeta(i,j) = inival
433 END DO
434
435# ifdef DIAGNOSTICS_BIO
436# if defined BIO_FENNEL
437 IF (ndbio2d.gt.0) THEN
438 DO idiag=1,ndbio2d
439 DO i=imin,imax
440 diags(ng) % DiaBio2d(i,j,idiag) = inival
441 END DO
442 END DO
443 END IF
444 IF (ndbio3d.gt.0) THEN
445 DO idiag=1,ndbio3d
446 DO k=1,n(ng)
447 DO i=imin,imax
448 diags(ng) % DiaBio3d(i,j,k,idiag) = inival
449 END DO
450 END DO
451 END DO
452 END IF
453# elif defined ECOSIM
454 IF (ndbio3d.gt.0) THEN
455 DO idiag=1,ndbio3d
456 DO k=1,ndbands
457 DO i=imin,imax
458 diags(ng) % DiaBio3d(i,j,k,idiag) = inival
459 END DO
460 END DO
461 END DO
462 END IF
463 IF (ndbio4d.gt.0) THEN
464 DO idiag=1,ndbio4d
465 DO iband=1,ndbands
466 DO k=1,n(ng)
467 DO i=imin,imax
468 diags(ng) % DiaBio4d(i,j,k,iband,idiag) = inival
469 END DO
470 END DO
471 END DO
472 END DO
473 END IF
474# endif
475# endif
476# ifdef DIAGNOSTICS_TS
477 DO idiag=1,ndt
478 DO itrc=1,nt(ng)
479 DO k=1,n(ng)
480 DO i=imin,imax
481 diags(ng) % DiaTrc(i,j,k,itrc,idiag) = inival
482 diags(ng) % DiaTwrk(i,j,k,itrc,idiag) = inival
483 END DO
484 END DO
485 END DO
486 END DO
487# endif
488# ifdef DIAGNOSTICS_UV
489 DO idiag=1,ndm2d
490 DO i=imin,imax
491 diags(ng) % DiaU2d(i,j,idiag) = inival
492 diags(ng) % DiaV2d(i,j,idiag) = inival
493 diags(ng) % DiaU2wrk(i,j,idiag) = inival
494 diags(ng) % DiaV2wrk(i,j,idiag) = inival
495# ifdef SOLVE3D
496 diags(ng) % DiaU2int(i,j,idiag) = inival
497 diags(ng) % DiaV2int(i,j,idiag) = inival
498# endif
499 END DO
500 END DO
501 DO idiag=1,ndm2d-1
502 DO i=imin,imax
503 diags(ng) % DiaRUbar(i,j,1,idiag) = inival
504 diags(ng) % DiaRUbar(i,j,2,idiag) = inival
505 diags(ng) % DiaRVbar(i,j,1,idiag) = inival
506 diags(ng) % DiaRVbar(i,j,2,idiag) = inival
507# ifdef SOLVE3D
508 diags(ng) % DiaRUfrc(i,j,1,idiag) = inival
509 diags(ng) % DiaRUfrc(i,j,2,idiag) = inival
510 diags(ng) % DiaRUfrc(i,j,3,idiag) = inival
511 diags(ng) % DiaRVfrc(i,j,1,idiag) = inival
512 diags(ng) % DiaRVfrc(i,j,2,idiag) = inival
513 diags(ng) % DiaRVfrc(i,j,3,idiag) = inival
514# endif
515 END DO
516 END DO
517# ifdef SOLVE3D
518 DO idiag=1,ndm3d
519 DO k=1,n(ng)
520 DO i=imin,imax
521 diags(ng) % DiaU3d(i,j,k,idiag) = inival
522 diags(ng) % DiaV3d(i,j,k,idiag) = inival
523 diags(ng) % DiaU3wrk(i,j,k,idiag) = inival
524 diags(ng) % DiaV3wrk(i,j,k,idiag) = inival
525 END DO
526 END DO
527 END DO
528 DO idiag=1,ndrhs
529 DO k=1,n(ng)
530 DO i=imin,imax
531 diags(ng) % DiaRU(i,j,k,1,idiag) = inival
532 diags(ng) % DiaRU(i,j,k,2,idiag) = inival
533 diags(ng) % DiaRV(i,j,k,1,idiag) = inival
534 diags(ng) % DiaRV(i,j,k,2,idiag) = inival
535 END DO
536 END DO
537 END DO
538# endif
539# endif
540 END DO
541!
542 RETURN
543 END SUBROUTINE initialize_diags
544#endif
545 END MODULE mod_diags
integer ndbands
Definition ecosim_mod.h:212
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
subroutine, public initialize_diags(ng, tile)
Definition mod_diags.F:364
subroutine, public allocate_diags(ng, lbi, ubi, lbj, ubj)
Definition mod_diags.F:105
subroutine, public deallocate_diags(ng)
Definition mod_diags.F:231
integer, parameter r8
Definition mod_kinds.F:28
integer ndrhs
Definition mod_param.F:590
integer ndbio2d
Definition mod_param.F:584
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
real(r8), dimension(:), allocatable dmem
Definition mod_param.F:137
integer ndbio4d
Definition mod_param.F:586
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer ndm3d
Definition mod_param.F:579
integer ngrids
Definition mod_param.F:113
integer ndt
Definition mod_param.F:574
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer ndbio3d
Definition mod_param.F:585
integer ndm2d
Definition mod_param.F:578