MOM6
RGC_tracer.F90
Go to the documentation of this file.
1 !> This module contains the routines used to set up a
2 !! dynamically passive tracer.
3 !! Set up and use passive tracers requires the following:
4 !! (1) register_RGC_tracer
5 !! (2) apply diffusion, physics/chemistry and advect the tracer
6 
7 !********+*********+*********+*********+*********+*********+*********+**
8 !* *
9 !* By Elizabeth Yankovsky, June 2019 *
10 !*********+*********+*********+*********+*********+*********+***********
11 
12 module rgc_tracer
13 
14 ! This file is part of MOM6. See LICENSE.md for the license.
15 
16 use mom_diag_mediator, only : diag_ctrl
17 use mom_error_handler, only : mom_error, fatal, warning
19 use mom_forcing_type, only : forcing
20 use mom_hor_index, only : hor_index_type
21 use mom_grid, only : ocean_grid_type
22 use mom_io, only : file_exists, read_data, slasher, vardesc, var_desc, query_vardesc
23 use mom_restart, only : mom_restart_cs
26 use mom_time_manager, only : time_type
30 use mom_variables, only : surface
33 
34 implicit none ; private
35 
36 #include <MOM_memory.h>
37 
38 !< Publicly available functions
41 
42 integer, parameter :: ntr = 1 !< The number of tracers in this module.
43 
44 !> tracer control structure
45 type, public :: rgc_tracer_cs ; private
46  logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
47  character(len = 200) :: tracer_ic_file !< The full path to the IC file, or " " to initialize internally.
48  type(time_type), pointer :: time !< A pointer to the ocean model's clock.
49  type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the tracer registry.
50  real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this package.
51  real, pointer :: tr_aux(:,:,:,:) => null() !< The masked tracer concentration.
52  real :: land_val(ntr) = -1.0 !< The value of tr used where land is masked out.
53  real :: lenlat !< the latitudinal or y-direction length of the domain.
54  real :: lenlon !< the longitudinal or x-direction length of the domain.
55  real :: csl !< The length of the continental shelf (x dir, km)
56  real :: lensponge !< the length of the sponge layer.
57  logical :: mask_tracers !< If true, tracers are masked out in massless layers.
58  logical :: use_sponge !< If true, sponges may be applied somewhere in the domain.
59  type(diag_ctrl), pointer :: diag !< A structure that is used to regulate the timing of diagnostic output.
60  type(vardesc) :: tr_desc(ntr) !< Descriptions and metadata for the tracers.
61 end type rgc_tracer_cs
62 
63 contains
64 
65 
66 !> This subroutine is used to register tracer fields
67 function register_rgc_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
68  type(hor_index_type), intent(in) :: hi !< A horizontal index type structure.
69  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
70  type(param_file_type), intent(in) :: param_file !<A structure indicating the open file to parse
71  !! for model parameter values.
72  type(rgc_tracer_cs), pointer :: cs !< A pointer that is set to point to the control
73  !! structure for this module (in/out).
74  type(tracer_registry_type), pointer :: tr_reg !< A pointer to the tracer registry.
75  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure.
76 
77  character(len=80) :: name, longname
78 ! This include declares and sets the variable "version".
79 #include "version_variable.h"
80  character(len=40) :: mdl = "RGC_tracer" ! This module's name.
81  character(len=200) :: inputdir
82  real, pointer :: tr_ptr(:,:,:) => null()
83  logical :: register_rgc_tracer
84  integer :: isd, ied, jsd, jed, nz, m
85  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
86 
87  if (associated(cs)) then
88  call mom_error(warning, "RGC_register_tracer called with an "// &
89  "associated control structure.")
90  return
91  endif
92  allocate(cs)
93 
94  ! Read all relevant parameters and write them to the model log.
95  call log_version(param_file, mdl, version, "")
96  call get_param(param_file, mdl, "RGC_TRACER_IC_FILE", cs%tracer_IC_file, &
97  "The name of a file from which to read the initial \n"//&
98  "conditions for the RGC tracers, or blank to initialize \n"//&
99  "them internally.", default=" ")
100  if (len_trim(cs%tracer_IC_file) >= 1) then
101  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
102  inputdir = slasher(inputdir)
103  cs%tracer_IC_file = trim(inputdir)//trim(cs%tracer_IC_file)
104  call log_param(param_file, mdl, "INPUTDIR/RGC_TRACER_IC_FILE", &
105  cs%tracer_IC_file)
106  endif
107  call get_param(param_file, mdl, "SPONGE", cs%use_sponge, &
108  "If true, sponges may be applied anywhere in the domain. \n"//&
109  "The exact location and properties of those sponges are \n"//&
110  "specified from MOM_initialization.F90.", default=.false.)
111 
112  call get_param(param_file, mdl, "LENLAT", cs%lenlat, &
113  "The latitudinal or y-direction length of the domain", &
114  fail_if_missing=.true., do_not_log=.true.)
115 
116  call get_param(param_file, mdl, "LENLON", cs%lenlon, &
117  "The longitudinal or x-direction length of the domain", &
118  fail_if_missing=.true., do_not_log=.true.)
119 
120  call get_param(param_file, mdl, "CONT_SHELF_LENGTH", cs%CSL, &
121  "The length of the continental shelf (x dir, km).", &
122  default=15.0)
123 
124  call get_param(param_file, mdl, "LENSPONGE", cs%lensponge, &
125  "The length of the sponge layer (km).", &
126  default=10.0)
127 
128  allocate(cs%tr(isd:ied,jsd:jed,nz,ntr)) ; cs%tr(:,:,:,:) = 0.0
129  if (cs%mask_tracers) then
130  allocate(cs%tr_aux(isd:ied,jsd:jed,nz,ntr)) ; cs%tr_aux(:,:,:,:) = 0.0
131  endif
132 
133  do m=1,ntr
134  if (m < 10) then ; write(name,'("tr_RGC",I1.1)') m
135  else ; write(name,'("tr_RGC",I2.2)') m ; endif
136  write(longname,'("Concentration of RGC Tracer ",I2.2)') m
137  cs%tr_desc(m) = var_desc(name, units="kg kg-1", longname=longname, caller=mdl)
138 
139  ! This is needed to force the compiler not to do a copy in the registration calls.
140  tr_ptr => cs%tr(:,:,:,m)
141  ! Register the tracer for horizontal advection & diffusion.
142  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
143  name=name, longname=longname, units="kg kg-1", &
144  registry_diags=.true., flux_units="kg/s", &
145  restart_cs=restart_cs)
146  enddo
147 
148  cs%tr_Reg => tr_reg
149  register_rgc_tracer = .true.
150 end function register_rgc_tracer
151 
152 !> Initializes the NTR tracer fields in tr(:,:,:,:)
153 !! and it sets up the tracer output.
154 subroutine initialize_rgc_tracer(restart, day, G, GV, h, diag, OBC, CS, &
155  layer_CSp, sponge_CSp)
156 
157  type(ocean_grid_type), intent(in) :: g !< Grid structure.
158  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
159  logical, intent(in) :: restart !< .true. if the fields have already
160  !! been read from a restart file.
161  type(time_type), target, intent(in) :: day !< Time of the start of the run.
162  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
163  intent(in) :: h !< Layer thickness, in m or kg m-2.
164  type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
165  type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
166  !! whether, where, and what open boundary
167  !! conditions are used. This is not being used for now.
168  type(rgc_tracer_cs), pointer :: cs !< The control structure returned by a previous
169  !! call to RGC_register_tracer.
170  type(sponge_cs), pointer :: layer_csp !< A pointer to the control structure
171  type(ale_sponge_cs), pointer :: sponge_csp !< A pointer to the control structure for the
172  !! sponges, if they are in use. Otherwise this may be unassociated.
173 
174  real, allocatable :: temp(:,:,:)
175  real, pointer, dimension(:,:,:) :: &
176  obc_tr1_u => null(), & ! These arrays should be allocated and set to
177  obc_tr1_v => null() ! specify the values of tracer 1 that should come
178  ! in through u- and v- points through the open
179  ! boundary conditions, in the same units as tr.
180  character(len=16) :: name ! A variable's name in a NetCDF file.
181  character(len=72) :: longname ! The long name of that variable.
182  character(len=48) :: units ! The dimensions of the variable.
183  character(len=48) :: flux_units ! The units for tracer fluxes, usually
184  ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
185  real, pointer :: tr_ptr(:,:,:) => null()
186  real :: h_neglect ! A thickness that is so small it is usually lost
187  ! in roundoff and can be neglected [H ~> m or kg-2].
188  real :: e(szk_(g)+1), e_top, e_bot, d_tr ! Heights [Z ~> m].
189  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
190  integer :: isdb, iedb, jsdb, jedb
191  integer :: nzdata
192 
193  if (.not.associated(cs)) return
194  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
195  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
196  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
197  h_neglect = gv%H_subroundoff
198 
199  cs%Time => day
200  cs%diag => diag
201 
202  if (.not.restart) then
203  if (len_trim(cs%tracer_IC_file) >= 1) then
204  ! Read the tracer concentrations from a netcdf file.
205  if (.not.file_exists(cs%tracer_IC_file, g%Domain)) &
206  call mom_error(fatal, "RGC_initialize_tracer: Unable to open "// &
207  cs%tracer_IC_file)
208  do m=1,ntr
209  call query_vardesc(cs%tr_desc(m), name, caller="initialize_RGC_tracer")
210  call read_data(cs%tracer_IC_file, trim(name), &
211  cs%tr(:,:,:,m), domain=g%Domain%mpp_domain)
212  enddo
213  else
214  do m=1,ntr
215  do k=1,nz ; do j=js,je ; do i=is,ie
216  cs%tr(i,j,k,m) = 0.0
217  enddo ; enddo ; enddo
218  enddo
219  m=1
220  do j=js,je ; do i=is,ie
221  !set tracer to 1.0 in the surface of the continental shelf
222  if (g%geoLonT(i,j) <= (cs%CSL)) then
223  cs%tr(i,j,1,m) = 1.0 !first layer
224  endif
225  enddo ; enddo
226 
227  endif
228  endif ! restart
229 
230  if ( cs%use_sponge ) then
231 ! If sponges are used, this damps values to zero in the offshore boundary.
232 ! For any tracers that are not damped in the sponge, the call
233 ! to set_up_sponge_field can simply be omitted.
234  if (associated(sponge_csp)) then !ALE mode
235  nzdata = get_ale_sponge_nz_data(sponge_csp)
236  if (nzdata>0) then
237  allocate(temp(g%isd:g%ied,g%jsd:g%jed,nzdata))
238  do k=1,nzdata ; do j=js,je ; do i=is,ie
239  if (g%geoLonT(i,j) >= (cs%lenlon - cs%lensponge) .AND. g%geoLonT(i,j) <= cs%lenlon) then
240  temp(i,j,k) = 0.0
241  endif
242  enddo ; enddo; enddo
243  do m=1,1
244  ! This is needed to force the compiler not to do a copy in the sponge calls.
245  tr_ptr => cs%tr(:,:,:,m)
246  call set_up_ale_sponge_field(temp, g, tr_ptr, sponge_csp)
247  enddo
248  deallocate(temp)
249  endif
250 
251  elseif (associated(layer_csp)) then !layer mode
252  if (nz>0) then
253  allocate(temp(g%isd:g%ied,g%jsd:g%jed,nz))
254  do k=1,nz ; do j=js,je ; do i=is,ie
255  if (g%geoLonT(i,j) >= (cs%lenlon - cs%lensponge) .AND. g%geoLonT(i,j) <= cs%lenlon) then
256  temp(i,j,k) = 0.0
257  endif
258  enddo ; enddo; enddo
259  do m=1,1
260  tr_ptr => cs%tr(:,:,:,m)
261  call set_up_sponge_field(temp, tr_ptr, g, nz, layer_csp)
262  enddo
263  deallocate(temp)
264  endif
265  else
266  call mom_error(fatal, "RGC_initialize_tracer: "// &
267  "The pointer to sponge_CSp must be associated if SPONGE is defined.")
268  endif !selecting mode/calling error if no pointer
269  endif !using sponge
270 
271 end subroutine initialize_rgc_tracer
272 
273 !> This subroutine applies diapycnal diffusion and any other column
274 !! tracer physics or chemistry to the tracers from this file.
275 !! This is a simple example of a set of advected passive tracers.
276 subroutine rgc_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, &
277  evap_CFL_limit, minimum_forcing_depth)
278  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
279  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
280  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
281  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
282  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
283  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
284  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
285  intent(in) :: ea !< an array to which the amount of fluid entrained
286  !! from the layer above during this call will be
287  !! added [H ~> m or kg m-2].
288  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
289  intent(in) :: eb !< an array to which the amount of fluid entrained
290  !! from the layer below during this call will be
291  !! added [H ~> m or kg m-2].
292  type(forcing), intent(in) :: fluxes !< A structure containing pointers to any possible
293  !! forcing fields. Unused fields have NULL ptrs.
294  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s].
295  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
296  type(rgc_tracer_cs), pointer :: cs !< The control structure returned by a previous call.
297  real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can be
298  !! fluxed out of the top layer in a timestep [nondim].
299  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which fluxes
300  !! can be applied [H ~> m or kg m-2].
301 
302 ! The arguments to this subroutine are redundant in that
303 ! h_new[k] = h_old[k] + ea[k] - eb[k-1] + eb[k] - ea[k+1]
304 
305  real :: b1(szi_(g)) ! b1 and c1 are variables used by the
306  real :: c1(szi_(g),szk_(g)) ! tridiagonal solver.
307  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
308  real :: in_flux(szi_(g),szj_(g),2) ! total amount of tracer to be injected
309 
310  integer :: i, j, k, is, ie, js, je, nz, m
311  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
312 
313  if (.not.associated(cs)) return
314 
315  in_flux(:,:,:) = 0.0
316  m=1
317  do j=js,je ; do i=is,ie
318  !set tracer to 1.0 in the surface of the continental shelf
319  if (g%geoLonT(i,j) <= (cs%CSL)) then
320  cs%tr(i,j,1,m) = 1.0 !first layer
321  endif
322  enddo ; enddo
323 
324  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
325  do m=1,ntr
326  do k=1,nz ;do j=js,je ; do i=is,ie
327  h_work(i,j,k) = h_old(i,j,k)
328  enddo ; enddo ; enddo;
329  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m) , dt, fluxes, h_work, &
330  evap_cfl_limit, minimum_forcing_depth, in_flux(:,:,m))
331 
332  call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
333  enddo
334  else
335  do m=1,ntr
336  call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
337  enddo
338  endif
339 
340 end subroutine rgc_tracer_column_physics
341 
342 subroutine rgc_tracer_end(CS)
343  type(rgc_tracer_cs), pointer :: cs !< The control structure returned by a previous call to RGC_register_tracer.
344  integer :: m
345 
346  if (associated(cs)) then
347  if (associated(cs%tr)) deallocate(cs%tr)
348  deallocate(cs)
349  endif
350 end subroutine rgc_tracer_end
351 
352 end module rgc_tracer
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
rgc_tracer::initialize_rgc_tracer
subroutine, public initialize_rgc_tracer(restart, day, G, GV, h, diag, OBC, CS, layer_CSp, sponge_CSp)
Initializes the NTR tracer fields in tr(:,:,:,:) and it sets up the tracer output.
Definition: RGC_tracer.F90:156
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
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
rgc_tracer
This module contains the routines used to set up a dynamically passive tracer. Set up and use passive...
Definition: RGC_tracer.F90:12
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
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_ale_sponge::ale_sponge_cs
ALE sponge control structure.
Definition: MOM_ALE_sponge.F90:84
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_ale_sponge::get_ale_sponge_nz_data
integer function, public get_ale_sponge_nz_data(CS)
Return the number of layers in the data with a fixed ALE sponge, or 0 if there are no sponge columns ...
Definition: MOM_ALE_sponge.F90:326
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_ale_sponge
This module contains the routines used to apply sponge layers when using the ALE mode.
Definition: MOM_ALE_sponge.F90:11
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_restart::mom_restart_cs
A restart registry and the control structure for restarts.
Definition: MOM_restart.F90:72
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_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
rgc_tracer::ntr
integer, parameter ntr
The number of tracers in this module.
Definition: RGC_tracer.F90:42
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_ale_sponge::set_up_ale_sponge_field
Store the reference profile at h points for a variable.
Definition: MOM_ALE_sponge.F90:33
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
rgc_tracer::rgc_tracer_column_physics
subroutine, public rgc_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: RGC_tracer.F90:278
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_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
rgc_tracer::rgc_tracer_cs
tracer control structure
Definition: RGC_tracer.F90:45
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
rgc_tracer::register_rgc_tracer
logical function, public register_rgc_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
This subroutine is used to register tracer fields.
Definition: RGC_tracer.F90:68
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
rgc_tracer::rgc_tracer_end
subroutine, public rgc_tracer_end(CS)
Definition: RGC_tracer.F90:343
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