MOM6
MOM_tracer_registry.F90
Go to the documentation of this file.
1 !> This module contains the tracer_registry_type and the subroutines
2 !! that handle registration of tracers and related subroutines.
3 !! The primary subroutine, register_tracer, is called to indicate the
4 !! tracers advected and diffused.
6 
7 ! This file is part of MOM6. See LICENSE.md for the license.
8 
9 ! use MOM_diag_mediator, only : diag_ctrl
10 use mom_coms, only : reproducing_sum
11 use mom_debugging, only : hchksum
12 use mom_diag_mediator, only : diag_ctrl, register_diag_field, post_data, safe_alloc_ptr
15 use mom_error_handler, only : mom_error, fatal, warning, mom_mesg, is_root_pe
17 use mom_hor_index, only : hor_index_type
18 use mom_grid, only : ocean_grid_type
22 use mom_time_manager, only : time_type
25 
26 implicit none ; private
27 
28 #include <MOM_memory.h>
29 
30 public register_tracer
35 public tracer_name_lookup
36 
37 !> The tracer type
38 type, public :: tracer_type
39 
40  real, dimension(:,:,:), pointer :: t => null() !< tracer concentration array [conc]
41 ! real :: OBC_inflow_conc= 0.0 !< tracer concentration for generic inflows [conc]
42 ! real, dimension(:,:,:), pointer :: OBC_in_u => NULL() !< structured values for flow into the domain
43 ! !! specified in OBCs through u-face of cell
44 ! real, dimension(:,:,:), pointer :: OBC_in_v => NULL() !< structured values for flow into the domain
45 ! !! specified in OBCs through v-face of cell
46 
47  real, dimension(:,:,:), pointer :: ad_x => null() !< diagnostic array for x-advective tracer flux
48  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
49  real, dimension(:,:,:), pointer :: ad_y => null() !< diagnostic array for y-advective tracer flux
50  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
51  real, dimension(:,:), pointer :: ad2d_x => null() !< diagnostic vertical sum x-advective tracer flux
52  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
53  real, dimension(:,:), pointer :: ad2d_y => null() !< diagnostic vertical sum y-advective tracer flux
54  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
55 
56  real, dimension(:,:,:), pointer :: df_x => null() !< diagnostic array for x-diffusive tracer flux
57  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
58  real, dimension(:,:,:), pointer :: df_y => null() !< diagnostic array for y-diffusive tracer flux
59  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
60  real, dimension(:,:,:), pointer :: lbd_dfx => null() !< diagnostic array for x-diffusive tracer flux
61  !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
62  real, dimension(:,:,:), pointer :: lbd_dfy => null() !< diagnostic array for y-diffusive tracer flux
63  !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
64  real, dimension(:,:), pointer :: lbd_dfx_2d => null() !< diagnostic array for x-diffusive tracer flux
65  !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
66  real, dimension(:,:), pointer :: lbd_dfy_2d => null() !< diagnostic array for y-diffusive tracer flux
67  !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
68  real, dimension(:,:), pointer :: lbd_bulk_df_x => null() !< diagnostic array for x-diffusive tracer flux
69  !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
70  real, dimension(:,:), pointer :: lbd_bulk_df_y => null() !< diagnostic array for y-diffusive tracer flux
71  !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
72  real, dimension(:,:), pointer :: df2d_x => null() !< diagnostic vertical sum x-diffusive flux
73  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
74  real, dimension(:,:), pointer :: df2d_y => null() !< diagnostic vertical sum y-diffusive flux
75  !! [conc H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
76 ! real, dimension(:,:), pointer :: df2d_conc_x => NULL() !< diagnostic vertical sum x-diffusive content flux
77 ! !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
78 ! real, dimension(:,:), pointer :: df2d_conc_y => NULL() !< diagnostic vertical sum y-diffusive content flux
79 ! !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1]
80 
81  real, dimension(:,:,:), pointer :: advection_xy => null() !< convergence of lateral advective tracer fluxes
82  !! [conc H T-1 ~> conc m s-1 or conc kg m-2 s-1]
83 ! real, dimension(:,:,:), pointer :: diff_cont_xy => NULL() !< convergence of lateral diffusive tracer fluxes
84 ! !! [conc H s-1 ~> conc m s-1 or conc kg m-2 s-1]
85 ! real, dimension(:,:,:), pointer :: diff_conc_xy => NULL() !< convergence of lateral diffusive tracer fluxes
86 ! !! expressed as a change in concentration [conc s-1]
87  real, dimension(:,:,:), pointer :: t_prev => null() !< tracer concentration array at a previous
88  !! timestep used for diagnostics [conc]
89  real, dimension(:,:,:), pointer :: trxh_prev => null() !< layer integrated tracer concentration array
90  !! at a previous timestep used for diagnostics
91 
92  character(len=32) :: name !< tracer name used for diagnostics and error messages
93  character(len=64) :: units !< Physical dimensions of the tracer concentration
94  character(len=240) :: longname !< Long name of the variable
95 ! type(vardesc), pointer :: vd => NULL() !< metadata describing the tracer
96  logical :: registry_diags = .false. !< If true, use the registry to set up the
97  !! diagnostics associated with this tracer.
98  character(len=64) :: cmor_name !< CMOR name of this tracer
99  character(len=64) :: cmor_units !< CMOR physical dimensions of the tracer
100  character(len=240) :: cmor_longname !< CMOR long name of the tracer
101  character(len=32) :: flux_nameroot = "" !< Short tracer name snippet used construct the
102  !! names of flux diagnostics.
103  character(len=64) :: flux_longname = "" !< A word or phrase used construct the long
104  !! names of flux diagnostics.
105  real :: flux_scale= 1.0 !< A scaling factor used to convert the fluxes
106  !! of this tracer to its desired units.
107  character(len=48) :: flux_units = "" !< The units for fluxes of this variable.
108  character(len=48) :: conv_units = "" !< The units for the flux convergence of this tracer.
109  real :: conv_scale = 1.0 !< A scaling factor used to convert the flux
110  !! convergence of this tracer to its desired units.
111  character(len=48) :: cmor_tendprefix = "" !< The CMOR variable prefix for tendencies of this
112  !! tracer, required because CMOR does not follow any
113  !! discernable pattern for these names.
114  integer :: ind_tr_squared = -1 !< The tracer registry index for the square of this tracer
115 
116  !### THESE CAPABILITIES HAVE NOT YET BEEN IMPLEMENTED.
117  ! logical :: advect_tr = .true. !< If true, this tracer should be advected
118  ! logical :: hordiff_tr = .true. !< If true, this tracer should experience epineutral diffusion
119  logical :: remap_tr = .true. !< If true, this tracer should be vertically remapped
120 
121  integer :: diag_form = 1 !< An integer indicating which template is to be used to label diagnostics.
122  !>@{ Diagnostic IDs
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
134  !!@}
135 end type tracer_type
136 
137 !> Type to carry basic tracer information
138 type, public :: tracer_registry_type
139  integer :: ntr = 0 !< number of registered tracers
140  type(tracer_type) :: tr(max_fields_) !< array of registered tracers
141 ! type(diag_ctrl), pointer :: diag !< structure to regulate timing of diagnostics
142  logical :: locked = .false. !< New tracers may be registered if locked=.false.
143  !! When locked=.true., no more tracers can be registered,
144  !! at which point common diagnostics can be set up
145  !! for the registered tracers
146 end type tracer_registry_type
147 
148 contains
149 
150 !> This subroutine registers a tracer to be advected and laterally diffused.
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)
158  type(hor_index_type), intent(in) :: hi !< horizontal index type
159  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
160  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
161  real, dimension(SZI_(HI),SZJ_(HI),SZK_(GV)), &
162  target :: tr_ptr !< target or pointer to the tracer array
163  type(param_file_type), intent(in) :: param_file !< file to parse for model parameter values
164  character(len=*), optional, intent(in) :: name !< Short tracer name
165  character(len=*), optional, intent(in) :: longname !< The long tracer name
166  character(len=*), optional, intent(in) :: units !< The units of this tracer
167  character(len=*), optional, intent(in) :: cmor_name !< CMOR name
168  character(len=*), optional, intent(in) :: cmor_units !< CMOR physical dimensions of variable
169  character(len=*), optional, intent(in) :: cmor_longname !< CMOR long name
170  type(vardesc), optional, intent(in) :: tr_desc !< A structure with metadata about the tracer
171 
172  real, optional, intent(in) :: obc_inflow !< the tracer for all inflows via OBC for which OBC_in_u
173  !! or OBC_in_v are not specified (units of tracer CONC)
174  real, dimension(:,:,:), optional, pointer :: obc_in_u !< tracer at inflows through u-faces of
175  !! tracer cells (units of tracer CONC)
176  real, dimension(:,:,:), optional, pointer :: obc_in_v !< tracer at inflows through v-faces of
177  !! tracer cells (units of tracer CONC)
178 
179  ! The following are probably not necessary if registry_diags is present and true.
180  real, dimension(:,:,:), optional, pointer :: ad_x !< diagnostic x-advective flux
181  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
182  real, dimension(:,:,:), optional, pointer :: ad_y !< diagnostic y-advective flux
183  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
184  real, dimension(:,:,:), optional, pointer :: df_x !< diagnostic x-diffusive flux
185  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
186  real, dimension(:,:,:), optional, pointer :: df_y !< diagnostic y-diffusive flux
187  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
188  real, dimension(:,:), optional, pointer :: ad_2d_x !< vert sum of diagnostic x-advect flux
189  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
190  real, dimension(:,:), optional, pointer :: ad_2d_y !< vert sum of diagnostic y-advect flux
191  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
192  real, dimension(:,:), optional, pointer :: df_2d_x !< vert sum of diagnostic x-diffuse flux
193  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
194  real, dimension(:,:), optional, pointer :: df_2d_y !< vert sum of diagnostic y-diffuse flux
195  !! [conc H L2 T-1 ~> CONC m3 s-1 or CONC kg s-1]
196 
197  real, dimension(:,:,:), optional, pointer :: advection_xy !< convergence of lateral advective tracer fluxes
198  logical, optional, intent(in) :: registry_diags !< If present and true, use the registry for
199  !! the diagnostics of this tracer.
200  character(len=*), optional, intent(in) :: flux_nameroot !< Short tracer name snippet used construct the
201  !! names of flux diagnostics.
202  character(len=*), optional, intent(in) :: flux_longname !< A word or phrase used construct the long
203  !! names of flux diagnostics.
204  character(len=*), optional, intent(in) :: flux_units !< The units for the fluxes of this tracer.
205  real, optional, intent(in) :: flux_scale !< A scaling factor used to convert the fluxes
206  !! of this tracer to its desired units.
207  character(len=*), optional, intent(in) :: convergence_units !< The units for the flux convergence of
208  !! this tracer.
209  real, optional, intent(in) :: convergence_scale !< A scaling factor used to convert the flux
210  !! convergence of this tracer to its desired units.
211  character(len=*), optional, intent(in) :: cmor_tendprefix !< The CMOR name for the layer-integrated
212  !! tendencies of this tracer.
213  integer, optional, intent(in) :: diag_form !< An integer (1 or 2, 1 by default) indicating the
214  !! character string template to use in
215  !! labeling diagnostics
216  type(mom_restart_cs), optional, pointer :: restart_cs !< A pointer to the restart control structure
217  !! this tracer will be registered for
218  !! restarts if this argument is present
219  logical, optional, intent(in) :: mandatory !< If true, this tracer must be read
220  !! from a restart file.
221 
222  logical :: mand
223  type(tracer_type), pointer :: tr=>null()
224  character(len=256) :: mesg ! Message for error messages.
225 
226  if (.not. associated(reg)) call tracer_registry_init(param_file, reg)
227 
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)
232  endif
233  reg%ntr = reg%ntr + 1
234 
235  tr => reg%Tr(reg%ntr)
236 
237  if (present(name)) then
238  tr%name = name
239  tr%longname = name ; if (present(longname)) tr%longname = longname
240  tr%units = "Conc" ; if (present(units)) tr%units = units
241 
242  tr%cmor_name = ""
243  if (present(cmor_name)) tr%cmor_name = cmor_name
244 
245  tr%cmor_units = tr%units
246  if (present(cmor_units)) tr%cmor_units = cmor_units
247 
248  tr%cmor_longname = ""
249  if (present(cmor_longname)) tr%cmor_longname = cmor_longname
250 
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
254  call query_vardesc(tr_desc, name=tr%name, units=tr%units, &
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
258  else
259  call mom_error(fatal,"MOM register_tracer: Either name or "//&
260  "tr_desc must be present when registering a tracer.")
261  endif
262 
263  if (reg%locked) call mom_error(fatal, &
264  "MOM register_tracer was called for variable "//trim(tr%name)//&
265  " with a locked tracer registry.")
266 
267  tr%flux_nameroot = tr%name
268  if (present(flux_nameroot)) then
269  if (len_trim(flux_nameroot) > 0) tr%flux_nameroot = flux_nameroot
270  endif
271 
272  tr%flux_longname = tr%longname
273  if (present(flux_longname)) then
274  if (len_trim(flux_longname) > 0) tr%flux_longname = flux_longname
275  endif
276 
277  tr%flux_units = ""
278  if (present(flux_units)) tr%flux_units = flux_units
279 
280  tr%flux_scale = 1.0
281  if (present(flux_scale)) tr%flux_scale = flux_scale
282 
283  tr%conv_units = ""
284  if (present(convergence_units)) tr%conv_units = convergence_units
285 
286  tr%cmor_tendprefix = ""
287  if (present(cmor_tendprefix)) tr%cmor_tendprefix = cmor_tendprefix
288 
289  tr%conv_scale = 1.0
290  if (present(convergence_scale)) then
291  tr%conv_scale = convergence_scale
292  elseif (present(flux_scale)) then
293  tr%conv_scale = flux_scale
294  endif
295 
296  tr%diag_form = 1
297  if (present(diag_form)) tr%diag_form = diag_form
298 
299  tr%t => tr_ptr
300 
301  if (present(registry_diags)) tr%registry_diags = registry_diags
302 
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
307 ! if (present(OBC_inflow)) Tr%OBC_inflow_conc = OBC_inflow
308 ! if (present(OBC_in_u)) then ; if (associated(OBC_in_u)) &
309 ! Tr%OBC_in_u => OBC_in_u ; endif
310 ! if (present(OBC_in_v)) then ; if (associated(OBC_in_v)) &
311 ! Tr%OBC_in_v => OBC_in_v ; 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
315 
316  if (present(advection_xy)) then ; if (associated(advection_xy)) tr%advection_xy => advection_xy ; endif
317 
318  if (present(restart_cs)) then ; if (associated(restart_cs)) then
319  ! Register this tracer to be read from and written to restart files.
320  mand = .true. ; if (present(mandatory)) mand = mandatory
321 
322  call register_restart_field(tr_ptr, tr%name, mand, restart_cs, &
323  longname=tr%longname, units=tr%units)
324  endif ; endif
325 
326 end subroutine register_tracer
327 
328 
329 !> This subroutine locks the tracer registry to prevent the addition of more
330 !! tracers. After locked=.true., can then register common diagnostics.
331 subroutine lock_tracer_registry(Reg)
332  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
333 
334  if (.not. associated(reg)) call mom_error(warning, &
335  "lock_tracer_registry called with an unassociated registry.")
336 
337  reg%locked = .true.
338 
339 end subroutine lock_tracer_registry
340 
341 !> register_tracer_diagnostics does a set of register_diag_field calls for any previously
342 !! registered in a tracer registry with a value of registry_diags set to .true.
343 subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, US, use_ALE)
344  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
345  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
346  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
347  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
348  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
349  intent(in) :: h !< Layer thicknesses
350  type(time_type), intent(in) :: time !< current model time
351  type(diag_ctrl), intent(in) :: diag !< structure to regulate diagnostic output
352  logical, intent(in) :: use_ale !< If true active diagnostics that only
353  !! apply to ALE configurations
354 
355  character(len=24) :: name ! A variable's name in a NetCDF file.
356  character(len=24) :: shortnm ! A shortened version of a variable's name for
357  ! creating additional diagnostics.
358  character(len=72) :: longname ! The long name of that tracer variable.
359  character(len=72) :: flux_longname ! The tracer name in the long names of fluxes.
360  character(len=48) :: units ! The dimensions of the tracer.
361  character(len=48) :: flux_units ! The units for fluxes, either
362  ! [units] m3 s-1 or [units] kg s-1.
363  character(len=48) :: conv_units ! The units for flux convergences, either
364  ! [units] m2 s-1 or [units] kg s-1.
365  character(len=48) :: unit2 ! The dimensions of the tracer squared
366  character(len=72) :: cmorname ! The CMOR name of this tracer.
367  character(len=120) :: cmor_longname ! The CMOR long name of that variable.
368  character(len=120) :: var_lname ! A temporary longname for a diagnostic.
369  character(len=120) :: cmor_var_lname ! The temporary CMOR long name for a diagnostic
370  character(len=72) :: cmor_varname ! The temporary CMOR name for a diagnostic
371  type(tracer_type), pointer :: tr=>null()
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
377 
378  if (.not. associated(reg)) call mom_error(fatal, "register_tracer_diagnostics: "//&
379  "register_tracer must be called before register_tracer_diagnostics")
380 
381  ntr_in = reg%ntr
382 
383  do m=1,ntr_in ; if (reg%Tr(m)%registry_diags) then
384  tr => reg%Tr(m)
385 ! call query_vardesc(Tr%vd, name, units=units, longname=longname, &
386 ! cmor_field_name=cmorname, cmor_longname=cmor_longname, &
387 ! caller="register_tracer_diagnostics")
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
393 
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
397 
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
401 
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))
405  else
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, &
409  cmor_standard_name=cmor_long_std(cmor_longname))
410  endif
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)
436  else
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, &
446  y_cell_method='sum')
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, &
450  x_cell_method='sum')
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, &
454  y_cell_method='sum')
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, &
458  x_cell_method='sum')
459  endif
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)
466 
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, &
479  y_cell_method='sum')
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, &
484  x_cell_method='sum')
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, &
489  y_cell_method='sum')
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, &
494  x_cell_method='sum')
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')
503 
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)
512 
513  tr%id_adv_xy = register_diag_field('ocean_model', trim(shortnm)//"_advection_xy", &
514  diag%axesTL, time, &
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", &
519  diag%axesT1, time, &
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)
524 
525  tr%id_tendency = register_diag_field('ocean_model', trim(shortnm)//'_tendency', &
526  diag%axesTL, time, &
527  'Net time tendency for '//trim(lowercase(longname)), trim(units)//' s-1', conversion=us%s_to_T)
528 
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
534  endif
535 
536  ! Neutral/Lateral diffusion convergence tendencies
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.)
541 
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')
546 
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.)
550 
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')
555  else
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.)
564 
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')
573 
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.)
578 
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')
583  endif
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)
587 
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)
591 
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)
600  else
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)
614  endif
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
620  endif
621 
622  ! Vertical regridding/remapping tendencies
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)
628 
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)
633 
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)
639 
640  endif
641 
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
648  ! Set up a new tracer for this tracer squared
649  m2 = reg%ntr+1
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
657  ! Augment the total number of tracers, including the squared tracers.
658  reg%ntr = reg%ntr + 1
659  endif
660  endif
661 
662  endif ; enddo
663 
664 end subroutine register_tracer_diagnostics
665 
666 subroutine preale_tracer_diagnostics(Reg, G, GV)
667  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
668  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
669  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
670 
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
673 
674  do m=1,reg%ntr ; if (reg%Tr(m)%ind_tr_squared > 0) then
675  m2 = reg%Tr(m)%ind_tr_squared
676  ! Update squared quantities
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
680  endif ; enddo
681 
682 end subroutine preale_tracer_diagnostics
683 
684 subroutine postale_tracer_diagnostics(Reg, G, GV, diag, dt)
685  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
686  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
687  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
688  type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
689  real, intent(in) :: dt !< total time interval for these diagnostics [T ~> s]
690 
691  real :: work(szi_(g),szj_(g),szk_(g))
692  real :: idt ! The inverse of the time step [T-1 ~> s-1]
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
695 
696  ! The "if" is to avoid NaNs if the diagnostic is called for a zero length interval
697  idt = 0.0 ; if (dt /= 0.0) idt = 1.0 / dt
698 
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))
702  ! Update squared quantities
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)
707  endif ; enddo
708 
709 end subroutine postale_tracer_diagnostics
710 
711 !> post_tracer_diagnostics does post_data calls for any diagnostics that are
712 !! being handled via the tracer registry.
713 subroutine post_tracer_diagnostics(Reg, h, diag_prev, diag, G, GV, dt)
714  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
715  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
716  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
717  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
718  intent(in) :: h !< Layer thicknesses
719  type(diag_grid_storage), intent(in) :: diag_prev !< Contains diagnostic grids from previous timestep
720  type(diag_ctrl), intent(inout) :: diag !< structure to regulate diagnostic output
721  real, intent(in) :: dt !< total time step for tracer updates [T ~> s]
722 
723  real :: work3d(szi_(g),szj_(g),szk_(gv))
724  real :: work2d(szi_(g),szj_(g))
725  real :: idt ! The inverse of the time step [T-1 ~> s-1]
726  type(tracer_type), pointer :: tr=>null()
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
729 
730  idt = 0.; if (dt/=0.) idt = 1.0 / dt ! The "if" is in case the diagnostic is called for a zero length interval
731 
732  ! Tendency diagnostics need to be posted on the grid from the last call to this routine
733  call diag_save_grids(diag)
734  call diag_copy_storage_to_diag(diag, diag_prev)
735  do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
736  tr => reg%Tr(m)
737  if (tr%id_tendency > 0) then
738  work3d(:,:,:) = 0.0
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)
744  endif
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
752  work2d(:,:) = 0.0
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)
757  endif
758  endif
759  endif ; enddo
760  call diag_restore_grids(diag)
761 
762 end subroutine post_tracer_diagnostics
763 
764 !> Post the advective and diffusive tendencies
765 subroutine post_tracer_transport_diagnostics(G, GV, Reg, h_diag, diag)
766  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
767  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
768  type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
769  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
770  intent(in) :: h_diag !< Layer thicknesses on which to post fields
771  type(diag_ctrl), intent(in) :: diag !< structure to regulate diagnostic output
772 
773  integer :: i, j, k, is, ie, js, je, nz, m
774  real :: work2d(szi_(g),szj_(g))
775  type(tracer_type), pointer :: tr=>null()
776 
777  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
778 
779  do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
780  tr => reg%Tr(m)
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
792  work2d(:,:) = 0.0
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)
797  endif
798  endif ; enddo
799 
801 
802 !> This subroutine writes out chksums for tracers.
803 subroutine mom_tracer_chksum(mesg, Tr, ntr, G)
804  character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
805  type(tracer_type), intent(in) :: tr(:) !< array of all of registered tracers
806  integer, intent(in) :: ntr !< number of registered tracers
807  type(ocean_grid_type), intent(in) :: g !< ocean grid structure
808 
809  integer :: is, ie, js, je, nz
810  integer :: i, j, k, m
811 
812  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
813  do m=1,ntr
814  call hchksum(tr(m)%t, mesg//trim(tr(m)%name), g%HI)
815  enddo
816 
817 end subroutine mom_tracer_chksum
818 
819 !> Calculates and prints the global inventory of all tracers in the registry.
820 subroutine mom_tracer_chkinv(mesg, G, h, Tr, ntr)
821  character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
822  type(ocean_grid_type), intent(in) :: g !< ocean grid structure
823  type(tracer_type), dimension(:), intent(in) :: tr !< array of all of registered tracers
824  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses
825  integer, intent(in) :: ntr !< number of registered tracers
826 
827  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: tr_inv !< Tracer inventory
828  real :: total_inv
829  integer :: is, ie, js, je, nz
830  integer :: i, j, k, m
831 
832  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
833  do m=1,ntr
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
839  enddo
840 
841 end subroutine mom_tracer_chkinv
842 
843 !> Find a tracer in the tracer registry by name.
844 subroutine tracer_name_lookup(Reg, tr_ptr, name)
845  type(tracer_registry_type), pointer :: reg !< pointer to tracer registry
846  type(tracer_type), pointer :: tr_ptr !< target or pointer to the tracer array
847  character(len=32), intent(in) :: name !< tracer name
848 
849  integer n
850  do n=1,reg%ntr
851  if (lowercase(reg%Tr(n)%name) == lowercase(name)) tr_ptr => reg%Tr(n)
852  enddo
853 
854 end subroutine tracer_name_lookup
855 
856 !> Initialize the tracer registry.
857 subroutine tracer_registry_init(param_file, Reg)
858  type(param_file_type), intent(in) :: param_file !< open file to parse for model parameters
859  type(tracer_registry_type), pointer :: reg !< pointer to tracer registry
860 
861  integer, save :: init_calls = 0
862 
863 ! This include declares and sets the variable "version".
864 #include "version_variable.h"
865  character(len=40) :: mdl = "MOM_tracer_registry" ! This module's name.
866  character(len=256) :: mesg ! Message for error messages.
867 
868  if (.not.associated(reg)) then ; allocate(reg)
869  else ; return ; endif
870 
871  ! Read all relevant parameters and write them to the model log.
872  call log_version(param_file, mdl, version, "")
873 
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
878  if (is_root_pe()) call mom_error(warning,"MOM_tracer"//mesg)
879  endif
880 
881 end subroutine tracer_registry_init
882 
883 
884 !> This routine closes the tracer registry module.
885 subroutine tracer_registry_end(Reg)
886  type(tracer_registry_type), pointer :: reg !< The tracer registry that will be deallocated
887  if (associated(reg)) deallocate(reg)
888 end subroutine tracer_registry_end
889 
890 end module mom_tracer_registry
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_diag_mediator::diag_restore_grids
subroutine, public diag_restore_grids(diag)
Restore the diagnostic grids from the temporary structure within diag.
Definition: MOM_diag_mediator.F90:3564
mom_tracer_registry::register_tracer
subroutine, public register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, cmor_name, cmor_units, cmor_longname, tr_desc, OBC_inflow, OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, ad_2d_x, ad_2d_y, df_2d_x, df_2d_y, advection_xy, registry_diags, flux_nameroot, flux_longname, flux_units, flux_scale, convergence_units, convergence_scale, cmor_tendprefix, diag_form, restart_CS, mandatory)
This subroutine registers a tracer to be advected and laterally diffused.
Definition: MOM_tracer_registry.F90:158
mom_io::cmor_long_std
character(len=len(longname)) function, public cmor_long_std(longname)
This function returns the CMOR standard name given a CMOR longname, based on the standard pattern of ...
Definition: MOM_io.F90:683
mom_tracer_registry::post_tracer_transport_diagnostics
subroutine, public post_tracer_transport_diagnostics(G, GV, Reg, h_diag, diag)
Post the advective and diffusive tendencies.
Definition: MOM_tracer_registry.F90:766
mom_io::query_vardesc
subroutine, public query_vardesc(vd, name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
This routine queries vardesc.
Definition: MOM_io.F90:699
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_tracer_registry::tracer_type
The tracer type.
Definition: MOM_tracer_registry.F90:38
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_diag_mediator::register_diag_field
integer function, public register_diag_field(module_name, field_name, axes_in, init_time, long_name, units, missing_value, range, mask_variant, standard_name, verbose, do_not_log, err_msg, interp_method, tile_count, cmor_field_name, cmor_long_name, cmor_units, cmor_standard_name, cell_methods, x_cell_method, y_cell_method, v_cell_method, conversion, v_extensive)
Returns the "diag_mediator" handle for a group (native, CMOR, z-coord, ...) of diagnostics derived fr...
Definition: MOM_diag_mediator.F90:1878
mom_string_functions
Handy functions for manipulating strings.
Definition: MOM_string_functions.F90:2
mom_diag_mediator::diag_save_grids
subroutine, public diag_save_grids(diag)
Save the current diagnostic grids in the temporary structure within diag.
Definition: MOM_diag_mediator.F90:3548
mom_tracer_registry::tracer_registry_end
subroutine, public tracer_registry_end(Reg)
This routine closes the tracer registry module.
Definition: MOM_tracer_registry.F90:886
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_tracer_registry
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
Definition: MOM_tracer_registry.F90:5
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_restart::mom_restart_cs
A restart registry and the control structure for restarts.
Definition: MOM_restart.F90:72
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_coms
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Definition: MOM_coms.F90:3
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_tracer_registry::mom_tracer_chkinv
subroutine, public mom_tracer_chkinv(mesg, G, h, Tr, ntr)
Calculates and prints the global inventory of all tracers in the registry.
Definition: MOM_tracer_registry.F90:821
mom_tracer_registry::register_tracer_diagnostics
subroutine, public register_tracer_diagnostics(Reg, h, Time, diag, G, GV, US, use_ALE)
register_tracer_diagnostics does a set of register_diag_field calls for any previously registered in ...
Definition: MOM_tracer_registry.F90:344
mom_diag_mediator::diag_grid_storage
Stores all the remapping grids and the model's native space thicknesses.
Definition: MOM_diag_mediator.F90:146
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_restart
The MOM6 facility for reading and writing restart files, and querying what has been read.
Definition: MOM_restart.F90:2
mom_tracer_registry::post_tracer_diagnostics
subroutine, public post_tracer_diagnostics(Reg, h, diag_prev, diag, G, GV, dt)
post_tracer_diagnostics does post_data calls for any diagnostics that are being handled via the trace...
Definition: MOM_tracer_registry.F90:714
mom_tracer_registry::preale_tracer_diagnostics
subroutine, public preale_tracer_diagnostics(Reg, G, GV)
Definition: MOM_tracer_registry.F90:667
mom_tracer_registry::postale_tracer_diagnostics
subroutine, public postale_tracer_diagnostics(Reg, G, GV, diag, dt)
Definition: MOM_tracer_registry.F90:685
mom_tracer_registry::tracer_registry_init
subroutine, public tracer_registry_init(param_file, Reg)
Initialize the tracer registry.
Definition: MOM_tracer_registry.F90:858
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_tracer_registry::tracer_name_lookup
subroutine, public tracer_name_lookup(Reg, tr_ptr, name)
Find a tracer in the tracer registry by name.
Definition: MOM_tracer_registry.F90:845
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains.
Definition: MOM_hor_index.F90:15
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_restart::register_restart_field
Register fields for restarts.
Definition: MOM_restart.F90:107
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_io::vardesc
Type for describing a variable, typically a tracer.
Definition: MOM_io.F90:53
mom_tracer_registry::mom_tracer_chksum
subroutine, public mom_tracer_chksum(mesg, Tr, ntr, G)
This subroutine writes out chksums for tracers.
Definition: MOM_tracer_registry.F90:804
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_string_functions::lowercase
character(len=len(input_string)) function, public lowercase(input_string)
Return a string in which all uppercase letters have been replaced by their lowercase counterparts.
Definition: MOM_string_functions.F90:24
mom_coms::reproducing_sum
Find an accurate and order-invariant sum of distributed 2d or 3d fields.
Definition: MOM_coms.F90:54
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_tracer_registry::lock_tracer_registry
subroutine, public lock_tracer_registry(Reg)
This subroutine locks the tracer registry to prevent the addition of more tracers....
Definition: MOM_tracer_registry.F90:332
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diag_mediator::diag_copy_storage_to_diag
subroutine, public diag_copy_storage_to_diag(diag, grid_storage)
Copy from the stored diagnostic arrays to the main diagnostic grids.
Definition: MOM_diag_mediator.F90:3530
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239