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_inflow)) Tr%OBC_inflow_conc = OBC_inflow
251! if (present(OBC_in_u)) then ; if (associated(OBC_in_u)) Tr%OBC_in_u => OBC_in_u ; endif
252! if (present(OBC_in_v)) then ; if (associated(OBC_in_v)) Tr%OBC_in_v => OBC_in_v ; endif
253 if (present(ad_2d_x)) then ; if (associated(ad_2d_x)) tr%ad2d_x => ad_2d_x ; endif
254 if (present(ad_2d_y)) then ; if (associated(ad_2d_y)) tr%ad2d_y => ad_2d_y ; endif
255 if (present(df_2d_x)) then ; if (associated(df_2d_x)) tr%df2d_x => df_2d_x ; 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] m2 s-1 or [units] kg 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", trim(flux_units), &
383 v_extensive=.true., &
384 x_cell_method='sum', 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, &
491 x_cell_method='sum', y_cell_method='sum', v_extensive=.true.)
492
493 tr%id_dfxy_cont_2d = register_diag_field("ocean_model", &
494 trim(shortnm)//'_dfxy_cont_tendency_2d', &
495 diag%axesT1, time, "Depth integrated neutral diffusion tracer content "//&
496 "tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
497 x_cell_method='sum', y_cell_method='sum')
498
499 tr%id_hbdxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_hbdxy_cont_tendency', &
500 diag%axesTL, time, "Horizontal boundary diffusion tracer content tendency for "//&
501 trim(shortnm), &
502 conv_units, conversion=tr%conv_scale*us%s_to_T, &
503 x_cell_method='sum', y_cell_method='sum', v_extensive=.true.)
504
505 tr%id_hbdxy_cont_2d = register_diag_field("ocean_model", &
506 trim(shortnm)//'_hbdxy_cont_tendency_2d', &
507 diag%axesT1, time, "Depth integrated horizontal boundary diffusion tracer content "//&
508 "tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
509 x_cell_method='sum', y_cell_method='sum')
510 else
511 cmor_var_lname = 'Tendency of '//trim(lowercase(cmor_longname))//' expressed as '//&
512 trim(lowercase(flux_longname))//&
513 ' content due to parameterized mesoscale neutral diffusion'
514 tr%id_dfxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_dfxy_cont_tendency', &
515 diag%axesTL, time, "Neutral diffusion tracer content tendency for "//trim(shortnm), &
516 conv_units, conversion=tr%conv_scale*us%s_to_T, &
517 cmor_field_name=trim(tr%cmor_tendprefix)//'pmdiff', &
518 cmor_long_name=trim(cmor_var_lname), &
519 cmor_standard_name=trim(cmor_long_std(cmor_var_lname)), &
520 x_cell_method='sum', y_cell_method='sum', v_extensive=.true.)
521
522 cmor_var_lname = 'Tendency of '//trim(lowercase(cmor_longname))//' expressed as '//&
523 trim(lowercase(flux_longname))//&
524 ' content due to parameterized mesoscale neutral diffusion'
525 tr%id_dfxy_cont_2d = register_diag_field("ocean_model", &
526 trim(shortnm)//'_dfxy_cont_tendency_2d', &
527 diag%axesT1, time, "Depth integrated neutral diffusion tracer "//&
528 "content tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
529 cmor_field_name=trim(tr%cmor_tendprefix)//'pmdiff_2d', &
530 cmor_long_name=trim(cmor_var_lname), &
531 cmor_standard_name=trim(cmor_long_std(cmor_var_lname)), &
532 x_cell_method='sum', y_cell_method='sum')
533
534 tr%id_hbdxy_cont = register_diag_field("ocean_model", trim(shortnm)//'_hbdxy_cont_tendency', &
535 diag%axesTL, time, &
536 "Horizontal boundary diffusion tracer content tendency for "//trim(shortnm), &
537 conv_units, conversion=tr%conv_scale*us%s_to_T, &
538 x_cell_method='sum', y_cell_method='sum', v_extensive=.true.)
539
540 tr%id_hbdxy_cont_2d = register_diag_field("ocean_model", &
541 trim(shortnm)//'_hbdxy_cont_tendency_2d', &
542 diag%axesT1, time, "Depth integrated horizontal boundary diffusion of tracer "//&
543 "content tendency for "//trim(shortnm), conv_units, conversion=tr%conv_scale*us%s_to_T, &
544 x_cell_method='sum', y_cell_method='sum')
545 endif
546 tr%id_dfxy_conc = register_diag_field("ocean_model", trim(shortnm)//'_dfxy_conc_tendency', &
547 diag%axesTL, time, "Neutral diffusion tracer concentration tendency for "//trim(shortnm), &
548 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
549
550 tr%id_hbdxy_conc = register_diag_field("ocean_model", trim(shortnm)//'_hbdxy_conc_tendency', &
551 diag%axesTL, time, &
552 "Horizontal diffusion tracer concentration tendency for "//trim(shortnm), &
553 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
554
555 var_lname = "Net time tendency for "//lowercase(flux_longname)
556 if (len_trim(tr%cmor_tendprefix) == 0) then
557 tr%id_trxh_tendency = register_diag_field('ocean_model', trim(shortnm)//'h_tendency', &
558 diag%axesTL, time, var_lname, conv_units, conversion=tr%conv_scale*us%s_to_T, &
559 v_extensive=.true.)
560 tr%id_trxh_tendency_2d = register_diag_field('ocean_model', trim(shortnm)//'h_tendency_2d', &
561 diag%axesT1, time, "Vertical sum of "//trim(lowercase(var_lname)), &
562 conv_units, conversion=tr%conv_scale*us%s_to_T)
563 else
564 cmor_var_lname = "Tendency of "//trim(cmor_longname)//" Expressed as "//&
565 trim(flux_longname)//" Content"
566 tr%id_trxh_tendency = register_diag_field('ocean_model', trim(shortnm)//'h_tendency', &
567 diag%axesTL, time, var_lname, conv_units, conversion=tr%conv_scale*us%s_to_T, &
568 cmor_field_name=trim(tr%cmor_tendprefix)//"tend", &
569 cmor_standard_name=cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname, &
570 v_extensive=.true.)
571 cmor_var_lname = trim(cmor_var_lname)//" Vertical Sum"
572 tr%id_trxh_tendency_2d = register_diag_field('ocean_model', trim(shortnm)//'h_tendency_2d', &
573 diag%axesT1, time, "Vertical sum of "//trim(lowercase(var_lname)), &
574 conv_units, conversion=tr%conv_scale*us%s_to_T, &
575 cmor_field_name=trim(tr%cmor_tendprefix)//"tend_2d", &
576 cmor_standard_name=cmor_long_std(cmor_var_lname), cmor_long_name=cmor_var_lname)
577 endif
578 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0)) then
579 call safe_alloc_ptr(tr%Trxh_prev,isd,ied,jsd,jed,nz)
580 do k=1,nz ; do j=js,je ; do i=is,ie
581 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
582 enddo ; enddo ; enddo
583 endif
584
585 ! Vertical regridding/remapping tendencies
586 if (use_ale .and. tr%remap_tr) then
587 var_lname = "Vertical remapping tracer concentration tendency for "//trim(reg%Tr(m)%name)
588 tr%id_remap_conc= register_diag_field('ocean_model', &
589 trim(tr%flux_nameroot)//'_tendency_vert_remap', diag%axesTL, time, var_lname, &
590 trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
591
592 var_lname = "Vertical remapping tracer content tendency for "//trim(reg%Tr(m)%flux_longname)
593 tr%id_remap_cont = register_diag_field('ocean_model', &
594 trim(tr%flux_nameroot)//'h_tendency_vert_remap', &
595 diag%axesTL, time, var_lname, conv_units, v_extensive=.true., conversion=tr%conv_scale*us%s_to_T)
596
597 var_lname = "Vertical sum of vertical remapping tracer content tendency for "//&
598 trim(reg%Tr(m)%flux_longname)
599 tr%id_remap_cont_2d = register_diag_field('ocean_model', &
600 trim(tr%flux_nameroot)//'h_tendency_vert_remap_2d', &
601 diag%axesT1, time, var_lname, conv_units, conversion=tr%conv_scale*us%s_to_T)
602
603 endif
604
605 if (use_ale .and. (reg%ntr<max_fields_) .and. tr%remap_tr) then
606 unit2 = trim(units)//"2"
607 if (index(units(1:len_trim(units))," ") > 0) unit2 = "("//trim(units)//")2"
608 tr%id_tr_vardec = register_diag_field('ocean_model', trim(shortnm)//"_vardec", diag%axesTL, &
609 time, "ALE variance decay for "//lowercase(longname), &
610 trim(unit2)//" s-1", conversion=tr%conc_scale**2*us%s_to_T)
611 if (tr%id_tr_vardec > 0) then
612 ! Set up a new tracer for this tracer squared
613 m2 = reg%ntr+1
614 tr%ind_tr_squared = m2
615 call safe_alloc_ptr(reg%Tr(m2)%t,isd,ied,jsd,jed,nz) ; reg%Tr(m2)%t(:,:,:) = 0.0
616 reg%Tr(m2)%name = trim(shortnm)//"2"
617 reg%Tr(m2)%longname = "Squared "//trim(longname)
618 reg%Tr(m2)%units = unit2
619 reg%Tr(m2)%registry_diags = .false.
620 reg%Tr(m2)%ind_tr_squared = -1
621 ! Augment the total number of tracers, including the squared tracers.
622 reg%ntr = reg%ntr + 1
623 endif
624 endif
625
626 ! KPP nonlocal term diagnostics
627 if (use_kpp) then
628 tr%id_net_surfflux = register_diag_field('ocean_model', tr%net_surfflux_name, diag%axesT1, time, &
629 tr%net_surfflux_longname, trim(units)//' m s-1', conversion=tr%conc_scale*gv%H_to_m*us%s_to_T)
630 tr%id_NLT_tendency = register_diag_field('ocean_model', "KPP_NLT_d"//trim(shortnm)//"dt", &
631 diag%axesTL, time, &
632 trim(longname)//' tendency due to non-local transport of '//trim(lowercase(flux_longname))//&
633 ', as calculated by [CVMix] KPP', trim(units)//' s-1', conversion=tr%conc_scale*us%s_to_T)
634 if (tr%conv_scale == 0.001*gv%H_to_kg_m2) then
635 conversion = gv%H_to_kg_m2
636 else
637 conversion = tr%conv_scale
638 endif
639 ! We actually want conversion=Tr%conv_scale for all tracers, but introducing the local variable
640 ! 'conversion' and setting it to GV%H_to_kg_m2 instead of 0.001*GV%H_to_kg_m2 for salt tracers
641 ! keeps changes introduced by this refactoring limited to round-off level; as it turns out,
642 ! there is a bug in the code and the NLT budget term for salinity is off by a factor of 10^3
643 ! so introducing the 0.001 here will fix that bug.
644 tr%id_NLT_budget = register_diag_field('ocean_model', tr%NLT_budget_name, &
645 diag%axesTL, time, &
646 trim(flux_longname)//&
647 ' content change due to non-local transport, as calculated by [CVMix] KPP', &
648 conv_units, conversion=conversion*us%s_to_T, v_extensive=.true.)
649 endif
650
651 endif ; enddo
652
653end subroutine register_tracer_diagnostics
654
655subroutine preale_tracer_diagnostics(Reg, G, GV)
656 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
657 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
658 type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
659
660 integer :: i, j, k, is, ie, js, je, nz, m, m2
661 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
662
663 do m=1,reg%ntr ; if (reg%Tr(m)%ind_tr_squared > 0) then
664 m2 = reg%Tr(m)%ind_tr_squared
665 ! Update squared quantities
666 do k=1,nz ; do j=js,je ; do i=is,ie
667 reg%Tr(m2)%T(i,j,k) = reg%Tr(m)%T(i,j,k)**2
668 enddo ; enddo ; enddo
669 endif ; enddo
670
671end subroutine preale_tracer_diagnostics
672
673subroutine postale_tracer_diagnostics(Reg, G, GV, diag, dt)
674 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
675 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
676 type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
677 type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
678 real, intent(in) :: dt !< total time interval for these diagnostics [T ~> s]
679
680 real :: work(szi_(g),szj_(g),szk_(gv)) ! Variance decay [CU2 T-1 ~> conc2 s-1]
681 real :: idt ! The inverse of the time step [T-1 ~> s-1]
682 integer :: i, j, k, is, ie, js, je, nz, m, m2
683 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
684
685 ! The "if" is to avoid NaNs if the diagnostic is called for a zero length interval
686 idt = 0.0 ; if (dt /= 0.0) idt = 1.0 / dt
687
688 do m=1,reg%ntr ; if (reg%Tr(m)%id_tr_vardec > 0) then
689 m2 = reg%Tr(m)%ind_tr_squared
690 if (m2 < 1) call mom_error(fatal, "Bad value of Tr%ind_tr_squared for "//trim(reg%Tr(m)%name))
691 ! Update squared quantities
692 do k=1,nz ; do j=js,je ; do i=is,ie
693 work(i,j,k) = (reg%Tr(m2)%T(i,j,k) - reg%Tr(m)%T(i,j,k)**2) * idt
694 enddo ; enddo ; enddo
695 call post_data(reg%Tr(m)%id_tr_vardec, work, diag)
696 endif ; enddo
697
698end subroutine postale_tracer_diagnostics
699
700!> Post tracer diganostics when that should only be posted when MOM's state
701!! is self-consistent (also referred to as 'synchronized')
702subroutine post_tracer_diagnostics_at_sync(Reg, h, diag_prev, diag, G, GV, dt)
703 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
704 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
705 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
706 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
707 intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
708 type(diag_grid_storage), intent(in) :: diag_prev !< Contains diagnostic grids from previous timestep
709 type(diag_ctrl), intent(inout) :: diag !< structure to regulate diagnostic output
710 real, intent(in) :: dt !< total time step for tracer updates [T ~> s]
711
712 real :: work3d(szi_(g),szj_(g),szk_(gv)) ! The time tendency of a diagnostic [CU T-1 ~> conc s-1]
713 real :: work2d(szi_(g),szj_(g)) ! The vertically integrated time tendency of a diagnostic
714 ! in [CU H T-1 ~> conc m s-1 or conc kg m-2 s-1]
715 real :: idt ! The inverse of the time step [T-1 ~> s-1]
716 type(tracer_type), pointer :: tr=>null()
717 integer :: i, j, k, is, ie, js, je, nz, m
718 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
719
720 idt = 0. ; if (dt/=0.) idt = 1.0 / dt ! The "if" is in case the diagnostic is called for a zero length interval
721
722 ! Tendency diagnostics need to be posted on the grid from the last call to this routine
723 call diag_save_grids(diag)
724 call diag_copy_storage_to_diag(diag, diag_prev)
725 do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
726 tr => reg%Tr(m)
727 if (tr%id_tr > 0) call post_data(tr%id_tr, tr%t, diag)
728 if (tr%id_tendency > 0) then
729 work3d(:,:,:) = 0.0
730 do k=1,nz ; do j=js,je ; do i=is,ie
731 work3d(i,j,k) = (tr%t(i,j,k) - tr%t_prev(i,j,k))*idt
732 tr%t_prev(i,j,k) = tr%t(i,j,k)
733 enddo ; enddo ; enddo
734 call post_data(tr%id_tendency, work3d, diag, alt_h=diag_prev%h_state)
735 endif
736 if ((tr%id_trxh_tendency > 0) .or. (tr%id_trxh_tendency_2d > 0)) then
737 do k=1,nz ; do j=js,je ; do i=is,ie
738 work3d(i,j,k) = (tr%t(i,j,k)*h(i,j,k) - tr%Trxh_prev(i,j,k)) * idt
739 tr%Trxh_prev(i,j,k) = tr%t(i,j,k) * h(i,j,k)
740 enddo ; enddo ; enddo
741 if (tr%id_trxh_tendency > 0) call post_data(tr%id_trxh_tendency, work3d, diag, &
742 alt_h=diag_prev%h_state)
743 if (tr%id_trxh_tendency_2d > 0) then
744 work2d(:,:) = 0.0
745 do k=1,nz ; do j=js,je ; do i=is,ie
746 work2d(i,j) = work2d(i,j) + work3d(i,j,k)
747 enddo ; enddo ; enddo
748 call post_data(tr%id_trxh_tendency_2d, work2d, diag)
749 endif
750 endif
751 endif ; enddo
752 call diag_restore_grids(diag)
753
755
756!> Post the advective and diffusive tendencies
757subroutine post_tracer_transport_diagnostics(G, GV, Reg, h_diag, diag)
758 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
759 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
760 type(tracer_registry_type), pointer :: reg !< pointer to the tracer registry
761 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
762 intent(in) :: h_diag !< Layer thicknesses on which to post fields [H ~> m or kg m-2]
763 type(diag_ctrl), intent(in) :: diag !< structure to regulate diagnostic output
764
765 integer :: i, j, k, is, ie, js, je, nz, m, khi
766 real :: work2d(szi_(g),szj_(g)) ! The vertically integrated convergence of lateral advective
767 ! tracer fluxes [CU H T-1 ~> conc m s-1 or conc kg m-2 s-1]
768 real :: frac_under_100m(szi_(g),szj_(g),szk_(gv)) ! weights used to compute 100m vertical integrals [nondim]
769 real :: ztop(szi_(g),szj_(g)) ! position of the top interface [H ~> m or kg m-2]
770 real :: zbot(szi_(g),szj_(g)) ! position of the bottom interface [H ~> m or kg m-2]
771 type(tracer_type), pointer :: tr=>null()
772
773 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
774
775 ! If any tracers are posting 100m vertical integrals, compute weights
776 frac_under_100m(:,:,:) = 0.0
777 ! khi will be the largest layer index corresponding where ztop < 100m and ztop >= 100m
778 ! in any column (we can reduce computation of 100m integrals by only looping through khi
779 ! rather than GV%ke)
780 khi = 0
781 do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
782 tr => reg%Tr(m)
783 if (tr%id_zint_100m > 0) then
784 zbot(:,:) = 0.0
785 do k=1, nz
786 do j=js,je ; do i=is,ie
787 ztop(i,j) = zbot(i,j)
788 zbot(i,j) = ztop(i,j) + h_diag(i,j,k)*gv%H_to_m
789 if (zbot(i,j) <= 100.0) then
790 frac_under_100m(i,j,k) = 1.0
791 elseif (ztop(i,j) < 100.0) then
792 frac_under_100m(i,j,k) = (100.0 - ztop(i,j)) / (zbot(i,j) - ztop(i,j))
793 else
794 frac_under_100m(i,j,k) = 0.0
795 endif
796 ! frac_under_100m(i,j,k) = max(0, min(1.0, (100.0 - ztop(i,j)) / (zbot(i,j) - ztop(i,j))))
797 enddo ; enddo
798 if (any(frac_under_100m(:,:,k) > 0)) khi = k
799 enddo
800 exit
801 endif
802 endif ; enddo
803
804 do m=1,reg%ntr ; if (reg%Tr(m)%registry_diags) then
805 tr => reg%Tr(m)
806 if (tr%id_tr_post_horzn> 0) call post_data(tr%id_tr_post_horzn, tr%t, diag)
807 if (tr%id_adx > 0) call post_data(tr%id_adx, tr%ad_x, diag, alt_h=h_diag)
808 if (tr%id_ady > 0) call post_data(tr%id_ady, tr%ad_y, diag, alt_h=h_diag)
809 if (tr%id_dfx > 0) call post_data(tr%id_dfx, tr%df_x, diag, alt_h=h_diag)
810 if (tr%id_dfy > 0) call post_data(tr%id_dfy, tr%df_y, diag, alt_h=h_diag)
811 if (tr%id_adx_2d > 0) call post_data(tr%id_adx_2d, tr%ad2d_x, diag)
812 if (tr%id_ady_2d > 0) call post_data(tr%id_ady_2d, tr%ad2d_y, diag)
813 if (tr%id_dfx_2d > 0) call post_data(tr%id_dfx_2d, tr%df2d_x, diag)
814 if (tr%id_dfy_2d > 0) call post_data(tr%id_dfy_2d, tr%df2d_y, diag)
815 if (tr%id_adv_xy > 0) call post_data(tr%id_adv_xy, tr%advection_xy, diag, alt_h=h_diag)
816 if (tr%id_adv_xy_2d > 0) then
817 work2d(:,:) = 0.0
818 do k=1,nz ; do j=js,je ; do i=is,ie
819 work2d(i,j) = work2d(i,j) + tr%advection_xy(i,j,k)
820 enddo ; enddo ; enddo
821 call post_data(tr%id_adv_xy_2d, work2d, diag)
822 endif
823
824 ! A few diagnostics introduce with MARBL driver
825 ! Compute full-depth vertical integral
826 if (tr%id_zint > 0) then
827 work2d(:,:) = 0.0
828 do k=1,nz ; do j=js,je ; do i=is,ie
829 work2d(i,j) = work2d(i,j) + (h_diag(i,j,k)*gv%H_to_m)*tr%t(i,j,k)
830 enddo ; enddo ; enddo
831 call post_data(tr%id_zint, work2d, diag)
832 endif
833
834 ! Compute 100m vertical integral
835 if (tr%id_zint_100m > 0) then
836 work2d(:,:) = 0.0
837 do k=1,khi ; do j=js,je ; do i=is,ie
838 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))
839 enddo ; enddo ; enddo
840 call post_data(tr%id_zint_100m, work2d, diag)
841 endif
842
843 ! Surface values of tracers
844 if (tr%id_SURF > 0) call post_data(tr%id_SURF, tr%t(:,:,1), diag)
845 endif ; enddo
846
848
849!> This subroutine writes out chksums for the first ntr registered tracers.
850subroutine tracer_array_chksum(mesg, Tr, ntr, G)
851 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
852 type(tracer_type), intent(in) :: Tr(:) !< array of all of registered tracers
853 integer, intent(in) :: ntr !< number of registered tracers
854 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
855
856 integer :: m
857
858 do m=1,ntr
859 call hchksum(tr(m)%t, mesg//trim(tr(m)%name), g%HI, unscale=tr(m)%conc_scale)
860 enddo
861
862end subroutine tracer_array_chksum
863
864!> This subroutine writes out chksums for all the registered tracers.
865subroutine tracer_reg_chksum(mesg, Reg, G)
866 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
867 type(tracer_registry_type), pointer :: Reg !< pointer to the tracer registry
868 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
869
870 integer :: m
871
872 if (.not.associated(reg)) return
873
874 do m=1,reg%ntr
875 call hchksum(reg%Tr(m)%t, mesg//trim(reg%Tr(m)%name), g%HI, unscale=reg%Tr(m)%conc_scale)
876 enddo
877
878end subroutine tracer_reg_chksum
879
880!> Calculates and prints the global inventory of the first ntr tracers in the registry.
881subroutine tracer_array_chkinv(mesg, G, GV, h, Tr, ntr)
882 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
883 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
884 type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure
885 type(tracer_type), dimension(:), intent(in) :: Tr !< array of all of registered tracers
886 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
887 integer, intent(in) :: ntr !< number of registered tracers
888
889 ! Local variables
890 real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1] or cell
891 ! masses to kg [kg H-1 L-2 ~> 1], depending on whether the Boussinesq approximation is used
892 real :: tr_inv(SZI_(G),SZJ_(G),SZK_(GV)) ! Volumetric or mass-based tracer inventory in
893 ! each cell [conc m3] or [conc kg]
894 real :: total_inv ! The total amount of tracer [conc m3] or [conc kg]
895 integer :: is, ie, js, je, nz
896 integer :: i, j, k, m
897
898 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
899 vol_scale = gv%H_to_MKS*g%US%L_to_m**2
900 do m=1,ntr
901 do k=1,nz ; do j=js,je ; do i=is,ie
902 tr_inv(i,j,k) = tr(m)%conc_scale*tr(m)%t(i,j,k) * &
903 (vol_scale * h(i,j,k) * g%areaT(i,j)*g%mask2dT(i,j))
904 enddo ; enddo ; enddo
905 total_inv = reproducing_sum(tr_inv, is+(1-g%isd), ie+(1-g%isd), js+(1-g%jsd), je+(1-g%jsd))
906 if (is_root_pe()) write(0,'(A,1X,A5,1X,ES25.16,1X,A)') &
907 "h-point: inventory", tr(m)%name, total_inv, mesg
908 enddo
909
910end subroutine tracer_array_chkinv
911
912
913!> Calculates and prints the global inventory of all tracers in the registry.
914subroutine tracer_reg_chkinv(mesg, G, GV, h, Reg)
915 character(len=*), intent(in) :: mesg !< message that appears on the chksum lines
916 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
917 type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure
918 type(tracer_registry_type), pointer :: Reg !< pointer to the tracer registry
919 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
920
921 ! Local variables
922 real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1] or cell
923 ! masses to kg [kg H-1 L-2 ~> 1], depending on whether the Boussinesq approximation is used
924 real :: tr_inv(SZI_(G),SZJ_(G),SZK_(GV)) ! Volumetric or mass-based tracer inventory in
925 ! each cell [conc m3] or [conc kg]
926 real :: total_inv ! The total amount of tracer [conc m3] or [conc kg]
927 integer :: is, ie, js, je, nz
928 integer :: i, j, k, m
929
930 if (.not.associated(reg)) return
931
932 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
933 vol_scale = gv%H_to_MKS*g%US%L_to_m**2
934 do m=1,reg%ntr
935 do k=1,nz ; do j=js,je ; do i=is,ie
936 tr_inv(i,j,k) = reg%Tr(m)%conc_scale*reg%Tr(m)%t(i,j,k) * &
937 (vol_scale * h(i,j,k) * g%areaT(i,j)*g%mask2dT(i,j))
938 enddo ; enddo ; enddo
939 total_inv = reproducing_sum(tr_inv, is+(1-g%isd), ie+(1-g%isd), js+(1-g%jsd), je+(1-g%jsd))
940 if (is_root_pe()) write(0,'(A,1X,A5,1X,ES25.16,1X,A)') &
941 "h-point: inventory", reg%Tr(m)%name, total_inv, mesg
942 enddo
943
944end subroutine tracer_reg_chkinv
945
946
947!> Find a tracer in the tracer registry by name.
948subroutine tracer_name_lookup(Reg, n, tr_ptr, name)
949 type(tracer_registry_type), pointer :: reg !< pointer to tracer registry
950 type(tracer_type), pointer :: tr_ptr !< target or pointer to the tracer array
951 character(len=32), intent(in) :: name !< tracer name
952 integer, intent(out) :: n !< index to tracer registery
953
954 do n=1,reg%ntr
955 if (lowercase(reg%Tr(n)%name) == lowercase(name)) then
956 tr_ptr => reg%Tr(n)
957 return
958 endif
959 enddo
960
961 call mom_error(fatal,"MOM cannot find registered tracer: "//name)
962
963end subroutine tracer_name_lookup
964
965!> Initialize the tracer registry.
966subroutine tracer_registry_init(param_file, Reg)
967 type(param_file_type), intent(in) :: param_file !< open file to parse for model parameters
968 type(tracer_registry_type), pointer :: reg !< pointer to tracer registry
969
970 integer, save :: init_calls = 0
971
972! This include declares and sets the variable "version".
973#include "version_variable.h"
974 character(len=40) :: mdl = "MOM_tracer_registry" ! This module's name.
975 character(len=256) :: mesg ! Message for error messages.
976
977 if (.not.associated(reg)) then ; allocate(reg)
978 else ; return ; endif
979
980 ! Read all relevant parameters and write them to the model log.
981 call log_version(param_file, mdl, version, "", all_default=.true.)
982
983 init_calls = init_calls + 1
984 if (init_calls > 1) then
985 write(mesg,'("tracer_registry_init called ",I0, &
986 &" times with different registry pointers.")') init_calls
987 if (is_root_pe()) call mom_error(warning,"MOM_tracer "//mesg)
988 endif
989
990end subroutine tracer_registry_init
991
992
993!> This routine closes the tracer registry module.
994subroutine tracer_registry_end(Reg)
995 type(tracer_registry_type), pointer :: reg !< The tracer registry that will be deallocated
996 if (associated(reg)) deallocate(reg)
997end subroutine tracer_registry_end
998
999end module mom_tracer_registry