MOM6
advection_test_tracer.F90
Go to the documentation of this file.
1 !> This tracer package is used to test advection schemes
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
20 use mom_variables, only : surface
22 
23 use coupler_types_mod, only : coupler_type_set_data, ind_csurf
25 
26 implicit none ; private
27 
28 #include <MOM_memory.h>
29 
33 
34 integer, parameter :: ntr = 11 !< The number of tracers in this module.
35 
36 !> The control structure for the advect_test_tracer module
37 type, public :: advection_test_tracer_cs ; private
38  integer :: ntr = ntr !< Number of tracers in this module
39  logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
40  character(len=200) :: tracer_ic_file !< The full path to the IC file, or " " to initialize internally.
41  type(time_type), pointer :: time => null() !< A pointer to the ocean model's clock.
42  type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the MOM tracer registry
43  real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this subroutine, in g m-3?
44  real :: land_val(ntr) = -1.0 !< The value of tr used where land is masked out.
45  logical :: use_sponge !< If true, sponges may be applied somewhere in the domain.
46  logical :: tracers_may_reinit !< If true, the tracers may be set up via the initialization code if
47  !! they are not found in the restart files. Otherwise it is a fatal error
48  !! if the tracers are not found in the restart files of a restarted run.
49  real :: x_origin !< Parameters describing the test functions
50  real :: x_width !< Parameters describing the test functions
51  real :: y_origin !< Parameters describing the test functions
52  real :: y_width !< Parameters describing the test functions
53 
54  integer, dimension(NTR) :: ind_tr !< Indices returned by aof_set_coupler_flux if it is used and
55  !! the surface tracer concentrations are to be provided to the coupler.
56 
57  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
58  !! regulate the timing of diagnostic output.
59  type(mom_restart_cs), pointer :: restart_csp => null() !< A pointer to the restart control structure.
60 
61  type(vardesc) :: tr_desc(ntr) !< Descriptions and metadata for the tracers
63 
64 contains
65 
66 !> Register tracer fields and subroutines to be used with MOM.
67 function register_advection_test_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 to parse for run-time parameters
71  type(advection_test_tracer_cs), pointer :: cs !< The control structure returned by a previous
72  !! call to register_advection_test_tracer.
73  type(tracer_registry_type), pointer :: tr_reg !< A pointer that is set to point to the control
74  !! structure for the tracer advection and
75  !! diffusion module
76  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure
77 
78  ! Local variables
79  character(len=80) :: name, longname
80 ! This include declares and sets the variable "version".
81 #include "version_variable.h"
82  character(len=40) :: mdl = "advection_test_tracer" ! This module's name.
83  character(len=200) :: inputdir
84  character(len=48) :: flux_units ! The units for tracer fluxes, usually
85  ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
86  real, pointer :: tr_ptr(:,:,:) => null()
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_advection_test_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 
101  call get_param(param_file, mdl, "ADVECTION_TEST_X_ORIGIN", cs%x_origin, &
102  "The x-coorindate of the center of the test-functions.", default=0.)
103  call get_param(param_file, mdl, "ADVECTION_TEST_Y_ORIGIN", cs%y_origin, &
104  "The y-coorindate of the center of the test-functions.", default=0.)
105  call get_param(param_file, mdl, "ADVECTION_TEST_X_WIDTH", cs%x_width, &
106  "The x-width of the test-functions.", default=0.)
107  call get_param(param_file, mdl, "ADVECTION_TEST_Y_WIDTH", cs%y_width, &
108  "The y-width of the test-functions.", default=0.)
109  call get_param(param_file, mdl, "ADVECTION_TEST_TRACER_IC_FILE", cs%tracer_IC_file, &
110  "The name of a file from which to read the initial "//&
111  "conditions for the tracers, or blank to initialize "//&
112  "them internally.", default=" ")
113 
114  if (len_trim(cs%tracer_IC_file) >= 1) then
115  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
116  cs%tracer_IC_file = trim(slasher(inputdir))//trim(cs%tracer_IC_file)
117  call log_param(param_file, mdl, "INPUTDIR/ADVECTION_TEST_TRACER_IC_FILE", &
118  cs%tracer_IC_file)
119  endif
120  call get_param(param_file, mdl, "SPONGE", cs%use_sponge, &
121  "If true, sponges may be applied anywhere in the domain. "//&
122  "The exact location and properties of those sponges are "//&
123  "specified from MOM_initialization.F90.", default=.false.)
124 
125  call get_param(param_file, mdl, "TRACERS_MAY_REINIT", cs%tracers_may_reinit, &
126  "If true, tracers may go through the initialization code "//&
127  "if they are not found in the restart files. Otherwise "//&
128  "it is a fatal error if the tracers are not found in the "//&
129  "restart files of a restarted run.", default=.false.)
130 
131 
132  allocate(cs%tr(isd:ied,jsd:jed,nz,ntr)) ; cs%tr(:,:,:,:) = 0.0
133 
134  do m=1,ntr
135  if (m < 10) then ; write(name,'("tr",I1.1)') m
136  else ; write(name,'("tr",I2.2)') m ; endif
137  write(longname,'("Concentration of Tracer ",I2.2)') m
138  cs%tr_desc(m) = var_desc(name, units="kg kg-1", longname=longname, caller=mdl)
139  if (gv%Boussinesq) then ; flux_units = "kg kg-1 m3 s-1"
140  else ; flux_units = "kg s-1" ; endif
141 
142 
143  ! This is needed to force the compiler not to do a copy in the registration
144  ! calls. Curses on the designers and implementers of Fortran90.
145  tr_ptr => cs%tr(:,:,:,m)
146  ! Register the tracer for horizontal advection, diffusion, and restarts.
147  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
148  name=name, longname=longname, units="kg kg-1", &
149  registry_diags=.true., flux_units=flux_units, &
150  restart_cs=restart_cs, mandatory=.not.cs%tracers_may_reinit)
151 
152  ! Set coupled_tracers to be true (hard-coded above) to provide the surface
153  ! values to the coupler (if any). This is meta-code and its arguments will
154  ! currently (deliberately) give fatal errors if it is used.
155  if (cs%coupled_tracers) &
156  cs%ind_tr(m) = aof_set_coupler_flux(trim(name)//'_flux', &
157  flux_type=' ', implementation=' ', caller="register_advection_test_tracer")
158  enddo
159 
160  cs%tr_Reg => tr_reg
161  cs%restart_CSp => restart_cs
164 
165 !> Initializes the NTR tracer fields in tr(:,:,:,:) and it sets up the tracer output.
166 subroutine initialize_advection_test_tracer(restart, day, G, GV, h,diag, OBC, CS, &
167  sponge_CSp)
168  logical, intent(in) :: restart !< .true. if the fields have already
169  !! been read from a restart file.
170  type(time_type), target, intent(in) :: day !< Time of the start of the run.
171  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
172  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
173  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
174  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
175  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate
176  !! diagnostic output.
177  type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
178  !! whether, where, and what open boundary
179  !! conditions are used.
180  type(advection_test_tracer_cs), pointer :: cs !< The control structure returned by a previous
181  !! call to register_advection_test_tracer.
182  type(sponge_cs), pointer :: sponge_csp !< Pointer to the control structure for the sponges.
183 
184  ! Local variables
185  real, allocatable :: temp(:,:,:)
186  real, pointer, dimension(:,:,:) :: &
187  obc_tr1_u => null(), & ! These arrays should be allocated and set to
188  obc_tr1_v => null() ! specify the values of tracer 1 that should come
189  ! in through u- and v- points through the open
190  ! boundary conditions, in the same units as tr.
191  character(len=16) :: name ! A variable's name in a NetCDF file.
192  character(len=72) :: longname ! The long name of that variable.
193  character(len=48) :: units ! The dimensions of the variable.
194  character(len=48) :: flux_units ! The units for tracer fluxes, usually
195  ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
196  real, pointer :: tr_ptr(:,:,:) => null()
197  real :: h_neglect ! A thickness that is so small it is usually lost
198  ! in roundoff and can be neglected [H ~> m or kg m-2].
199  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
200  integer :: isdb, iedb, jsdb, jedb
201  real :: tmpx, tmpy, locx, locy
202 
203  if (.not.associated(cs)) return
204  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
205  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
206  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
207  h_neglect = gv%H_subroundoff
208 
209  cs%diag => diag
210  cs%ntr = ntr
211  do m=1,ntr
212  call query_vardesc(cs%tr_desc(m), name=name, &
213  caller="initialize_advection_test_tracer")
214  if ((.not.restart) .or. (cs%tracers_may_reinit .and. .not. &
215  query_initialized(cs%tr(:,:,:,m), name, cs%restart_CSp))) then
216  do k=1,nz ; do j=js,je ; do i=is,ie
217  cs%tr(i,j,k,m) = 0.0
218  enddo ; enddo ; enddo
219  k=1 ! Square wave
220  do j=js,je ; do i=is,ie
221  if (abs(g%geoLonT(i,j)-cs%x_origin)<0.5*cs%x_width .and. &
222  abs(g%geoLatT(i,j)-cs%y_origin)<0.5*cs%y_width) cs%tr(i,j,k,m) = 1.0
223  enddo ; enddo
224  k=2 ! Triangle wave
225  do j=js,je ; do i=is,ie
226  locx=abs(g%geoLonT(i,j)-cs%x_origin)/cs%x_width
227  locy=abs(g%geoLatT(i,j)-cs%y_origin)/cs%y_width
228  cs%tr(i,j,k,m) = max(0.0, 1.0-locx)*max(0.0, 1.0-locy)
229  enddo ; enddo
230  k=3 ! Cosine bell
231  do j=js,je ; do i=is,ie
232  locx=min(1.0, abs(g%geoLonT(i,j)-cs%x_origin)/cs%x_width)*(acos(0.0)*2.)
233  locy=min(1.0, abs(g%geoLatT(i,j)-cs%y_origin)/cs%y_width)*(acos(0.0)*2.)
234  cs%tr(i,j,k,m) = (1.0+cos(locx))*(1.0+cos(locy))*0.25
235  enddo ; enddo
236  k=4 ! Cylinder
237  do j=js,je ; do i=is,ie
238  locx=abs(g%geoLonT(i,j)-cs%x_origin)/cs%x_width
239  locy=abs(g%geoLatT(i,j)-cs%y_origin)/cs%y_width
240  if (locx**2+locy**2<=1.0) cs%tr(i,j,k,m) = 1.0
241  enddo ; enddo
242  k=5 ! Cut cylinder
243  do j=js,je ; do i=is,ie
244  locx=(g%geoLonT(i,j)-cs%x_origin)/cs%x_width
245  locy=(g%geoLatT(i,j)-cs%y_origin)/cs%y_width
246  if (locx**2+locy**2<=1.0) cs%tr(i,j,k,m) = 1.0
247  if (locx>0.0.and.abs(locy)<0.2) cs%tr(i,j,k,m) = 0.0
248  enddo ; enddo
249  endif ! restart
250  enddo
251 
252 
254 
255 
256 !> Applies diapycnal diffusion and any other column tracer physics or chemistry to the tracers
257 !! from this package. This is a simple example of a set of advected passive tracers.
258 subroutine advection_test_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, &
259  evap_CFL_limit, minimum_forcing_depth)
260  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
261  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
262  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
263  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
264  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
265  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
266  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
267  intent(in) :: ea !< an array to which the amount of fluid entrained
268  !! from the layer above during this call will be
269  !! added [H ~> m or kg m-2].
270  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
271  intent(in) :: eb !< an array to which the amount of fluid entrained
272  !! from the layer below during this call will be
273  !! added [H ~> m or kg m-2].
274  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
275  !! and tracer forcing fields. Unused fields have NULL ptrs.
276  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
277  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
278  type(advection_test_tracer_cs), pointer :: cs !< The control structure returned by a previous
279  !! call to register_advection_test_tracer.
280  real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can
281  !! be fluxed out of the top layer in a timestep [nondim]
282  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
283  !! fluxes can be applied [H ~> m or kg m-2]
284 ! This subroutine applies diapycnal diffusion and any other column
285 ! tracer physics or chemistry to the tracers from this file.
286 ! This is a simple example of a set of advected passive tracers.
287 
288 ! The arguments to this subroutine are redundant in that
289 ! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
290  ! Local variables
291  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
292  real :: b1(szi_(g)) ! b1 and c1 are variables used by the
293  real :: c1(szi_(g),szk_(g)) ! tridiagonal solver.
294  integer :: i, j, k, is, ie, js, je, nz, m
295  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
296 
297  if (.not.associated(cs)) return
298 
299  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
300  do m=1,cs%ntr
301  do k=1,nz ;do j=js,je ; do i=is,ie
302  h_work(i,j,k) = h_old(i,j,k)
303  enddo ; enddo ; enddo
304  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
305  evap_cfl_limit, minimum_forcing_depth)
306  call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
307  enddo
308  else
309  do m=1,ntr
310  call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
311  enddo
312  endif
313 
315 
316 !> This subroutine extracts the surface fields from this tracer package that
317 !! are to be shared with the atmosphere in coupled configurations.
318 !! This particular tracer package does not report anything back to the coupler.
319 subroutine advection_test_tracer_surface_state(state, h, G, CS)
320  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
321  type(surface), intent(inout) :: state !< A structure containing fields that
322  !! describe the surface state of the ocean.
323  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
324  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
325  type(advection_test_tracer_cs), pointer :: cs !< The control structure returned by a previous
326  !! call to register_advection_test_tracer.
327 
328  ! This particular tracer package does not report anything back to the coupler.
329  ! The code that is here is just a rough guide for packages that would.
330 
331  integer :: m, is, ie, js, je, isd, ied, jsd, jed
332  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
333  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
334 
335  if (.not.associated(cs)) return
336 
337  if (cs%coupled_tracers) then
338  do m=1,cs%ntr
339  ! This call loads the surface values into the appropriate array in the
340  ! coupler-type structure.
341  call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
342  state%tr_fields, idim=(/isd, is, ie, ied/), &
343  jdim=(/jsd, js, je, jed/) )
344  enddo
345  endif
346 
348 
349 !> Calculate the mass-weighted integral of all tracer stocks, returning the number of stocks it has calculated.
350 !! If the stock_index is present, only the stock corresponding to that coded index is returned.
351 function advection_test_stock(h, stocks, G, GV, CS, names, units, stock_index)
352  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
353  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
354  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each
355  !! tracer, in kg times concentration units [kg conc].
356  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
357  type(advection_test_tracer_cs), pointer :: cs !< The control structure returned by a previous
358  !! call to register_advection_test_tracer.
359  character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated.
360  character(len=*), dimension(:), intent(out) :: units !< the units of the stocks calculated.
361  integer, optional, intent(in) :: stock_index !< the coded index of a specific stock being sought.
362  integer :: advection_test_stock !< the number of stocks calculated here.
363 
364  integer :: i, j, k, is, ie, js, je, nz, m
365  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
366 
368  if (.not.associated(cs)) return
369  if (cs%ntr < 1) return
370 
371  if (present(stock_index)) then ; if (stock_index > 0) then
372  ! Check whether this stock is available from this routine.
373 
374  ! No stocks from this routine are being checked yet. Return 0.
375  return
376  endif ; endif
377 
378  do m=1,cs%ntr
379  call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="advection_test_stock")
380  stocks(m) = 0.0
381  do k=1,nz ; do j=js,je ; do i=is,ie
382  stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
383  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
384  enddo ; enddo ; enddo
385  stocks(m) = gv%H_to_kg_m2 * stocks(m)
386  enddo
387  advection_test_stock = cs%ntr
388 
389 end function advection_test_stock
390 
391 !> Deallocate memory associated with this module
392 subroutine advection_test_tracer_end(CS)
393  type(advection_test_tracer_cs), pointer :: cs !< The control structure returned by a previous
394  !! call to register_advection_test_tracer.
395  integer :: m
396 
397  if (associated(cs)) then
398  if (associated(cs%tr)) deallocate(cs%tr)
399  deallocate(cs)
400  endif
401 end subroutine advection_test_tracer_end
402 
403 end module advection_test_tracer
advection_test_tracer::advection_test_tracer_cs
The control structure for the advect_test_tracer module.
Definition: advection_test_tracer.F90:37
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
advection_test_tracer::advection_test_tracer_surface_state
subroutine, public advection_test_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: advection_test_tracer.F90:320
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
advection_test_tracer::initialize_advection_test_tracer
subroutine, public initialize_advection_test_tracer(restart, day, G, GV, h, diag, OBC, CS, sponge_CSp)
Initializes the NTR tracer fields in tr(:,:,:,:) and it sets up the tracer output.
Definition: advection_test_tracer.F90:168
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
advection_test_tracer::advection_test_tracer_column_physics
subroutine, public advection_test_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, evap_CFL_limit, minimum_forcing_depth)
Applies diapycnal diffusion and any other column tracer physics or chemistry to the tracers from this...
Definition: advection_test_tracer.F90:260
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
advection_test_tracer::ntr
integer, parameter ntr
The number of tracers in this module.
Definition: advection_test_tracer.F90:34
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
advection_test_tracer::register_advection_test_tracer
logical function, public register_advection_test_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
Register tracer fields and subroutines to be used with MOM.
Definition: advection_test_tracer.F90:68
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
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
advection_test_tracer::advection_test_tracer_end
subroutine, public advection_test_tracer_end(CS)
Deallocate memory associated with this module.
Definition: advection_test_tracer.F90:393
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
advection_test_tracer
This tracer package is used to test advection schemes.
Definition: advection_test_tracer.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
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
advection_test_tracer::advection_test_stock
integer function, public advection_test_stock(h, stocks, G, GV, CS, names, units, stock_index)
Calculate the mass-weighted integral of all tracer stocks, returning the number of stocks it has calc...
Definition: advection_test_tracer.F90:352
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
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