Go to the documentation of this file.
19 implicit none ;
private
21 #include <MOM_memory.h>
27 integer :: continuity_scheme
42 subroutine continuity(u, v, hin, h, uh, vh, dt, G, GV, US, CS, uhbt, vhbt, OBC, &
43 visc_rem_u, visc_rem_v, u_cor, v_cor, BT_cont)
46 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
48 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
50 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
52 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
54 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
57 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
60 real,
intent(in) :: dt
63 real,
dimension(SZIB_(G),SZJ_(G)), &
64 optional,
intent(in) :: uhbt
66 real,
dimension(SZI_(G),SZJB_(G)), &
67 optional,
intent(in) :: vhbt
70 optional,
pointer :: obc
71 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
72 optional,
intent(in) :: visc_rem_u
76 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
77 optional,
intent(in) :: visc_rem_v
81 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
82 optional,
intent(out) :: u_cor
84 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
85 optional,
intent(out) :: v_cor
88 optional,
pointer :: bt_cont
91 if (
present(visc_rem_u) .neqv.
present(visc_rem_v))
call mom_error(fatal, &
92 "MOM_continuity: Either both visc_rem_u and visc_rem_v or neither"// &
93 " one must be present in call to continuity.")
94 if (
present(u_cor) .neqv.
present(v_cor))
call mom_error(fatal, &
95 "MOM_continuity: Either both u_cor and v_cor or neither"// &
96 " one must be present in call to continuity.")
99 call continuity_ppm(u, v, hin, h, uh, vh, dt, g, gv, us, cs%PPM_CSp, uhbt, vhbt, obc, &
100 visc_rem_u, visc_rem_v, u_cor, v_cor, bt_cont=bt_cont)
102 call mom_error(fatal,
"continuity: Unrecognized value of continuity_scheme")
109 type(time_type),
target,
intent(in) :: time
114 type(
diag_ctrl),
target,
intent(inout) :: diag
117 #include "version_variable.h"
118 character(len=40) :: mdl =
"MOM_continuity"
119 character(len=20) :: tmpstr
121 if (
associated(cs))
then
122 call mom_error(warning,
"continuity_init called with associated control structure.")
129 call get_param(param_file, mdl,
"CONTINUITY_SCHEME", tmpstr, &
130 "CONTINUITY_SCHEME selects the discretization for the "//&
131 "continuity solver. The only valid value currently is: \n"//&
132 "\t PPM - use a positive-definite (or monotonic) \n"//&
133 "\t piecewise parabolic reconstruction solver.", &
136 tmpstr =
uppercase(tmpstr) ; cs%continuity_scheme = 0
137 select case (trim(tmpstr))
140 call mom_mesg(
'continuity_init: CONTINUITY_SCHEME ="'//trim(tmpstr)//
'"', 0)
141 call mom_mesg(
"continuity_init: The only valid value is currently "// &
143 call mom_error(fatal,
"continuity_init: Unrecognized setting "// &
144 "#define CONTINUITY_SCHEME "//trim(tmpstr)//
" found in input file.")
148 call continuity_ppm_init(time, g, gv, us, param_file, diag, cs%PPM_CSp)
162 stencil = continuity_ppm_stencil(cs%PPM_CSp)
Provides a transparent vertical ocean grid type and supporting routines.
Control structure for mom_continuity.
subroutine, public continuity_ppm_init(Time, G, GV, US, param_file, diag, CS)
Initializes continuity_ppm_cs.
An overloaded interface to log version information about modules.
integer function, public continuity_stencil(CS)
continuity_stencil returns the continuity solver stencil size
Handy functions for manipulating strings.
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
integer, parameter ppm_scheme
Enumerated constant to select PPM.
subroutine, public continuity_init(Time, G, GV, US, param_file, diag, CS)
Initializes continuity_cs.
A structure that can be parsed to read and document run-time parameters.
Container for information about the summed layer transports and how they will vary as the barotropic ...
An overloaded interface to read and log the values of various types of parameters.
character(len=len(input_string)) function, public uppercase(input_string)
Return a string in which all uppercase letters have been replaced by their lowercase counterparts.
Describes various unit conversion factors.
subroutine, public continuity(u, v, hin, h, uh, vh, dt, G, GV, US, CS, uhbt, vhbt, OBC, visc_rem_u, visc_rem_v, u_cor, v_cor, BT_cont)
Time steps the layer thicknesses, using a monotonically limited, directionally split PPM scheme,...
Solve the layer continuity equation using the PPM method for layer fluxes.
subroutine, public continuity_ppm_end(CS)
Destructor for continuity_ppm_cs.
subroutine, public continuity_end(CS)
Destructor for continuity_cs.
Describes the vertical ocean grid, including unit conversion factors.
Provides transparent structures with groups of MOM6 variables and supporting routines.
Control structure for mom_continuity_ppm.
Controls where open boundary conditions are applied.
The MOM6 facility to parse input files for runtime parameters.
Provides the ocean grid type.
character(len=20), parameter ppm_string
String to select PPM.
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Solve the layer continuity equation.
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...
Routines for error handling and I/O management.
Ocean grid type. See mom_grid for details.
integer function, public continuity_ppm_stencil(CS)
continuity_PPM_stencil returns the continuity solver stencil size
subroutine, public continuity_ppm(u, v, hin, h, uh, vh, dt, G, GV, US, CS, uhbt, vhbt, OBC, visc_rem_u, visc_rem_v, u_cor, v_cor, BT_cont)
Time steps the layer thicknesses, using a monotonically limit, directionally split PPM scheme,...