14 use mom_diag_mediator,
only : diag_copy_storage_to_diag, diag_save_grids, diag_restore_grids
25 implicit none ;
private
27 #include <MOM_memory.h>
29 public register_tracer
30 public mom_tracer_chksum, mom_tracer_chkinv
31 public register_tracer_diagnostics, post_tracer_diagnostics, post_tracer_transport_diagnostics
32 public preale_tracer_diagnostics, postale_tracer_diagnostics
33 public tracer_registry_init, lock_tracer_registry, tracer_registry_end
34 public tracer_name_lookup
39 real,
dimension(:,:,:),
pointer :: t => null()
46 real,
dimension(:,:,:),
pointer :: ad_x => null()
48 real,
dimension(:,:,:),
pointer :: ad_y => null()
50 real,
dimension(:,:),
pointer :: ad2d_x => null()
52 real,
dimension(:,:),
pointer :: ad2d_y => null()
55 real,
dimension(:,:,:),
pointer :: df_x => null()
57 real,
dimension(:,:,:),
pointer :: df_y => null()
59 real,
dimension(:,:),
pointer :: df2d_x => null()
61 real,
dimension(:,:),
pointer :: df2d_y => null()
63 real,
dimension(:,:),
pointer :: df2d_conc_x => null()
65 real,
dimension(:,:),
pointer :: df2d_conc_y => null()
68 real,
dimension(:,:,:),
pointer :: advection_xy => null()
70 real,
dimension(:,:,:),
pointer :: diff_cont_xy => null()
72 real,
dimension(:,:,:),
pointer :: diff_conc_xy => null()
74 real,
dimension(:,:,:),
pointer :: t_prev => null()
76 real,
dimension(:,:,:),
pointer :: trxh_prev => null()
79 character(len=32) :: name
80 character(len=64) :: units
81 character(len=240) :: longname
83 logical :: registry_diags = .false.
85 character(len=64) :: cmor_name
86 character(len=64) :: cmor_units
87 character(len=240) :: cmor_longname
88 character(len=32) :: flux_nameroot =
""
90 character(len=64) :: flux_longname =
""
92 real :: flux_scale= 1.0
94 character(len=48) :: flux_units =
""
95 character(len=48) :: conv_units =
""
96 real :: conv_scale = 1.0
98 character(len=48) :: cmor_tendprefix =
""
101 integer :: ind_tr_squared = -1
104 logical :: advect_tr = .true.
105 logical :: hordiff_tr = .true.
106 logical :: remap_tr = .true.
108 integer :: diag_form = 1
110 integer :: id_tr = -1
111 integer :: id_adx = -1, id_ady = -1, id_dfx = -1, id_dfy = -1
112 integer :: id_adx_2d = -1, id_ady_2d = -1, id_dfx_2d = -1, id_dfy_2d = -1
113 integer :: id_adv_xy = -1, id_adv_xy_2d = -1
114 integer :: id_dfxy_cont = -1, id_dfxy_cont_2d = -1, id_dfxy_conc = -1
115 integer :: id_remap_conc = -1, id_remap_cont = -1, id_remap_cont_2d = -1
116 integer :: id_tendency = -1, id_trxh_tendency = -1, id_trxh_tendency_2d = -1
117 integer :: id_tr_vardec = -1
126 logical :: locked = .false.
135 subroutine register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, &
136 cmor_name, cmor_units, cmor_longname, tr_desc, &
137 OBC_inflow, OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, &
138 ad_2d_x, ad_2d_y, df_2d_x, df_2d_y, advection_xy, registry_diags, &
139 flux_nameroot, flux_longname, flux_units, flux_scale, &
140 convergence_units, convergence_scale, cmor_tendprefix, diag_form, &
141 restart_CS, mandatory)
145 real,
dimension(SZI_(HI),SZJ_(HI),SZK_(GV)), &
148 character(len=*),
optional,
intent(in) :: name
149 character(len=*),
optional,
intent(in) :: longname
150 character(len=*),
optional,
intent(in) :: units
151 character(len=*),
optional,
intent(in) :: cmor_name
152 character(len=*),
optional,
intent(in) :: cmor_units
153 character(len=*),
optional,
intent(in) :: cmor_longname
154 type(
vardesc),
optional,
intent(in) :: tr_desc
156 real,
optional,
intent(in) :: obc_inflow
158 real,
dimension(:,:,:),
optional,
pointer :: obc_in_u
160 real,
dimension(:,:,:),
optional,
pointer :: obc_in_v
164 real,
dimension(:,:,:),
optional,
pointer :: ad_x
165 real,
dimension(:,:,:),
optional,
pointer :: ad_y
166 real,
dimension(:,:,:),
optional,
pointer :: df_x
167 real,
dimension(:,:,:),
optional,
pointer :: df_y
168 real,
dimension(:,:),
optional,
pointer :: ad_2d_x
170 real,
dimension(:,:),
optional,
pointer :: ad_2d_y
172 real,
dimension(:,:),
optional,
pointer :: df_2d_x
174 real,
dimension(:,:),
optional,
pointer :: df_2d_y
177 real,
dimension(:,:,:),
optional,
pointer :: advection_xy
178 logical,
optional,
intent(in) :: registry_diags
180 character(len=*),
optional,
intent(in) :: flux_nameroot
182 character(len=*),
optional,
intent(in) :: flux_longname
184 character(len=*),
optional,
intent(in) :: flux_units
185 real,
optional,
intent(in) :: flux_scale
187 character(len=*),
optional,
intent(in) :: convergence_units
189 real,
optional,
intent(in) :: convergence_scale
191 character(len=*),
optional,
intent(in) :: cmor_tendprefix
193 integer,
optional,
intent(in) :: diag_form
199 logical,
optional,
intent(in) :: mandatory
204 character(len=256) :: mesg
206 if (.not.
associated(reg))
call tracer_registry_init(param_file, reg)
208 if (reg%ntr>=max_fields_)
then
209 write(mesg,
'("Increase MAX_FIELDS_ in MOM_memory.h to at least ",I3," to allow for &
210 &all the tracers being registered via register_tracer.")') reg%ntr+1
211 call mom_error(fatal,
"MOM register_tracer: "//mesg)
213 reg%ntr = reg%ntr + 1
215 tr => reg%Tr(reg%ntr)
217 if (
present(name))
then
219 tr%longname = name ;
if (
present(longname)) tr%longname = longname
220 tr%units =
"Conc" ;
if (
present(units)) tr%units = units
223 if (
present(cmor_name)) tr%cmor_name = cmor_name
225 tr%cmor_units = tr%units
226 if (
present(cmor_units)) tr%cmor_units = cmor_units
228 tr%cmor_longname =
""
229 if (
present(cmor_longname)) tr%cmor_longname = cmor_longname
231 if (
present(tr_desc))
call mom_error(warning,
"MOM register_tracer: "//&
232 "It is a bad idea to use both name and tr_desc when registring "//trim(name))
233 elseif (
present(tr_desc))
then
234 call query_vardesc(tr_desc, name=tr%name, units=tr%units, &
235 longname=tr%longname, cmor_field_name=tr%cmor_name, &
236 cmor_longname=tr%cmor_longname, caller=
"register_tracer")
237 tr%cmor_units = tr%units
239 call mom_error(fatal,
"MOM register_tracer: Either name or "//&
240 "tr_desc must be present when registering a tracer.")
243 if (reg%locked)
call mom_error(fatal, &
244 "MOM register_tracer was called for variable "//trim(tr%name)//&
245 " with a locked tracer registry.")
247 tr%flux_nameroot = tr%name
248 if (
present(flux_nameroot))
then
249 if (len_trim(flux_nameroot) > 0) tr%flux_nameroot = flux_nameroot
252 tr%flux_longname = tr%longname
253 if (
present(flux_longname))
then
254 if (len_trim(flux_longname) > 0) tr%flux_longname = flux_longname
258 if (
present(flux_units)) tr%flux_units = flux_units
261 if (
present(flux_scale)) tr%flux_scale = flux_scale
264 if (
present(convergence_units)) tr%conv_units = convergence_units
266 tr%cmor_tendprefix =
""
267 if (
present(cmor_tendprefix)) tr%cmor_tendprefix = cmor_tendprefix
270 if (
present(convergence_scale))
then
271 tr%conv_scale = convergence_scale
272 elseif (
present(flux_scale))
then
273 tr%conv_scale = flux_scale
277 if (
present(diag_form)) tr%diag_form = diag_form
281 if (
present(registry_diags)) tr%registry_diags = registry_diags
283 if (
present(ad_x))
then ;
if (
associated(ad_x)) tr%ad_x => ad_x ;
endif
284 if (
present(ad_y))
then ;
if (
associated(ad_y)) tr%ad_y => ad_y ;
endif
285 if (
present(df_x))
then ;
if (
associated(df_x)) tr%df_x => df_x ;
endif
286 if (
present(df_y))
then ;
if (
associated(df_y)) tr%df_y => df_y ;
endif
292 if (
present(ad_2d_x))
then ;
if (
associated(ad_2d_x)) tr%ad2d_x => ad_2d_x ;
endif
293 if (
present(ad_2d_y))
then ;
if (
associated(ad_2d_y)) tr%ad2d_y => ad_2d_y ;
endif
294 if (
present(df_2d_x))
then ;
if (
associated(df_2d_x)) tr%df2d_x => df_2d_x ;
endif
296 if (
present(advection_xy))
then ;
if (
associated(advection_xy)) tr%advection_xy => advection_xy ;
endif
298 if (
present(restart_cs))
then ;
if (
associated(restart_cs))
then
300 mand = .true. ;
if (
present(mandatory)) mand = mandatory
303 longname=tr%longname, units=tr%units)
306 end subroutine register_tracer
311 subroutine lock_tracer_registry(Reg)
314 if (.not.
associated(reg))
call mom_error(warning, &
315 "lock_tracer_registry called with an unassociated registry.")
319 end subroutine lock_tracer_registry
323 subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, use_ALE)
327 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
329 type(time_type),
intent(in) :: time
331 logical,
intent(in) :: use_ale
334 character(len=24) :: name
335 character(len=24) :: shortnm
337 character(len=72) :: longname
338 character(len=72) :: flux_longname
339 character(len=48) :: units
340 character(len=48) :: flux_units
342 character(len=48) :: conv_units
344 character(len=48) :: unit2
345 character(len=72) :: cmorname
346 character(len=120) :: cmor_longname
347 character(len=120) :: var_lname
348 character(len=120) :: cmor_var_lname
349 character(len=72) :: cmor_varname
351 integer :: i, j, k, is, ie, js, je, nz, m, m2, ntr_in
352 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
353 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
354 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
355 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
357 if (.not.
associated(reg))
call mom_error(fatal,
"register_tracer_diagnostics: "//&
358 "register_tracer must be called before register_tracer_diagnostics")
362 do m=1,ntr_in ;
if (reg%Tr(m)%registry_diags)
then
367 name = tr%name ; units=adjustl(tr%units) ; longname = tr%longname
368 cmorname = tr%cmor_name ; cmor_longname = tr%cmor_longname
369 shortnm = tr%flux_nameroot
370 flux_longname = tr%flux_longname
371 if (len_trim(cmor_longname) == 0) cmor_longname = longname
373 if (len_trim(tr%flux_units) > 0)
then ; flux_units = tr%flux_units
374 elseif (gv%Boussinesq)
then ; flux_units = trim(units)//
" m3 s-1"
375 else ; flux_units = trim(units)//
" kg s-1" ;
endif
377 if (len_trim(tr%conv_units) > 0)
then ; conv_units = tr%conv_units
378 elseif (gv%Boussinesq)
then ; conv_units = trim(units)//
" m s-1"
379 else ; conv_units = trim(units)//
" kg m-2 s-1" ;
endif
381 if (len_trim(cmorname) == 0)
then
382 tr%id_tr = register_diag_field(
"ocean_model", trim(name), diag%axesTL, &
383 time, trim(longname), trim(units))
385 tr%id_tr = register_diag_field(
"ocean_model", trim(name), diag%axesTL, &
386 time, trim(longname), trim(units), cmor_field_name=cmorname, &
387 cmor_long_name=cmor_longname, cmor_units=tr%cmor_units, &
388 cmor_standard_name=cmor_long_std(cmor_longname))
390 if (tr%diag_form == 1)
then
391 tr%id_adx = register_diag_field(
"ocean_model", trim(shortnm)//
"_adx", &
392 diag%axesCuL, time, trim(flux_longname)//
" advective zonal flux" , &
393 trim(flux_units), v_extensive = .true., y_cell_method =
'sum')
394 tr%id_ady = register_diag_field(
"ocean_model", trim(shortnm)//
"_ady", &
395 diag%axesCvL, time, trim(flux_longname)//
" advective meridional flux" , &
396 trim(flux_units), v_extensive = .true., x_cell_method =
'sum')
397 tr%id_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_dfx", &
398 diag%axesCuL, time, trim(flux_longname)//
" diffusive zonal flux" , &
399 trim(flux_units), v_extensive = .true., y_cell_method =
'sum')
400 tr%id_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_dfy", &
401 diag%axesCvL, time, trim(flux_longname)//
" diffusive zonal flux" , &
402 trim(flux_units), v_extensive = .true., x_cell_method =
'sum')
404 tr%id_adx = register_diag_field(
"ocean_model", trim(shortnm)//
"_adx", &
405 diag%axesCuL, time,
"Advective (by residual mean) Zonal Flux of "//trim(flux_longname), &
406 flux_units, v_extensive=.true., conversion=tr%flux_scale, y_cell_method =
'sum')
407 tr%id_ady = register_diag_field(
"ocean_model", trim(shortnm)//
"_ady", &
408 diag%axesCvL, time,
"Advective (by residual mean) Meridional Flux of "//trim(flux_longname), &
409 flux_units, v_extensive=.true., conversion=tr%flux_scale, x_cell_method =
'sum')
410 tr%id_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffx", &
411 diag%axesCuL, time,
"Diffusive Zonal Flux of "//trim(flux_longname), &
412 flux_units, v_extensive=.true., conversion=tr%flux_scale, y_cell_method =
'sum')
413 tr%id_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffy", &
414 diag%axesCvL, time,
"Diffusive Meridional Flux of "//trim(flux_longname), &
415 flux_units, v_extensive=.true., conversion=tr%flux_scale, x_cell_method =
'sum')
417 if (tr%id_adx > 0)
call safe_alloc_ptr(tr%ad_x,isdb,iedb,jsd,jed,nz)
418 if (tr%id_ady > 0)
call safe_alloc_ptr(tr%ad_y,isd,ied,jsdb,jedb,nz)
419 if (tr%id_dfx > 0)
call safe_alloc_ptr(tr%df_x,isdb,iedb,jsd,jed,nz)
420 if (tr%id_dfy > 0)
call safe_alloc_ptr(tr%df_y,isd,ied,jsdb,jedb,nz)
422 tr%id_adx_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_adx_2d", &
423 diag%axesCu1, time, &
424 "Vertically Integrated Advective Zonal Flux of "//trim(flux_longname), &
425 flux_units, conversion=tr%flux_scale, y_cell_method =
'sum')
426 tr%id_ady_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_ady_2d", &
427 diag%axesCv1, time, &
428 "Vertically Integrated Advective Meridional Flux of "//trim(flux_longname), &
429 flux_units, conversion=tr%flux_scale, x_cell_method =
'sum')
430 tr%id_dfx_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffx_2d", &
431 diag%axesCu1, time, &
432 "Vertically Integrated Diffusive Zonal Flux of "//trim(flux_longname), &
433 flux_units, conversion=tr%flux_scale, y_cell_method =
'sum')
434 tr%id_dfy_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffy_2d", &
435 diag%axesCv1, time, &
436 "Vertically Integrated Diffusive Meridional Flux of "//trim(flux_longname), &
437 flux_units, conversion=tr%flux_scale, x_cell_method =
'sum')
439 if (tr%id_adx_2d > 0)
call safe_alloc_ptr(tr%ad2d_x,isdb,iedb,jsd,jed)
440 if (tr%id_ady_2d > 0)
call safe_alloc_ptr(tr%ad2d_y,isd,ied,jsdb,jedb)
441 if (tr%id_dfx_2d > 0)
call safe_alloc_ptr(tr%df2d_x,isdb,iedb,jsd,jed)
442 if (tr%id_dfy_2d > 0)
call safe_alloc_ptr(tr%df2d_y,isd,ied,jsdb,jedb)
444 tr%id_adv_xy = register_diag_field(
'ocean_model', trim(shortnm)//
"_advection_xy", &
446 'Horizontal convergence of residual mean advective fluxes of '//&
447 trim(lowercase(flux_longname)), conv_units, v_extensive=.true., &
448 conversion=tr%conv_scale)
449 tr%id_adv_xy_2d = register_diag_field(
'ocean_model', trim(shortnm)//
"_advection_xy_2d", &
451 'Vertical sum of horizontal convergence of residual mean advective fluxes of '//&
452 trim(lowercase(flux_longname)), conv_units, conversion=tr%conv_scale)
453 if ((tr%id_adv_xy > 0) .or. (tr%id_adv_xy_2d > 0)) &
454 call safe_alloc_ptr(tr%advection_xy,isd,ied,jsd,jed,nz)
456 tr%id_tendency = register_diag_field(
'ocean_model', trim(shortnm)//
'_tendency', &
458 'Net time tendency for '//trim(lowercase(longname)), trim(units)//
' s-1')
460 if (tr%id_tendency > 0)
then
461 call safe_alloc_ptr(tr%t_prev,isd,ied,jsd,jed,nz)
462 do k=1,nz ;
do j=js,je ;
do i=is,ie
463 tr%t_prev(i,j,k) = tr%t(i,j,k)
464 enddo ;
enddo ;
enddo
468 if (tr%diag_form == 1)
then
469 tr%id_dfxy_cont = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency', &
470 diag%axesTL, time,
"Lateral or neutral diffusion tracer content tendency for "//trim(shortnm), &
471 conv_units, conversion=tr%conv_scale, x_cell_method=
'sum', y_cell_method=
'sum', v_extensive=.true.)
473 tr%id_dfxy_cont_2d = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency_2d', &
474 diag%axesT1, time,
"Depth integrated lateral or neutral diffusion tracer concentration "//&
475 "tendency for "//trim(shortnm), conv_units, conversion = tr%conv_scale, &
476 x_cell_method =
'sum', y_cell_method =
'sum')
478 cmor_var_lname =
'Tendency of '//trim(lowercase(cmor_longname))//&
479 ' expressed as '//trim(lowercase(flux_longname))//&
480 ' content due to parameterized mesoscale diffusion'
481 tr%id_dfxy_cont = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency', &
482 diag%axesTL, time,
"Lateral or neutral diffusion tracer concentration tendency for "//trim(shortnm), &
483 conv_units, conversion = tr%conv_scale, cmor_field_name = trim(tr%cmor_tendprefix)//
'pmdiff', &
484 cmor_long_name = trim(cmor_var_lname), cmor_standard_name = trim(cmor_long_std(cmor_var_lname)), &
485 x_cell_method =
'sum', y_cell_method =
'sum', v_extensive = .true.)
487 cmor_var_lname =
'Tendency of '//trim(lowercase(cmor_longname))//
' expressed as '//&
488 trim(lowercase(flux_longname))//
' content due to parameterized mesoscale diffusion'
489 tr%id_dfxy_cont_2d = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency_2d', &
490 diag%axesT1, time,
"Depth integrated lateral or neutral diffusion tracer "//&
491 "concentration tendency for "//trim(shortnm), conv_units, &
492 conversion=tr%conv_scale, cmor_field_name=trim(tr%cmor_tendprefix)//
'pmdiff_2d', &
493 cmor_long_name=trim(cmor_var_lname), cmor_standard_name=trim(cmor_long_std(cmor_var_lname)), &
494 x_cell_method=
'sum', y_cell_method=
'sum')
496 tr%id_dfxy_conc = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_conc_tendency', &
497 diag%axesTL, time,
"Lateral (neutral) tracer concentration tendency for "//trim(shortnm), &
500 var_lname =
"Net time tendency for "//lowercase(flux_longname)
501 if (len_trim(tr%cmor_tendprefix) == 0)
then
502 tr%id_trxh_tendency = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency', &
503 diag%axesTL, time, var_lname, conv_units, v_extensive=.true.)
504 tr%id_trxh_tendency_2d = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency_2d', &
505 diag%axesT1, time,
"Vertical sum of "//trim(lowercase(var_lname)), conv_units)
507 cmor_var_lname =
"Tendency of "//trim(cmor_longname)//
" Expressed as "//&
508 trim(flux_longname)//
" Content"
509 tr%id_trxh_tendency = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency', &
510 diag%axesTL, time, var_lname, conv_units, &
511 cmor_field_name=trim(tr%cmor_tendprefix)//
"tend", &
512 cmor_standard_name=cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname, &
513 v_extensive=.true., conversion=tr%conv_scale)
514 cmor_var_lname = trim(cmor_var_lname)//
" Vertical Sum"
515 tr%id_trxh_tendency_2d = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency_2d', &
516 diag%axesT1, time,
"Vertical sum of "//trim(lowercase(var_lname)), conv_units, &
517 cmor_field_name=trim(tr%cmor_tendprefix)//
"tend_2d", &
518 cmor_standard_name=cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname, &
519 conversion=tr%conv_scale)
521 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0))
then
522 call safe_alloc_ptr(tr%Trxh_prev,isd,ied,jsd,jed,nz)
523 do k=1,nz ;
do j=js,je ;
do i=is,ie
524 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
525 enddo ;
enddo ;
enddo
529 if (use_ale .and. tr%remap_tr)
then
530 var_lname =
"Vertical remapping tracer concentration tendency for "//trim(reg%Tr(m)%name)
531 tr%id_remap_conc= register_diag_field(
'ocean_model', &
532 trim(tr%flux_nameroot)//
'_tendency_vert_remap', diag%axesTL, time, var_lname, &
535 var_lname =
"Vertical remapping tracer content tendency for "//trim(reg%Tr(m)%flux_longname)
536 tr%id_remap_cont = register_diag_field(
'ocean_model', &
537 trim(tr%flux_nameroot)//
'h_tendency_vert_remap', &
538 diag%axesTL, time, var_lname, flux_units, v_extensive=.true., conversion = tr%conv_scale)
540 var_lname =
"Vertical sum of vertical remapping tracer content tendency for "//&
541 trim(reg%Tr(m)%flux_longname)
542 tr%id_remap_cont_2d = register_diag_field(
'ocean_model', &
543 trim(tr%flux_nameroot)//
'h_tendency_vert_remap_2d', &
544 diag%axesT1, time, var_lname, flux_units, conversion = tr%conv_scale)
548 if (use_ale .and. (reg%ntr<max_fields_) .and. tr%remap_tr)
then
549 unit2 = trim(units)//
"2"
550 if (index(units(1:len_trim(units)),
" ") > 0) unit2 =
"("//trim(units)//
")2"
551 tr%id_tr_vardec = register_diag_field(
'ocean_model', trim(shortnm)//
"_vardec", diag%axesTL, &
552 time,
"ALE variance decay for "//lowercase(longname), trim(unit2)//
" s-1")
553 if (tr%id_tr_vardec > 0)
then
556 tr%ind_tr_squared = m2
557 call safe_alloc_ptr(reg%Tr(m2)%t,isd,ied,jsd,jed,nz) ; reg%Tr(m2)%t(:,:,:) = 0.0
558 reg%Tr(m2)%name = trim(shortnm)//
"2"
559 reg%Tr(m2)%longname =
"Squared "//trim(longname)
560 reg%Tr(m2)%units = unit2
561 reg%Tr(m2)%registry_diags = .false.
562 reg%Tr(m2)%ind_tr_squared = -1
564 reg%ntr = reg%ntr + 1
570 end subroutine register_tracer_diagnostics
572 subroutine preale_tracer_diagnostics(Reg, G, GV)
577 integer :: i, j, k, is, ie, js, je, nz, m, m2
578 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
580 do m=1,reg%ntr ;
if (reg%Tr(m)%ind_tr_squared > 0)
then
581 m2 = reg%Tr(m)%ind_tr_squared
583 do k=1,nz ;
do j=js,je ;
do i=is,ie
584 reg%Tr(m2)%T(i,j,k) = reg%Tr(m)%T(i,j,k)**2
585 enddo ;
enddo ;
enddo
588 end subroutine preale_tracer_diagnostics
590 subroutine postale_tracer_diagnostics(Reg, G, GV, diag, dt)
595 real,
intent(in) :: dt
597 real :: work(szi_(g),szj_(g),szk_(g))
599 integer :: i, j, k, is, ie, js, je, nz, m, m2
600 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
603 idt = 0.0 ;
if (dt /= 0.0) idt = 1.0 / dt
605 do m=1,reg%ntr ;
if (reg%Tr(m)%id_tr_vardec > 0)
then
606 m2 = reg%Tr(m)%ind_tr_squared
607 if (m2 < 1)
call mom_error(fatal,
"Bad value of Tr%ind_tr_squared for "//trim(reg%Tr(m)%name))
609 do k=1,nz ;
do j=js,je ;
do i=is,ie
610 work(i,j,k) = (reg%Tr(m2)%T(i,j,k) - reg%Tr(m)%T(i,j,k)**2) * idt
611 enddo ;
enddo ;
enddo
612 call post_data(reg%Tr(m)%id_tr_vardec, work, diag)
615 end subroutine postale_tracer_diagnostics
619 subroutine post_tracer_diagnostics(Reg, h, diag_prev, diag, G, GV, dt)
623 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
627 real,
intent(in) :: dt
629 real :: work3d(szi_(g),szj_(g),szk_(gv))
630 real :: work2d(szi_(g),szj_(g))
633 integer :: i, j, k, is, ie, js, je, nz, m
634 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
636 idt = 0.;
if (dt/=0.) idt = 1.0 / dt
639 call diag_save_grids(diag)
640 call diag_copy_storage_to_diag(diag, diag_prev)
641 do m=1,reg%ntr ;
if (reg%Tr(m)%registry_diags)
then
643 if (tr%id_tendency > 0)
then
645 do k=1,nz ;
do j=js,je ;
do i=is,ie
646 work3d(i,j,k) = (tr%t(i,j,k) - tr%t_prev(i,j,k))*idt
647 tr%t_prev(i,j,k) = tr%t(i,j,k)
648 enddo ;
enddo ;
enddo
649 call post_data(tr%id_tendency, work3d, diag, alt_h = diag_prev%h_state)
651 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0))
then
652 do k=1,nz ;
do j=js,je ;
do i=is,ie
653 work3d(i,j,k) = (tr%t(i,j,k)*h(i,j,k) - tr%Trxh_prev(i,j,k)) * idt
654 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
655 enddo ;
enddo ;
enddo
656 if (tr%id_trxh_tendency > 0)
call post_data(tr%id_trxh_tendency, work3d, diag, alt_h = diag_prev%h_state)
657 if (tr%id_trxh_tendency_2d > 0)
then
659 do k=1,nz ;
do j=js,je ;
do i=is,ie
660 work2d(i,j) = work2d(i,j) + work3d(i,j,k)
661 enddo ;
enddo ;
enddo
662 call post_data(tr%id_trxh_tendency_2d, work2d, diag)
666 call diag_restore_grids(diag)
668 end subroutine post_tracer_diagnostics
671 subroutine post_tracer_transport_diagnostics(G, GV, Reg, h_diag, diag)
675 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
679 integer :: i, j, k, is, ie, js, je, nz, m
680 real :: work2d(szi_(g),szj_(g))
683 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
685 do m=1,reg%ntr ;
if (reg%Tr(m)%registry_diags)
then
687 if (tr%id_tr > 0)
call post_data(tr%id_tr, tr%t, diag)
688 if (tr%id_adx > 0)
call post_data(tr%id_adx, tr%ad_x, diag, alt_h = h_diag)
689 if (tr%id_ady > 0)
call post_data(tr%id_ady, tr%ad_y, diag, alt_h = h_diag)
690 if (tr%id_dfx > 0)
call post_data(tr%id_dfx, tr%df_x, diag, alt_h = h_diag)
691 if (tr%id_dfy > 0)
call post_data(tr%id_dfy, tr%df_y, diag, alt_h = h_diag)
692 if (tr%id_adx_2d > 0)
call post_data(tr%id_adx_2d, tr%ad2d_x, diag)
693 if (tr%id_ady_2d > 0)
call post_data(tr%id_ady_2d, tr%ad2d_y, diag)
694 if (tr%id_dfx_2d > 0)
call post_data(tr%id_dfx_2d, tr%df2d_x, diag)
695 if (tr%id_dfy_2d > 0)
call post_data(tr%id_dfy_2d, tr%df2d_y, diag)
696 if (tr%id_adv_xy > 0)
call post_data(tr%id_adv_xy, tr%advection_xy, diag, alt_h = h_diag)
697 if (tr%id_adv_xy_2d > 0)
then
699 do k=1,nz ;
do j=js,je ;
do i=is,ie
700 work2d(i,j) = work2d(i,j) + tr%advection_xy(i,j,k)
701 enddo ;
enddo ;
enddo
702 call post_data(tr%id_adv_xy_2d, work2d, diag)
706 end subroutine post_tracer_transport_diagnostics
709 subroutine mom_tracer_chksum(mesg, Tr, ntr, G)
710 character(len=*),
intent(in) :: mesg
712 integer,
intent(in) :: ntr
715 integer :: is, ie, js, je, nz
716 integer :: i, j, k, m
718 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
720 call hchksum(tr(m)%t, mesg//trim(tr(m)%name), g%HI)
723 end subroutine mom_tracer_chksum
726 subroutine mom_tracer_chkinv(mesg, G, h, Tr, ntr)
727 character(len=*),
intent(in) :: mesg
730 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
731 integer,
intent(in) :: ntr
733 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: tr_inv
735 integer :: is, ie, js, je, nz
736 integer :: i, j, k, m
738 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
740 do k=1,nz ;
do j=js,je ;
do i=is,ie
741 tr_inv(i,j,k) = tr(m)%t(i,j,k)*h(i,j,k)*g%areaT(i,j)*g%mask2dT(i,j)
742 enddo ;
enddo ;
enddo
743 total_inv =
reproducing_sum(tr_inv, is+(1-g%isd), ie+(1-g%isd), js+(1-g%jsd), je+(1-g%jsd))
744 if (is_root_pe())
write(0,
'(A,1X,A5,1X,ES25.16,1X,A)')
"h-point: inventory", tr(m)%name, total_inv, mesg
747 end subroutine mom_tracer_chkinv
750 subroutine tracer_name_lookup(Reg, tr_ptr, name)
753 character(len=32),
intent(in) :: name
757 if (lowercase(reg%Tr(n)%name) == lowercase(name)) tr_ptr => reg%Tr(n)
760 end subroutine tracer_name_lookup
763 subroutine tracer_registry_init(param_file, Reg)
767 integer,
save :: init_calls = 0
770 #include "version_variable.h"
771 character(len=40) :: mdl =
"MOM_tracer_registry"
772 character(len=256) :: mesg
774 if (.not.
associated(reg))
then ;
allocate(reg)
775 else ;
return ;
endif
780 init_calls = init_calls + 1
781 if (init_calls > 1)
then
782 write(mesg,
'("tracer_registry_init called ",I3, &
783 &" times with different registry pointers.")') init_calls
784 if (is_root_pe())
call mom_error(warning,
"MOM_tracer"//mesg)
787 end subroutine tracer_registry_init
791 subroutine tracer_registry_end(Reg)
793 if (
associated(reg))
deallocate(reg)
794 end subroutine tracer_registry_end