MOM6
mom_verticalgrid Module Reference

Detailed Description

Provides a transparent vertical ocean grid type and supporting routines.

Data Types

type  verticalgrid_type
 Describes the vertical ocean grid, including unit conversion factors. More...
 

Functions/Subroutines

subroutine, public verticalgridinit (param_file, GV, US)
 Allocates and initializes the ocean model vertical grid structure. More...
 
subroutine, public fix_restart_scaling (GV)
 Set the scaling factors for restart files to the scaling factors for this run. More...
 
character(len=48) function, public get_thickness_units (GV)
 Returns the model's thickness units, usually m or kg/m^2. More...
 
character(len=48) function, public get_flux_units (GV)
 Returns the model's thickness flux units, usually m^3/s or kg/s. More...
 
character(len=48) function, public get_tr_flux_units (GV, tr_units, tr_vol_conc_units, tr_mass_conc_units)
 Returns the model's tracer flux units. More...
 
subroutine, public setverticalgridaxes (Rlay, GV, scale)
 This sets the coordinate data for the "layer mode" of the isopycnal model. More...
 
subroutine, public verticalgridend (GV)
 Deallocates the model's vertical grid structure. More...
 

Function/Subroutine Documentation

◆ fix_restart_scaling()

subroutine, public mom_verticalgrid::fix_restart_scaling ( type(verticalgrid_type), intent(inout)  GV)

Set the scaling factors for restart files to the scaling factors for this run.

Parameters
[in,out]gvThe ocean's vertical grid structure

Definition at line 183 of file MOM_verticalGrid.F90.

183  type(verticalGrid_type), intent(inout) :: GV !< The ocean's vertical grid structure
184 
185  gv%m_to_H_restart = gv%m_to_H

◆ get_flux_units()

character(len=48) function, public mom_verticalgrid::get_flux_units ( type(verticalgrid_type), intent(in)  GV)

Returns the model's thickness flux units, usually m^3/s or kg/s.

Returns
The thickness flux units
Parameters
[in]gvThe ocean's vertical grid structure

Definition at line 205 of file MOM_verticalGrid.F90.

205  character(len=48) :: get_flux_units !< The thickness flux units
206  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
207  ! This subroutine returns the appropriate units for thickness fluxes,
208  ! depending on whether the model is Boussinesq or not and the scaling for
209  ! the vertical thickness.
210 
211  if (gv%Boussinesq) then
212  get_flux_units = "m3 s-1"
213  else
214  get_flux_units = "kg s-1"
215  endif

Referenced by mom_dynamics_split_rk2::initialize_dyn_split_rk2(), mom_dynamics_unsplit::initialize_dyn_unsplit(), mom_dynamics_unsplit_rk2::initialize_dyn_unsplit_rk2(), mom_dynamics_split_rk2::register_restarts_dyn_split_rk2(), mom_dynamics_unsplit::register_restarts_dyn_unsplit(), mom_dynamics_unsplit_rk2::register_restarts_dyn_unsplit_rk2(), and mom::set_restart_fields().

Here is the caller graph for this function:

◆ get_thickness_units()

character(len=48) function, public mom_verticalgrid::get_thickness_units ( type(verticalgrid_type), intent(in)  GV)

Returns the model's thickness units, usually m or kg/m^2.

Returns
The vertical thickness units
Parameters
[in]gvThe ocean's vertical grid structure

Definition at line 190 of file MOM_verticalGrid.F90.

190  character(len=48) :: get_thickness_units !< The vertical thickness units
191  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
192  ! This subroutine returns the appropriate units for thicknesses,
193  ! depending on whether the model is Boussinesq or not and the scaling for
194  ! the vertical thickness.
195 
196  if (gv%Boussinesq) then
197  get_thickness_units = "m"
198  else
199  get_thickness_units = "kg m-2"
200  endif

Referenced by mom_ale::ale_register_diags(), mom_diabatic_driver::diabatic_driver_init(), mom_geothermal::geothermal_init(), mom::register_diags(), mom_diagnostics::register_transport_diags(), and mom::set_restart_fields().

Here is the caller graph for this function:

◆ get_tr_flux_units()

character(len=48) function, public mom_verticalgrid::get_tr_flux_units ( type(verticalgrid_type), intent(in)  GV,
character(len=*), intent(in), optional  tr_units,
character(len=*), intent(in), optional  tr_vol_conc_units,
character(len=*), intent(in), optional  tr_mass_conc_units 
)

Returns the model's tracer flux units.

Returns
The model's flux units for a tracer.
Parameters
[in]gvThe ocean's vertical grid structure.
[in]tr_unitsUnits for a tracer, for example Celsius or PSU.
[in]tr_vol_conc_unitsThe concentration units per unit volume, for example if the units are umol m-3, tr_vol_conc_units would be umol.
[in]tr_mass_conc_unitsThe concentration units per unit mass of sea water, for example if the units are mol kg-1, tr_vol_conc_units would be mol.

Definition at line 220 of file MOM_verticalGrid.F90.

220  character(len=48) :: get_tr_flux_units !< The model's flux units
221  !! for a tracer.
222  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical
223  !! grid structure.
224  character(len=*), optional, intent(in) :: tr_units !< Units for a tracer, for example
225  !! Celsius or PSU.
226  character(len=*), optional, intent(in) :: tr_vol_conc_units !< The concentration units per unit
227  !! volume, for example if the units are
228  !! umol m-3, tr_vol_conc_units would
229  !! be umol.
230  character(len=*), optional, intent(in) :: tr_mass_conc_units !< The concentration units per unit
231  !! mass of sea water, for example if
232  !! the units are mol kg-1,
233  !! tr_vol_conc_units would be mol.
234 
235  ! This subroutine returns the appropriate units for thicknesses and fluxes,
236  ! depending on whether the model is Boussinesq or not and the scaling for
237  ! the vertical thickness.
238  integer :: cnt
239 
240  cnt = 0
241  if (present(tr_units)) cnt = cnt+1
242  if (present(tr_vol_conc_units)) cnt = cnt+1
243  if (present(tr_mass_conc_units)) cnt = cnt+1
244 
245  if (cnt == 0) call mom_error(fatal, "get_tr_flux_units: One of the three "//&
246  "arguments tr_units, tr_vol_conc_units, or tr_mass_conc_units "//&
247  "must be present.")
248  if (cnt > 1) call mom_error(fatal, "get_tr_flux_units: Only one of "//&
249  "tr_units, tr_vol_conc_units, and tr_mass_conc_units may be present.")
250  if (present(tr_units)) then
251  if (gv%Boussinesq) then
252  get_tr_flux_units = trim(tr_units)//" m3 s-1"
253  else
254  get_tr_flux_units = trim(tr_units)//" kg s-1"
255  endif
256  endif
257  if (present(tr_vol_conc_units)) then
258  if (gv%Boussinesq) then
259  get_tr_flux_units = trim(tr_vol_conc_units)//" s-1"
260  else
261  get_tr_flux_units = trim(tr_vol_conc_units)//" m-3 kg s-1"
262  endif
263  endif
264  if (present(tr_mass_conc_units)) then
265  if (gv%Boussinesq) then
266  get_tr_flux_units = trim(tr_mass_conc_units)//" kg-1 m3 s-1"
267  else
268  get_tr_flux_units = trim(tr_mass_conc_units)//" s-1"
269  endif
270  endif
271 

References mom_error_handler::mom_error().

Referenced by mom::initialize_mom().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ setverticalgridaxes()

subroutine, public mom_verticalgrid::setverticalgridaxes ( real, dimension(gv%ke), intent(in)  Rlay,
type(verticalgrid_type), intent(inout)  GV,
real, intent(in)  scale 
)

This sets the coordinate data for the "layer mode" of the isopycnal model.

Parameters
[in,out]gvThe container for vertical grid data
[in]rlayThe layer target density [R ~> kg m-3]
[in]scaleA unit scaling factor for Rlay

Definition at line 276 of file MOM_verticalGrid.F90.

276  type(verticalGrid_type), intent(inout) :: GV !< The container for vertical grid data
277  real, dimension(GV%ke), intent(in) :: Rlay !< The layer target density [R ~> kg m-3]
278  real, intent(in) :: scale !< A unit scaling factor for Rlay
279  ! Local variables
280  integer :: k, nk
281 
282  nk = gv%ke
283 
284  gv%zAxisLongName = 'Target Potential Density'
285  gv%zAxisUnits = 'kg m-3'
286  do k=1,nk ; gv%sLayer(k) = scale*rlay(k) ; enddo
287  if (nk > 1) then
288  gv%sInterface(1) = scale * (1.5*rlay(1) - 0.5*rlay(2))
289  do k=2,nk ; gv%sInterface(k) = scale * 0.5*( rlay(k-1) + rlay(k) ) ; enddo
290  gv%sInterface(nk+1) = scale * (1.5*rlay(nk) - 0.5*rlay(nk-1))
291  else
292  gv%sInterface(1) = 0.0 ; gv%sInterface(nk+1) = 2.0*scale*rlay(nk)
293  endif
294 

Referenced by mom_coord_initialization::mom_initialize_coord().

Here is the caller graph for this function:

◆ verticalgridend()

subroutine, public mom_verticalgrid::verticalgridend ( type(verticalgrid_type), pointer  GV)

Deallocates the model's vertical grid structure.

Parameters
gvThe ocean's vertical grid structure

Definition at line 299 of file MOM_verticalGrid.F90.

299  type(verticalGrid_type), pointer :: GV !< The ocean's vertical grid structure
300 
301  deallocate( gv%g_prime, gv%Rlay )
302  deallocate( gv%sInterface , gv%sLayer )
303  deallocate( gv )
304 

◆ verticalgridinit()

subroutine, public mom_verticalgrid::verticalgridinit ( type(param_file_type), intent(in)  param_file,
type(verticalgrid_type), pointer  GV,
type(unit_scale_type), intent(in)  US 
)

Allocates and initializes the ocean model vertical grid structure.

Parameters
[in]param_fileParameter file handle/type
gvThe container for vertical grid data
[in]usA dimensional unit scaling type

Definition at line 76 of file MOM_verticalGrid.F90.

