MOM6
dye_example.F90
Go to the documentation of this file.
1 !> A tracer package for using dyes to diagnose regional flows.
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
12 use mom_io, only : file_exists, read_data, slasher, vardesc, var_desc, query_vardesc
16 use mom_time_manager, only : time_type
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 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
36 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
37 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
38 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
39 
40 !> The control structure for the regional dyes tracer package
41 type, public :: dye_tracer_cs ; private
42  integer :: ntr !< The number of tracers that are actually used.
43  logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
44  real, allocatable, dimension(:) :: dye_source_minlon !< Minimum longitude of region dye will be injected.
45  real, allocatable, dimension(:) :: dye_source_maxlon !< Maximum longitude of region dye will be injected.
46  real, allocatable, dimension(:) :: dye_source_minlat !< Minimum latitude of region dye will be injected.
47  real, allocatable, dimension(:) :: dye_source_maxlat !< Maximum latitude of region dye will be injected.
48  real, allocatable, dimension(:) :: dye_source_mindepth !< Minimum depth of region dye will be injected [Z ~> m].
49  real, allocatable, dimension(:) :: dye_source_maxdepth !< Maximum depth of region dye will be injected [Z ~> m].
50  type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the tracer registry
51  real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this subroutine, in g m-3?
52 
53  integer, allocatable, dimension(:) :: ind_tr !< Indices returned by aof_set_coupler_flux if it is used and the
54  !! surface tracer concentrations are to be provided to the coupler.
55 
56  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
57  !! regulate the timing of diagnostic output.
58  type(mom_restart_cs), pointer :: restart_csp => null() !< A pointer to the restart control structure
59 
60  type(vardesc), allocatable :: tr_desc(:) !< Descriptions and metadata for the tracers
61  logical :: tracers_may_reinit = .false. !< If true the tracers may be initialized if not found in a restart file
62 end type dye_tracer_cs
63 
64 contains
65 
66 !> This subroutine is used to register tracer fields and subroutines
67 !! to be used with MOM.
68 function register_dye_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS)
69  type(hor_index_type), intent(in) :: hi !< A horizontal index type structure.
70  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
71  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
72  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
73  type(dye_tracer_cs), pointer :: cs !< A pointer that is set to point to the control
74  !! structure for this module
75  type(tracer_registry_type), pointer :: tr_reg !< A pointer that is set to point to the control
76  !! structure for the tracer advection and diffusion module.
77  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure.
78 
79 ! Local variables
80 ! This include declares and sets the variable "version".
81 #include "version_variable.h"
82  character(len=40) :: mdl = "regional_dyes" ! This module's name.
83  character(len=200) :: inputdir ! The directory where the input files are.
84  character(len=48) :: var_name ! The variable's name.
85  character(len=48) :: desc_name ! The variable's descriptor.
86  real, pointer :: tr_ptr(:,:,:) => null()
87  logical :: register_dye_tracer
88  integer :: isd, ied, jsd, jed, nz, m
89  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
90 
91  if (associated(cs)) then
92  call mom_error(warning, "register_dye_tracer called with an "// &
93  "associated control structure.")
94  return
95  endif
96  allocate(cs)
97 
98  ! Read all relevant parameters and write them to the model log.
99  call log_version(param_file, mdl, version, "")
100  call get_param(param_file, mdl, "NUM_DYE_TRACERS", cs%ntr, &
101  "The number of dye tracers in this run. Each tracer "//&
102  "should have a separate region.", default=0)
103  allocate(cs%dye_source_minlon(cs%ntr), &
104  cs%dye_source_maxlon(cs%ntr), &
105  cs%dye_source_minlat(cs%ntr), &
106  cs%dye_source_maxlat(cs%ntr), &
107  cs%dye_source_mindepth(cs%ntr), &
108  cs%dye_source_maxdepth(cs%ntr))
109  allocate(cs%ind_tr(cs%ntr))
110  allocate(cs%tr_desc(cs%ntr))
111 
112  cs%dye_source_minlon(:) = -1.e30
113  call get_param(param_file, mdl, "DYE_SOURCE_MINLON", cs%dye_source_minlon, &
114  "This is the starting longitude at which we start injecting dyes.", &
115  fail_if_missing=.true.)
116  if (minval(cs%dye_source_minlon(:)) < -1.e29) &
117  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINLON ")
118 
119  cs%dye_source_maxlon(:) = -1.e30
120  call get_param(param_file, mdl, "DYE_SOURCE_MAXLON", cs%dye_source_maxlon, &
121  "This is the ending longitude at which we finish injecting dyes.", &
122  fail_if_missing=.true.)
123  if (minval(cs%dye_source_maxlon(:)) < -1.e29) &
124  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXLON ")
125 
126  cs%dye_source_minlat(:) = -1.e30
127  call get_param(param_file, mdl, "DYE_SOURCE_MINLAT", cs%dye_source_minlat, &
128  "This is the starting latitude at which we start injecting dyes.", &
129  fail_if_missing=.true.)
130  if (minval(cs%dye_source_minlat(:)) < -1.e29) &
131  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINLAT ")
132 
133  cs%dye_source_maxlat(:) = -1.e30
134  call get_param(param_file, mdl, "DYE_SOURCE_MAXLAT", cs%dye_source_maxlat, &
135  "This is the ending latitude at which we finish injecting dyes.", &
136  fail_if_missing=.true.)
137  if (minval(cs%dye_source_maxlat(:)) < -1.e29) &
138  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXLAT ")
139 
140  cs%dye_source_mindepth(:) = -1.e30
141  call get_param(param_file, mdl, "DYE_SOURCE_MINDEPTH", cs%dye_source_mindepth, &
142  "This is the minimum depth at which we inject dyes.", &
143  units="m", scale=us%m_to_Z, fail_if_missing=.true.)
144  if (minval(cs%dye_source_mindepth(:)) < -1.e29*us%m_to_Z) &
145  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINDEPTH")
146 
147  cs%dye_source_maxdepth(:) = -1.e30
148  call get_param(param_file, mdl, "DYE_SOURCE_MAXDEPTH", cs%dye_source_maxdepth, &
149  "This is the maximum depth at which we inject dyes.", &
150  units="m", scale=us%m_to_Z, fail_if_missing=.true.)
151  if (minval(cs%dye_source_maxdepth(:)) < -1.e29*us%m_to_Z) &
152  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXDEPTH ")
153 
154  allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr)) ; cs%tr(:,:,:,:) = 0.0
155 
156  do m = 1, cs%ntr
157  write(var_name(:),'(A,I3.3)') "dye",m
158  write(desc_name(:),'(A,I3.3)') "Dye Tracer ",m
159  cs%tr_desc(m) = var_desc(trim(var_name), "conc", trim(desc_name), caller=mdl)
160 
161  ! This is needed to force the compiler not to do a copy in the registration
162  ! calls. Curses on the designers and implementers of Fortran90.
163  tr_ptr => cs%tr(:,:,:,m)
164  call query_vardesc(cs%tr_desc(m), name=var_name, &
165  caller="register_dye_tracer")
166  ! Register the tracer for horizontal advection, diffusion, and restarts.
167  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
168  tr_desc=cs%tr_desc(m), registry_diags=.true., &
169  restart_cs=restart_cs, mandatory=.not.cs%tracers_may_reinit)
170 
171  ! Set coupled_tracers to be true (hard-coded above) to provide the surface
172  ! values to the coupler (if any). This is meta-code and its arguments will
173  ! currently (deliberately) give fatal errors if it is used.
174  if (cs%coupled_tracers) &
175  cs%ind_tr(m) = aof_set_coupler_flux(trim(var_name)//'_flux', &
176  flux_type=' ', implementation=' ', caller="register_dye_tracer")
177  enddo
178 
179  cs%tr_Reg => tr_reg
180  cs%restart_CSp => restart_cs
181  register_dye_tracer = .true.
182 end function register_dye_tracer
183 
184 !> This subroutine initializes the CS%ntr tracer fields in tr(:,:,:,:)
185 !! and it sets up the tracer output.
186 subroutine initialize_dye_tracer(restart, day, G, GV, h, diag, OBC, CS, sponge_CSp)
187  logical, intent(in) :: restart !< .true. if the fields have already been
188  !! read from a restart file.
189  type(time_type), target, intent(in) :: day !< Time of the start of the run.
190  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
191  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
192  real, dimension(NIMEM_,NJMEM_,NKMEM_), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
193  type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
194  type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
195  !! whether, where, and what open boundary
196  !! conditions are used.
197  type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
198  !! call to register_dye_tracer.
199  type(sponge_cs), pointer :: sponge_csp !< A pointer to the control structure
200  !! for the sponges, if they are in use.
201 
202 ! Local variables
203  character(len=24) :: name ! A variable's name in a NetCDF file.
204  character(len=72) :: longname ! The long name of that variable.
205  character(len=48) :: units ! The dimensions of the variable.
206  character(len=48) :: flux_units ! The units for age tracer fluxes, either
207  ! years m3 s-1 or years kg s-1.
208  logical :: ok
209  integer :: i, j, k, m
210  real :: z_bot, z_center
211 
212  if (.not.associated(cs)) return
213  if (cs%ntr < 1) return
214 
215  cs%diag => diag
216 
217  ! Establish location of source
218  do m= 1, cs%ntr
219  do j=g%jsd,g%jed ; do i=g%isd,g%ied
220  ! A dye is set dependent on the center of the cell being inside the rectangular box.
221  if (cs%dye_source_minlon(m)<g%geoLonT(i,j) .and. &
222  cs%dye_source_maxlon(m)>=g%geoLonT(i,j) .and. &
223  cs%dye_source_minlat(m)<g%geoLatT(i,j) .and. &
224  cs%dye_source_maxlat(m)>=g%geoLatT(i,j) .and. &
225  g%mask2dT(i,j) > 0.0 ) then
226  z_bot = -g%bathyT(i,j)
227  do k = gv%ke, 1, -1
228  z_center = z_bot + 0.5*h(i,j,k)*gv%H_to_Z
229  if ( z_center > -cs%dye_source_maxdepth(m) .and. &
230  z_center < -cs%dye_source_mindepth(m) ) then
231  cs%tr(i,j,k,m) = 1.0
232  endif
233  z_bot = z_bot + h(i,j,k)*gv%H_to_Z
234  enddo
235  endif
236  enddo ; enddo
237  enddo
238 
239 end subroutine initialize_dye_tracer
240 
241 !> This subroutine applies diapycnal diffusion and any other column
242 !! tracer physics or chemistry to the tracers from this file.
243 !! This is a simple example of a set of advected passive tracers.
244 !! The arguments to this subroutine are redundant in that
245 !! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
246 subroutine dye_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, &
247  evap_CFL_limit, minimum_forcing_depth)
248  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
249  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
250  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
251  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
252  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
253  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
254  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
255  intent(in) :: ea !< an array to which the amount of fluid entrained
256  !! from the layer above during this call will be
257  !! added [H ~> m or kg m-2].
258  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
259  intent(in) :: eb !< an array to which the amount of fluid entrained
260  !! from the layer below during this call will be
261  !! added [H ~> m or kg m-2].
262  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
263  !! and tracer forcing fields. Unused fields have NULL ptrs.
264  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
265  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
266  type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
267  !! call to register_dye_tracer.
268  real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can
269  !! be fluxed out of the top layer in a timestep [nondim]
270  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
271  !! fluxes can be applied [H ~> m or kg m-2]
272 
273 ! Local variables
274  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
275  real :: sfc_val ! The surface value for the tracers.
276  real :: isecs_per_year ! The number of seconds in a year.
277  real :: year ! The time in years.
278  integer :: secs, days ! Integer components of the time type.
279  integer :: i, j, k, is, ie, js, je, nz, m
280  real :: z_bot, z_center
281 
282  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
283 
284  if (.not.associated(cs)) return
285  if (cs%ntr < 1) return
286 
287  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
288  do m=1,cs%ntr
289  do k=1,nz ;do j=js,je ; do i=is,ie
290  h_work(i,j,k) = h_old(i,j,k)
291  enddo ; enddo ; enddo
292  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
293  evap_cfl_limit, minimum_forcing_depth)
294  call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
295  enddo
296  else
297  do m=1,cs%ntr
298  call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
299  enddo
300  endif
301 
302  do m=1,cs%ntr
303  do j=g%jsd,g%jed ; do i=g%isd,g%ied
304  ! A dye is set dependent on the center of the cell being inside the rectangular box.
305  if (cs%dye_source_minlon(m)<g%geoLonT(i,j) .and. &
306  cs%dye_source_maxlon(m)>=g%geoLonT(i,j) .and. &
307  cs%dye_source_minlat(m)<g%geoLatT(i,j) .and. &
308  cs%dye_source_maxlat(m)>=g%geoLatT(i,j) .and. &
309  g%mask2dT(i,j) > 0.0 ) then
310  z_bot = -g%bathyT(i,j)
311  do k=nz,1,-1
312  z_center = z_bot + 0.5*h_new(i,j,k)*gv%H_to_Z
313  if ( z_center > -cs%dye_source_maxdepth(m) .and. &
314  z_center < -cs%dye_source_mindepth(m) ) then
315  cs%tr(i,j,k,m) = 1.0
316  endif
317  z_bot = z_bot + h_new(i,j,k)*gv%H_to_Z
318  enddo
319  endif
320  enddo ; enddo
321  enddo
322 
323 end subroutine dye_tracer_column_physics
324 
325 !> This function calculates the mass-weighted integral of all tracer stocks,
326 !! returning the number of stocks it has calculated. If the stock_index
327 !! is present, only the stock corresponding to that coded index is returned.
328 function dye_stock(h, stocks, G, GV, CS, names, units, stock_index)
329  real, dimension(NIMEM_,NJMEM_,NKMEM_), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
330  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of
331  !! each tracer, in kg times concentration units [kg conc].
332  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
333  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
334  type(dye_tracer_cs), pointer :: cs !< The control structure returned by a
335  !! previous call to register_dye_tracer.
336  character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated.
337  character(len=*), dimension(:), intent(out) :: units !< the units of the stocks calculated.
338  integer, optional, intent(in) :: stock_index !< the coded index of a specific stock
339  !! being sought.
340  integer :: dye_stock !< Return value: the number of stocks
341  !! calculated here.
342 
343 ! Local variables
344  integer :: i, j, k, is, ie, js, je, nz, m
345  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
346 
347  dye_stock = 0
348  if (.not.associated(cs)) return
349  if (cs%ntr < 1) return
350 
351  if (present(stock_index)) then ; if (stock_index > 0) then
352  ! Check whether this stock is available from this routine.
353 
354  ! No stocks from this routine are being checked yet. Return 0.
355  return
356  endif ; endif
357 
358  do m=1,cs%ntr
359  call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="dye_stock")
360  units(m) = trim(units(m))//" kg"
361  stocks(m) = 0.0
362  do k=1,nz ; do j=js,je ; do i=is,ie
363  stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
364  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
365  enddo ; enddo ; enddo
366  stocks(m) = gv%H_to_kg_m2 * stocks(m)
367  enddo
368  dye_stock = cs%ntr
369 
370 end function dye_stock
371 
372 !> This subroutine extracts the surface fields from this tracer package that
373 !! are to be shared with the atmosphere in coupled configurations.
374 !! This particular tracer package does not report anything back to the coupler.
375 subroutine dye_tracer_surface_state(state, h, G, CS)
376  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
377  type(surface), intent(inout) :: state !< A structure containing fields that
378  !! describe the surface state of the ocean.
379  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
380  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
381  type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
382  !! call to register_dye_tracer.
383 
384  ! This particular tracer package does not report anything back to the coupler.
385  ! The code that is here is just a rough guide for packages that would.
386 
387  integer :: m, is, ie, js, je, isd, ied, jsd, jed
388  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
389  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
390 
391  if (.not.associated(cs)) return
392 
393  if (cs%coupled_tracers) then
394  do m=1,cs%ntr
395  ! This call loads the surface values into the appropriate array in the
396  ! coupler-type structure.
397  call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
398  state%tr_fields, idim=(/isd, is, ie, ied/), &
399  jdim=(/jsd, js, je, jed/) )
400  enddo
401  endif
402 
403 end subroutine dye_tracer_surface_state
404 
405 !> Clean up any allocated memory after the run.
406 subroutine regional_dyes_end(CS)
407  type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
408  !! call to register_dye_tracer.
409  integer :: m
410 
411  if (associated(cs)) then
412  if (associated(cs%tr)) deallocate(cs%tr)
413  deallocate(cs)
414  endif
415 end subroutine regional_dyes_end
416 
417 !> \namespace regional_dyes
418 !!
419 !! This file contains an example of the code that is needed to set
420 !! up and use a set (in this case two) of dynamically passive tracers
421 !! for diagnostic purposes. The tracers here are dye tracers which
422 !! are set to 1 within the geographical region specified. The depth
423 !! which a tracer is set is determined by calculating the depth from
424 !! the seafloor upwards through the column.
425 
426 end module regional_dyes
regional_dyes::dye_tracer_column_physics
subroutine, public dye_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, evap_CFL_limit, minimum_forcing_depth)
This subroutine applies diapycnal diffusion and any other column tracer physics or chemistry to the t...
Definition: dye_example.F90:248
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
regional_dyes
A tracer package for using dyes to diagnose regional flows.
Definition: dye_example.F90:2
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
regional_dyes::initialize_dye_tracer
subroutine, public initialize_dye_tracer(restart, day, G, GV, h, diag, OBC, CS, sponge_CSp)
This subroutine initializes the CSntr tracer fields in tr(:,:,:,:) and it sets up the tracer output.
Definition: dye_example.F90:187
mom_io::query_vardesc
subroutine, public query_vardesc(vd, name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
This routine queries vardesc.
Definition: MOM_io.F90:699
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_tracer_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
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
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
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
regional_dyes::regional_dyes_end
subroutine, public regional_dyes_end(CS)
Clean up any allocated memory after the run.
Definition: dye_example.F90:407
regional_dyes::register_dye_tracer
logical function, public register_dye_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS)
This subroutine is used to register tracer fields and subroutines to be used with MOM.
Definition: dye_example.F90:69
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
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
regional_dyes::dye_stock
integer function, public dye_stock(h, stocks, G, GV, CS, names, units, stock_index)
This function calculates the mass-weighted integral of all tracer stocks, returning the number of sto...
Definition: dye_example.F90:329
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
regional_dyes::dye_tracer_cs
The control structure for the regional dyes tracer package.
Definition: dye_example.F90:41
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
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
regional_dyes::dye_tracer_surface_state
subroutine, public dye_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: dye_example.F90:376
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
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