26 implicit none ;
private
28 #include <MOM_memory.h>
40 real,
dimension(:,:,:),
pointer :: t => null()
47 real,
dimension(:,:,:),
pointer :: ad_x => null()
49 real,
dimension(:,:,:),
pointer :: ad_y => null()
51 real,
dimension(:,:),
pointer :: ad2d_x => null()
53 real,
dimension(:,:),
pointer :: ad2d_y => null()
56 real,
dimension(:,:,:),
pointer :: df_x => null()
58 real,
dimension(:,:,:),
pointer :: df_y => null()
60 real,
dimension(:,:,:),
pointer :: lbd_dfx => null()
62 real,
dimension(:,:,:),
pointer :: lbd_dfy => null()
64 real,
dimension(:,:),
pointer :: lbd_dfx_2d => null()
66 real,
dimension(:,:),
pointer :: lbd_dfy_2d => null()
68 real,
dimension(:,:),
pointer :: lbd_bulk_df_x => null()
70 real,
dimension(:,:),
pointer :: lbd_bulk_df_y => null()
72 real,
dimension(:,:),
pointer :: df2d_x => null()
74 real,
dimension(:,:),
pointer :: df2d_y => null()
81 real,
dimension(:,:,:),
pointer :: advection_xy => null()
87 real,
dimension(:,:,:),
pointer :: t_prev => null()
89 real,
dimension(:,:,:),
pointer :: trxh_prev => null()
92 character(len=32) :: name
93 character(len=64) :: units
94 character(len=240) :: longname
96 logical :: registry_diags = .false.
98 character(len=64) :: cmor_name
99 character(len=64) :: cmor_units
100 character(len=240) :: cmor_longname
101 character(len=32) :: flux_nameroot =
""
103 character(len=64) :: flux_longname =
""
105 real :: flux_scale= 1.0
107 character(len=48) :: flux_units =
""
108 character(len=48) :: conv_units =
""
109 real :: conv_scale = 1.0
111 character(len=48) :: cmor_tendprefix =
""
114 integer :: ind_tr_squared = -1
119 logical :: remap_tr = .true.
121 integer :: diag_form = 1
123 integer :: id_tr = -1
124 integer :: id_adx = -1, id_ady = -1, id_dfx = -1, id_dfy = -1
125 integer :: id_lbd_bulk_dfx = -1, id_lbd_bulk_dfy = -1, id_lbd_dfx = -1, id_lbd_dfy = -1
126 integer :: id_lbd_dfx_2d = -1 , id_lbd_dfy_2d = -1
127 integer :: id_adx_2d = -1, id_ady_2d = -1, id_dfx_2d = -1, id_dfy_2d = -1
128 integer :: id_adv_xy = -1, id_adv_xy_2d = -1
129 integer :: id_dfxy_cont = -1, id_dfxy_cont_2d = -1, id_dfxy_conc = -1
130 integer :: id_lbdxy_cont = -1, id_lbdxy_cont_2d = -1, id_lbdxy_conc = -1
131 integer :: id_remap_conc = -1, id_remap_cont = -1, id_remap_cont_2d = -1
132 integer :: id_tendency = -1, id_trxh_tendency = -1, id_trxh_tendency_2d = -1
133 integer :: id_tr_vardec = -1
142 logical :: locked = .false.
151 subroutine register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, &
152 cmor_name, cmor_units, cmor_longname, tr_desc, &
153 OBC_inflow, OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, &
154 ad_2d_x, ad_2d_y, df_2d_x, df_2d_y, advection_xy, registry_diags, &
155 flux_nameroot, flux_longname, flux_units, flux_scale, &
156 convergence_units, convergence_scale, cmor_tendprefix, diag_form, &
157 restart_CS, mandatory)
161 real,
dimension(SZI_(HI),SZJ_(HI),SZK_(GV)), &
164 character(len=*),
optional,
intent(in) :: name
165 character(len=*),
optional,
intent(in) :: longname
166 character(len=*),
optional,
intent(in) :: units
167 character(len=*),
optional,
intent(in) :: cmor_name
168 character(len=*),
optional,
intent(in) :: cmor_units
169 character(len=*),
optional,
intent(in) :: cmor_longname
170 type(
vardesc),
optional,
intent(in) :: tr_desc
172 real,
optional,
intent(in) :: obc_inflow
174 real,
dimension(:,:,:),
optional,
pointer :: obc_in_u
176 real,
dimension(:,:,:),
optional,
pointer :: obc_in_v
180 real,
dimension(:,:,:),
optional,
pointer :: ad_x
182 real,
dimension(:,:,:),
optional,
pointer :: ad_y
184 real,
dimension(:,:,:),
optional,
pointer :: df_x
186 real,
dimension(:,:,:),
optional,
pointer :: df_y
188 real,
dimension(:,:),
optional,
pointer :: ad_2d_x
190 real,
dimension(:,:),
optional,
pointer :: ad_2d_y
192 real,
dimension(:,:),
optional,
pointer :: df_2d_x
194 real,
dimension(:,:),
optional,
pointer :: df_2d_y
197 real,
dimension(:,:,:),
optional,
pointer :: advection_xy
198 logical,
optional,
intent(in) :: registry_diags
200 character(len=*),
optional,
intent(in) :: flux_nameroot
202 character(len=*),
optional,
intent(in) :: flux_longname
204 character(len=*),
optional,
intent(in) :: flux_units
205 real,
optional,
intent(in) :: flux_scale
207 character(len=*),
optional,
intent(in) :: convergence_units
209 real,
optional,
intent(in) :: convergence_scale
211 character(len=*),
optional,
intent(in) :: cmor_tendprefix
213 integer,
optional,
intent(in) :: diag_form
219 logical,
optional,
intent(in) :: mandatory
224 character(len=256) :: mesg
228 if (reg%ntr>=max_fields_)
then
229 write(mesg,
'("Increase MAX_FIELDS_ in MOM_memory.h to at least ",I3," to allow for &
230 &all the tracers being registered via register_tracer.")') reg%ntr+1
231 call mom_error(fatal,
"MOM register_tracer: "//mesg)
233 reg%ntr = reg%ntr + 1
235 tr => reg%Tr(reg%ntr)
237 if (
present(name))
then
239 tr%longname = name ;
if (
present(longname)) tr%longname = longname
240 tr%units =
"Conc" ;
if (
present(units)) tr%units = units
243 if (
present(cmor_name)) tr%cmor_name = cmor_name
245 tr%cmor_units = tr%units
246 if (
present(cmor_units)) tr%cmor_units = cmor_units
248 tr%cmor_longname =
""
249 if (
present(cmor_longname)) tr%cmor_longname = cmor_longname
251 if (
present(tr_desc))
call mom_error(warning,
"MOM register_tracer: "//&
252 "It is a bad idea to use both name and tr_desc when registring "//trim(name))
253 elseif (
present(tr_desc))
then
255 longname=tr%longname, cmor_field_name=tr%cmor_name, &
256 cmor_longname=tr%cmor_longname, caller=
"register_tracer")
257 tr%cmor_units = tr%units
259 call mom_error(fatal,
"MOM register_tracer: Either name or "//&
260 "tr_desc must be present when registering a tracer.")
264 "MOM register_tracer was called for variable "//trim(tr%name)//&
265 " with a locked tracer registry.")
267 tr%flux_nameroot = tr%name
268 if (
present(flux_nameroot))
then
269 if (len_trim(flux_nameroot) > 0) tr%flux_nameroot = flux_nameroot
272 tr%flux_longname = tr%longname
273 if (
present(flux_longname))
then
274 if (len_trim(flux_longname) > 0) tr%flux_longname = flux_longname
278 if (
present(flux_units)) tr%flux_units = flux_units
281 if (
present(flux_scale)) tr%flux_scale = flux_scale
284 if (
present(convergence_units)) tr%conv_units = convergence_units
286 tr%cmor_tendprefix =
""
287 if (
present(cmor_tendprefix)) tr%cmor_tendprefix = cmor_tendprefix
290 if (
present(convergence_scale))
then
291 tr%conv_scale = convergence_scale
292 elseif (
present(flux_scale))
then
293 tr%conv_scale = flux_scale
297 if (
present(diag_form)) tr%diag_form = diag_form
301 if (
present(registry_diags)) tr%registry_diags = registry_diags
303 if (
present(ad_x))
then ;
if (
associated(ad_x)) tr%ad_x => ad_x ;
endif
304 if (
present(ad_y))
then ;
if (
associated(ad_y)) tr%ad_y => ad_y ;
endif
305 if (
present(df_x))
then ;
if (
associated(df_x)) tr%df_x => df_x ;
endif
306 if (
present(df_y))
then ;
if (
associated(df_y)) tr%df_y => df_y ;
endif
312 if (
present(ad_2d_x))
then ;
if (
associated(ad_2d_x)) tr%ad2d_x => ad_2d_x ;
endif
313 if (
present(ad_2d_y))
then ;
if (
associated(ad_2d_y)) tr%ad2d_y => ad_2d_y ;
endif
314 if (
present(df_2d_x))
then ;
if (
associated(df_2d_x)) tr%df2d_x => df_2d_x ;
endif
316 if (
present(advection_xy))
then ;
if (
associated(advection_xy)) tr%advection_xy => advection_xy ;
endif
318 if (
present(restart_cs))
then ;
if (
associated(restart_cs))
then
320 mand = .true. ;
if (
present(mandatory)) mand = mandatory
323 longname=tr%longname, units=tr%units)
334 if (.not.
associated(reg))
call mom_error(warning, &
335 "lock_tracer_registry called with an unassociated registry.")
348 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
350 type(time_type),
intent(in) :: time
352 logical,
intent(in) :: use_ale
355 character(len=24) :: name
356 character(len=24) :: shortnm
358 character(len=72) :: longname
359 character(len=72) :: flux_longname
360 character(len=48) :: units
361 character(len=48) :: flux_units
363 character(len=48) :: conv_units
365 character(len=48) :: unit2
366 character(len=72) :: cmorname
367 character(len=120) :: cmor_longname
368 character(len=120) :: var_lname
369 character(len=120) :: cmor_var_lname
370 character(len=72) :: cmor_varname
372 integer :: i, j, k, is, ie, js, je, nz, m, m2, ntr_in
373 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
374 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
375 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
376 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
378 if (.not.
associated(reg))
call mom_error(fatal,
"register_tracer_diagnostics: "//&
379 "register_tracer must be called before register_tracer_diagnostics")
383 do m=1,ntr_in ;
if (reg%Tr(m)%registry_diags)
then
388 name = tr%name ; units=adjustl(tr%units) ; longname = tr%longname
389 cmorname = tr%cmor_name ; cmor_longname = tr%cmor_longname
390 shortnm = tr%flux_nameroot
391 flux_longname = tr%flux_longname
392 if (len_trim(cmor_longname) == 0) cmor_longname = longname
394 if (len_trim(tr%flux_units) > 0)
then ; flux_units = tr%flux_units
395 elseif (gv%Boussinesq)
then ; flux_units = trim(units)//
" m3 s-1"
396 else ; flux_units = trim(units)//
" kg s-1" ;
endif
398 if (len_trim(tr%conv_units) > 0)
then ; conv_units = tr%conv_units
399 elseif (gv%Boussinesq)
then ; conv_units = trim(units)//
" m s-1"
400 else ; conv_units = trim(units)//
" kg m-2 s-1" ;
endif
402 if (len_trim(cmorname) == 0)
then
403 tr%id_tr = register_diag_field(
"ocean_model", trim(name), diag%axesTL, &
404 time, trim(longname), trim(units))
406 tr%id_tr = register_diag_field(
"ocean_model", trim(name), diag%axesTL, &
407 time, trim(longname), trim(units), cmor_field_name=cmorname, &
408 cmor_long_name=cmor_longname, cmor_units=tr%cmor_units, &
411 if (tr%diag_form == 1)
then
412 tr%id_adx = register_diag_field(
"ocean_model", trim(shortnm)//
"_adx", &
413 diag%axesCuL, time, trim(flux_longname)//
" advective zonal flux" , &
414 trim(flux_units), v_extensive = .true., y_cell_method =
'sum', &
415 conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T)
416 tr%id_ady = register_diag_field(
"ocean_model", trim(shortnm)//
"_ady", &
417 diag%axesCvL, time, trim(flux_longname)//
" advective meridional flux" , &
418 trim(flux_units), v_extensive = .true., x_cell_method =
'sum', &
419 conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T)
420 tr%id_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_dfx", &
421 diag%axesCuL, time, trim(flux_longname)//
" diffusive zonal flux" , &
422 trim(flux_units), v_extensive = .true., y_cell_method =
'sum', &
423 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
424 tr%id_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_dfy", &
425 diag%axesCvL, time, trim(flux_longname)//
" diffusive zonal flux" , &
426 trim(flux_units), v_extensive = .true., x_cell_method =
'sum', &
427 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
428 tr%id_lbd_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_diffx", &
429 diag%axesCuL, time, trim(flux_longname)//
" diffusive zonal flux from the lateral boundary diffusion "&
430 "scheme", trim(flux_units), v_extensive = .true., y_cell_method =
'sum', &
431 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
432 tr%id_lbd_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_diffy", &
433 diag%axesCvL, time, trim(flux_longname)//
" diffusive meridional flux from the lateral boundary diffusion"&
434 " scheme", trim(flux_units), v_extensive = .true., x_cell_method =
'sum', &
435 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
437 tr%id_adx = register_diag_field(
"ocean_model", trim(shortnm)//
"_adx", &
438 diag%axesCuL, time,
"Advective (by residual mean) Zonal Flux of "//trim(flux_longname), &
439 flux_units, v_extensive=.true., conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, y_cell_method =
'sum')
440 tr%id_ady = register_diag_field(
"ocean_model", trim(shortnm)//
"_ady", &
441 diag%axesCvL, time,
"Advective (by residual mean) Meridional Flux of "//trim(flux_longname), &
442 flux_units, v_extensive=.true., conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, x_cell_method =
'sum')
443 tr%id_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffx", &
444 diag%axesCuL, time,
"Diffusive Zonal Flux of "//trim(flux_longname), &
445 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
447 tr%id_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffy", &
448 diag%axesCvL, time,
"Diffusive Meridional Flux of "//trim(flux_longname), &
449 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
451 tr%id_lbd_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_diffx", &
452 diag%axesCuL, time,
"Lateral Boundary Diffusive Zonal Flux of "//trim(flux_longname), &
453 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
455 tr%id_lbd_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_diffy", &
456 diag%axesCvL, time,
"Lateral Boundary Diffusive Meridional Flux of "//trim(flux_longname), &
457 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
460 if (tr%id_adx > 0)
call safe_alloc_ptr(tr%ad_x,isdb,iedb,jsd,jed,nz)
461 if (tr%id_ady > 0)
call safe_alloc_ptr(tr%ad_y,isd,ied,jsdb,jedb,nz)
462 if (tr%id_dfx > 0)
call safe_alloc_ptr(tr%df_x,isdb,iedb,jsd,jed,nz)
463 if (tr%id_dfy > 0)
call safe_alloc_ptr(tr%df_y,isd,ied,jsdb,jedb,nz)
464 if (tr%id_lbd_dfx > 0)
call safe_alloc_ptr(tr%lbd_dfx,isdb,iedb,jsd,jed,nz)
465 if (tr%id_lbd_dfy > 0)
call safe_alloc_ptr(tr%lbd_dfy,isd,ied,jsdb,jedb,nz)
467 tr%id_adx_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_adx_2d", &
468 diag%axesCu1, time, &
469 "Vertically Integrated Advective Zonal Flux of "//trim(flux_longname), &
470 flux_units, conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, y_cell_method =
'sum')
471 tr%id_ady_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_ady_2d", &
472 diag%axesCv1, time, &
473 "Vertically Integrated Advective Meridional Flux of "//trim(flux_longname), &
474 flux_units, conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, x_cell_method =
'sum')
475 tr%id_dfx_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffx_2d", &
476 diag%axesCu1, time, &
477 "Vertically Integrated Diffusive Zonal Flux of "//trim(flux_longname), &
478 flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
480 tr%id_dfy_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_diffy_2d", &
481 diag%axesCv1, time, &
482 "Vertically Integrated Diffusive Meridional Flux of "//trim(flux_longname), &
483 flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
485 tr%id_lbd_bulk_dfx = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_bulk_diffx", &
486 diag%axesCu1, time, &
487 "Total Bulk Diffusive Zonal Flux of "//trim(flux_longname), &
488 flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
490 tr%id_lbd_bulk_dfy = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_bulk_diffy", &
491 diag%axesCv1, time, &
492 "Total Bulk Diffusive Meridional Flux of "//trim(flux_longname), &
493 flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
495 tr%id_lbd_dfx_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_diffx_2d", &
496 diag%axesCu1, time,
"Vertically-integrated zonal diffusive flux from the lateral boundary diffusion "//&
497 "scheme for "//trim(flux_longname), flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
498 y_cell_method =
'sum')
499 tr%id_lbd_dfy_2d = register_diag_field(
"ocean_model", trim(shortnm)//
"_lbd_diffy_2d", &
500 diag%axesCv1, time,
"Vertically-integrated meridional diffusive flux from the lateral boundary diffusion "//&
501 "scheme for "//trim(flux_longname), flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
502 x_cell_method =
'sum')
504 if (tr%id_adx_2d > 0)
call safe_alloc_ptr(tr%ad2d_x,isdb,iedb,jsd,jed)
505 if (tr%id_ady_2d > 0)
call safe_alloc_ptr(tr%ad2d_y,isd,ied,jsdb,jedb)
506 if (tr%id_dfx_2d > 0)
call safe_alloc_ptr(tr%df2d_x,isdb,iedb,jsd,jed)
507 if (tr%id_dfy_2d > 0)
call safe_alloc_ptr(tr%df2d_y,isd,ied,jsdb,jedb)
508 if (tr%id_lbd_bulk_dfx > 0)
call safe_alloc_ptr(tr%lbd_bulk_df_x,isdb,iedb,jsd,jed)
509 if (tr%id_lbd_bulk_dfy > 0)
call safe_alloc_ptr(tr%lbd_bulk_df_y,isd,ied,jsdb,jedb)
510 if (tr%id_lbd_dfx_2d > 0)
call safe_alloc_ptr(tr%lbd_dfx_2d,isdb,iedb,jsd,jed)
511 if (tr%id_lbd_dfy_2d > 0)
call safe_alloc_ptr(tr%lbd_dfy_2d,isd,ied,jsdb,jedb)
513 tr%id_adv_xy = register_diag_field(
'ocean_model', trim(shortnm)//
"_advection_xy", &
515 'Horizontal convergence of residual mean advective fluxes of '//&
516 trim(
lowercase(flux_longname)), conv_units, v_extensive=.true., &
517 conversion=tr%conv_scale*us%s_to_T)
518 tr%id_adv_xy_2d = register_diag_field(
'ocean_model', trim(shortnm)//
"_advection_xy_2d", &
520 'Vertical sum of horizontal convergence of residual mean advective fluxes of '//&
521 trim(
lowercase(flux_longname)), conv_units, conversion=tr%conv_scale*us%s_to_T)
522 if ((tr%id_adv_xy > 0) .or. (tr%id_adv_xy_2d > 0)) &
523 call safe_alloc_ptr(tr%advection_xy,isd,ied,jsd,jed,nz)
525 tr%id_tendency = register_diag_field(
'ocean_model', trim(shortnm)//
'_tendency', &
527 'Net time tendency for '//trim(
lowercase(longname)), trim(units)//
' s-1', conversion=us%s_to_T)
529 if (tr%id_tendency > 0)
then
530 call safe_alloc_ptr(tr%t_prev,isd,ied,jsd,jed,nz)
531 do k=1,nz ;
do j=js,je ;
do i=is,ie
532 tr%t_prev(i,j,k) = tr%t(i,j,k)
533 enddo ;
enddo ;
enddo
537 if (tr%diag_form == 1)
then
538 tr%id_dfxy_cont = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency', &
539 diag%axesTL, time,
"Neutral diffusion tracer content tendency for "//trim(shortnm), &
540 conv_units, conversion=tr%conv_scale*us%s_to_T, x_cell_method=
'sum', y_cell_method=
'sum', v_extensive=.true.)
542 tr%id_dfxy_cont_2d = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency_2d', &
543 diag%axesT1, time,
"Depth integrated neutral diffusion tracer content "//&
544 "tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
545 x_cell_method=
'sum', y_cell_method=
'sum')
547 tr%id_lbdxy_cont = register_diag_field(
"ocean_model", trim(shortnm)//
'_lbdxy_cont_tendency', &
548 diag%axesTL, time,
"Lateral diffusion tracer content tendency for "//trim(shortnm), &
549 conv_units, conversion=tr%conv_scale*us%s_to_T, x_cell_method=
'sum', y_cell_method=
'sum', v_extensive=.true.)
551 tr%id_lbdxy_cont_2d = register_diag_field(
"ocean_model", trim(shortnm)//
'_lbdxy_cont_tendency_2d', &
552 diag%axesT1, time,
"Depth integrated lateral diffusion tracer content "//&
553 "tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
554 x_cell_method=
'sum', y_cell_method=
'sum')
556 cmor_var_lname =
'Tendency of '//trim(
lowercase(cmor_longname))//&
557 ' expressed as '//trim(
lowercase(flux_longname))//&
558 ' content due to parameterized mesoscale neutral diffusion'
559 tr%id_dfxy_cont = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency', &
560 diag%axesTL, time,
"Neutral diffusion tracer content tendency for "//trim(shortnm), &
561 conv_units, conversion=tr%conv_scale*us%s_to_T, cmor_field_name = trim(tr%cmor_tendprefix)//
'pmdiff', &
562 cmor_long_name = trim(cmor_var_lname), cmor_standard_name = trim(
cmor_long_std(cmor_var_lname)), &
563 x_cell_method =
'sum', y_cell_method =
'sum', v_extensive = .true.)
565 cmor_var_lname =
'Tendency of '//trim(
lowercase(cmor_longname))//
' expressed as '//&
566 trim(
lowercase(flux_longname))//
' content due to parameterized mesoscale neutral diffusion'
567 tr%id_dfxy_cont_2d = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_cont_tendency_2d', &
568 diag%axesT1, time,
"Depth integrated neutral diffusion tracer "//&
569 "content tendency for "//trim(shortnm), conv_units, &
570 conversion=tr%conv_scale*us%s_to_T, cmor_field_name=trim(tr%cmor_tendprefix)//
'pmdiff_2d', &
571 cmor_long_name=trim(cmor_var_lname), cmor_standard_name=trim(
cmor_long_std(cmor_var_lname)), &
572 x_cell_method=
'sum', y_cell_method=
'sum')
574 tr%id_lbdxy_cont = register_diag_field(
"ocean_model", trim(shortnm)//
'_lbdxy_cont_tendency', &
575 diag%axesTL, time,
"Lateral diffusion tracer content tendency for "//trim(shortnm), &
576 conv_units, conversion=tr%conv_scale*us%s_to_T, &
577 x_cell_method =
'sum', y_cell_method =
'sum', v_extensive = .true.)
579 tr%id_lbdxy_cont_2d = register_diag_field(
"ocean_model", trim(shortnm)//
'_lbdxy_cont_tendency_2d', &
580 diag%axesT1, time,
"Depth integrated lateral diffusion tracer "//&
581 "content tendency for "//trim(shortnm), conv_units, &
582 conversion=tr%conv_scale*us%s_to_T, x_cell_method=
'sum', y_cell_method=
'sum')
584 tr%id_dfxy_conc = register_diag_field(
"ocean_model", trim(shortnm)//
'_dfxy_conc_tendency', &
585 diag%axesTL, time,
"Neutral diffusion tracer concentration tendency for "//trim(shortnm), &
586 trim(units)//
' s-1', conversion=us%s_to_T)
588 tr%id_lbdxy_conc = register_diag_field(
"ocean_model", trim(shortnm)//
'_lbdxy_conc_tendency', &
589 diag%axesTL, time,
"Lateral diffusion tracer concentration tendency for "//trim(shortnm), &
590 trim(units)//
' s-1', conversion=us%s_to_T)
592 var_lname =
"Net time tendency for "//
lowercase(flux_longname)
593 if (len_trim(tr%cmor_tendprefix) == 0)
then
594 tr%id_trxh_tendency = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency', &
595 diag%axesTL, time, var_lname, conv_units, v_extensive=.true., &
596 conversion=tr%conv_scale*us%s_to_T)
597 tr%id_trxh_tendency_2d = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency_2d', &
598 diag%axesT1, time,
"Vertical sum of "//trim(
lowercase(var_lname)), conv_units, &
599 conversion=tr%conv_scale*us%s_to_T)
601 cmor_var_lname =
"Tendency of "//trim(cmor_longname)//
" Expressed as "//&
602 trim(flux_longname)//
" Content"
603 tr%id_trxh_tendency = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency', &
604 diag%axesTL, time, var_lname, conv_units, &
605 cmor_field_name=trim(tr%cmor_tendprefix)//
"tend", &
606 cmor_standard_name=
cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname, &
607 v_extensive=.true., conversion=tr%conv_scale*us%s_to_T)
608 cmor_var_lname = trim(cmor_var_lname)//
" Vertical Sum"
609 tr%id_trxh_tendency_2d = register_diag_field(
'ocean_model', trim(shortnm)//
'h_tendency_2d', &
610 diag%axesT1, time,
"Vertical sum of "//trim(
lowercase(var_lname)), conv_units, &
611 cmor_field_name=trim(tr%cmor_tendprefix)//
"tend_2d", &
612 cmor_standard_name=
cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname, &
613 conversion=tr%conv_scale*us%s_to_T)
615 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0))
then
616 call safe_alloc_ptr(tr%Trxh_prev,isd,ied,jsd,jed,nz)
617 do k=1,nz ;
do j=js,je ;
do i=is,ie
618 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
619 enddo ;
enddo ;
enddo
623 if (use_ale .and. tr%remap_tr)
then
624 var_lname =
"Vertical remapping tracer concentration tendency for "//trim(reg%Tr(m)%name)
625 tr%id_remap_conc= register_diag_field(
'ocean_model', &
626 trim(tr%flux_nameroot)//
'_tendency_vert_remap', diag%axesTL, time, var_lname, &
627 trim(units)//
' s-1', conversion=us%s_to_T)
629 var_lname =
"Vertical remapping tracer content tendency for "//trim(reg%Tr(m)%flux_longname)
630 tr%id_remap_cont = register_diag_field(
'ocean_model', &
631 trim(tr%flux_nameroot)//
'h_tendency_vert_remap', &
632 diag%axesTL, time, var_lname, flux_units, v_extensive=.true., conversion=tr%conv_scale*us%s_to_T)
634 var_lname =
"Vertical sum of vertical remapping tracer content tendency for "//&
635 trim(reg%Tr(m)%flux_longname)
636 tr%id_remap_cont_2d = register_diag_field(
'ocean_model', &
637 trim(tr%flux_nameroot)//
'h_tendency_vert_remap_2d', &
638 diag%axesT1, time, var_lname, flux_units, conversion=tr%conv_scale*us%s_to_T)
642 if (use_ale .and. (reg%ntr<max_fields_) .and. tr%remap_tr)
then
643 unit2 = trim(units)//
"2"
644 if (index(units(1:len_trim(units)),
" ") > 0) unit2 =
"("//trim(units)//
")2"
645 tr%id_tr_vardec = register_diag_field(
'ocean_model', trim(shortnm)//
"_vardec", diag%axesTL, &
646 time,
"ALE variance decay for "//
lowercase(longname), trim(unit2)//
" s-1", conversion=us%s_to_T)
647 if (tr%id_tr_vardec > 0)
then
650 tr%ind_tr_squared = m2
651 call safe_alloc_ptr(reg%Tr(m2)%t,isd,ied,jsd,jed,nz) ; reg%Tr(m2)%t(:,:,:) = 0.0
652 reg%Tr(m2)%name = trim(shortnm)//
"2"
653 reg%Tr(m2)%longname =
"Squared "//trim(longname)
654 reg%Tr(m2)%units = unit2
655 reg%Tr(m2)%registry_diags = .false.
656 reg%Tr(m2)%ind_tr_squared = -1
658 reg%ntr = reg%ntr + 1
671 integer :: i, j, k, is, ie, js, je, nz, m, m2
672 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
674 do m=1,reg%ntr ;
if (reg%Tr(m)%ind_tr_squared > 0)
then
675 m2 = reg%Tr(m)%ind_tr_squared
677 do k=1,nz ;
do j=js,je ;
do i=is,ie
678 reg%Tr(m2)%T(i,j,k) = reg%Tr(m)%T(i,j,k)**2
679 enddo ;
enddo ;
enddo
689 real,
intent(in) :: dt
691 real :: work(szi_(g),szj_(g),szk_(g))
693 integer :: i, j, k, is, ie, js, je, nz, m, m2
694 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
697 idt = 0.0 ;
if (dt /= 0.0) idt = 1.0 / dt
699 do m=1,reg%ntr ;
if (reg%Tr(m)%id_tr_vardec > 0)
then
700 m2 = reg%Tr(m)%ind_tr_squared
701 if (m2 < 1)
call mom_error(fatal,
"Bad value of Tr%ind_tr_squared for "//trim(reg%Tr(m)%name))
703 do k=1,nz ;
do j=js,je ;
do i=is,ie
704 work(i,j,k) = (reg%Tr(m2)%T(i,j,k) - reg%Tr(m)%T(i,j,k)**2) * idt
705 enddo ;
enddo ;
enddo
706 call post_data(reg%Tr(m)%id_tr_vardec, work, diag)
717 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
721 real,
intent(in) :: dt
723 real :: work3d(szi_(g),szj_(g),szk_(gv))
724 real :: work2d(szi_(g),szj_(g))
727 integer :: i, j, k, is, ie, js, je, nz, m
728 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
730 idt = 0.;
if (dt/=0.) idt = 1.0 / dt
735 do m=1,reg%ntr ;
if (reg%Tr(m)%registry_diags)
then
737 if (tr%id_tendency > 0)
then
739 do k=1,nz ;
do j=js,je ;
do i=is,ie
740 work3d(i,j,k) = (tr%t(i,j,k) - tr%t_prev(i,j,k))*idt
741 tr%t_prev(i,j,k) = tr%t(i,j,k)
742 enddo ;
enddo ;
enddo
743 call post_data(tr%id_tendency, work3d, diag, alt_h=diag_prev%h_state)
745 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0))
then
746 do k=1,nz ;
do j=js,je ;
do i=is,ie
747 work3d(i,j,k) = (tr%t(i,j,k)*h(i,j,k) - tr%Trxh_prev(i,j,k)) * idt
748 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
749 enddo ;
enddo ;
enddo
750 if (tr%id_trxh_tendency > 0)
call post_data(tr%id_trxh_tendency, work3d, diag, alt_h=diag_prev%h_state)
751 if (tr%id_trxh_tendency_2d > 0)
then
753 do k=1,nz ;
do j=js,je ;
do i=is,ie
754 work2d(i,j) = work2d(i,j) + work3d(i,j,k)
755 enddo ;
enddo ;
enddo
756 call post_data(tr%id_trxh_tendency_2d, work2d, diag)
769 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
773 integer :: i, j, k, is, ie, js, je, nz, m
774 real :: work2d(szi_(g),szj_(g))
777 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
779 do m=1,reg%ntr ;
if (reg%Tr(m)%registry_diags)
then
781 if (tr%id_tr > 0)
call post_data(tr%id_tr, tr%t, diag)
782 if (tr%id_adx > 0)
call post_data(tr%id_adx, tr%ad_x, diag, alt_h=h_diag)
783 if (tr%id_ady > 0)
call post_data(tr%id_ady, tr%ad_y, diag, alt_h=h_diag)
784 if (tr%id_dfx > 0)
call post_data(tr%id_dfx, tr%df_x, diag, alt_h=h_diag)
785 if (tr%id_dfy > 0)
call post_data(tr%id_dfy, tr%df_y, diag, alt_h=h_diag)
786 if (tr%id_adx_2d > 0)
call post_data(tr%id_adx_2d, tr%ad2d_x, diag)
787 if (tr%id_ady_2d > 0)
call post_data(tr%id_ady_2d, tr%ad2d_y, diag)
788 if (tr%id_dfx_2d > 0)
call post_data(tr%id_dfx_2d, tr%df2d_x, diag)
789 if (tr%id_dfy_2d > 0)
call post_data(tr%id_dfy_2d, tr%df2d_y, diag)
790 if (tr%id_adv_xy > 0)
call post_data(tr%id_adv_xy, tr%advection_xy, diag, alt_h=h_diag)
791 if (tr%id_adv_xy_2d > 0)
then
793 do k=1,nz ;
do j=js,je ;
do i=is,ie
794 work2d(i,j) = work2d(i,j) + tr%advection_xy(i,j,k)
795 enddo ;
enddo ;
enddo
796 call post_data(tr%id_adv_xy_2d, work2d, diag)
804 character(len=*),
intent(in) :: mesg
806 integer,
intent(in) :: ntr
809 integer :: is, ie, js, je, nz
810 integer :: i, j, k, m
812 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
814 call hchksum(tr(m)%t, mesg//trim(tr(m)%name), g%HI)
821 character(len=*),
intent(in) :: mesg
824 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
825 integer,
intent(in) :: ntr
827 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: tr_inv
829 integer :: is, ie, js, je, nz
830 integer :: i, j, k, m
832 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
834 do k=1,nz ;
do j=js,je ;
do i=is,ie
835 tr_inv(i,j,k) = tr(m)%t(i,j,k)*h(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)*g%mask2dT(i,j)
836 enddo ;
enddo ;
enddo
837 total_inv =
reproducing_sum(tr_inv, is+(1-g%isd), ie+(1-g%isd), js+(1-g%jsd), je+(1-g%jsd))
838 if (
is_root_pe())
write(0,
'(A,1X,A5,1X,ES25.16,1X,A)')
"h-point: inventory", tr(m)%name, total_inv, mesg
847 character(len=32),
intent(in) :: name
861 integer,
save :: init_calls = 0
864 #include "version_variable.h"
865 character(len=40) :: mdl =
"MOM_tracer_registry"
866 character(len=256) :: mesg
868 if (.not.
associated(reg))
then ;
allocate(reg)
869 else ;
return ;
endif
874 init_calls = init_calls + 1
875 if (init_calls > 1)
then
876 write(mesg,
'("tracer_registry_init called ",I3, &
877 &" times with different registry pointers.")') init_calls
887 if (
associated(reg))
deallocate(reg)