mom_tracer_advect module reference¶
This module contains the subroutines that advect tracers along coordinate surfaces.
Data Types¶
Control structure for this module. |
Functions/Subroutines¶
This routine time steps the tracer concentration using a monotonic, conservative, weakly diffusive scheme. |
|
This subroutine does 1-d flux-form advection in the zonal direction using a monotonic piecewise linear scheme. |
|
This subroutine does 1-d flux-form advection using a monotonic piecewise linear scheme. |
|
Initialize lateral tracer advection module. |
|
Close the tracer advection module. |
Detailed Description¶
This program contains the subroutines that advect tracers horizontally (i.e. along layers).
section_mom_advect_intro¶
advect_tracer advects tracer concentrations using a combination of the modified flux advection scheme from Easter (Mon. Wea. Rev., 1993) with tracer distributions given by the monotonic modified van Leer scheme proposed by Lin et al. (Mon. Wea. Rev., 1994). This scheme conserves the total amount of tracer while avoiding spurious maxima and minima of the tracer concentration. If a higher order accuracy scheme is needed, suggest monotonic piecewise parabolic method, as described in Carpenter et al. (MWR, 1990).
advect_tracer has 4 arguments, described below. This subroutine determines the volume of a layer in a grid cell at the previous instance when the tracer concentration was changed, so it is essential that the volume fluxes should be correct. It is also important that the tracer advection occurs before each calculation of the diabatic forcing.
Type Documentation¶
-
type
mom_tracer_advect/
tracer_advect_cs
¶ Control structure for this module.
- Type fields:
%
dt
[real] :: The baroclinic dynamics time step [T ~> s].%
diag
[type( diag_ctrl ),pointer] :: A structure that is used to regulate the timing of diagnostic output.%
debug
[logical] :: If true, write verbose checksums for debugging purposes.%
useppm
[logical] :: If true, use PPM instead of PLM.%
usehuynh
[logical] :: If true, use the Huynh scheme for PPM interface values.%
pass_uhr_vhr_t_hprev
[type(group_pass_type)] :: A structure used for group passes.
Function/Subroutine Documentation¶
-
subroutine
mom_tracer_advect/
advect_tracer
(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first_in, vol_prev, max_iter_in, update_vol_prev, uhr_out, vhr_out)¶ This routine time steps the tracer concentration using a monotonic, conservative, weakly diffusive scheme.
- Parameters:
g :: [inout] ocean grid structure
gv :: [in] ocean vertical grid structure
h_end :: [in] Layer thickness after advection [H ~> m or kg m-2]
uhtr :: [in] Accumulated volume or mass flux through the
vhtr :: [in] Accumulated volume or mass flux through the
obc :: specifies whether, where, and what OBCs are used
dt :: [in] time increment [T ~> s]
us :: [in] A dimensional unit scaling type
cs :: control structure for module
reg :: pointer to tracer registry
x_first_in :: [in] If present, indicate whether to update first in the x- or y-direction.
vol_prev :: [inout] Cell volume before advection [H L2 ~> m3 or kg].
max_iter_in :: [in] The maximum number of iterations
update_vol_prev :: [in] If present and true, update vol_prev to return its value after the tracer have been updated.
uhr_out :: [out] Remaining accumulated volume or mass fluxes
vhr_out :: [out] Remaining accumulated volume or mass fluxes
- Call to:
advect_x
advect_y
id_clock_advect
id_clock_pass
id_clock_sync
mom_error_handler::mom_error
- Called from:
mom_offline_main::offline_advection_ale
mom_offline_main::offline_advection_layer
mom_offline_main::offline_redistribute_residual
-
subroutine
mom_tracer_advect/
advect_x
(Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, is, ie, js, je, k, G, GV, US, usePPM, useHuynh)¶ This subroutine does 1-d flux-form advection in the zonal direction using a monotonic piecewise linear scheme.
- Parameters:
g :: [inout] The ocean’s grid structure
gv :: [in] The ocean’s vertical grid structure
ntr :: [in] The number of tracers
tr :: [inout] The array of registered tracers to work on
hprev :: [inout] cell volume at the end of previous tracer change [H L2 ~> m3 or kg]
uhr :: [inout] accumulated volume/mass flux through the zonal face [H L2 ~> m3 or kg]
uh_neglect :: [in] A tiny zonal mass flux that can be neglected [H L2 ~> m3 or kg]
obc :: specifies whether, where, and what OBCs are used
domore_u :: [inout] If true, there is more advection to be done in this u-row
idt :: [in] The inverse of dt [T-1 ~> s-1]
is :: [in] The starting tracer i-index to work on
ie :: [in] The ending tracer i-index to work on
js :: [in] The starting tracer j-index to work on
je :: [in] The ending tracer j-index to work on
k :: [in] The k-level to work on
us :: [in] A dimensional unit scaling type
useppm :: [in] If true, use PPM instead of PLM
usehuynh :: [in] If true, use the Huynh scheme for PPM interface values
- Called from:
-
subroutine
mom_tracer_advect/
advect_y
(Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, is, ie, js, je, k, G, GV, US, usePPM, useHuynh)¶ This subroutine does 1-d flux-form advection using a monotonic piecewise linear scheme.
- Parameters:
g :: [inout] The ocean’s grid structure
gv :: [in] The ocean’s vertical grid structure
ntr :: [in] The number of tracers
tr :: [inout] The array of registered tracers to work on
hprev :: [inout] cell volume at the end of previous tracer change [H L2 ~> m3 or kg]
vhr :: [inout] accumulated volume/mass flux through the meridional face [H L2 ~> m3 or kg]
vh_neglect :: [inout] A tiny meridional mass flux that can be neglected [H L2 ~> m3 or kg]
obc :: specifies whether, where, and what OBCs are used
domore_v :: [inout] If true, there is more advection to be done in this v-row
idt :: [in] The inverse of dt [T-1 ~> s-1]
is :: [in] The starting tracer i-index to work on
ie :: [in] The ending tracer i-index to work on
js :: [in] The starting tracer j-index to work on
je :: [in] The ending tracer j-index to work on
k :: [in] The k-level to work on
us :: [in] A dimensional unit scaling type
useppm :: [in] If true, use PPM instead of PLM
usehuynh :: [in] If true, use the Huynh scheme for PPM interface values
- Called from:
-
subroutine
mom_tracer_advect/
tracer_advect_init
(Time, G, US, param_file, diag, CS)¶ Initialize lateral tracer advection module.
- Parameters:
time :: [in] current model time
g :: [in] ocean grid structure
us :: [in] A dimensional unit scaling type
param_file :: [in] open file to parse for model parameters
diag :: [inout] regulates diagnostic output
cs :: module control structure
- Call to:
id_clock_advect
id_clock_pass
id_clock_sync
mom_error_handler::mom_error
-
subroutine
mom_tracer_advect/
tracer_advect_end
(CS)¶ Close the tracer advection module.
- Parameters:
cs :: module control structure
- Called from: