MOM6
mom_coord_initialization Module Reference

Detailed Description

Initializes fixed aspects of the related to its vertical coordinate.

Functions/Subroutines

subroutine, public mom_initialize_coord (GV, US, PF, write_geom, output_dir, tv, max_depth)
 MOM_initialize_coord sets up time-invariant quantities related to MOM6's vertical coordinate. More...
 
subroutine set_coord_from_gprime (Rlay, g_prime, GV, US, param_file)
 Sets the layer densities (Rlay) and the interface reduced gravities (g). More...
 
subroutine set_coord_from_layer_density (Rlay, g_prime, GV, US, param_file)
 Sets the layer densities (Rlay) and the interface reduced gravities (g). More...
 
subroutine set_coord_from_ts_ref (Rlay, g_prime, GV, US, param_file, eqn_of_state, P_Ref)
 Sets the layer densities (Rlay) and the interface reduced gravities (g) from a profile of g'. More...
 
subroutine set_coord_from_ts_profile (Rlay, g_prime, GV, US, param_file, eqn_of_state, P_Ref)
 Sets the layer densities (Rlay) and the interface reduced gravities (g) from a T-S profile. More...
 
subroutine set_coord_from_ts_range (Rlay, g_prime, GV, US, param_file, eqn_of_state, P_Ref)
 Sets the layer densities (Rlay) and the interface reduced gravities (g) from a linear T-S profile. More...
 
subroutine set_coord_from_file (Rlay, g_prime, GV, US, param_file)
 
subroutine set_coord_linear (Rlay, g_prime, GV, US, param_file)
 Sets the layer densities (Rlay) and the interface reduced gravities (g) according to a linear profile starting at a reference surface layer density and spanning a range of densities to the bottom defined by the parameter RLAY_RANGE (defaulting to 2.0 if not defined) More...
 
subroutine set_coord_to_none (Rlay, g_prime, GV, US, param_file)
 Sets Rlay to Rho0 and g_prime to zero except for the free surface. This is for use only in ALE mode where Rlay should not be used and g_prime(1) alone might be used. More...
 
subroutine write_vertgrid_file (GV, US, param_file, directory)
 Writes out a file containing any available data related to the vertical grid used by the MOM ocean model. More...
 

Variables

character(len=40) mdl = "MOM_coord_initialization"
 This module's name. More...
 

Function/Subroutine Documentation

◆ mom_initialize_coord()

subroutine, public mom_coord_initialization::mom_initialize_coord ( type(verticalgrid_type), intent(inout)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  PF,
logical, intent(in)  write_geom,
character(len=*), intent(in)  output_dir,
type(thermo_var_ptrs), intent(inout)  tv,
real, intent(in)  max_depth 
)

MOM_initialize_coord sets up time-invariant quantities related to MOM6's vertical coordinate.

Parameters
[in,out]gvOcean vertical grid structure.
[in]usA dimensional unit scaling type
[in]pfA structure indicating the open file to parse for model parameter values.
[in]write_geomIf true, write grid geometry files.
[in]output_dirThe directory into which to write files.
[in,out]tvThe thermodynamic variable structure.
[in]max_depthThe ocean's maximum depth [Z ~> m].

Definition at line 40 of file MOM_coord_initialization.F90.

40  type(verticalGrid_type), intent(inout) :: GV !< Ocean vertical grid structure.
41  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
42  type(param_file_type), intent(in) :: PF !< A structure indicating the open file
43  !! to parse for model parameter values.
44  logical, intent(in) :: write_geom !< If true, write grid geometry files.
45  character(len=*), intent(in) :: output_dir !< The directory into which to write files.
46  type(thermo_var_ptrs), intent(inout) :: tv !< The thermodynamic variable structure.
47  real, intent(in) :: max_depth !< The ocean's maximum depth [Z ~> m].
48  ! Local
49  character(len=200) :: config
50  logical :: debug
51 ! This include declares and sets the variable "version".
52 #include "version_variable.h"
53  integer :: nz
54  type(EOS_type), pointer :: eos => null()
55 
56  if (associated(tv%eqn_of_state)) eos => tv%eqn_of_state
57 
58  nz = gv%ke
59 
60  call calltree_enter("MOM_initialize_coord(), MOM_coord_initialization.F90")
61  call log_version(pf, mdl, version, "")
62  call get_param(pf, mdl, "DEBUG", debug, default=.false.)
63 
64 ! Set-up the layer densities, GV%Rlay, and reduced gravities, GV%g_prime.
65  call get_param(pf, mdl, "COORD_CONFIG", config, &
66  "This specifies how layers are to be defined: \n"//&
67  " \t ALE or none - used to avoid defining layers in ALE mode \n"//&
68  " \t file - read coordinate information from the file \n"//&
69  " \t\t specified by (COORD_FILE).\n"//&
70  " \t BFB - Custom coords for buoyancy-forced basin case \n"//&
71  " \t\t based on SST_S, T_BOT and DRHO_DT.\n"//&
72  " \t linear - linear based on interfaces not layers \n"//&
73  " \t layer_ref - linear based on layer densities \n"//&
74  " \t ts_ref - use reference temperature and salinity \n"//&
75  " \t ts_range - use range of temperature and salinity \n"//&
76  " \t\t (T_REF and S_REF) to determine surface density \n"//&
77  " \t\t and GINT calculate internal densities. \n"//&
78  " \t gprime - use reference density (RHO_0) for surface \n"//&
79  " \t\t density and GINT calculate internal densities. \n"//&
80  " \t ts_profile - use temperature and salinity profiles \n"//&
81  " \t\t (read from COORD_FILE) to set layer densities. \n"//&
82  " \t USER - call a user modified routine.", &
83  fail_if_missing=.true.)
84  select case ( trim(config) )
85  case ("gprime")
86  call set_coord_from_gprime(gv%Rlay, gv%g_prime, gv, us, pf)
87  case ("layer_ref")
88  call set_coord_from_layer_density(gv%Rlay, gv%g_prime, gv, us, pf)
89  case ("linear")
90  call set_coord_linear(gv%Rlay, gv%g_prime, gv, us, pf)
91  case ("ts_ref")
92  call set_coord_from_ts_ref(gv%Rlay, gv%g_prime, gv, us, pf, eos, tv%P_Ref)
93  case ("ts_profile")
94  call set_coord_from_ts_profile(gv%Rlay, gv%g_prime, gv, us, pf, eos, tv%P_Ref)
95  case ("ts_range")
96  call set_coord_from_ts_range(gv%Rlay, gv%g_prime, gv, us, pf, eos, tv%P_Ref)
97  case ("file")
98  call set_coord_from_file(gv%Rlay, gv%g_prime, gv, us, pf)
99  case ("USER")
100  call user_set_coord(gv%Rlay, gv%g_prime, gv, us, pf, eos)
101  case ("BFB")
102  call bfb_set_coord(gv%Rlay, gv%g_prime, gv, us, pf, eos)
103  case ("none", "ALE")
104  call set_coord_to_none(gv%Rlay, gv%g_prime, gv, us, pf)
105  case default ; call mom_error(fatal,"MOM_initialize_coord: "// &
106  "Unrecognized coordinate setup"//trim(config))
107  end select
108  if (debug) call chksum(us%R_to_kg_m3*gv%Rlay(:), "MOM_initialize_coord: Rlay ", 1, nz)
109  if (debug) call chksum(us%m_to_Z*us%L_to_m**2*us%s_to_T**2*gv%g_prime(:), "MOM_initialize_coord: g_prime ", 1, nz)
110  call setverticalgridaxes( gv%Rlay, gv, scale=us%R_to_kg_m3 )
111 
112 ! Copy the maximum depth across from the input argument
113  gv%max_depth = max_depth
114 
115 ! Write out all of the grid data used by this run.
116  if (write_geom) call write_vertgrid_file(gv, us, pf, output_dir)
117 
118  call calltree_leave('MOM_initialize_coord()')
119 

References bfb_initialization::bfb_set_coord(), mom_error_handler::calltree_enter(), mom_error_handler::calltree_leave(), mdl, set_coord_from_file(), set_coord_from_gprime(), set_coord_from_layer_density(), set_coord_from_ts_profile(), set_coord_from_ts_range(), set_coord_from_ts_ref(), set_coord_linear(), set_coord_to_none(), mom_verticalgrid::setverticalgridaxes(), user_initialization::user_set_coord(), and write_vertgrid_file().

Referenced by mom_oda_driver_mod::init_oda(), and mom::initialize_mom().

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

◆ set_coord_from_file()

subroutine mom_coord_initialization::set_coord_from_file ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file 
)
private
Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters

Definition at line 385 of file MOM_coord_initialization.F90.

385  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
386  !! (potential density) [R ~> kg m-3].
387  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
388  !! [L2 Z-1 T-2 ~> m s-2].
389  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
390  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
391  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
392  ! Local variables
393  real :: g_fs ! Reduced gravity across the free surface [m s-2].
394  integer :: k, nz
395  character(len=40) :: mdl = "set_coord_from_file" ! This subroutine's name.
396  character(len=40) :: coord_var
397  character(len=200) :: filename,coord_file,inputdir ! Strings for file/path
398  nz = gv%ke
399 
400  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
401 
402  call get_param(param_file, mdl, "GFS", g_fs, &
403  "The reduced gravity at the free surface.", units="m s-2", &
404  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
405  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
406  inputdir = slasher(inputdir)
407  call get_param(param_file, mdl, "COORD_FILE", coord_file, &
408  "The file from which the coordinate densities are read.", &
409  fail_if_missing=.true.)
410  call get_param(param_file, mdl, "COORD_VAR", coord_var, &
411  "The variable in COORD_FILE that is to be used for the "//&
412  "coordinate densities.", default="Layer")
413  filename = trim(inputdir)//trim(coord_file)
414  call log_param(param_file, mdl, "INPUTDIR/COORD_FILE", filename)
415  if (.not.file_exists(filename)) call mom_error(fatal, &
416  " set_coord_from_file: Unable to open "//trim(filename))
417 
418  call read_axis_data(filename, coord_var, rlay)
419  do k=1,nz ; rlay(k) = us%kg_m3_to_R*rlay(k) ; enddo
420  g_prime(1) = g_fs
421  do k=2,nz ; g_prime(k) = (gv%g_Earth/(gv%Rho0)) * (rlay(k) - rlay(k-1)) ; enddo
422  do k=1,nz ; if (g_prime(k) <= 0.0) then
423  call mom_error(fatal, "MOM_initialization set_coord_from_file: "//&
424  "Zero or negative g_primes read from variable "//"Layer"//" in file "//&
425  trim(filename))
426  endif ; enddo
427 
428  call calltree_leave(trim(mdl)//'()')

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_from_gprime()

subroutine mom_coord_initialization::set_coord_from_gprime ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file 
)
private

Sets the layer densities (Rlay) and the interface reduced gravities (g).

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters

Definition at line 126 of file MOM_coord_initialization.F90.

126  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
127  !! (potential density) [R ~> kg m-3].
128  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
129  !! [L2 Z-1 T-2 ~> m s-2].
130  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
131  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
132  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
133  ! Local variables
134  real :: g_int ! Reduced gravities across the internal interfaces [m s-2].
135  real :: g_fs ! Reduced gravity across the free surface [m s-2].
136  character(len=40) :: mdl = "set_coord_from_gprime" ! This subroutine's name.
137  integer :: k, nz
138  nz = gv%ke
139 
140  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
141 
142  call get_param(param_file, mdl, "GFS" , g_fs, &
143  "The reduced gravity at the free surface.", units="m s-2", &
144  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
145  call get_param(param_file, mdl, "GINT", g_int, &
146  "The reduced gravity across internal interfaces.", &
147  units="m s-2", fail_if_missing=.true., scale=us%m_s_to_L_T**2*us%Z_to_m)
148 
149  g_prime(1) = g_fs
150  do k=2,nz ; g_prime(k) = g_int ; enddo
151  rlay(1) = gv%Rho0
152  do k=2,nz ; rlay(k) = rlay(k-1) + g_prime(k)*(gv%Rho0/gv%g_Earth) ; enddo
153 
154  call calltree_leave(trim(mdl)//'()')
155 

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_from_layer_density()

subroutine mom_coord_initialization::set_coord_from_layer_density ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file 
)
private

Sets the layer densities (Rlay) and the interface reduced gravities (g).

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters

Definition at line 160 of file MOM_coord_initialization.F90.

160  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
161  !! (potential density) [R ~> kg m-3].
162  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
163  !! [L2 Z-1 T-2 ~> m s-2].
164  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
165  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
166  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
167  ! Local variables
168  real :: g_fs ! Reduced gravity across the free surface [m s-2].
169  real :: Rlay_Ref! The surface layer's target density [kg m-3].
170  real :: RLay_range ! The range of densities [kg m-3].
171  character(len=40) :: mdl = "set_coord_from_layer_density" ! This subroutine's name.
172  integer :: k, nz
173  nz = gv%ke
174 
175  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
176 
177  call get_param(param_file, mdl, "GFS", g_fs, &
178  "The reduced gravity at the free surface.", units="m s-2", &
179  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
180  call get_param(param_file, mdl, "LIGHTEST_DENSITY", rlay_ref, &
181  "The reference potential density used for layer 1.", &
182  units="kg m-3", default=us%R_to_kg_m3*gv%Rho0, scale=us%kg_m3_to_R)
183  call get_param(param_file, mdl, "DENSITY_RANGE", rlay_range, &
184  "The range of reference potential densities in the layers.", &
185  units="kg m-3", default=2.0, scale=us%kg_m3_to_R)
186 
187  g_prime(1) = g_fs
188  rlay(1) = rlay_ref
189  do k=2,nz
190  rlay(k) = rlay(k-1) + rlay_range/(real(nz-1))
191  enddo
192 ! These statements set the interface reduced gravities. !
193  do k=2,nz
194  g_prime(k) = (gv%g_Earth/(gv%Rho0)) * (rlay(k) - rlay(k-1))
195  enddo
196 
197  call calltree_leave(trim(mdl)//'()')

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_from_ts_profile()

subroutine mom_coord_initialization::set_coord_from_ts_profile ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(eos_type), pointer  eqn_of_state,
real, intent(in)  P_Ref 
)
private

Sets the layer densities (Rlay) and the interface reduced gravities (g) from a T-S profile.

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters
eqn_of_stateinteger that selects equation of state.
[in]p_refThe coordinate-density reference pressure [Pa].

Definition at line 254 of file MOM_coord_initialization.F90.

254  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
255  !! (potential density) [R ~> kg m-3].
256  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
257  !! [L2 Z-1 T-2 ~> m s-2].
258  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
259  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
260  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time
261  !! parameters
262  type(EOS_type), pointer :: eqn_of_state !< integer that selects equation of state.
263  real, intent(in) :: P_Ref !< The coordinate-density reference pressure [Pa].
264  ! Local variables
265  real, dimension(GV%ke) :: T0, S0, Pref
266  real :: g_fs ! Reduced gravity across the free surface [m s-2].
267  integer :: k, nz
268  character(len=40) :: mdl = "set_coord_from_TS_profile" ! This subroutine's name.
269  character(len=200) :: filename, coord_file, inputdir ! Strings for file/path
270  nz = gv%ke
271 
272  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
273 
274  call get_param(param_file, mdl, "GFS", g_fs, &
275  "The reduced gravity at the free surface.", units="m s-2", &
276  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
277  call get_param(param_file, mdl, "COORD_FILE", coord_file, &
278  "The file from which the coordinate temperatures and "//&
279  "salinities are read.", fail_if_missing=.true.)
280 
281  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
282  filename = trim(slasher(inputdir))//trim(coord_file)
283  call log_param(param_file, mdl, "INPUTDIR/COORD_FILE", filename)
284 
285  call mom_read_data(filename,"PTEMP",t0(:))
286  call mom_read_data(filename,"SALT",s0(:))
287 
288  if (.not.file_exists(filename)) call mom_error(fatal, &
289  " set_coord_from_TS_profile: Unable to open " //trim(filename))
290 ! These statements set the interface reduced gravities. !
291  g_prime(1) = g_fs
292  do k=1,nz ; pref(k) = p_ref ; enddo
293  call calculate_density(t0, s0, pref, rlay, 1, nz, eqn_of_state, scale=us%kg_m3_to_R)
294  do k=2,nz; g_prime(k) = (gv%g_Earth/(gv%Rho0)) * (rlay(k) - rlay(k-1)) ; enddo
295 
296  call calltree_leave(trim(mdl)//'()')

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_from_ts_range()

subroutine mom_coord_initialization::set_coord_from_ts_range ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(eos_type), pointer  eqn_of_state,
real, intent(in)  P_Ref 
)
private

Sets the layer densities (Rlay) and the interface reduced gravities (g) from a linear T-S profile.

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters
eqn_of_stateinteger that selects equation of state
[in]p_refThe coordinate-density reference pressure [Pa]

Definition at line 302 of file MOM_coord_initialization.F90.

302  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
303  !! (potential density) [R ~> kg m-3].
304  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
305  !! [L2 Z-1 T-2 ~> m s-2].
306  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
307  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
308  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time
309  !! parameters
310  type(EOS_type), pointer :: eqn_of_state !< integer that selects equation of state
311  real, intent(in) :: P_Ref !< The coordinate-density reference pressure [Pa]
312 
313  ! Local variables
314  real, dimension(GV%ke) :: T0, S0, Pref
315  real :: S_Ref, S_Light, S_Dense ! Salinity range parameters [ppt].
316  real :: T_Ref, T_Light, T_Dense ! Temperature range parameters [decC].
317  real :: res_rat ! The ratio of density space resolution in the denser part
318  ! of the range to that in the lighter part of the range.
319  ! Setting this greater than 1 increases the resolution for
320  ! the denser water.
321  real :: g_fs ! Reduced gravity across the free surface [m s-2].
322  real :: a1, frac_dense, k_frac
323  integer :: k, nz, k_light
324  character(len=40) :: mdl = "set_coord_from_TS_range" ! This subroutine's name.
325  character(len=200) :: filename, coord_file, inputdir ! Strings for file/path
326  nz = gv%ke
327 
328  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
329 
330  call get_param(param_file, mdl, "T_REF", t_ref, &
331  "The default initial temperatures.", units="degC", default=10.0)
332  call get_param(param_file, mdl, "TS_RANGE_T_LIGHT", t_light, &
333  "The initial temperature of the lightest layer when "//&
334  "COORD_CONFIG is set to ts_range.", units="degC", default=t_ref)
335  call get_param(param_file, mdl, "TS_RANGE_T_DENSE", t_dense, &
336  "The initial temperature of the densest layer when "//&
337  "COORD_CONFIG is set to ts_range.", units="degC", default=t_ref)
338 
339  call get_param(param_file, mdl, "S_REF", s_ref, &
340  "The default initial salinities.", units="PSU", default=35.0)
341  call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", s_light, &
342  "The initial lightest salinities when COORD_CONFIG "//&
343  "is set to ts_range.", default = s_ref, units="PSU")
344  call get_param(param_file, mdl, "TS_RANGE_S_DENSE", s_dense, &
345  "The initial densest salinities when COORD_CONFIG "//&
346  "is set to ts_range.", default = s_ref, units="PSU")
347 
348  call get_param(param_file, mdl, "TS_RANGE_RESOLN_RATIO", res_rat, &
349  "The ratio of density space resolution in the densest "//&
350  "part of the range to that in the lightest part of the "//&
351  "range when COORD_CONFIG is set to ts_range. Values "//&
352  "greater than 1 increase the resolution of the denser water.",&
353  default=1.0, units="nondim")
354 
355  call get_param(param_file, mdl, "GFS", g_fs, &
356  "The reduced gravity at the free surface.", units="m s-2", &
357  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
358 
359  k_light = gv%nk_rho_varies + 1
360 
361  ! Set T0(k) to range from T_LIGHT to T_DENSE, and simliarly for S0(k).
362  t0(k_light) = t_light ; s0(k_light) = s_light
363  a1 = 2.0 * res_rat / (1.0 + res_rat)
364  do k=k_light+1,nz
365  k_frac = real(k-k_light)/real(nz-k_light)
366  frac_dense = a1 * k_frac + (1.0 - a1) * k_frac**2
367  t0(k) = frac_dense * (t_dense - t_light) + t_light
368  s0(k) = frac_dense * (s_dense - s_light) + s_light
369  enddo
370 
371  g_prime(1) = g_fs
372  do k=1,nz ; pref(k) = p_ref ; enddo
373  call calculate_density(t0, s0, pref, rlay, k_light, nz-k_light+1, eqn_of_state, scale=us%kg_m3_to_R)
374  ! Extrapolate target densities for the variable density mixed and buffer layers.
375  do k=k_light-1,1,-1
376  rlay(k) = 2.0*rlay(k+1) - rlay(k+2)
377  enddo
378  do k=2,nz ; g_prime(k) = (gv%g_Earth/(gv%Rho0)) * (rlay(k) - rlay(k-1)) ; enddo
379 
380  call calltree_leave(trim(mdl)//'()')

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_from_ts_ref()

subroutine mom_coord_initialization::set_coord_from_ts_ref ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(eos_type), pointer  eqn_of_state,
real, intent(in)  P_Ref 
)
private

Sets the layer densities (Rlay) and the interface reduced gravities (g) from a profile of g'.

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters
eqn_of_stateinteger selecting the equation of state.
[in]p_refThe coordinate-density reference pressure [Pa].

Definition at line 203 of file MOM_coord_initialization.F90.

203  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
204  !! (potential density) [R ~> kg m-3].
205  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
206  !! [L2 Z-1 T-2 ~> m s-2].
207  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
208  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
209  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time
210  !! parameters
211  type(EOS_type), pointer :: eqn_of_state !< integer selecting the equation of state.
212  real, intent(in) :: P_Ref !< The coordinate-density reference pressure [Pa].
213  ! Local variables
214  real :: T_ref ! Reference temperature
215  real :: S_ref ! Reference salinity
216  real :: g_int ! Reduced gravities across the internal interfaces [m s-2].
217  real :: g_fs ! Reduced gravity across the free surface [m s-2].
218  character(len=40) :: mdl = "set_coord_from_TS_ref" ! This subroutine's name.
219  integer :: k, nz
220  nz = gv%ke
221 
222  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
223 
224  call get_param(param_file, mdl, "T_REF", t_ref, &
225  "The initial temperature of the lightest layer.", units="degC", &
226  fail_if_missing=.true.)
227  call get_param(param_file, mdl, "S_REF", s_ref, &
228  "The initial salinities.", units="PSU", default=35.0)
229  call get_param(param_file, mdl, "GFS", g_fs, &
230  "The reduced gravity at the free surface.", units="m s-2", &
231  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
232  call get_param(param_file, mdl, "GINT", g_int, &
233  "The reduced gravity across internal interfaces.", &
234  units="m s-2", fail_if_missing=.true., scale=us%m_s_to_L_T**2*us%Z_to_m)
235  !
236 ! These statements set the interface reduced gravities. !
237  g_prime(1) = g_fs
238  do k=2,nz ; g_prime(k) = g_int ; enddo
239 
240 ! The uppermost layer's density is set here. Subsequent layers' !
241 ! densities are determined from this value and the g values. !
242 ! T0 = 28.228 ; S0 = 34.5848 ; Pref = P_Ref
243  call calculate_density(t_ref, s_ref, p_ref, rlay(1), eqn_of_state, scale=us%kg_m3_to_R)
244 
245 ! These statements set the layer densities. !
246  do k=2,nz ; rlay(k) = rlay(k-1) + g_prime(k)*(gv%Rho0/gv%g_Earth) ; enddo
247 
248  call calltree_leave(trim(mdl)//'()')

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_linear()

subroutine mom_coord_initialization::set_coord_linear ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file 
)
private

Sets the layer densities (Rlay) and the interface reduced gravities (g) according to a linear profile starting at a reference surface layer density and spanning a range of densities to the bottom defined by the parameter RLAY_RANGE (defaulting to 2.0 if not defined)

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters

Definition at line 437 of file MOM_coord_initialization.F90.

437  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
438  !! (potential density) [R ~> kg m-3].
439  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces
440  !! [L2 Z-1 T-2 ~> m s-2].
441  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
442  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
443  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
444  ! Local variables
445  character(len=40) :: mdl = "set_coord_linear" ! This subroutine
446  real :: Rlay_ref, Rlay_range, g_fs
447  integer :: k, nz
448  nz = gv%ke
449 
450  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
451 
452  call get_param(param_file, mdl, "LIGHTEST_DENSITY", rlay_ref, &
453  "The reference potential density used for the surface interface.", &
454  units="kg m-3", default=us%R_to_kg_m3*gv%Rho0, scale=us%kg_m3_to_R)
455  call get_param(param_file, mdl, "DENSITY_RANGE", rlay_range, &
456  "The range of reference potential densities across all interfaces.", &
457  units="kg m-3", default=2.0, scale=us%kg_m3_to_R)
458  call get_param(param_file, mdl, "GFS", g_fs, &
459  "The reduced gravity at the free surface.", units="m s-2", &
460  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
461 
462  ! This following sets the target layer densities such that a the
463  ! surface interface has density Rlay_ref and the bottom
464  ! is Rlay_range larger
465  do k=1,nz
466  rlay(k) = rlay_ref + rlay_range*((real(k)-0.5)/real(nz))
467  enddo
468  ! These statements set the interface reduced gravities.
469  g_prime(1) = g_fs
470  do k=2,nz
471  g_prime(k) = (gv%g_Earth/(gv%Rho0)) * (rlay(k) - rlay(k-1))
472  enddo
473 
474  call calltree_leave(trim(mdl)//'()')

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ set_coord_to_none()

subroutine mom_coord_initialization::set_coord_to_none ( real, dimension(:), intent(out)  Rlay,
real, dimension(:), intent(out)  g_prime,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file 
)
private

Sets Rlay to Rho0 and g_prime to zero except for the free surface. This is for use only in ALE mode where Rlay should not be used and g_prime(1) alone might be used.

Parameters
[out]rlayThe layers' target coordinate values (potential density) [R ~> kg m-3].
[out]g_primeThe reduced gravity across the interfaces, [L2 Z-1 T-2 ~> m s-2].
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters

Definition at line 481 of file MOM_coord_initialization.F90.

481  real, dimension(:), intent(out) :: Rlay !< The layers' target coordinate values
482  !! (potential density) [R ~> kg m-3].
483  real, dimension(:), intent(out) :: g_prime !< The reduced gravity across the interfaces,
484  !! [L2 Z-1 T-2 ~> m s-2].
485  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
486  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
487  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
488  ! Local variables
489  real :: g_fs ! Reduced gravity across the free surface [m s-2].
490  character(len=40) :: mdl = "set_coord_to_none" ! This subroutine's name.
491  integer :: k, nz
492  nz = gv%ke
493 
494  call calltree_enter(trim(mdl)//"(), MOM_coord_initialization.F90")
495 
496  call get_param(param_file, mdl, "GFS" , g_fs, &
497  "The reduced gravity at the free surface.", units="m s-2", &
498  default=gv%mks_g_Earth, scale=us%m_s_to_L_T**2*us%Z_to_m)
499 
500  g_prime(1) = g_fs
501  do k=2,nz ; g_prime(k) = 0. ; enddo
502  rlay(1) = gv%Rho0
503  do k=2,nz ; rlay(k) = rlay(k-1) + g_prime(k)*(gv%Rho0/gv%g_Earth) ; enddo
504 
505  call calltree_leave(trim(mdl)//'()')
506 

References mom_error_handler::calltree_enter(), and mom_error_handler::calltree_leave().

Referenced by mom_initialize_coord().

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

◆ write_vertgrid_file()

subroutine mom_coord_initialization::write_vertgrid_file ( type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
character(len=*), intent(in)  directory 
)
private

Writes out a file containing any available data related to the vertical grid used by the MOM ocean model.

Parameters
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters
[in]directoryThe directory into which to place the file.

Definition at line 512 of file MOM_coord_initialization.F90.

512  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
513  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
514  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
515  character(len=*), intent(in) :: directory !< The directory into which to place the file.
516  ! Local variables
517  character(len=240) :: filepath
518  type(vardesc) :: vars(2)
519  type(fieldtype) :: fields(2)
520  integer :: unit
521 
522  filepath = trim(directory) // trim("Vertical_coordinate")
523 
524  vars(1) = var_desc("R","kilogram meter-3","Target Potential Density",'1','L','1')
525  vars(2) = var_desc("g","meter second-2","Reduced gravity",'1','L','1')
526 
527  call create_file(unit, trim(filepath), vars, 2, fields, single_file, gv=gv)
528 
529  call write_field(unit, fields(1), us%R_to_kg_m3*gv%Rlay(:))
530  call write_field(unit, fields(2), us%L_T_to_m_s**2*us%m_to_Z*gv%g_prime(:))
531 
532  call close_file(unit)
533 

References mom_io::var_desc().

Referenced by mom_initialize_coord().

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

Variable Documentation

◆ mdl

character(len=40) mom_coord_initialization::mdl = "MOM_coord_initialization"
private

This module's name.

Definition at line 33 of file MOM_coord_initialization.F90.

33 character(len=40) :: mdl = "MOM_coord_initialization" !< This module's name.

Referenced by mom_initialize_coord().