76  type(param_file_type), intent(in) :: param_file !< Parameter file handle/type
77  type(verticalGrid_type), pointer :: GV !< The container for vertical grid data
78  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
79  ! This routine initializes the verticalGrid_type structure (GV).
80  ! All memory is allocated but not necessarily set to meaningful values until later.
81 
82  ! Local variables
83  integer :: nk, H_power
84  real :: H_rescale_factor
85  ! This include declares and sets the variable "version".
86 # include "version_variable.h"
87  character(len=16) :: mdl = 'MOM_verticalGrid'
88 
89  if (associated(gv)) call mom_error(fatal, &
90  'verticalGridInit: called with an associated GV pointer.')
91  allocate(gv)
92 
93  ! Read all relevant parameters and write them to the model log.
94  call log_version(param_file, mdl, version, &
95  "Parameters providing information about the vertical grid.")
96  call get_param(param_file, mdl, "G_EARTH", gv%mks_g_Earth, &
97  "The gravitational acceleration of the Earth.", &
98  units="m s-2", default = 9.80)
99  call get_param(param_file, mdl, "RHO_0", gv%Rho0, &
100  "The mean ocean density used with BOUSSINESQ true to "//&
101  "calculate accelerations and the mass for conservation "//&
102  "properties, or with BOUSSINSEQ false to convert some "//&
103  "parameters from vertical units of m to kg m-2.", &
104  units="kg m-3", default=1035.0, scale=us%kg_m3_to_R)
105  call get_param(param_file, mdl, "BOUSSINESQ", gv%Boussinesq, &
106  "If true, make the Boussinesq approximation.", default=.true.)
107  call get_param(param_file, mdl, "ANGSTROM", gv%Angstrom_m, &
108  "The minimum layer thickness, usually one-Angstrom.", &
109  units="m", default=1.0e-10)
110  call get_param(param_file, mdl, "H_RESCALE_POWER", h_power, &
111  "An integer power of 2 that is used to rescale the model's "//&
112  "intenal units of thickness. Valid values range from -300 to 300.", &
113  units="nondim", default=0, debuggingparam=.true.)
114  if (abs(h_power) > 300) call mom_error(fatal, "verticalGridInit: "//&
115  "H_RESCALE_POWER is outside of the valid range of -300 to 300.")
116  h_rescale_factor = 1.0
117  if (h_power /= 0) h_rescale_factor = 2.0**h_power
118  if (.not.gv%Boussinesq) then
119  call get_param(param_file, mdl, "H_TO_KG_M2", gv%H_to_kg_m2,&
120  "A constant that translates thicknesses from the model's "//&
121  "internal units of thickness to kg m-2.", units="kg m-2 H-1", &
122  default=1.0)
123  gv%H_to_kg_m2 = gv%H_to_kg_m2 * h_rescale_factor
124  else
125  call get_param(param_file, mdl, "H_TO_M", gv%H_to_m, &
126  "A constant that translates the model's internal "//&
127  "units of thickness into m.", units="m H-1", default=1.0)
128  gv%H_to_m = gv%H_to_m * h_rescale_factor
129  endif
130  gv%g_Earth = us%m_to_L**2*us%Z_to_m*us%T_to_s**2 * gv%mks_g_Earth
131 #ifdef STATIC_MEMORY_
132  ! Here NK_ is a macro, while nk is a variable.
133  call get_param(param_file, mdl, "NK", nk, &
134  "The number of model layers.", units="nondim", &
135  static_value=nk_)
136  if (nk /= nk_) call mom_error(fatal, "verticalGridInit: " // &
137  "Mismatched number of layers NK_ between MOM_memory.h and param_file")
138 
139 #else
140  call get_param(param_file, mdl, "NK", nk, &
141  "The number of model layers.", units="nondim", fail_if_missing=.true.)
142 #endif
143  gv%ke = nk
144 
145  if (gv%Boussinesq) then
146  gv%H_to_kg_m2 = us%R_to_kg_m3*gv%Rho0 * gv%H_to_m
147  gv%kg_m2_to_H = 1.0 / gv%H_to_kg_m2
148  gv%m_to_H = 1.0 / gv%H_to_m
149  gv%Angstrom_H = gv%m_to_H * gv%Angstrom_m
150  gv%H_to_MKS = gv%H_to_m
151  else
152  gv%kg_m2_to_H = 1.0 / gv%H_to_kg_m2
153  gv%m_to_H = us%R_to_kg_m3*gv%Rho0 * gv%kg_m2_to_H
154  gv%H_to_m = gv%H_to_kg_m2 / (us%R_to_kg_m3*gv%Rho0)
155  gv%Angstrom_H = gv%Angstrom_m*1000.0*gv%kg_m2_to_H
156  gv%H_to_MKS = gv%H_to_kg_m2
157  endif
158  gv%H_subroundoff = 1e-20 * max(gv%Angstrom_H,gv%m_to_H*1e-17)
159  gv%H_to_Pa = gv%mks_g_Earth * gv%H_to_kg_m2
160 
161  gv%H_to_Z = gv%H_to_m * us%m_to_Z
162  gv%Z_to_H = us%Z_to_m * gv%m_to_H
163  gv%Angstrom_Z = us%m_to_Z * gv%Angstrom_m
164 
165  gv%H_to_RZ = gv%H_to_kg_m2 * us%kg_m3_to_R * us%m_to_Z
166  gv%RZ_to_H = gv%kg_m2_to_H * us%R_to_kg_m3 * us%Z_to_m
167 
168 ! Log derivative values.
169  call log_param(param_file, mdl, "M to THICKNESS", gv%m_to_H*h_rescale_factor)
170  call log_param(param_file, mdl, "M to THICKNESS rescaled by 2^-n", gv%m_to_H)
171  call log_param(param_file, mdl, "THICKNESS to M rescaled by 2^n", gv%H_to_m)
172 
173  allocate( gv%sInterface(nk+1) )
174  allocate( gv%sLayer(nk) )
175  allocate( gv%g_prime(nk+1) ) ; gv%g_prime(:) = 0.0
176  ! The extent of Rlay should be changed to nk?
177  allocate( gv%Rlay(nk+1) ) ; gv%Rlay(:) = 0.0
178 

References mom_error_handler::mom_error().

Referenced by mom_oda_driver_mod::init_oda().

Here is the call graph for this function:
Here is the caller graph for this function: