MOM_tracer_registry.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> This module contains subroutines that handle registration of tracers
6!! and related subroutines. The primary subroutine, register_tracer, is
7!! called to indicate the tracers advected and diffused.
8!! It also makes public the types defined in MOM_tracer_types.
10
11! use MOM_diag_mediator, only : diag_ctrl
12use mom_coms, only : reproducing_sum
13use mom_debugging, only : hchksum
14use mom_diag_mediator, only : diag_ctrl, register_diag_field, post_data, safe_alloc_ptr
17use mom_error_handler, only : mom_error, fatal, warning, mom_mesg, is_root_pe
18use mom_file_parser, only : get_param, log_version, param_file_type
20use mom_grid, only : ocean_grid_type
22use mom_restart, only : register_restart_field, mom_restart_cs
24use mom_time_manager, only : time_type
28
29implicit none ; private
30
31#include <MOM_memory.h>
32
33public register_tracer
34public mom_tracer_chksum, mom_tracer_chkinv
41
42!> Write out checksums for registered tracers
43interface mom_tracer_chksum
45end interface mom_tracer_chksum
46
47!> Calculate and print the global inventories of registered tracers
48interface mom_tracer_chkinv
50end interface mom_tracer_chkinv
51
52contains
53
54!> This subroutine registers a tracer to be advected and laterally diffused.
55subroutine register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, &
56 cmor_name, cmor_units, cmor_longname, net_surfflux_name, &
57 NLT_budget_name, net_surfflux_longname, tr_desc, OBC_inflow, &
58 OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, ad_2d_x, ad_2d_y, &
59 df_2d_x, df_2d_y, advection_xy, registry_diags, &
60 conc_scale, flux_nameroot, flux_longname, flux_units, flux_scale, &
61 convergence_units, convergence_scale, cmor_tendprefix, diag_form, &
62 restart_CS, mandatory, underflow_conc, Tr_out, advect_scheme)
63 type(hor_index_type), intent(in) :: hi !< horizontal index type
64 type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
65 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
66 real, dimension(SZI_(HI),SZJ_(HI),SZK_(GV)), &
67 target :: tr_ptr !< target or pointer to the tracer array [CU ~> conc]
68 type(param_file_type), intent(in) :: param_file !< file to parse for model parameter values
69 character(len=*), optional, intent(in) :: name !< Short tracer name
70 character(len=*), optional, intent(in) :: longname !< The long tracer name
71 character(len=*), optional, intent(in) :: units !< The units of this tracer
72 character(len=*), optional, intent(in) :: cmor_name !< CMOR name
73 character(len=*), optional, intent(in) :: cmor_units !< CMOR physical dimensions of variable
74 character(len=*), optional, intent(in) :: cmor_longname !< CMOR long name
75 character(len=*), optional, intent(in) :: net_surfflux_name !< Name for net_surfflux diag
76 character(len=*), optional, intent(in) :: nlt_budget_name !< Name for NLT_budget diag
77 character(len=*), optional, intent(in) :: net_surfflux_longname !< Long name for net_surfflux diag
78 type(vardesc), optional, intent(in) :: tr_desc !< A structure with metadata about the tracer
79
80 real, optional, intent(in) :: obc_inflow !< the tracer for all inflows via OBC for which OBC_in_u
81 !! or OBC_in_v are not specified [CU ~> conc]
82 real, dimension(:,:,:), optional, pointer :: obc_in_u !< tracer at inflows through u-faces of
83 !! tracer cells [CU ~> conc]
84 real, dimension(:,:,:), optional, pointer :: obc_in_v !< tracer at inflows through v-faces of
85 !! tracer cells [CU ~> conc]
86
87 ! The following are probably not necessary if registry_diags is present and true.
88 real, dimension(:,:,:), optional, pointer :: ad_x !< diagnostic x-advective flux
89 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
90 real, dimension(:,:,:), optional, pointer :: ad_y !< diagnostic y-advective flux
91 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
92 real, dimension(:,:,:), optional, pointer :: df_x !< diagnostic x-diffusive flux
93 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
94 real, dimension(:,:,:), optional, pointer :: df_y !< diagnostic y-diffusive flux
95 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
96 real, dimension(:,:), optional, pointer :: ad_2d_x !< vert sum of diagnostic x-advect flux
97 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
98 real, dimension(:,:), optional, pointer :: ad_2d_y !< vert sum of diagnostic y-advect flux
99 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
100 real, dimension(:,:), optional, pointer :: df_2d_x !< vert sum of diagnostic x-diffuse flux
101 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
102 real, dimension(:,:), optional, pointer :: df_2d_y !< vert sum of diagnostic y-diffuse flux
103 !! [CU H L2 T-1 ~> conc m3 s-1 or conc kg s-1]
104
105 real, dimension(:,:,:), optional, pointer :: advection_xy !< convergence of lateral advective tracer fluxes
106 !! [CU H T-1 ~> conc m s-1 or conc kg m-2 s-1]
107 logical, optional, intent(in) :: registry_diags !< If present and true, use the registry for
108 !! the diagnostics of this tracer.
109 real, optional, intent(in) :: conc_scale !< A scaling factor used to convert the concentration
110 !! of this tracer to its desired units [conc CU-1 ~> 1]
111 character(len=*), optional, intent(in) :: flux_nameroot !< Short tracer name snippet used construct the
112 !! names of flux diagnostics.
113 character(len=*), optional, intent(in) :: flux_longname !< A word or phrase used construct the long
114 !! names of flux diagnostics.
115 character(len=*), optional, intent(in) :: flux_units !< The units for the fluxes of this tracer.
116 real, optional, intent(in) :: flux_scale !< A scaling factor used to convert the fluxes
117 !! of this tracer to its desired units
118 !! [conc m CU-1 H-1 ~> 1] or [conc kg m-2 CU-1 H-1 ~> 1]
119 character(len=*), optional, intent(in) :: convergence_units !< The units for the flux convergence of
120 !! this tracer.
121 real, optional, intent(in) :: convergence_scale !< A scaling factor used to convert the flux
122 !! convergence of this tracer to its desired units.
123 !! [conc m CU-1 H-1 ~> 1] or [conc kg m-2 CU-1 H-1 ~> 1]
124 character(len=*), optional, intent(in) :: cmor_tendprefix !< The CMOR name for the layer-integrated
125 !! tendencies of this tracer.
126 integer, optional, intent(in) :: diag_form !< An integer (1 or 2, 1 by default) indicating the
127 !! character string template to use in
128 !! labeling diagnostics
129 type(mom_restart_cs), optional, intent(inout) :: restart_cs !< MOM restart control struct
130 logical, optional, intent(in) :: mandatory !< If true, this tracer must be read
131 !! from a restart file.
132 real, optional, intent(in) :: underflow_conc !< A tiny concentration, below which the tracer
133 !! concentration underflows to 0 [CU ~> conc].
134 type(tracer_type), optional, pointer :: tr_out !< If present, returns pointer into registry
135
136 integer, optional, intent(in) :: advect_scheme !< Advection scheme for this tracer, the default is -1
137 !! indicating to use the scheme from MOM_tracer_advect
138
139 logical :: mand
140 type(tracer_type), pointer :: tr=>null()
141 character(len=256) :: mesg ! Message for error messages.
142
143 if (.not. associated(reg)) call tracer_registry_init(param_file, reg)
144
145 if (reg%ntr>=max_fields_) then
146 write(mesg,'("Increase MAX_FIELDS_ in MOM_memory.h to at least ",I0," to allow for &
147 &all the tracers being registered via register_tracer.")') reg%ntr+1
148 call mom_error(fatal,"MOM register_tracer: "//mesg)
149 endif
150 reg%ntr = reg%ntr + 1
151
152 tr => reg%Tr(reg%ntr)
153 if (present(tr_out)) tr_out => reg%Tr(reg%ntr)
154
155 if (present(name)) then
156 tr%name = name
157 tr%longname = name ; if (present(longname)) tr%longname = longname
158 tr%units = "Conc" ; if (present(units)) tr%units = units
159
160 tr%cmor_name = ""
161 if (present(cmor_name)) tr%cmor_name = cmor_name
162
163 tr%cmor_units = tr%units
164 if (present(cmor_units)) tr%cmor_units = cmor_units
165
166 tr%cmor_longname = ""
167 if (present(cmor_longname)) tr%cmor_longname = cmor_longname
168
169 if (present(tr_desc)) call mom_error(warning, "MOM register_tracer: "//&
170 "It is a bad idea to use both name and tr_desc when registring "//trim(name))
171 elseif (present(tr_desc)) then
172 call query_vardesc(tr_desc, name=tr%name, units=tr%units, &
173 longname=tr%longname, cmor_field_name=tr%cmor_name, &
174 cmor_longname=tr%cmor_longname, caller="register_tracer")
175 tr%cmor_units = tr%units
176 else
177 call mom_error(fatal,"MOM register_tracer: Either name or "//&
178 "tr_desc must be present when registering a tracer.")
179 endif
180
181 if (reg%locked) call mom_error(fatal, &
182 "MOM register_tracer was called for variable "//trim(tr%name)//&
183 " with a locked tracer registry.")
184
185 tr%conc_scale = 1.0
186 if (present(conc_scale)) tr%conc_scale = conc_scale
187
188 tr%conc_underflow = 0.0
189 if (present(underflow_conc)) tr%conc_underflow = underflow_conc
190
191 tr%flux_nameroot = tr%name
192 if (present(flux_nameroot)) then
193 if (len_trim(flux_nameroot) > 0) tr%flux_nameroot = flux_nameroot
194 endif
195
196 tr%flux_longname = tr%longname
197 if (present(flux_longname)) then
198 if (len_trim(flux_longname) > 0) tr%flux_longname = flux_longname
199 endif
200
201 tr%net_surfflux_name = "KPP_net"//trim(tr%name)
202 if (present(net_surfflux_name)) then
203 tr%net_surfflux_name = net_surfflux_name
204 endif
205
206 tr%NLT_budget_name = 'KPP_NLT_'//trim(tr%flux_nameroot)//'_budget'
207 if (present(nlt_budget_name)) then
208 tr%NLT_budget_name = nlt_budget_name
209 endif
210
211 tr%net_surfflux_longname = 'Effective net surface '//trim(lowercase(tr%flux_longname))//&
212 ' flux, as used by [CVMix] KPP'
213 if (present(net_surfflux_longname)) then
214 tr%net_surfflux_longname = net_surfflux_longname
215 endif
216
217 tr%flux_units = ""
218 if (present(flux_units)) tr%flux_units = flux_units
219
220 tr%flux_scale = gv%H_to_MKS*tr%conc_scale
221 if (present(flux_scale)) tr%flux_scale = flux_scale
222
223 tr%conv_units = ""
224 if (present(convergence_units)) tr%conv_units = convergence_units
225
226 tr%cmor_tendprefix = ""
227 if (present(cmor_tendprefix)) tr%cmor_tendprefix = cmor_tendprefix
228
229 tr%conv_scale = gv%H_to_MKS*tr%conc_scale
230 if (present(convergence_scale)) then
231 tr%conv_scale = convergence_scale
232 elseif (present(flux_scale)) then
233 tr%conv_scale = flux_scale
234 endif
235
236 tr%diag_form = 1
237 if (present(diag_form)) tr%diag_form = diag_form
238
239 tr%advect_scheme = -1
240 if (present(advect_scheme)) tr%advect_scheme = advect_scheme
241
242 tr%t => tr_ptr
243
244 if (present(registry_diags)) tr%registry_diags = registry_diags
245
246 if (present(ad_x)) then ; if (associated(ad_x)) tr%ad_x => ad_x ; endif
247 if (present(ad_y)) then ; if (associated(ad_y)) tr%ad_y => ad_y ; endif
248 if (present(df_x)) then ; if (associated(df_x)) tr%df_x => df_x ; endif
249 if (present(df_y)) then ; if (associated(df_y)) tr%df_y => df_y ; endif
250! if (present(OBC_in_u)) then ; if (associated(OBC_in_u)) Tr%OBC_in_u => OBC_in_u ; endif
251! if (present(OBC_in_v)) then ; if (associated(OBC_in_v)) Tr%OBC_in_v => OBC_in_v ; endif
252 if (present(ad_2d_x)) then ; if (associated(ad_2d_x)) tr%ad2d_x => ad_2d_x ; endif
253 if (present(ad_2d_y)) then ; if (associated(ad_2d_y)) tr%ad2d_y => ad_2d_y ; endif
254 if (present(df_2d_x)) then ; if (associated(df_2d_x)) tr%df2d_x => df_2d_x ; endif
255 if (present(df_2d_y)) then ; if (associated(df_2d_y)) tr%df2d_y => df_2d_y ; endif
256
257 if (present(advection_xy)) then
258 if (associated(advection_xy)) tr%advection_xy => advection_xy
259 endif
260
261 if (present(restart_cs)) then
262 ! Register this tracer to be read from and written to restart files.
263 mand = .true. ; if (present(mandatory)) mand = mandatory
264
265 call register_restart_field(tr_ptr, tr%name, mand, restart_cs, &
266 longname=tr%longname, units=tr%units, conversion=conc_scale)
267 endif
268end subroutine register_tracer
269
270
271!> This subroutine locks the tracer registry to prevent the addition of more
272!! tracers. After locked=.true., can then register common diagnostics.
273subroutine lock_tracer_registry(Reg)
274 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
275
276 if (.not. associated(reg)) call mom_error(warning, &
277 "lock_tracer_registry called with an unassociated registry.")
278
279 reg%locked = .true.
280
281end subroutine lock_tracer_registry
282
283!> register_tracer_diagnostics does a set of register_diag_field calls for any previously
284!! registered in a tracer registry with a value of registry_diags set to .true.
285subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, US, use_ALE, use_KPP)
286 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
287 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
288 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
289 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
290 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
291 intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
292 type(time_type), intent(in) :: time !< current model time
293 type(diag_ctrl), intent(in) :: diag !< structure to regulate diagnostic output
294 logical, intent(in) :: use_ale !< If true active diagnostics that only
295 !! apply to ALE configurations
296 logical, intent(in) :: use_kpp !< If true active diagnostics that only
297 !! apply to CVMix KPP mixings
298
299 character(len=24) :: name ! A variable's name in a NetCDF file.
300 character(len=24) :: shortnm ! A shortened version of a variable's name for
301 ! creating additional diagnostics.
302 character(len=72) :: longname ! The long name of that tracer variable.
303 character(len=72) :: flux_longname ! The tracer name in the long names of fluxes.
304 character(len=48) :: units ! The dimensions of the tracer.
305 character(len=48) :: flux_units ! The units for fluxes, either
306 ! [units] m3 s-1 or [units] kg s-1.
307 character(len=48) :: conv_units ! The units for flux convergences, either
308 ! [units] m s-1 or [units] kg m-2 s-1.
309 character(len=48) :: unit2 ! The dimensions of the tracer squared
310 character(len=72) :: cmorname ! The CMOR name of this tracer.
311 character(len=120) :: cmor_longname ! The CMOR long name of that variable.
312 character(len=120) :: var_lname ! A temporary longname for a diagnostic.
313 character(len=120) :: cmor_var_lname ! The temporary CMOR long name for a diagnostic
314 real :: conversion ! Temporary term while we address a bug [conc m CU-1 H-1 ~> 1] or [conc kg m-2 CU-1 H-1 ~> 1]
315 type(tracer_type), pointer :: tr=>null()
316 integer :: i, j, k, is, ie, js, je, nz, m, m2, ntr_in
317 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
318 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
319 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
320 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
321
322 if (.not. associated(reg)) call mom_error(fatal, "register_tracer_diagnostics: "//&
323 "register_tracer must be called before register_tracer_diagnostics")
324
325 ntr_in = reg%ntr
326
327 do m=1,ntr_in ; if (reg%Tr(m)%registry_diags) then
328 tr => reg%Tr(m)
329! call query_vardesc(Tr%vd, name, units=units, longname=longname, &
330! cmor_field_name=cmorname, cmor_longname=cmor_longname, &
331! caller="register_tracer_diagnostics")
332 name = tr%name ; units=adjustl(tr%units) ; longname = tr%longname
333 cmorname = tr%cmor_name ; cmor_longname = tr%cmor_longname
334 shortnm = tr%flux_nameroot
335 flux_longname = tr%flux_longname
336 if (len_trim(cmor_longname) == 0) cmor_longname = longname
337
338 if (len_trim(tr%flux_units) > 0) then ; flux_units = tr%flux_units
339 elseif (gv%Boussinesq) then ; flux_units = trim(units)//" m3 s-1"
340 else ; flux_units = trim(units)//" kg s-1" ; endif
341
342 if (len_trim(tr%conv_units) > 0) then ; conv_units = tr%conv_units
343 elseif (gv%Boussinesq) then ; conv_units = trim(units)//" m s-1"
344 else ; conv_units = trim(units)//" kg m-2 s-1" ; endif
345
346 if (len_trim(cmorname) == 0) then
347 tr%id_tr = register_diag_field("ocean_model", trim(name), diag%axesTL, &
348 time, trim(longname), trim(units), conversion=tr%conc_scale)
349 else
350 tr%id_tr = register_diag_field("ocean_model", trim(name), diag%axesTL, &
351 time, trim(longname), trim(units), conversion=tr%conc_scale, &
352 cmor_field_name=cmorname, cmor_long_name=cmor_longname, &
353 cmor_units=tr%cmor_units, cmor_standard_name=cmor_long_std(cmor_longname))
354 endif
355 tr%id_tr_post_horzn = register_diag_field("ocean_model", &
356 trim(name)//"_post_horzn", diag%axesTL, time, &
357 trim(longname)//" after horizontal transport (advection/diffusion) has occurred", &
358 trim(units), conversion=tr%conc_scale)
359 if (tr%diag_form == 1) then
360 tr%id_adx = register_diag_field("ocean_model", trim(shortnm)//"_adx", &
361 diag%axesCuL, time, trim(flux_longname)//" advective zonal flux" , &
362 trim(flux_units), v_extensive=.true., y_cell_method='sum', &
363 conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T)
364 tr%id_ady = register_diag_field("ocean_model", trim(shortnm)//"_ady", &
365 diag%axesCvL, time, trim(flux_longname)//" advective meridional flux" , &
366 trim(flux_units), v_extensive=.true., x_cell_method='sum', &
367 conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T)
368 tr%id_dfx = register_diag_field("ocean_model", trim(shortnm)//"_dfx", &
369 diag%axesCuL, time, trim(flux_longname)//" diffusive zonal flux" , &
370 trim(flux_units), v_extensive=.true., y_cell_method='sum', &
371 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
372 tr%id_dfy = register_diag_field("ocean_model", trim(shortnm)//"_dfy", &
373 diag%axesCvL, time, trim(flux_longname)//" diffusive meridional flux" , &
374 trim(flux_units), v_extensive=.true., x_cell_method='sum', &
375 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
376 tr%id_hbd_dfx = register_diag_field("ocean_model", trim(shortnm)//"_hbd_diffx", &
377 diag%axesCuL, time, trim(flux_longname)//" diffusive zonal flux " //&
378 "from the horizontal boundary diffusion scheme", trim(flux_units), v_extensive=.true., &
379 y_cell_method='sum', conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
380 tr%id_hbd_dfy = register_diag_field("ocean_model", trim(shortnm)//"_hbd_diffy", &
381 diag%axesCvL, time, trim(flux_longname)//" diffusive meridional " //&
382 "flux from the horizontal boundary diffusion scheme", &
383 trim(flux_units), v_extensive=.true., x_cell_method='sum', &
384 conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T)
385 else
386 tr%id_adx = register_diag_field("ocean_model", trim(shortnm)//"_adx", &
387 diag%axesCuL, time, "Advective (by residual mean) Zonal Flux of "//trim(flux_longname), &
388 flux_units, v_extensive=.true., conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, y_cell_method='sum')
389 tr%id_ady = register_diag_field("ocean_model", trim(shortnm)//"_ady", &
390 diag%axesCvL, time, "Advective (by residual mean) Meridional Flux of "//trim(flux_longname), &
391 flux_units, v_extensive=.true., conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, x_cell_method='sum')
392 tr%id_dfx = register_diag_field("ocean_model", trim(shortnm)//"_diffx", &
393 diag%axesCuL, time, "Diffusive Zonal Flux of "//trim(flux_longname), &
394 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
395 y_cell_method='sum')
396 tr%id_dfy = register_diag_field("ocean_model", trim(shortnm)//"_diffy", &
397 diag%axesCvL, time, "Diffusive Meridional Flux of "//trim(flux_longname), &
398 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
399 x_cell_method='sum')
400 tr%id_hbd_dfx = register_diag_field("ocean_model", trim(shortnm)//"_hbd_diffx", &
401 diag%axesCuL, time, &
402 "Horizontal Boundary Diffusive Zonal Flux of "//trim(flux_longname), &
403 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
404 y_cell_method='sum')
405 tr%id_hbd_dfy = register_diag_field("ocean_model", trim(shortnm)//"_hbd_diffy", &
406 diag%axesCvL, time, &
407 "Horizontal Boundary Diffusive Meridional Flux of "//trim(flux_longname), &
408 flux_units, v_extensive=.true., conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
409 x_cell_method='sum')
410 endif
411 tr%id_zint = register_diag_field("ocean_model", trim(shortnm)//"_zint", &
412 diag%axesT1, time, &
413 "Thickness-weighted integral of " // trim(longname), &
414 trim(units) // " m")
415 tr%id_zint_100m = register_diag_field("ocean_model", trim(shortnm)//"_zint_100m", &
416 diag%axesT1, time, &
417 "Thickness-weighted integral of "// trim(longname) // " over top 100m", &
418 trim(units) // " m")
419 tr%id_surf = register_diag_field("ocean_model", trim(shortnm)//"_SURF", &
420 diag%axesT1, time, "Surface values of "// trim(longname), trim(units))
421 if (tr%id_adx > 0) call safe_alloc_ptr(tr%ad_x,isdb,iedb,jsd,jed,nz)
422 if (tr%id_ady > 0) call safe_alloc_ptr(tr%ad_y,isd,ied,jsdb,jedb,nz)
423 if (tr%id_dfx > 0) call safe_alloc_ptr(tr%df_x,isdb,iedb,jsd,jed,nz)
424 if (tr%id_dfy > 0) call safe_alloc_ptr(tr%df_y,isd,ied,jsdb,jedb,nz)
425 if (tr%id_hbd_dfx > 0) call safe_alloc_ptr(tr%hbd_dfx,isdb,iedb,jsd,jed,nz)
426 if (tr%id_hbd_dfy > 0) call safe_alloc_ptr(tr%hbd_dfy,isd,ied,jsdb,jedb,nz)
427
428 tr%id_adx_2d = register_diag_field("ocean_model", trim(shortnm)//"_adx_2d", &
429 diag%axesCu1, time, &
430 "Vertically Integrated Advective Zonal Flux of "//trim(flux_longname), &
431 flux_units, conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, y_cell_method='sum')
432 tr%id_ady_2d = register_diag_field("ocean_model", trim(shortnm)//"_ady_2d", &
433 diag%axesCv1, time, &
434 "Vertically Integrated Advective Meridional Flux of "//trim(flux_longname), &
435 flux_units, conversion=tr%flux_scale*(us%L_to_m**2)*us%s_to_T, x_cell_method='sum')
436 tr%id_dfx_2d = register_diag_field("ocean_model", trim(shortnm)//"_diffx_2d", &
437 diag%axesCu1, time, &
438 "Vertically Integrated Diffusive Zonal Flux of "//trim(flux_longname), &
439 flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
440 y_cell_method='sum')
441 tr%id_dfy_2d = register_diag_field("ocean_model", trim(shortnm)//"_diffy_2d", &
442 diag%axesCv1, time, &
443 "Vertically Integrated Diffusive Meridional Flux of "//trim(flux_longname), &
444 flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
445 x_cell_method='sum')
446 tr%id_hbd_dfx_2d = register_diag_field("ocean_model", trim(shortnm)//"_hbd_diffx_2d", &
447 diag%axesCu1, time, "Vertically-integrated zonal diffusive flux from the horizontal boundary diffusion "//&
448 "scheme for "//trim(flux_longname), flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
449 y_cell_method='sum')
450 tr%id_hbd_dfy_2d = register_diag_field("ocean_model", trim(shortnm)//"_hbd_diffy_2d", &
451 diag%axesCv1, time, "Vertically-integrated meridional diffusive flux from the horizontal boundary diffusion "//&
452 "scheme for "//trim(flux_longname), flux_units, conversion=(us%L_to_m**2)*tr%flux_scale*us%s_to_T, &
453 x_cell_method='sum')
454
455 if (tr%id_adx_2d > 0) call safe_alloc_ptr(tr%ad2d_x,isdb,iedb,jsd,jed)
456 if (tr%id_ady_2d > 0) call safe_alloc_ptr(tr%ad2d_y,isd,ied,jsdb,jedb)
457 if (tr%id_dfx_2d > 0) call safe_alloc_ptr(tr%df2d_x,isdb,iedb,jsd,jed)
458 if (tr%id_dfy_2d > 0) call safe_alloc_ptr(tr%df2d_y,isd,ied,jsdb,jedb)
459 if (tr%id_hbd_dfx_2d > 0) call safe_alloc_ptr(tr%hbd_dfx_2d,isdb,iedb,jsd,jed)
460 if (tr%id_hbd_dfy_2d > 0) call safe_alloc_ptr(tr%hbd_dfy_2d,isd,ied,jsdb,jedb)
461
462 tr%id_adv_xy = register_diag_field('ocean_model', trim(shortnm)//"_advection_xy", &
463 diag%axesTL, time, &
464 'Horizontal convergence of residual mean advective fluxes of '//&
465 trim(lowercase(flux_longname)), &
466 conv_units, v_extensive=.true., conversion=tr%conv_scale*us%s_to_T)
467 tr%id_adv_xy_2d = register_diag_field('ocean_model', trim(shortnm)//"_advection_xy_2d", &
468 diag%axesT1, time, &
469 'Vertical sum of horizontal convergence of residual mean advective fluxes of '//&
470 trim(lowercase(flux_longname)), conv_units, conversion=tr%conv_scale*us%s_to_T)
471 if ((tr%id_adv_xy > 0) .or. (tr%id_adv_xy_2d > 0)) &
472 call safe_alloc_ptr(tr%advection_xy,isd,ied,jsd,jed,nz)
473
474 tr%id_tendency = register_diag_field('ocean_model', trim(shortnm)//'_tendency', &
475 diag%axesTL, time, &
476 'Net time tendency for '//trim(lowercase(longname)), &
477 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
478
479 if (tr%id_tendency > 0) then
480 call safe_alloc_ptr(tr%t_prev,isd,ied,jsd,jed,nz)
481 do k=1,nz ; do j=js,je ; do i=is,ie
482 tr%t_prev(i,j,k) = tr%t(i,j,k)
483 enddo ; enddo ; enddo
484 endif
485
486 ! Neutral/Horizontal diffusion convergence tendencies
487 if (tr%diag_form == 1) then
488 tr%id_dfxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_dfxy_cont_tendency', &
489 diag%axesTL, time, "Neutral diffusion tracer content tendency for "//trim(shortnm), &
490 conv_units, conversion=tr%conv_scale*us%s_to_T, v_extensive=.true.)
491
492 tr%id_dfxy_cont_2d = register_diag_field("ocean_model", &
493 trim(shortnm)//'_dfxy_cont_tendency_2d', &
494 diag%axesT1, time, "Depth integrated neutral diffusion tracer content "//&
495 "tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T)
496
497 tr%id_hbdxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_hbdxy_cont_tendency', &
498 diag%axesTL, time, "Horizontal boundary diffusion tracer content tendency for "//&
499 trim(shortnm), &
500 conv_units, conversion=tr%conv_scale*us%s_to_T, v_extensive=.true.)
501
502 tr%id_hbdxy_cont_2d = register_diag_field("ocean_model", &
503 trim(shortnm)//'_hbdxy_cont_tendency_2d', &
504 diag%axesT1, time, "Depth integrated horizontal boundary diffusion tracer content "//&
505 "tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T)
506 else
507 cmor_var_lname = 'Tendency of '//trim(lowercase(cmor_longname))//' expressed as '//&
508 trim(lowercase(flux_longname))//&
509 ' content due to parameterized mesoscale neutral diffusion'
510 tr%id_dfxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_dfxy_cont_tendency', &
511 diag%axesTL, time, "Neutral diffusion tracer content tendency for "//trim(shortnm), &
512 conv_units, conversion=tr%conv_scale*us%s_to_T, v_extensive=.true., &
513 cmor_field_name=trim(tr%cmor_tendprefix)//'pmdiff', &
514 cmor_long_name=trim(cmor_var_lname), &
515 cmor_standard_name=trim(cmor_long_std(cmor_var_lname)))
516
517 cmor_var_lname = 'Tendency of '//trim(lowercase(cmor_longname))//' expressed as '//&
518 trim(lowercase(flux_longname))//&
519 ' content due to parameterized mesoscale neutral diffusion'
520 tr%id_dfxy_cont_2d = register_diag_field("ocean_model", &
521 trim(shortnm)//'_dfxy_cont_tendency_2d', &
522 diag%axesT1, time, "Depth integrated neutral diffusion tracer "//&
523 "content tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
524 cmor_field_name=trim(tr%cmor_tendprefix)//'pmdiff_2d', &
525 cmor_long_name=trim(cmor_var_lname), &
526 cmor_standard_name=trim(cmor_long_std(cmor_var_lname)))
527
528 tr%id_hbdxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_hbdxy_cont_tendency', &
529 diag%axesTL, time, &
530 "Horizontal boundary diffusion tracer content tendency for "//trim(shortnm), &
531 conv_units, conversion=tr%conv_scale*us%s_to_T, v_extensive=.true.)
532
533 tr%id_hbdxy_cont_2d = register_diag_field("ocean_model", &
534 trim(shortnm)//'_hbdxy_cont_tendency_2d', &
535 diag%axesT1, time, "Depth integrated horizontal boundary diffusion of tracer "//&
536 "content tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T)
537 endif
538 tr%id_dfxy_conc = register_diag_field("ocean_model", trim(shortnm)//'_dfxy_conc_tendency', &
539 diag%axesTL, time, "Neutral diffusion tracer concentration tendency for "//trim(shortnm), &
540 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
541
542 tr%id_hbdxy_conc = register_diag_field("ocean_model", trim(shortnm)//'_hbdxy_conc_tendency', &
543 diag%axesTL, time, &
544 "Horizontal diffusion tracer concentration tendency for "//trim(shortnm), &
545 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
546
547 var_lname = "Net time tendency for "//lowercase(flux_longname)
548 if (len_trim(tr%cmor_tendprefix) == 0) then
549 tr%id_trxh_tendency = register_diag_field('ocean_model', trim(shortnm)//'h_tendency', &
550 diag%axesTL, time, var_lname, conv_units, conversion=tr%conv_scale*us%s_to_T, &
551 v_extensive=.true.)
552 tr%id_trxh_tendency_2d = register_diag_field('ocean_model', trim(shortnm)//'h_tendency_2d', &
553 diag%axesT1, time, "Vertical sum of "//trim(lowercase(var_lname)), &
554 conv_units, conversion=tr%conv_scale*us%s_to_T)
555 else
556 cmor_var_lname = "Tendency of "//trim(cmor_longname)//" Expressed as "//&
557 trim(flux_longname)//" Content"
558 tr%id_trxh_tendency = register_diag_field('ocean_model', trim(shortnm)//'h_tendency', &
559 diag%axesTL, time, var_lname, conv_units, conversion=tr%conv_scale*us%s_to_T, &
560 cmor_field_name=trim(tr%cmor_tendprefix)//"tend", &
561 cmor_standard_name=cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname, &
562 v_extensive=.true.)
563 cmor_var_lname = trim(cmor_var_lname)//" Vertical Sum"
564 tr%id_trxh_tendency_2d = register_diag_field('ocean_model', trim(shortnm)//'h_tendency_2d', &
565 diag%axesT1, time, "Vertical sum of "//trim(lowercase(var_lname)), &
566 conv_units, conversion=tr%conv_scale*us%s_to_T, &
567 cmor_field_name=trim(tr%cmor_tendprefix)//"tend_2d", &
568 cmor_standard_name=cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname)
569 endif
570 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0)) then
571 call safe_alloc_ptr(tr%Trxh_prev,isd,ied,jsd,jed,nz)
572 do k=1,nz ; do j=js,je ; do i=is,ie
573 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
574 enddo ; enddo ; enddo
575 endif
576
577 ! Vertical regridding/remapping tendencies
578 if (use_ale .and. tr%remap_tr) then
579 var_lname = "Vertical remapping tracer concentration tendency for "//trim(reg%Tr(m)%name)
580 tr%id_remap_conc= register_diag_field('ocean_model', &
581 trim(tr%flux_nameroot)//'_tendency_vert_remap', diag%axesTL, time, var_lname, &
582 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
583
584 var_lname = "Vertical remapping tracer content tendency for "//trim(reg%Tr(m)%flux_longname)
585 tr%id_remap_cont = register_diag_field('ocean_model', &
586 trim(tr%flux_nameroot)//'h_tendency_vert_remap', &
587 diag%axesTL, time, var_lname, conv_units, v_extensive=.true., conversion=tr%conv_scale*us%s_to_T)
588
589 var_lname = "Vertical sum of vertical remapping tracer content tendency for "//&
590 trim(reg%Tr(m)%flux_longname)
591 tr%id_remap_cont_2d = register_diag_field('ocean_model', &
592 trim(tr%flux_nameroot)//'h_tendency_vert_remap_2d', &
593 diag%axesT1, time, var_lname, conv_units, conversion=tr%conv_scale*us%s_to_T)
594
595 endif
596
597 if (use_ale .and. (reg%ntr<max_fields_) .and. tr%remap_tr) then
598 unit2 = trim(units)//"2"
599 if (index(units(1:len_trim(units))," ") > 0) unit2 = "("//trim(units)//")2"
600 tr%id_tr_vardec = register_diag_field('ocean_model', trim(shortnm)//"_vardec", diag%axesTL, &
601 time, "ALE variance decay for "//lowercase(longname), &
602 trim(unit2)//" s-1", conversion=tr%conc_scale**2*us%s_to_T)
603 if (tr%id_tr_vardec > 0) then
604 ! Set up a new tracer for this tracer squared
605 m2 = reg%ntr+1
606 tr%ind_tr_squared = m2
607 call safe_alloc_ptr(reg%Tr(m2)%t,isd,ied,jsd,jed,nz) ; reg%Tr(m2)%t(:,:,:) = 0.0
608 reg%Tr(m2)%name = trim(shortnm)//"2"
609 reg%Tr(m2)%longname = "Squared "//trim(longname)
610 reg%Tr(m2)%units = unit2
611 reg%Tr(m2)%registry_diags = .false.
612 reg%Tr(m2)%ind_tr_squared = -1
613 ! Augment the total number of tracers, including the squared tracers.
614 reg%ntr = reg%ntr + 1
615 endif
616 endif
617
618 ! KPP nonlocal term diagnostics
619 if (use_kpp) then
620 tr%id_net_surfflux = register_diag_field('ocean_model', tr%net_surfflux_name, diag%axesT1, time, &
621 tr%net_surfflux_longname, trim(units)//' m s-1', conversion=tr%conc_scale*gv%H_to_m*us%s_to_T)
622 tr%id_NLT_tendency = register_diag_field('ocean_model', "KPP_NLT_d"//trim(shortnm)//"dt", &
623 diag%axesTL, time, &
624 trim(longname)//' tendency due to non-local transport of '//trim(lowercase(flux_longname))//&
625 ', as calculated by [CVMix] KPP', trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
626 if (tr%conv_scale == 0.001*gv%H_to_kg_m2) then
627 conversion = gv%H_to_kg_m2
628 else
629 conversion = tr%conv_scale
630 endif
631 ! We actually want conversion=Tr%conv_scale for all tracers, but introducing the local variable
632 ! 'conversion' and setting it to GV%H_to_kg_m2 instead of 0.001*GV%H_to_kg_m2 for salt tracers
633 ! keeps changes introduced by this refactoring limited to round-off level; as it turns out,
634 ! there is a bug in the code and the NLT budget term for salinity is off by a factor of 10^3
635 ! so introducing the 0.001 here will fix that bug.
636 tr%id_NLT_budget = register_diag_field('ocean_model', tr%NLT_budget_name, &
637 diag%axesTL, time, &
638 trim(flux_longname)//&
639 ' content change due to non-local transport, as calculated by [CVMix] KPP', &
640 conv_units, conversion=conversion*us%s_to_T, v_extensive=.true.)
641 endif
642
643 endif ; enddo
644
645end subroutine register_tracer_diagnostics
646
647subroutine preale_tracer_diagnostics(Reg, G, GV)
648 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
649 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
650 type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
651
652 integer :: i, j, k, is, ie, js, je, nz, m, m2
653 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
654
655 do m=1,reg%ntr ; if (reg%Tr(m)%ind_tr_squared > 0) then
656 m2 = reg%Tr(m)%ind_tr_squared
657 ! Update squared quantities
658 do k=1,nz ; do j=js,je ; do i=is,ie
659 reg%Tr(m2)%T(i,j,k) = reg%Tr(m)%T(i,j,k)**2
660 enddo ; enddo ; enddo
661 endif ; enddo
662
663end subroutine preale_tracer_diagnostics
664
665subroutine postale_tracer_diagnostics(Reg, G, GV, diag, dt)
666 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
667 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
668 type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
669 type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
670 real, intent(in) :: dt !< total time interval for these diagnostics [T ~> s]
671
672 real :: work(szi_(g),szj_(g),szk_(gv)) ! Variance decay [CU2 T-1 ~> conc2 s-1]
673 real :: idt ! The inverse of the time step [T-1 ~> s-1]
674 integer :: i, j, k, is, ie, js, je, nz, m, m2
675 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
676
677 ! The "if" is to avoid NaNs if the diagnostic is called for a zero length interval
678 idt = 0.0 ; if (dt /= 0.0) idt = 1.0 / dt
679
680 do m=1,reg%ntr ; if (reg%Tr(m)%id_tr_vardec > 0) then
681 m2 = reg%Tr(m)%ind_tr_squared
682 if (m2 < 1) call mom_error(fatal, "Bad value of Tr%ind_tr_squared for "//trim(reg%Tr(m)%name))
683 ! Update squared quantities
684 do k=1,nz ; do j=js,je ; do i=is,ie
685 work(i,j,k) = (reg%Tr(m2)%T(i,j,k) - reg%Tr(m)%T(i,j,k)**2) * idt
686 enddo ; enddo ; enddo
687 call post_data(reg%Tr(m)%id_tr_vardec, work, diag)
688 endif ; enddo
689
690end subroutine postale_tracer_diagnostics
691
692!> Post tracer diganostics when that should only be posted when MOM's state
693!! is self-consistent (also referred to as 'synchronized')
694subroutine post_tracer_diagnostics_at_sync(Reg, h, diag_prev, diag, G, GV, dt)
695 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
696 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
697 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
698 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
699 intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
700 type(diag_grid_storage), intent(in) :: diag_prev !< Contains diagnostic grids from previous timestep
701 type(diag_ctrl), intent(inout) :: diag !< structure to regulate diagnostic output
702 real, intent(in) :: dt !< total time step for tracer updates [T ~> s]
703
704 real :: work3d(szi_(g),szj_(g),szk_(gv)) ! The time tendency of a diagnostic [CU T-1 ~> conc s-1]
705 real :: work2d(szi_(g),szj_(g)) ! The vertically integrated time tendency of a diagnostic
706 ! in [CU H T-1 ~> conc m s-1 or conc kg m-2 s-1]
707 real :: idt ! The inverse of the time step [T-1 ~> s-1]
708 type(tracer_type), pointer :: tr=>null()
709 integer :: i, j, k, is, ie, js, je, nz, m
710 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
711
712 idt = 0. ; if (dt/=0.) idt = 1.0 / dt ! The "if" is in case the diagnostic is called for a zero length interval
713
714 ! Tendency diagnostics need to be posted on the grid from the last call to this routine
715 call diag_save_grids(diag)
716 call diag_copy_storage_to_diag(diag, diag_prev)
717 do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
718 tr => reg%Tr(m)
719 if (tr%id_tr > 0) call post_data(tr%id_tr, tr%t, diag)
720 if (tr%id_tendency > 0) then
721 work3d(:,:,:) = 0.0
722 do k=1,nz ; do j=js,je ; do i=is,ie
723 work3d(i,j,k) = (tr%t(i,j,k) - tr%t_prev(i,j,k))*idt
724 tr%t_prev(i,j,k) = tr%t(i,j,k)
725 enddo ; enddo ; enddo
726 call post_data(tr%id_tendency, work3d, diag, alt_h=diag_prev%h_state)
727 endif
728 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0)) then
729 do k=1,nz ; do j=js,je ; do i=is,ie
730 work3d(i,j,k) = (tr%t(i,j,k)*h(i,j,k) - tr%Trxh_prev(i,j,k)) * idt
731 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
732 enddo ; enddo ; enddo
733 if (tr%id_trxh_tendency > 0) call post_data(tr%id_trxh_tendency, work3d, diag, &
734 alt_h=diag_prev%h_state)
735 if (tr%id_trxh_tendency_2d > 0) then
736 work2d(:,:) = 0.0
737 do k=1,nz ; do j=js,je ; do i=is,ie
738 work2d(i,j) = work2d(i,j) + work3d(i,j,k)
739 enddo ; enddo ; enddo
740 call post_data(tr%id_trxh_tendency_2d, work2d, diag)
741 endif
742 endif
743 endif ; enddo
744 call diag_restore_grids(diag)
745
747
748!> Post the advective and diffusive tendencies
749subroutine post_tracer_transport_diagnostics(G, GV, Reg, h_diag, diag)
750 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
751 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
752 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
753 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
754 intent(in) :: h_diag !< Layer thicknesses on which to post fields [H ~> m or kg m-2]
755 type(diag_ctrl), intent(in) :: diag !< structure to regulate diagnostic output
756
757 integer :: i, j, k, is, ie, js, je, nz, m, khi
758 real :: work2d(szi_(g),szj_(g)) ! The vertically integrated convergence of lateral advective
759 ! tracer fluxes [CU H T-1 ~> conc m s-1 or conc kg m-2 s-1]
760 real :: frac_under_100m(szi_(g),szj_(g),szk_(gv)) ! weights used to compute 100m vertical integrals [nondim]
761 real :: ztop(szi_(g),szj_(g)) ! position of the top interface [H ~> m or kg m-2]
762 real :: zbot(szi_(g),szj_(g)) ! position of the bottom interface [H ~> m or kg m-2]
763 type(tracer_type), pointer :: tr=>null()
764
765 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
766
767 ! If any tracers are posting 100m vertical integrals, compute weights
768 frac_under_100m(:,:,:) = 0.0
769 ! khi will be the largest layer index corresponding where ztop < 100m and ztop >= 100m
770 ! in any column (we can reduce computation of 100m integrals by only looping through khi
771 ! rather than GV%ke)
772 khi = 0
773 do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
774 tr => reg%Tr(m)
775 if (tr%id_zint_100m > 0) then
776 zbot(:,:) = 0.0
777 do k=1, nz
778 do j=js,je ; do i=is,ie
779 ztop(i,j) = zbot(i,j)
780 zbot(i,j) = ztop(i,j) + h_diag(i,j,k)*gv%H_to_m
781 if (zbot(i,j) <= 100.0) then
782 frac_under_100m(i,j,k) = 1.0
783 elseif (ztop(i,j) < 100.0) then
784 frac_under_100m(i,j,k) = (100.0 - ztop(i,j)) / (zbot(i,j) - ztop(i,j))
785 else
786 frac_under_100m(i,j,k) = 0.0
787 endif
788 ! frac_under_100m(i,j,k) = max(0, min(1.0, (100.0 - ztop(i,j)) / (zbot(i,j) - ztop(i,j))))
789 enddo ; enddo
790 if (any(frac_under_100m(:,:,k) > 0)) khi = k
791 enddo
792 exit
793 endif
794 endif ; enddo
795
796 do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
797 tr => reg%Tr(m)
798 if (tr%id_tr_post_horzn> 0) call post_data(tr%id_tr_post_horzn, tr%t, diag)
799 if (tr%id_adx > 0) call post_data(tr%id_adx, tr%ad_x, diag, alt_h=h_diag)
800 if (tr%id_ady > 0) call post_data(tr%id_ady, tr%ad_y, diag, alt_h=h_diag)
801 if (tr%id_dfx > 0) call post_data(tr%id_dfx, tr%df_x, diag, alt_h=h_diag)
802 if (tr%id_dfy > 0) call post_data(tr%id_dfy, tr%df_y, diag, alt_h=h_diag)
803 if (tr%id_adx_2d > 0) call post_data(tr%id_adx_2d, tr%ad2d_x, diag)
804 if (tr%id_ady_2d > 0) call post_data(tr%id_ady_2d, tr%ad2d_y, diag)
805 if (tr%id_dfx_2d > 0) call post_data(tr%id_dfx_2d, tr%df2d_x, diag)
806 if (tr%id_dfy_2d > 0) call post_data(tr%id_dfy_2d, tr%df2d_y, diag)
807 if (tr%id_adv_xy > 0) call post_data(tr%id_adv_xy, tr%advection_xy, diag, alt_h=h_diag)
808 if (tr%id_adv_xy_2d > 0) then
809 work2d(:,:) = 0.0
810 do k=1,nz ; do j=js,je ; do i=is,ie
811 work2d(i,j) = work2d(i,j) + tr%advection_xy(i,j,k)
812 enddo ; enddo ; enddo
813 call post_data(tr%id_adv_xy_2d, work2d, diag)
814 endif
815
816 ! A few diagnostics introduce with MARBL driver
817 ! Compute full-depth vertical integral
818 if (tr%id_zint > 0) then
819 work2d(:,:) = 0.0
820 do k=1,nz ; do j=js,je ; do i=is,ie
821 work2d(i,j) = work2d(i,j) + (h_diag(i,j,k)*gv%H_to_m)*tr%t(i,j,k)
822 enddo ; enddo ; enddo
823 call post_data(tr%id_zint, work2d, diag)
824 endif
825
826 ! Compute 100m vertical integral
827 if (tr%id_zint_100m > 0) then
828 work2d(:,:) = 0.0
829 do k=1,khi ; do j=js,je ; do i=is,ie
830 work2d(i,j) = work2d(i,j) + frac_under_100m(i,j,k)*((h_diag(i,j,k)*gv%H_to_m)*tr%t(i,j,k))
831 enddo ; enddo ; enddo
832 call post_data(tr%id_zint_100m, work2d, diag)
833 endif
834
835 ! Surface values of tracers
836 if (tr%id_SURF > 0) call post_data(tr%id_SURF, tr%t(:,:,1), diag)
837 endif ; enddo
838
840
841!> This subroutine writes out chksums for the first ntr registered tracers.
842subroutine tracer_array_chksum(mesg, Tr, ntr, G)
843 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
844 type(tracer_type), intent(in) :: Tr(:) !< array of all of registered tracers
845 integer, intent(in) :: ntr !< number of registered tracers
846 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
847
848 integer :: m
849
850 do m=1,ntr
851 call hchksum(tr(m)%t, mesg//trim(tr(m)%name), g%HI, unscale=tr(m)%conc_scale)
852 enddo
853
854end subroutine tracer_array_chksum
855
856!> This subroutine writes out chksums for all the registered tracers.
857subroutine tracer_reg_chksum(mesg, Reg, G)
858 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
859 type(tracer_registry_type), pointer :: Reg !< pointer to the tracer registry
860 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
861
862 integer :: m
863
864 if (.not.associated(reg)) return
865
866 do m=1,reg%ntr
867 call hchksum(reg%Tr(m)%t, mesg//trim(reg%Tr(m)%name), g%HI, unscale=reg%Tr(m)%conc_scale)
868 enddo
869
870end subroutine tracer_reg_chksum
871
872!> Calculates and prints the global inventory of the first ntr tracers in the registry.
873subroutine tracer_array_chkinv(mesg, G, GV, h, Tr, ntr)
874 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
875 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
876 type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure
877 type(tracer_type), dimension(:), intent(in) :: Tr !< array of all of registered tracers
878 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
879 integer, intent(in) :: ntr !< number of registered tracers
880
881 ! Local variables
882 real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1] or cell
883 ! masses to kg [kg H-1 L-2 ~> 1], depending on whether the Boussinesq approximation is used
884 real :: tr_inv(SZI_(G),SZJ_(G),SZK_(GV)) ! Volumetric or mass-based tracer inventory in
885 ! each cell [conc m3] or [conc kg]
886 real :: total_inv ! The total amount of tracer [conc m3] or [conc kg]
887 integer :: is, ie, js, je, nz
888 integer :: i, j, k, m
889
890 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
891 vol_scale = gv%H_to_MKS*g%US%L_to_m**2
892 do m=1,ntr
893 do k=1,nz ; do j=js,je ; do i=is,ie
894 tr_inv(i,j,k) = tr(m)%conc_scale*tr(m)%t(i,j,k) * &
895 (vol_scale * h(i,j,k) * g%areaT(i,j)*g%mask2dT(i,j))
896 enddo ; enddo ; enddo
897 total_inv = reproducing_sum(tr_inv, is+(1-g%isd), ie+(1-g%isd), js+(1-g%jsd), je+(1-g%jsd))
898 if (is_root_pe()) write(0,'(A,1X,A5,1X,ES25.16,1X,A)') &
899 "h-point: inventory", tr(m)%name, total_inv, mesg
900 enddo
901
902end subroutine tracer_array_chkinv
903
904
905!> Calculates and prints the global inventory of all tracers in the registry.
906subroutine tracer_reg_chkinv(mesg, G, GV, h, Reg)
907 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
908 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
909 type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure
910 type(tracer_registry_type), pointer :: Reg !< pointer to the tracer registry
911 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
912
913 ! Local variables
914 real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1] or cell
915 ! masses to kg [kg H-1 L-2 ~> 1], depending on whether the Boussinesq approximation is used
916 real :: tr_inv(SZI_(G),SZJ_(G),SZK_(GV)) ! Volumetric or mass-based tracer inventory in
917 ! each cell [conc m3] or [conc kg]
918 real :: total_inv ! The total amount of tracer [conc m3] or [conc kg]
919 integer :: is, ie, js, je, nz
920 integer :: i, j, k, m
921
922 if (.not.associated(reg)) return
923
924 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
925 vol_scale = gv%H_to_MKS*g%US%L_to_m**2
926 do m=1,reg%ntr
927 do k=1,nz ; do j=js,je ; do i=is,ie
928 tr_inv(i,j,k) = reg%Tr(m)%conc_scale*reg%Tr(m)%t(i,j,k) * &
929 (vol_scale * h(i,j,k) * g%areaT(i,j)*g%mask2dT(i,j))
930 enddo ; enddo ; enddo
931 total_inv = reproducing_sum(tr_inv, is+(1-g%isd), ie+(1-g%isd), js+(1-g%jsd), je+(1-g%jsd))
932 if (is_root_pe()) write(0,'(A,1X,A5,1X,ES25.16,1X,A)') &
933 "h-point: inventory", reg%Tr(m)%name, total_inv, mesg
934 enddo
935
936end subroutine tracer_reg_chkinv
937
938
939!> Find a tracer in the tracer registry by name.
940subroutine tracer_name_lookup(Reg, n, tr_ptr, name)
941 type(tracer_registry_type), pointer :: reg !< pointer to tracer registry
942 type(tracer_type), pointer :: tr_ptr !< target or pointer to the tracer array
943 character(len=32), intent(in) :: name !< tracer name
944 integer, intent(out) :: n !< index to tracer registery
945
946 do n=1,reg%ntr
947 if (lowercase(reg%Tr(n)%name) == lowercase(name)) then
948 tr_ptr => reg%Tr(n)
949 return
950 endif
951 enddo
952
953 call mom_error(fatal,"MOM cannot find registered tracer: "//name)
954
955end subroutine tracer_name_lookup
956
957!> Initialize the tracer registry.
958subroutine tracer_registry_init(param_file, Reg)
959 type(param_file_type), intent(in) :: param_file !< open file to parse for model parameters
960 type(tracer_registry_type), pointer :: reg !< pointer to tracer registry
961
962 integer, save :: init_calls = 0
963
964! This include declares and sets the variable "version".
965#include "version_variable.h"
966 character(len=40) :: mdl = "MOM_tracer_registry" ! This module's name.
967 character(len=256) :: mesg ! Message for error messages.
968
969 if (.not.associated(reg)) then ; allocate(reg)
970 else ; return ; endif
971
972 ! Read all relevant parameters and write them to the model log.
973 call log_version(param_file, mdl, version, "", all_default=.true.)
974
975 init_calls = init_calls + 1
976 if (init_calls > 1) then
977 write(mesg,'("tracer_registry_init called ",I0, &
978 &" times with different registry pointers.")') init_calls
979 if (is_root_pe()) call mom_error(warning,"MOM_tracer "//mesg)
980 endif
981
982end subroutine tracer_registry_init
983
984
985!> This routine closes the tracer registry module.
986subroutine tracer_registry_end(Reg)
987 type(tracer_registry_type), pointer :: reg !< The tracer registry that will be deallocated
988 if (associated(reg)) deallocate(reg)
989end subroutine tracer_registry_end
990
991end module mom_tracer_registry