MOM6
dumbbell_initialization.F90
Go to the documentation of this file.
1 !> Configures the model for the idealized dumbbell test case.
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_domains, only : sum_across_pes
10 use mom_get_input, only : directories
11 use mom_grid, only : ocean_grid_type
22 
23 implicit none ; private
24 
25 #include <MOM_memory.h>
26 
27 character(len=40) :: mdl = "dumbbell_initialization" !< This module's name.
28 
33 
34 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
35 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
36 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
37 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
38 
39 contains
40 
41 !> Initialization of topography.
42 subroutine dumbbell_initialize_topography( D, G, param_file, max_depth )
43  type(dyn_horgrid_type), intent(in) :: g !< The dynamic horizontal grid type
44  real, dimension(G%isd:G%ied,G%jsd:G%jed), &
45  intent(out) :: d !< Ocean bottom depth in the units of depth_max
46  type(param_file_type), intent(in) :: param_file !< Parameter file structure
47  real, intent(in) :: max_depth !< Maximum ocean depth in arbitrary units
48 
49  ! Local variables
50  integer :: i, j
51  real :: x, y, delta, dblen, dbfrac
52 
53  call get_param(param_file, mdl,"DUMBBELL_LEN",dblen, &
54  'Lateral Length scale for dumbbell.',&
55  units='k', default=600., do_not_log=.false.)
56  call get_param(param_file, mdl,"DUMBBELL_FRACTION",dbfrac, &
57  'Meridional fraction for narrow part of dumbbell.',&
58  units='nondim', default=0.5, do_not_log=.false.)
59 
60  if (g%x_axis_units == 'm') then
61  dblen=dblen*1.e3
62  endif
63 
64  do j=g%jsc,g%jec ; do i=g%isc,g%iec
65  ! Compute normalized zonal coordinates (x,y=0 at center of domain)
66  x = ( g%geoLonT(i,j) ) / dblen
67  y = ( g%geoLatT(i,j) ) / g%len_lat
68  d(i,j) = g%max_depth
69  if ((x>=-0.25 .and. x<=0.25) .and. (y <= -0.5*dbfrac .or. y >= 0.5*dbfrac)) then
70  d(i,j) = 0.0
71  endif
72  enddo ; enddo
73 
74 end subroutine dumbbell_initialize_topography
75 
76 !> Initializes the layer thicknesses to be uniform in the dumbbell test case
77 subroutine dumbbell_initialize_thickness ( h, G, GV, US, param_file, just_read_params)
78  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
79  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
80  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
81  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
82  intent(out) :: h !< The thickness that is being initialized [H ~> m or kg m-2].
83  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
84  !! to parse for model parameter values.
85  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
86  !! only read parameters without changing h.
87 
88  real :: e0(szk_(g)+1) ! The resting interface heights [Z ~> m], usually
89  ! negative because it is positive upward.
90  real :: eta1d(szk_(g)+1)! Interface height relative to the sea surface
91  ! positive upward [Z ~> m].
92  real :: min_thickness ! The minimum layer thicknesses [Z ~> m].
93  real :: s_surf, s_range, s_ref, s_light, s_dense ! Various salinities [ppt].
94  real :: eta_ic_quanta ! The granularity of quantization of intial interface heights [Z-1 ~> m-1].
95  ! This include declares and sets the variable "version".
96 # include "version_variable.h"
97  character(len=20) :: verticalcoordinate
98  logical :: just_read ! If true, just read parameters but set nothing.
99  integer :: i, j, k, is, ie, js, je, nz
100 
101  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
102 
103  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
104 
105  if (.not.just_read) &
106  call mom_mesg("MOM_initialization.F90, initialize_thickness_uniform: setting thickness")
107 
108  if (.not.just_read) call log_version(param_file, mdl, version, "")
109  call get_param(param_file, mdl,"MIN_THICKNESS", min_thickness, &
110  'Minimum thickness for layer',&
111  units='m', default=1.0e-3, do_not_log=just_read, scale=us%m_to_Z)
112  call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalcoordinate, &
113  default=default_coordinate_mode, do_not_log=just_read)
114 
115  ! WARNING: this routine specifies the interface heights so that the last layer
116  ! is vanished, even at maximum depth. In order to have a uniform
117  ! layer distribution, use this line of code within the loop:
118  ! e0(k) = -G%max_depth * real(k-1) / real(nz)
119  ! To obtain a thickness distribution where the last layer is
120  ! vanished and the other thicknesses uniformly distributed, use:
121  ! e0(k) = -G%max_depth * real(k-1) / real(nz-1)
122  !do k=1,nz+1
123  ! e0(k) = -G%max_depth * real(k-1) / real(nz)
124  !enddo
125 
126  select case ( coordinatemode(verticalcoordinate) )
127 
128  case ( regridding_layer, regridding_rho ) ! Initial thicknesses for isopycnal coordinates
129  call get_param(param_file, mdl,"INITIAL_SSS", s_surf, default=34., do_not_log=.true.)
130  call get_param(param_file, mdl,"INITIAL_S_RANGE", s_range, default=2., do_not_log=.true.)
131  call get_param(param_file, mdl, "S_REF", s_ref, default=35.0, do_not_log=.true.)
132  call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", s_light, default = s_ref, do_not_log=.true.)
133  call get_param(param_file, mdl, "TS_RANGE_S_DENSE", s_dense, default = s_ref, do_not_log=.true.)
134  call get_param(param_file, mdl, "INTERFACE_IC_QUANTA", eta_ic_quanta, &
135  "The granularity of initial interface height values "//&
136  "per meter, to avoid sensivity to order-of-arithmetic changes.", &
137  default=2048.0, units="m-1", scale=us%Z_to_m, do_not_log=just_read)
138  if (just_read) return ! All run-time parameters have been read, so return.
139 
140  do k=1,nz+1
141  ! Salinity of layer k is S_light + (k-1)/(nz-1) * (S_dense - S_light)
142  ! Salinity of interface K is S_light + (K-3/2)/(nz-1) * (S_dense - S_light)
143  ! Salinity at depth z should be S(z) = S_surf - S_range * z/max_depth
144  ! Equating: S_surf - S_range * z/max_depth = S_light + (K-3/2)/(nz-1) * (S_dense - S_light)
145  ! Equating: - S_range * z/max_depth = S_light - S_surf + (K-3/2)/(nz-1) * (S_dense - S_light)
146  ! Equating: z/max_depth = - ( S_light - S_surf + (K-3/2)/(nz-1) * (S_dense - S_light) ) / S_range
147  e0(k) = - g%max_depth * ( ( s_light - s_surf ) + ( s_dense - s_light ) * &
148  ( (real(k)-1.5) / real(nz-1) ) ) / s_range
149  ! Force round numbers ... the above expression has irrational factors ...
150  if (eta_ic_quanta > 0.0) &
151  e0(k) = nint(eta_ic_quanta*e0(k)) / eta_ic_quanta
152  e0(k) = min(real(1-k)*gv%Angstrom_Z, e0(k)) ! Bound by surface
153  e0(k) = max(-g%max_depth, e0(k)) ! Bound by bottom
154  enddo
155  do j=js,je ; do i=is,ie
156  eta1d(nz+1) = -g%bathyT(i,j)
157  do k=nz,1,-1
158  eta1d(k) = e0(k)
159  if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z)) then
160  eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
161  h(i,j,k) = gv%Angstrom_H
162  else
163  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
164  endif
165  enddo
166  enddo ; enddo
167 
168  case ( regridding_zstar ) ! Initial thicknesses for z coordinates
169  if (just_read) return ! All run-time parameters have been read, so return.
170  do j=js,je ; do i=is,ie
171  eta1d(nz+1) = -g%bathyT(i,j)
172  do k=nz,1,-1
173  eta1d(k) = -g%max_depth * real(k-1) / real(nz)
174  if (eta1d(k) < (eta1d(k+1) + min_thickness)) then
175  eta1d(k) = eta1d(k+1) + min_thickness
176  h(i,j,k) = gv%Z_to_H * min_thickness
177  else
178  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
179  endif
180  enddo
181  enddo ; enddo
182 
183  case ( regridding_sigma ) ! Initial thicknesses for sigma coordinates
184  if (just_read) return ! All run-time parameters have been read, so return.
185  do j=js,je ; do i=is,ie
186  h(i,j,:) = gv%Z_to_H * g%bathyT(i,j) / dfloat(nz)
187  enddo ; enddo
188 
189 end select
190 
191 end subroutine dumbbell_initialize_thickness
192 
193 !> Initial values for temperature and salinity for the dumbbell test case
194 subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, param_file, &
195  eqn_of_state, just_read_params)
196  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
197  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
198  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: t !< Potential temperature [degC]
199  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: s !< Salinity [ppt]
200  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(in) :: h !< Layer thickness [H ~> m or kg m-2]
201  type(param_file_type), intent(in) :: param_file !< Parameter file structure
202  type(eos_type), pointer :: eqn_of_state !< Equation of state structure
203  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
204  !! only read parameters without changing h.
205 
206  ! Local variables
207  integer :: i, j, k, is, ie, js, je, nz, k_light
208  real :: xi0, xi1, dxi, r, s_surf, t_surf, s_range, t_range
209  real :: x, y, dblen
210  real :: t_ref, t_light, t_dense, s_ref, s_light, s_dense, a1, frac_dense, k_frac, res_rat
211  logical :: just_read ! If true, just read parameters but set nothing.
212  character(len=20) :: verticalcoordinate, density_profile
213 
214  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
215 
216  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
217 
218  t_surf = 20.0
219 
220  call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalcoordinate, &
221  default=default_coordinate_mode, do_not_log=just_read)
222  call get_param(param_file, mdl,"INITIAL_DENSITY_PROFILE", density_profile, &
223  'Initial profile shape. Valid values are "linear", "parabolic" '// &
224  'and "exponential".', default='linear', do_not_log=just_read)
225  call get_param(param_file, mdl,"DUMBBELL_SREF", s_surf, &
226  'DUMBBELL REFERENCE SALINITY', units='1e-3', default=34., do_not_log=just_read)
227  call get_param(param_file, mdl,"DUMBBELL_S_RANGE", s_range, &
228  'DUMBBELL salinity range (right-left)', units='1e-3', &
229  default=2., do_not_log=just_read)
230  call get_param(param_file, mdl,"DUMBBELL_LEN",dblen, &
231  'Lateral Length scale for dumbbell ',&
232  units='k', default=600., do_not_log=just_read)
233 
234  if (g%x_axis_units == 'm') then
235  dblen=dblen*1.e3
236  endif
237 
238  do j=g%jsc,g%jec
239  do i=g%isc,g%iec
240  ! Compute normalized zonal coordinates (x,y=0 at center of domain)
241  x = ( g%geoLonT(i,j) ) / dblen
242  do k=1,nz
243  t(i,j,k)=t_surf
244  enddo
245  if (x>=0. ) then
246  do k=1,nz
247  s(i,j,k)=s_surf + 0.5*s_range
248  enddo
249  endif
250  if (x<0. ) then
251  do k=1,nz
252  s(i,j,k)=s_surf - 0.5*s_range
253  enddo
254  endif
255 
256  enddo
257  enddo
258 
260 
261 !> Initialize the restoring sponges for the dumbbell test case
262 subroutine dumbbell_initialize_sponges(G, GV, US, tv, param_file, use_ALE, CSp, ACSp)
263  type(ocean_grid_type), intent(in) :: g !< Horizontal grid control structure
264  type(verticalgrid_type), intent(in) :: gv !< Vertical grid control structure
265  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
266  type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables
267  type(param_file_type), intent(in) :: param_file !< Parameter file structure
268  logical, intent(in) :: use_ale !< ALE flag
269  type(sponge_cs), pointer :: csp !< Layered sponge control structure pointer
270  type(ale_sponge_cs), pointer :: acsp !< ALE sponge control structure pointer
271 
272  real :: sponge_time_scale
273 
274  real, dimension(SZI_(G),SZJ_(G)) :: idamp ! inverse damping timescale
275  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h, t, s ! sponge thicknesses, temp and salt
276  real, dimension(SZK_(GV)+1) :: e0, eta1d ! interface positions for ALE sponge
277 
278  integer :: i, j, k, nz
279  real :: x, zi, zmid, dist, min_thickness, dblen
280  real :: mld, s_ref, s_range, s_dense, t_ref, sill_height
281  call get_param(param_file, mdl,"DUMBBELL_LEN",dblen, &
282  'Lateral Length scale for dumbbell ',&
283  units='k', default=600., do_not_log=.true.)
284 
285  if (g%x_axis_units == 'm') then
286  dblen=dblen*1.e3
287  endif
288 
289  nz = gv%ke
290 
291  call get_param(param_file, mdl, "DUMBBELL_SPONGE_TIME_SCALE", sponge_time_scale, &
292  "The time scale in the reservoir for restoring. If zero, the sponge is disabled.", &
293  units="s", default=0.)
294  call get_param(param_file, mdl, "DUMBBELL_SREF", s_ref, do_not_log=.true.)
295  call get_param(param_file, mdl, "DUMBBELL_S_RANGE", s_range, do_not_log=.true.)
296  call get_param(param_file, mdl,"MIN_THICKNESS", min_thickness, &
297  'Minimum thickness for layer',&
298  units='m', default=1.0e-3, do_not_log=.true., scale=us%m_to_Z)
299 
300  ! no active sponges
301  if (sponge_time_scale <= 0.) return
302 
303  ! everywhere is initially unsponged
304  idamp(:,:) = 0.0
305 
306  do j = g%jsc, g%jec
307  do i = g%isc,g%iec
308  if (g%mask2dT(i,j) > 0.) then
309  ! nondimensional x position
310  x = (g%geoLonT(i,j) ) / dblen
311  if (x > 0.25 .or. x < -0.25) then
312  ! scale restoring by depth into sponge
313  idamp(i,j) = 1. / sponge_time_scale
314  endif
315  endif
316  enddo
317  enddo
318 
319  if (use_ale) then
320  ! construct a uniform grid for the sponge
321  do j=g%jsc,g%jec ; do i=g%isc,g%iec
322  eta1d(nz+1) = -g%bathyT(i,j)
323  do k=nz,1,-1
324  eta1d(k) = -g%max_depth * real(k-1) / real(nz)
325  if (eta1d(k) < (eta1d(k+1) + min_thickness)) then
326  eta1d(k) = eta1d(k+1) + min_thickness
327  h(i,j,k) = gv%Z_to_H * min_thickness
328  else
329  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
330  endif
331  enddo
332  enddo ; enddo
333 
334  call initialize_ale_sponge(idamp, g, param_file, acsp, h, nz)
335 
336  ! construct temperature and salinity for the sponge
337  ! start with initial condition
338  s(:,:,:) = 0.0
339 
340  do j=g%jsc,g%jec ; do i=g%isc,g%iec
341  ! Compute normalized zonal coordinates (x,y=0 at center of domain)
342  x = ( g%geoLonT(i,j) ) / dblen
343  if (x>=0.25 ) then
344  do k=1,nz
345  s(i,j,k)=s_ref + 0.5*s_range
346  enddo
347  endif
348  if (x<=-0.25 ) then
349  do k=1,nz
350  s(i,j,k)=s_ref - 0.5*s_range
351  enddo
352  endif
353  enddo ; enddo
354  endif
355 
356  if (associated(tv%S)) call set_up_ale_sponge_field(s, g, tv%S, acsp)
357 
358 end subroutine dumbbell_initialize_sponges
359 
360 end module dumbbell_initialization
dumbbell_initialization::mdl
character(len=40) mdl
This module's name.
Definition: dumbbell_initialization.F90:27
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_sponge::set_up_sponge_field
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
Definition: MOM_sponge.F90:214
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_get_input::directories
Container for paths and parameter file names.
Definition: MOM_get_input.F90:20
mom_ale_sponge::initialize_ale_sponge
Ddetermine the number of points which are within sponges in this computational domain.
Definition: MOM_ALE_sponge.F90:48
mom_dyn_horgrid
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
Definition: MOM_dyn_horgrid.F90:3
mom_ale_sponge::ale_sponge_cs
ALE sponge control structure.
Definition: MOM_ALE_sponge.F90:84
dumbbell_initialization::dumbbell_initialize_thickness
subroutine, public dumbbell_initialize_thickness(h, G, GV, US, param_file, just_read_params)
Initializes the layer thicknesses to be uniform in the dumbbell test case.
Definition: dumbbell_initialization.F90:78
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_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_tracer_registry
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
Definition: MOM_tracer_registry.F90:5
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_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
regrid_consts::regridding_layer
integer, parameter regridding_layer
Layer mode identifier.
Definition: regrid_consts.F90:13
mom_ale_sponge
This module contains the routines used to apply sponge layers when using the ALE mode.
Definition: MOM_ALE_sponge.F90:11
dumbbell_initialization::dumbbell_initialize_sponges
subroutine, public dumbbell_initialize_sponges(G, GV, US, tv, param_file, use_ALE, CSp, ACSp)
Initialize the restoring sponges for the dumbbell test case.
Definition: dumbbell_initialization.F90:263
mom_get_input
Reads the only Fortran name list needed to boot-strap the model.
Definition: MOM_get_input.F90:6
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
regrid_consts
Contains constants for interpreting input parameters that control regridding.
Definition: regrid_consts.F90:2
dumbbell_initialization
Configures the model for the idealized dumbbell test case.
Definition: dumbbell_initialization.F90:2
mom_eos::eos_type
A control structure for the equation of state.
Definition: MOM_EOS.F90:91
mom_ale_sponge::set_up_ale_sponge_field
Store the reference profile at h points for a variable.
Definition: MOM_ALE_sponge.F90:33
regrid_consts::regridding_sigma
integer, parameter regridding_sigma
Sigma coordinates identifier.
Definition: regrid_consts.F90:16
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
dumbbell_initialization::dumbbell_initialize_temperature_salinity
subroutine, public dumbbell_initialize_temperature_salinity(T, S, h, G, GV, param_file, eqn_of_state, just_read_params)
Initial values for temperature and salinity for the dumbbell test case.
Definition: dumbbell_initialization.F90:196
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_sponge
Implements sponge regions in isopycnal mode.
Definition: MOM_sponge.F90:2
regrid_consts::regridding_zstar
integer, parameter regridding_zstar
z* coordinates identifier
Definition: regrid_consts.F90:14
mom_eos::calculate_density_derivs
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Definition: MOM_EOS.F90:70
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_sponge::sponge_cs
This control structure holds memory and parameters for the MOM_sponge module.
Definition: MOM_sponge.F90:41
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
regrid_consts::coordinatemode
integer function coordinatemode(string)
Parse a string parameter specifying the coordinate mode and return the appropriate enumerated integer...
Definition: regrid_consts.F90:54
dumbbell_initialization::dumbbell_initialize_topography
subroutine, public dumbbell_initialize_topography(D, G, param_file, max_depth)
Initialization of topography.
Definition: dumbbell_initialization.F90:43
regrid_consts::default_coordinate_mode
character(len= *), parameter default_coordinate_mode
Default coordinate mode.
Definition: regrid_consts.F90:35
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
regrid_consts::regridding_rho
integer, parameter regridding_rho
Density coordinates identifier.
Definition: regrid_consts.F90:15
mom_sponge::initialize_sponge
subroutine, public initialize_sponge(Iresttime, int_height, G, param_file, CS, GV, Iresttime_i_mean, int_height_i_mean)
This subroutine determines the number of points which are within sponges in this computational domain...
Definition: MOM_sponge.F90:90
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_dyn_horgrid::dyn_horgrid_type
Describes the horizontal ocean grid with only dynamic memory arrays.
Definition: MOM_dyn_horgrid.F90:23
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:60