MOM6
MOM_continuity.F90
Go to the documentation of this file.
1 !> Solve the layer continuity equation.
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
9 use mom_diag_mediator, only : time_type, diag_ctrl
10 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
13 use mom_grid, only : ocean_grid_type
16 use mom_variables, only : bt_cont_type
18 
19 implicit none ; private
20 
21 #include <MOM_memory.h>
22 
24 
25 !> Control structure for mom_continuity
26 type, public :: continuity_cs ; private
27  integer :: continuity_scheme !< Selects the discretization for the continuity solver.
28  !! Valid values are:
29  !! - PPM - A directionally split piecewise parabolic reconstruction solver.
30  !! The default, PPM, seems most appropriate for use with our current
31  !! time-splitting strategies.
32  type(continuity_ppm_cs), pointer :: ppm_csp => null() !< Control structure for mom_continuity_ppm
33 end type continuity_cs
34 
35 integer, parameter :: ppm_scheme = 1 !< Enumerated constant to select PPM
36 character(len=20), parameter :: ppm_string = "PPM" !< String to select PPM
37 
38 contains
39 
40 !> Time steps the layer thicknesses, using a monotonically limited, directionally split PPM scheme,
41 !! based on Lin (1994).
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)
44  type(ocean_grid_type), intent(inout) :: g !< Ocean grid structure.
45  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure.
46  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
47  intent(in) :: u !< Zonal velocity [L T-1 ~> m s-1].
48  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
49  intent(in) :: v !< Meridional velocity [L T-1 ~> m s-1].
50  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
51  intent(in) :: hin !< Initial layer thickness [H ~> m or kg m-2].
52  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
53  intent(inout) :: h !< Final layer thickness [H ~> m or kg m-2].
54  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
55  intent(out) :: uh !< Volume flux through zonal faces =
56  !! u*h*dy [H L2 T-1 ~> m3 s-1 or kg s-1].
57  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
58  intent(out) :: vh !< Volume flux through meridional faces =
59  !! v*h*dx [H L2 T-1 ~> m3 s-1 or kg s-1].
60  real, intent(in) :: dt !< Time increment [T ~> s].
61  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
62  type(continuity_cs), pointer :: cs !< Control structure for mom_continuity.
63  real, dimension(SZIB_(G),SZJ_(G)), &
64  optional, intent(in) :: uhbt !< The vertically summed volume
65  !! flux through zonal faces [H L2 T-1 ~> m3 s-1 or kg s-1].
66  real, dimension(SZI_(G),SZJB_(G)), &
67  optional, intent(in) :: vhbt !< The vertically summed volume
68  !! flux through meridional faces [H L2 T-1 ~> m3 s-1 or kg s-1].
69  type(ocean_obc_type), &
70  optional, pointer :: obc !< Open boundaries control structure.
71  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
72  optional, intent(in) :: visc_rem_u !< Both the fraction of
73  !! zonal momentum that remains after a time-step of viscosity, and the fraction of a time-step's
74  !! worth of a barotropic acceleration that a layer experiences after viscosity is applied.
75  !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom).
76  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
77  optional, intent(in) :: visc_rem_v !< Both the fraction of
78  !! meridional momentum that remains after a time-step of viscosity, and the fraction of a time-step's
79  !! worth of a barotropic acceleration that a layer experiences after viscosity is applied.
80  !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom).
81  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
82  optional, intent(out) :: u_cor !< The zonal velocities that
83  !! give uhbt as the depth-integrated transport [L T-1 ~> m s-1].
84  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
85  optional, intent(out) :: v_cor !< The meridional velocities that
86  !! give vhbt as the depth-integrated transport [L T-1 ~> m s-1].
87  type(bt_cont_type), &
88  optional, pointer :: bt_cont !< A structure with elements
89  !! that describe the effective open face areas as a function of barotropic flow.
90 
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.")
97 
98  if (cs%continuity_scheme == ppm_scheme) then
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)
101  else
102  call mom_error(fatal, "continuity: Unrecognized value of continuity_scheme")
103  endif
104 
105 end subroutine continuity
106 
107 !> Initializes continuity_cs
108 subroutine continuity_init(Time, G, GV, US, param_file, diag, CS)
109  type(time_type), target, intent(in) :: time !< Current model time.
110  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure.
111  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure.
112  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
113  type(param_file_type), intent(in) :: param_file !< Parameter file handles.
114  type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure.
115  type(continuity_cs), pointer :: cs !< Control structure for mom_continuity.
116 ! This include declares and sets the variable "version".
117 #include "version_variable.h"
118  character(len=40) :: mdl = "MOM_continuity" ! This module's name.
119  character(len=20) :: tmpstr
120 
121  if (associated(cs)) then
122  call mom_error(warning, "continuity_init called with associated control structure.")
123  return
124  endif
125  allocate(cs)
126 
127  ! Read all relevant parameters and write them to the model log.
128  call log_version(param_file, mdl, version, "")
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.", &
134  default=ppm_string)
135 
136  tmpstr = uppercase(tmpstr) ; cs%continuity_scheme = 0
137  select case (trim(tmpstr))
138  case (ppm_string) ; cs%continuity_scheme = ppm_scheme
139  case default
140  call mom_mesg('continuity_init: CONTINUITY_SCHEME ="'//trim(tmpstr)//'"', 0)
141  call mom_mesg("continuity_init: The only valid value is currently "// &
142  trim(ppm_string), 0)
143  call mom_error(fatal, "continuity_init: Unrecognized setting "// &
144  "#define CONTINUITY_SCHEME "//trim(tmpstr)//" found in input file.")
145  end select
146 
147  if (cs%continuity_scheme == ppm_scheme) then
148  call continuity_ppm_init(time, g, gv, us, param_file, diag, cs%PPM_CSp)
149  endif
150 
151 end subroutine continuity_init
152 
153 
154 !> continuity_stencil returns the continuity solver stencil size
155 function continuity_stencil(CS) result(stencil)
156  type(continuity_cs), pointer :: cs !< Module's control structure.
157  integer :: stencil !< The continuity solver stencil size with the current settings.
158 
159  stencil = 1
160 
161  if (cs%continuity_scheme == ppm_scheme) then
162  stencil = continuity_ppm_stencil(cs%PPM_CSp)
163  endif
164 
165 end function continuity_stencil
166 
167 !> Destructor for continuity_cs.
168 subroutine continuity_end(CS)
169  type(continuity_cs), pointer :: cs !< Control structure for mom_continuity.
170 
171  if (cs%continuity_scheme == ppm_scheme) then
172  call continuity_ppm_end(cs%PPM_CSp)
173  endif
174 
175  deallocate(cs)
176 
177 end subroutine continuity_end
178 
179 end module mom_continuity
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_continuity::continuity_cs
Control structure for mom_continuity.
Definition: MOM_continuity.F90:26
mom_continuity_ppm::continuity_ppm_init
subroutine, public continuity_ppm_init(Time, G, GV, US, param_file, diag, CS)
Initializes continuity_ppm_cs.
Definition: MOM_continuity_PPM.F90:2181
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_continuity::continuity_stencil
integer function, public continuity_stencil(CS)
continuity_stencil returns the continuity solver stencil size
Definition: MOM_continuity.F90:156
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_string_functions
Handy functions for manipulating strings.
Definition: MOM_string_functions.F90:2
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_continuity::ppm_scheme
integer, parameter ppm_scheme
Enumerated constant to select PPM.
Definition: MOM_continuity.F90:35
mom_continuity::continuity_init
subroutine, public continuity_init(Time, G, GV, US, param_file, diag, CS)
Initializes continuity_cs.
Definition: MOM_continuity.F90:109
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_variables::bt_cont_type
Container for information about the summed layer transports and how they will vary as the barotropic ...
Definition: MOM_variables.F90:260
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_string_functions::uppercase
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.
Definition: MOM_string_functions.F90:42
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_continuity::continuity
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,...
Definition: MOM_continuity.F90:44
mom_continuity_ppm
Solve the layer continuity equation using the PPM method for layer fluxes.
Definition: MOM_continuity_PPM.F90:2
mom_continuity_ppm::continuity_ppm_end
subroutine, public continuity_ppm_end(CS)
Destructor for continuity_ppm_cs.
Definition: MOM_continuity_PPM.F90:2286
mom_continuity::continuity_end
subroutine, public continuity_end(CS)
Destructor for continuity_cs.
Definition: MOM_continuity.F90:169
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_continuity_ppm::continuity_ppm_cs
Control structure for mom_continuity_ppm.
Definition: MOM_continuity_PPM.F90:28
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_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_continuity::ppm_string
character(len=20), parameter ppm_string
String to select PPM.
Definition: MOM_continuity.F90:36
mom_open_boundary::ocean_obc_type
Open-boundary data.
Definition: MOM_open_boundary.F90:195
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_continuity
Solve the layer continuity equation.
Definition: MOM_continuity.F90:2
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
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_continuity_ppm::continuity_ppm_stencil
integer function, public continuity_ppm_stencil(CS)
continuity_PPM_stencil returns the continuity solver stencil size
Definition: MOM_continuity_PPM.F90:2277
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
mom_continuity_ppm::continuity_ppm
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,...
Definition: MOM_continuity_PPM.F90:78