MOM6
ideal_age_example.F90
Go to the documentation of this file.
1 !> A tracer package of ideal age tracers
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
7 use mom_error_handler, only : mom_error, fatal, warning
9 use mom_forcing_type, only : forcing
10 use mom_grid, only : ocean_grid_type
11 use mom_hor_index, only : hor_index_type
16 use mom_time_manager, only : time_type, time_type_to_real
21 use mom_variables, only : surface
23 
24 use coupler_types_mod, only : coupler_type_set_data, ind_csurf
26 
27 implicit none ; private
28 
29 #include <MOM_memory.h>
30 
34 
35 integer, parameter :: ntr_max = 3 !< the maximum number of tracers in this module.
36 
37 !> The control structure for the ideal_age_tracer package
38 type, public :: ideal_age_tracer_cs ; private
39  integer :: ntr !< The number of tracers that are actually used.
40  logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
41  integer :: nkml !< The number of layers in the mixed layer. The ideal
42  !1 age tracers are reset in the top nkml layers.
43  character(len=200) :: ic_file !< The file in which the age-tracer initial values
44  !! can be found, or an empty string for internal initialization.
45  logical :: z_ic_file !< If true, the IC_file is in Z-space. The default is false.
46  type(time_type), pointer :: time => null() !< A pointer to the ocean model's clock.
47  type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the tracer registry
48  real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this package, in g m-3?
49  real, dimension(NTR_MAX) :: ic_val = 0.0 !< The (uniform) initial condition value.
50  real, dimension(NTR_MAX) :: young_val = 0.0 !< The value assigned to tr at the surface.
51  real, dimension(NTR_MAX) :: land_val = -1.0 !< The value of tr used where land is masked out.
52  real, dimension(NTR_MAX) :: sfc_growth_rate !< The exponential growth rate for the surface value [year-1].
53  real, dimension(NTR_MAX) :: tracer_start_year !< The year in which tracers start aging, or at which the
54  !! surface value equals young_val, in years.
55  logical :: tracers_may_reinit !< If true, these tracers be set up via the initialization code if
56  !! they are not found in the restart files.
57  logical :: tracer_ages(ntr_max) !< Indicates whether each tracer ages.
58 
59  integer, dimension(NTR_MAX) :: ind_tr !< Indices returned by aof_set_coupler_flux if it is used and the
60  !! surface tracer concentrations are to be provided to the coupler.
61 
62  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
63  !! regulate the timing of diagnostic output.
64  type(mom_restart_cs), pointer :: restart_csp => null() !< A pointer to the restart controls structure
65 
66  type(vardesc) :: tr_desc(ntr_max) !< Descriptions and metadata for the tracers
67 end type ideal_age_tracer_cs
68 
69 contains
70 
71 !> Register the ideal age tracer fields to be used with MOM.
72 function register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
73  type(hor_index_type), intent(in) :: hi !< A horizontal index type structure
74  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
75  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
76  type(ideal_age_tracer_cs), pointer :: cs !< The control structure returned by a previous
77  !! call to register_ideal_age_tracer.
78  type(tracer_registry_type), pointer :: tr_reg !< A pointer that is set to point to the control
79  !! structure for the tracer advection and
80  !! diffusion module
81  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure
82 
83 ! This include declares and sets the variable "version".
84 #include "version_variable.h"
85  character(len=40) :: mdl = "ideal_age_example" ! This module's name.
86  character(len=200) :: inputdir ! The directory where the input files are.
87  character(len=48) :: var_name ! The variable's name.
88  real, pointer :: tr_ptr(:,:,:) => null()
89  logical :: register_ideal_age_tracer
90  logical :: do_ideal_age, do_vintage, do_ideal_age_dated
91  integer :: isd, ied, jsd, jed, nz, m
92  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
93 
94  if (associated(cs)) then
95  call mom_error(warning, "register_ideal_age_tracer called with an "// &
96  "associated control structure.")
97  return
98  endif
99  allocate(cs)
100 
101  ! Read all relevant parameters and write them to the model log.
102  call log_version(param_file, mdl, version, "")
103  call get_param(param_file, mdl, "DO_IDEAL_AGE", do_ideal_age, &
104  "If true, use an ideal age tracer that is set to 0 age "//&
105  "in the mixed layer and ages at unit rate in the interior.", &
106  default=.true.)
107  call get_param(param_file, mdl, "DO_IDEAL_VINTAGE", do_vintage, &
108  "If true, use an ideal vintage tracer that is set to an "//&
109  "exponentially increasing value in the mixed layer and "//&
110  "is conserved thereafter.", default=.false.)
111  call get_param(param_file, mdl, "DO_IDEAL_AGE_DATED", do_ideal_age_dated, &
112  "If true, use an ideal age tracer that is everywhere 0 "//&
113  "before IDEAL_AGE_DATED_START_YEAR, but the behaves like "//&
114  "the standard ideal age tracer - i.e. is set to 0 age in "//&
115  "the mixed layer and ages at unit rate in the interior.", &
116  default=.false.)
117 
118 
119  call get_param(param_file, mdl, "AGE_IC_FILE", cs%IC_file, &
120  "The file in which the age-tracer initial values can be "//&
121  "found, or an empty string for internal initialization.", &
122  default=" ")
123  if ((len_trim(cs%IC_file) > 0) .and. (scan(cs%IC_file,'/') == 0)) then
124  ! Add the directory if CS%IC_file is not already a complete path.
125  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
126  cs%IC_file = trim(slasher(inputdir))//trim(cs%IC_file)
127  call log_param(param_file, mdl, "INPUTDIR/AGE_IC_FILE", cs%IC_file)
128  endif
129  call get_param(param_file, mdl, "AGE_IC_FILE_IS_Z", cs%Z_IC_file, &
130  "If true, AGE_IC_FILE is in depth space, not layer space", &
131  default=.false.)
132  call get_param(param_file, mdl, "TRACERS_MAY_REINIT", cs%tracers_may_reinit, &
133  "If true, tracers may go through the initialization code "//&
134  "if they are not found in the restart files. Otherwise "//&
135  "it is a fatal error if the tracers are not found in the "//&
136  "restart files of a restarted run.", default=.false.)
137 
138  cs%ntr = 0
139  if (do_ideal_age) then
140  cs%ntr = cs%ntr + 1 ; m = cs%ntr
141  cs%tr_desc(m) = var_desc("age", "yr", "Ideal Age Tracer", cmor_field_name="agessc", caller=mdl)
142  cs%tracer_ages(m) = .true. ; cs%sfc_growth_rate(m) = 0.0
143  cs%IC_val(m) = 0.0 ; cs%young_val(m) = 0.0 ; cs%tracer_start_year(m) = 0.0
144  endif
145 
146  if (do_vintage) then
147  cs%ntr = cs%ntr + 1 ; m = cs%ntr
148  cs%tr_desc(m) = var_desc("vintage", "yr", "Exponential Vintage Tracer", &
149  caller=mdl)
150  cs%tracer_ages(m) = .false. ; cs%sfc_growth_rate(m) = 1.0/30.0
151  cs%IC_val(m) = 0.0 ; cs%young_val(m) = 1e-20 ; cs%tracer_start_year(m) = 0.0
152  call get_param(param_file, mdl, "IDEAL_VINTAGE_START_YEAR", cs%tracer_start_year(m), &
153  "The date at which the ideal vintage tracer starts.", &
154  units="years", default=0.0)
155  endif
156 
157  if (do_ideal_age_dated) then
158  cs%ntr = cs%ntr + 1 ; m = cs%ntr
159  cs%tr_desc(m) = var_desc("age_dated","yr","Ideal Age Tracer with a Start Date",&
160  caller=mdl)
161  cs%tracer_ages(m) = .true. ; cs%sfc_growth_rate(m) = 0.0
162  cs%IC_val(m) = 0.0 ; cs%young_val(m) = 0.0 ; cs%tracer_start_year(m) = 0.0
163  call get_param(param_file, mdl, "IDEAL_AGE_DATED_START_YEAR", cs%tracer_start_year(m), &
164  "The date at which the dated ideal age tracer starts.", &
165  units="years", default=0.0)
166  endif
167 
168  allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr)) ; cs%tr(:,:,:,:) = 0.0
169 
170  do m=1,cs%ntr
171  ! This is needed to force the compiler not to do a copy in the registration
172  ! calls. Curses on the designers and implementers of Fortran90.
173  tr_ptr => cs%tr(:,:,:,m)
174  call query_vardesc(cs%tr_desc(m), name=var_name, &
175  caller="register_ideal_age_tracer")
176  ! Register the tracer for horizontal advection, diffusion, and restarts.
177  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, tr_desc=cs%tr_desc(m), &
178  registry_diags=.true., restart_cs=restart_cs, &
179  mandatory=.not.cs%tracers_may_reinit, &
180  flux_scale=gv%H_to_m)
181 
182  ! Set coupled_tracers to be true (hard-coded above) to provide the surface
183  ! values to the coupler (if any). This is meta-code and its arguments will
184  ! currently (deliberately) give fatal errors if it is used.
185  if (cs%coupled_tracers) &
186  cs%ind_tr(m) = aof_set_coupler_flux(trim(var_name)//'_flux', &
187  flux_type=' ', implementation=' ', caller="register_ideal_age_tracer")
188  enddo
189 
190  cs%tr_Reg => tr_reg
191  cs%restart_CSp => restart_cs
193 end function register_ideal_age_tracer
194 
195 !> Sets the ideal age traces to their initial values and sets up the tracer output
196 subroutine initialize_ideal_age_tracer(restart, day, G, GV, US, h, diag, OBC, CS, &
197  sponge_CSp)
198  logical, intent(in) :: restart !< .true. if the fields have already
199  !! been read from a restart file.
200  type(time_type), target, intent(in) :: day !< Time of the start of the run.
201  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
202  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
203  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
204  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
205  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
206  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate
207  !! diagnostic output.
208  type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
209  !! whether, where, and what open boundary
210  !! conditions are used.
211  type(ideal_age_tracer_cs), pointer :: cs !< The control structure returned by a previous
212  !! call to register_ideal_age_tracer.
213  type(sponge_cs), pointer :: sponge_csp !< Pointer to the control structure for the sponges.
214 
215 ! This subroutine initializes the CS%ntr tracer fields in tr(:,:,:,:)
216 ! and it sets up the tracer output.
217 
218  ! Local variables
219  character(len=24) :: name ! A variable's name in a NetCDF file.
220  character(len=72) :: longname ! The long name of that variable.
221  character(len=48) :: units ! The dimensions of the variable.
222  character(len=48) :: flux_units ! The units for age tracer fluxes, either
223  ! years m3 s-1 or years kg s-1.
224  character(len=72) :: cmorname ! The CMOR name of that variable.
225  logical :: ok
226  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
227  integer :: isdb, iedb, jsdb, jedb
228 
229  if (.not.associated(cs)) return
230  if (cs%ntr < 1) return
231  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
232  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
233  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
234 
235  cs%Time => day
236  cs%diag => diag
237  cs%nkml = max(gv%nkml,1)
238 
239  do m=1,cs%ntr
240  call query_vardesc(cs%tr_desc(m), name=name, &
241  caller="initialize_ideal_age_tracer")
242  if ((.not.restart) .or. (cs%tracers_may_reinit .and. .not. &
243  query_initialized(cs%tr(:,:,:,m), name, cs%restart_CSp))) then
244 
245  if (len_trim(cs%IC_file) > 0) then
246  ! Read the tracer concentrations from a netcdf file.
247  if (.not.file_exists(cs%IC_file, g%Domain)) &
248  call mom_error(fatal, "initialize_ideal_age_tracer: "// &
249  "Unable to open "//cs%IC_file)
250 
251  if (cs%Z_IC_file) then
252  ok = tracer_z_init(cs%tr(:,:,:,m), h, cs%IC_file, name,&
253  g, us, -1e34, 0.0) ! CS%land_val(m))
254  if (.not.ok) then
255  ok = tracer_z_init(cs%tr(:,:,:,m), h, cs%IC_file, &
256  trim(name), g, us, -1e34, 0.0) ! CS%land_val(m))
257  if (.not.ok) call mom_error(fatal,"initialize_ideal_age_tracer: "//&
258  "Unable to read "//trim(name)//" from "//&
259  trim(cs%IC_file)//".")
260  endif
261  else
262  call mom_read_data(cs%IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
263  endif
264  else
265  do k=1,nz ; do j=js,je ; do i=is,ie
266  if (g%mask2dT(i,j) < 0.5) then
267  cs%tr(i,j,k,m) = cs%land_val(m)
268  else
269  cs%tr(i,j,k,m) = cs%IC_val(m)
270  endif
271  enddo ; enddo ; enddo
272  endif
273 
274  endif ! restart
275  enddo ! Tracer loop
276 
277  if (associated(obc)) then
278  ! Steal from updated DOME in the fullness of time.
279  endif
280 
281 end subroutine initialize_ideal_age_tracer
282 
283 !> Applies diapycnal diffusion, aging and regeneration at the surface to the ideal age tracers
284 subroutine ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, &
285  evap_CFL_limit, minimum_forcing_depth)
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  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
289  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
290  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
291  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
292  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
293  intent(in) :: ea !< an array to which the amount of fluid entrained
294  !! from the layer above during this call will be
295  !! added [H ~> m or kg m-2].
296  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
297  intent(in) :: eb !< an array to which the amount of fluid entrained
298  !! from the layer below during this call will be
299  !! added [H ~> m or kg m-2].
300  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
301  !! and tracer forcing fields. Unused fields have NULL ptrs.
302  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
303  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
304  type(ideal_age_tracer_cs), pointer :: cs !< The control structure returned by a previous
305  !! call to register_ideal_age_tracer.
306  real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can
307  !! be fluxed out of the top layer in a timestep [nondim]
308  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
309  !! fluxes can be applied [H ~> m or kg m-2]
310 ! This subroutine applies diapycnal diffusion and any other column
311 ! tracer physics or chemistry to the tracers from this file.
312 ! This is a simple example of a set of advected passive tracers.
313 
314 ! The arguments to this subroutine are redundant in that
315 ! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
316  ! Local variables
317  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
318  real :: sfc_val ! The surface value for the tracers.
319  real :: isecs_per_year ! The inverse of the amount of time in a year [T-1 ~> s-1]
320  real :: year ! The time in years.
321  integer :: i, j, k, is, ie, js, je, nz, m
322  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
323 
324  if (.not.associated(cs)) return
325  if (cs%ntr < 1) return
326 
327  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
328  do m=1,cs%ntr
329  do k=1,nz ;do j=js,je ; do i=is,ie
330  h_work(i,j,k) = h_old(i,j,k)
331  enddo ; enddo ; enddo
332  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
333  evap_cfl_limit, minimum_forcing_depth)
334  call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
335  enddo
336  else
337  do m=1,cs%ntr
338  call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
339  enddo
340  endif
341 
342  isecs_per_year = 1.0 / (365.0*86400.0*us%s_to_T)
343  ! Set the surface value of tracer 1 to increase exponentially
344  ! with a 30 year time scale.
345  year = us%s_to_T*time_type_to_real(cs%Time) * isecs_per_year
346 
347  do m=1,cs%ntr
348  if (cs%sfc_growth_rate(m) == 0.0) then
349  sfc_val = cs%young_val(m)
350  else
351  sfc_val = cs%young_val(m) * &
352  exp((year-cs%tracer_start_year(m)) * cs%sfc_growth_rate(m))
353  endif
354  do k=1,cs%nkml ; do j=js,je ; do i=is,ie
355  if (g%mask2dT(i,j) > 0.5) then
356  cs%tr(i,j,k,m) = sfc_val
357  else
358  cs%tr(i,j,k,m) = cs%land_val(m)
359  endif
360  enddo ; enddo ; enddo
361  enddo
362  do m=1,cs%ntr ; if (cs%tracer_ages(m) .and. &
363  (year>=cs%tracer_start_year(m))) then
364 !$OMP parallel do default(none) shared(is,ie,js,je,CS,nz,G,dt,Isecs_per_year,m)
365  do k=cs%nkml+1,nz ; do j=js,je ; do i=is,ie
366  cs%tr(i,j,k,m) = cs%tr(i,j,k,m) + g%mask2dT(i,j)*dt*isecs_per_year
367  enddo ; enddo ; enddo
368  endif ; enddo
369 
370 end subroutine ideal_age_tracer_column_physics
371 
372 !> Calculates the mass-weighted integral of all tracer stocks, returning the number of stocks it
373 !! has calculated. If stock_index is present, only the stock corresponding to that coded index is found.
374 function ideal_age_stock(h, stocks, G, GV, CS, names, units, stock_index)
375  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
376  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
377  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
378  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each
379  !! tracer, in kg times concentration units [kg conc].
380  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
381  type(ideal_age_tracer_cs), pointer :: cs !< The control structure returned by a previous
382  !! call to register_ideal_age_tracer.
383  character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated.
384  character(len=*), dimension(:), intent(out) :: units !< the units of the stocks calculated.
385  integer, optional, intent(in) :: stock_index !< the coded index of a specific stock
386  !! being sought.
387  integer :: ideal_age_stock !< The number of stocks calculated here.
388 ! This function calculates the mass-weighted integral of all tracer stocks,
389 ! returning the number of stocks it has calculated. If the stock_index
390 ! is present, only the stock corresponding to that coded index is returned.
391 
392  integer :: i, j, k, is, ie, js, je, nz, m
393  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
394 
395  ideal_age_stock = 0
396  if (.not.associated(cs)) return
397  if (cs%ntr < 1) return
398 
399  if (present(stock_index)) then ; if (stock_index > 0) then
400  ! Check whether this stock is available from this routine.
401 
402  ! No stocks from this routine are being checked yet. Return 0.
403  return
404  endif ; endif
405 
406  do m=1,cs%ntr
407  call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="ideal_age_stock")
408  units(m) = trim(units(m))//" kg"
409  stocks(m) = 0.0
410  do k=1,nz ; do j=js,je ; do i=is,ie
411  stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
412  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
413  enddo ; enddo ; enddo
414  stocks(m) = gv%H_to_kg_m2 * stocks(m)
415  enddo
416  ideal_age_stock = cs%ntr
417 
418 end function ideal_age_stock
419 
420 !> This subroutine extracts the surface fields from this tracer package that
421 !! are to be shared with the atmosphere in coupled configurations.
422 !! This particular tracer package does not report anything back to the coupler.
423 subroutine ideal_age_tracer_surface_state(state, h, G, CS)
424  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
425  type(surface), intent(inout) :: state !< A structure containing fields that
426  !! describe the surface state of the ocean.
427  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
428  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
429  type(ideal_age_tracer_cs), pointer :: cs !< The control structure returned by a previous
430  !! call to register_ideal_age_tracer.
431 
432  ! This particular tracer package does not report anything back to the coupler.
433  ! The code that is here is just a rough guide for packages that would.
434 
435  integer :: m, is, ie, js, je, isd, ied, jsd, jed
436  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
437  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
438 
439  if (.not.associated(cs)) return
440 
441  if (cs%coupled_tracers) then
442  do m=1,cs%ntr
443  ! This call loads the surface values into the appropriate array in the
444  ! coupler-type structure.
445  call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
446  state%tr_fields, idim=(/isd, is, ie, ied/), &
447  jdim=(/jsd, js, je, jed/) )
448  enddo
449  endif
450 
451 end subroutine ideal_age_tracer_surface_state
452 
453 !> Deallocate any memory associated with this tracer package
454 subroutine ideal_age_example_end(CS)
455  type(ideal_age_tracer_cs), pointer :: cs !< The control structure returned by a previous
456  !! call to register_ideal_age_tracer.
457 
458  integer :: m
459 
460  if (associated(cs)) then
461  if (associated(cs%tr)) deallocate(cs%tr)
462  deallocate(cs)
463  endif
464 end subroutine ideal_age_example_end
465 
466 !> \namespace ideal_age_example
467 !!
468 !! Originally by Robert Hallberg, 2002
469 !!
470 !! This file contains an example of the code that is needed to set
471 !! up and use a set (in this case two) of dynamically passive tracers
472 !! for diagnostic purposes. The tracers here are an ideal age tracer
473 !! that ages at a rate of 1/year once it is isolated from the surface,
474 !! and a vintage tracer, whose surface concentration grows exponen-
475 !! with time with a 30-year timescale (similar to CFCs).
476 
477 end module ideal_age_example
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_tracer_registry::register_tracer
subroutine, public register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, cmor_name, cmor_units, cmor_longname, tr_desc, OBC_inflow, OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, ad_2d_x, ad_2d_y, df_2d_x, df_2d_y, advection_xy, registry_diags, flux_nameroot, flux_longname, flux_units, flux_scale, convergence_units, convergence_scale, cmor_tendprefix, diag_form, restart_CS, mandatory)
This subroutine registers a tracer to be advected and laterally diffused.
Definition: MOM_tracer_registry.F90:158
mom_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
ideal_age_example::ideal_age_stock
integer function, public ideal_age_stock(h, stocks, G, GV, CS, names, units, stock_index)
Calculates the mass-weighted integral of all tracer stocks, returning the number of stocks it has cal...
Definition: ideal_age_example.F90:375
mom_io::var_desc
type(vardesc) function, public var_desc(name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
Returns a vardesc type whose elements have been filled with the provided fields. The argument name is...
Definition: MOM_io.F90:600
mom_io::query_vardesc
subroutine, public query_vardesc(vd, name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
This routine queries vardesc.
Definition: MOM_io.F90:699
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
ideal_age_example
A tracer package of ideal age tracers.
Definition: ideal_age_example.F90:2
mom_tracer_diabatic::applytracerboundaryfluxesinout
subroutine, public applytracerboundaryfluxesinout(G, GV, Tr, dt, fluxes, h, evap_CFL_limit, minimum_forcing_depth, in_flux_optional, out_flux_optional, update_h_opt)
This routine is modeled after applyBoundaryFluxesInOut in MOM_diabatic_aux.F90 NOTE: Please note that...
Definition: MOM_tracer_diabatic.F90:230
mom_sponge::set_up_sponge_field
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
Definition: MOM_sponge.F90:214
ideal_age_example::ideal_age_tracer_surface_state
subroutine, public ideal_age_tracer_surface_state(state, h, G, CS)
This subroutine extracts the surface fields from this tracer package that are to be shared with the a...
Definition: ideal_age_example.F90:424
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
atmos_ocean_fluxes_mod
A dummy version of atmos_ocean_fluxes_mod module for use when the vastly larger FMS package is not ne...
Definition: atmos_ocean_fluxes.F90:3
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_tracer_registry
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
Definition: MOM_tracer_registry.F90:5
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
ideal_age_example::ideal_age_example_end
subroutine, public ideal_age_example_end(CS)
Deallocate any memory associated with this tracer package.
Definition: ideal_age_example.F90:455
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_restart::mom_restart_cs
A restart registry and the control structure for restarts.
Definition: MOM_restart.F90:72
mom_tracer_z_init
Used to initialize tracers from a depth- (or z*-) space file.
Definition: MOM_tracer_Z_init.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_tracer_diabatic
This module contains routines that implement physical fluxes of tracers (e.g. due to surface fluxes o...
Definition: MOM_tracer_diabatic.F90:4
ideal_age_example::ideal_age_tracer_cs
The control structure for the ideal_age_tracer package.
Definition: ideal_age_example.F90:38
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_tracer_diabatic::tracer_vertdiff
subroutine, public tracer_vertdiff(h_old, ea, eb, dt, tr, G, GV, sfc_flux, btm_flux, btm_reservoir, sink_rate, convert_flux_in)
This subroutine solves a tridiagonal equation for the final tracer concentrations after the dual-entr...
Definition: MOM_tracer_diabatic.F90:27
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_restart
The MOM6 facility for reading and writing restart files, and querying what has been read.
Definition: MOM_restart.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
ideal_age_example::register_ideal_age_tracer
logical function, public register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
Register the ideal age tracer fields to be used with MOM.
Definition: ideal_age_example.F90:73
mom_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
mom_open_boundary
Controls where open boundary conditions are applied.
Definition: MOM_open_boundary.F90:2
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_sponge
Implements sponge regions in isopycnal mode.
Definition: MOM_sponge.F90:2
mom_tracer_z_init::tracer_z_init
logical function, public tracer_z_init(tr, h, filename, tr_name, G, US, missing_val, land_val)
This function initializes a tracer by reading a Z-space file, returning .true. if this appears to hav...
Definition: MOM_tracer_Z_init.F90:31
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains.
Definition: MOM_hor_index.F90:15
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_open_boundary::ocean_obc_type
Open-boundary data.
Definition: MOM_open_boundary.F90:195
mom_sponge::sponge_cs
This control structure holds memory and parameters for the MOM_sponge module.
Definition: MOM_sponge.F90:41
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_forcing_type::forcing
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Definition: MOM_forcing_type.F90:50
mom_io::vardesc
Type for describing a variable, typically a tracer.
Definition: MOM_io.F90:53
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
ideal_age_example::initialize_ideal_age_tracer
subroutine, public initialize_ideal_age_tracer(restart, day, G, GV, US, h, diag, OBC, CS, sponge_CSp)
Sets the ideal age traces to their initial values and sets up the tracer output.
Definition: ideal_age_example.F90:198
ideal_age_example::ntr_max
integer, parameter ntr_max
the maximum number of tracers in this module.
Definition: ideal_age_example.F90:35
atmos_ocean_fluxes_mod::aof_set_coupler_flux
integer function, public aof_set_coupler_flux(name, flux_type, implementation, atm_tr_index, param, flag, mol_wt, ice_restart_file, ocean_restart_file, units, caller, verbosity)
This subroutine duplicates an interface used by the FMS coupler, but only returns a value of -1....
Definition: atmos_ocean_fluxes.F90:18
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
mom_restart::query_initialized
Indicate whether a field has been read from a restart file.
Definition: MOM_restart.F90:116
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
ideal_age_example::ideal_age_tracer_column_physics
subroutine, public ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, evap_CFL_limit, minimum_forcing_depth)
Applies diapycnal diffusion, aging and regeneration at the surface to the ideal age tracers.
Definition: ideal_age_example.F90:286
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239