MOM6
MOM_EOS.F90
1 !> Provides subroutines for quantities specific to the equation of state
2 module mom_eos
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_eos_linear, only : calculate_density_linear, calculate_spec_vol_linear
7 use mom_eos_linear, only : calculate_density_derivs_linear
8 use mom_eos_linear, only : calculate_specvol_derivs_linear, int_density_dz_linear
9 use mom_eos_linear, only : calculate_density_second_derivs_linear
10 use mom_eos_linear, only : calculate_compress_linear, int_spec_vol_dp_linear
11 use mom_eos_wright, only : calculate_density_wright, calculate_spec_vol_wright
12 use mom_eos_wright, only : calculate_density_derivs_wright
13 use mom_eos_wright, only : calculate_specvol_derivs_wright, int_density_dz_wright
14 use mom_eos_wright, only : calculate_compress_wright, int_spec_vol_dp_wright
15 use mom_eos_wright, only : calculate_density_second_derivs_wright
16 use mom_eos_unesco, only : calculate_density_unesco, calculate_spec_vol_unesco
17 use mom_eos_unesco, only : calculate_density_derivs_unesco, calculate_density_unesco
18 use mom_eos_unesco, only : calculate_compress_unesco
19 use mom_eos_nemo, only : calculate_density_nemo
20 use mom_eos_nemo, only : calculate_density_derivs_nemo, calculate_density_nemo
21 use mom_eos_nemo, only : calculate_compress_nemo
22 use mom_eos_teos10, only : calculate_density_teos10, calculate_spec_vol_teos10
23 use mom_eos_teos10, only : calculate_density_derivs_teos10
24 use mom_eos_teos10, only : calculate_specvol_derivs_teos10
25 use mom_eos_teos10, only : calculate_density_second_derivs_teos10
26 use mom_eos_teos10, only : calculate_compress_teos10
27 use mom_eos_teos10, only : gsw_sp_from_sr, gsw_pt_from_ct
28 use mom_tfreeze, only : calculate_tfreeze_linear, calculate_tfreeze_millero
29 use mom_tfreeze, only : calculate_tfreeze_teos10
30 use mom_error_handler, only : mom_error, fatal, warning, mom_mesg
31 use mom_file_parser, only : get_param, log_version, param_file_type
32 use mom_hor_index, only : hor_index_type
33 use mom_string_functions, only : uppercase
34 use mom_unit_scaling, only : unit_scale_type
35 
36 implicit none ; private
37 
38 #include <MOM_memory.h>
39 
40 public eos_allocate
41 public eos_domain
42 public eos_end
43 public eos_init
44 public eos_manual_init
45 public eos_quadrature
46 public eos_use_linear
47 public analytic_int_density_dz
48 public analytic_int_specific_vol_dp
49 public calculate_compress
50 public calculate_density
51 public calculate_density_derivs
52 public calculate_density_second_derivs
53 public calculate_spec_vol
54 public calculate_specific_vol_derivs
55 public calculate_tfreeze
56 public convert_temp_salt_for_teos10
57 public extract_member_eos
58 public gsw_sp_from_sr
59 public gsw_pt_from_ct
60 public query_compressible
61 
62 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
63 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
64 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
65 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
66 
67 !> Calculates density of sea water from T, S and P
68 interface calculate_density
69  module procedure calculate_density_scalar, calculate_density_array, calculate_density_1d
70  module procedure calculate_stanley_density_scalar, calculate_stanley_density_array
71  module procedure calculate_stanley_density_1d
72 end interface calculate_density
73 
74 !> Calculates specific volume of sea water from T, S and P
75 interface calculate_spec_vol
76  module procedure calc_spec_vol_scalar, calculate_spec_vol_array, &
77  calc_spec_vol_1d
78 end interface calculate_spec_vol
79 
80 !> Calculate the derivatives of density with temperature and salinity from T, S, and P
81 interface calculate_density_derivs
82  module procedure calculate_density_derivs_scalar, calculate_density_derivs_array, &
83  calculate_density_derivs_1d
84 end interface calculate_density_derivs
85 
86 !> Calculate the derivatives of specific volume with temperature and salinity from T, S, and P
87 interface calculate_specific_vol_derivs
88  module procedure calculate_spec_vol_derivs_array, calc_spec_vol_derivs_1d
89 end interface calculate_specific_vol_derivs
90 
91 !> Calculates the second derivatives of density with various combinations of temperature,
92 !! salinity, and pressure from T, S and P
93 interface calculate_density_second_derivs
94  module procedure calculate_density_second_derivs_scalar, calculate_density_second_derivs_array
95 end interface calculate_density_second_derivs
96 
97 !> Calculates the freezing point of sea water from T, S and P
98 interface calculate_tfreeze
99  module procedure calculate_tfreeze_scalar, calculate_tfreeze_array
100 end interface calculate_tfreeze
101 
102 !> Calculates the compressibility of water from T, S, and P
103 interface calculate_compress
104  module procedure calculate_compress_scalar, calculate_compress_array
105 end interface calculate_compress
106 
107 !> A control structure for the equation of state
108 type, public :: eos_type ; private
109  integer :: form_of_eos = 0 !< The equation of state to use.
110  integer :: form_of_tfreeze = 0 !< The expression for the potential temperature
111  !! of the freezing point.
112  logical :: eos_quadrature !< If true, always use the generic (quadrature)
113  !! code for the integrals of density.
114  logical :: compressible = .true. !< If true, in situ density is a function of pressure.
115 ! The following parameters are used with the linear equation of state only.
116  real :: rho_t0_s0 !< The density at T=0, S=0 [kg m-3]
117  real :: drho_dt !< The partial derivative of density with temperature [kg m-3 degC-1]
118  real :: drho_ds !< The partial derivative of density with salinity [kg m-3 ppt-1]
119 ! The following parameters are use with the linear expression for the freezing
120 ! point only.
121  real :: tfr_s0_p0 !< The freezing potential temperature at S=0, P=0 [degC]
122  real :: dtfr_ds !< The derivative of freezing point with salinity [degC ppt-1]
123  real :: dtfr_dp !< The derivative of freezing point with pressure [degC Pa-1]
124 
125 ! Unit conversion factors (normally used for dimensional testing but could also allow for
126 ! change of units of arguments to functions)
127  real :: m_to_z = 1. !< A constant that translates distances in meters to the units of depth.
128  real :: kg_m3_to_r = 1. !< A constant that translates kilograms per meter cubed to the units of density.
129  real :: r_to_kg_m3 = 1. !< A constant that translates the units of density to kilograms per meter cubed.
130  real :: rl2_t2_to_pa = 1.!< Convert pressures from R L2 T-2 to Pa.
131  real :: l_t_to_m_s = 1. !< Convert lateral velocities from L T-1 to m s-1.
132 
133 ! logical :: test_EOS = .true. ! If true, test the equation of state
134 end type eos_type
135 
136 ! The named integers that might be stored in eqn_of_state_type%form_of_EOS.
137 integer, parameter, public :: eos_linear = 1 !< A named integer specifying an equation of state
138 integer, parameter, public :: eos_unesco = 2 !< A named integer specifying an equation of state
139 integer, parameter, public :: eos_wright = 3 !< A named integer specifying an equation of state
140 integer, parameter, public :: eos_teos10 = 4 !< A named integer specifying an equation of state
141 integer, parameter, public :: eos_nemo = 5 !< A named integer specifying an equation of state
142 
143 character*(10), parameter :: eos_linear_string = "LINEAR" !< A string for specifying the equation of state
144 character*(10), parameter :: eos_unesco_string = "UNESCO" !< A string for specifying the equation of state
145 character*(10), parameter :: eos_wright_string = "WRIGHT" !< A string for specifying the equation of state
146 character*(10), parameter :: eos_teos10_string = "TEOS10" !< A string for specifying the equation of state
147 character*(10), parameter :: eos_nemo_string = "NEMO" !< A string for specifying the equation of state
148 character*(10), parameter :: eos_default = eos_wright_string !< The default equation of state
149 
150 integer, parameter :: tfreeze_linear = 1 !< A named integer specifying a freezing point expression
151 integer, parameter :: tfreeze_millero = 2 !< A named integer specifying a freezing point expression
152 integer, parameter :: tfreeze_teos10 = 3 !< A named integer specifying a freezing point expression
153 character*(10), parameter :: tfreeze_linear_string = "LINEAR" !< A string for specifying the freezing point expression
154 character*(10), parameter :: tfreeze_millero_string = "MILLERO_78" !< A string for specifying
155  !! freezing point expression
156 character*(10), parameter :: tfreeze_teos10_string = "TEOS10" !< A string for specifying the freezing point expression
157 character*(10), parameter :: tfreeze_default = tfreeze_linear_string !< The default freezing point expression
158 
159 contains
160 
161 !> Calls the appropriate subroutine to calculate density of sea water for scalar inputs.
162 !! If rho_ref is present, the anomaly with respect to rho_ref is returned. The pressure and
163 !! density can be rescaled with the US. If both the US and scale arguments are present the density
164 !! scaling uses the product of the two scaling factors.
165 subroutine calculate_density_scalar(T, S, pressure, rho, EOS, rho_ref, scale)
166  real, intent(in) :: T !< Potential temperature referenced to the surface [degC]
167  real, intent(in) :: S !< Salinity [ppt]
168  real, intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
169  real, intent(out) :: rho !< Density (in-situ if pressure is local) [kg m-3] or [R ~> kg m-3]
170  type(EOS_type), pointer :: EOS !< Equation of state structure
171  real, optional, intent(in) :: rho_ref !< A reference density [kg m-3]
172  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density in
173  !! combination with scaling given by US [various]
174 
175  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
176  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
177 
178  if (.not.associated(eos)) call mom_error(fatal, &
179  "calculate_density_scalar called with an unassociated EOS_type EOS.")
180 
181  p_scale = eos%RL2_T2_to_Pa
182 
183  select case (eos%form_of_EOS)
184  case (eos_linear)
185  call calculate_density_linear(t, s, p_scale*pressure, rho, &
186  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS, rho_ref)
187  case (eos_unesco)
188  call calculate_density_unesco(t, s, p_scale*pressure, rho, rho_ref)
189  case (eos_wright)
190  call calculate_density_wright(t, s, p_scale*pressure, rho, rho_ref)
191  case (eos_teos10)
192  call calculate_density_teos10(t, s, p_scale*pressure, rho, rho_ref)
193  case (eos_nemo)
194  call calculate_density_nemo(t, s, p_scale*pressure, rho, rho_ref)
195  case default
196  call mom_error(fatal, "calculate_density_scalar: EOS is not valid.")
197  end select
198 
199  rho_scale = eos%kg_m3_to_R
200  if (present(scale)) rho_scale = rho_scale * scale
201  rho = rho_scale * rho
202 
203 end subroutine calculate_density_scalar
204 
205 !> Calls the appropriate subroutine to calculate density of sea water for scalar inputs
206 !! including the variance of T, S and covariance of T-S.
207 !! The calculation uses only the second order correction in a series as discussed
208 !! in Stanley et al., 2020.
209 !! If rho_ref is present, the anomaly with respect to rho_ref is returned. The
210 !! density can be rescaled using rho_ref.
211 subroutine calculate_stanley_density_scalar(T, S, pressure, Tvar, TScov, Svar, rho, EOS, rho_ref, scale)
212  real, intent(in) :: T !< Potential temperature referenced to the surface [degC]
213  real, intent(in) :: S !< Salinity [ppt]
214  real, intent(in) :: Tvar !< Variance of potential temperature referenced to the surface [degC2]
215  real, intent(in) :: TScov !< Covariance of potential temperature and salinity [degC ppt]
216  real, intent(in) :: Svar !< Variance of salinity [ppt2]
217  real, intent(in) :: pressure !< Pressure [Pa]
218  real, intent(out) :: rho !< Density (in-situ if pressure is local) [kg m-3] or [R ~> kg m-3]
219  type(EOS_type), pointer :: EOS !< Equation of state structure
220  real, optional, intent(in) :: rho_ref !< A reference density [kg m-3].
221  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
222  !! from kg m-3 to the desired units [R m3 kg-1]
223  ! Local variables
224  real :: d2RdTT, d2RdST, d2RdSS, d2RdSp, d2RdTp ! Second derivatives of density wrt T,S,p
225  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
226  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
227 
228  if (.not.associated(eos)) call mom_error(fatal, &
229  "calculate_stanley_density_scalar called with an unassociated EOS_type EOS.")
230 
231  p_scale = eos%RL2_T2_to_Pa
232 
233  select case (eos%form_of_EOS)
234  case (eos_linear)
235  call calculate_density_linear(t, s, p_scale*pressure, rho, &
236  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS, rho_ref)
237  call calculate_density_second_derivs_linear(t, s, pressure, d2rdss, d2rdst, &
238  d2rdtt, d2rdsp, d2rdtp)
239  case (eos_wright)
240  call calculate_density_wright(t, s, p_scale*pressure, rho, rho_ref)
241  call calculate_density_second_derivs_wright(t, s, pressure, d2rdss, d2rdst, &
242  d2rdtt, d2rdsp, d2rdtp)
243  case (eos_teos10)
244  call calculate_density_teos10(t, s, p_scale*pressure, rho, rho_ref)
245  call calculate_density_second_derivs_teos10(t, s, pressure, d2rdss, d2rdst, &
246  d2rdtt, d2rdsp, d2rdtp)
247  case default
248  call mom_error(fatal, "calculate_stanley_density_scalar: EOS is not valid.")
249  end select
250 
251  ! Equation 25 of Stanley et al., 2020.
252  rho = rho + ( 0.5 * d2rdtt * tvar + ( d2rdst * tscov + 0.5 * d2rdss * svar ) )
253 
254  rho_scale = eos%kg_m3_to_R
255  if (present(scale)) rho_scale = rho_scale * scale
256  rho = rho_scale * rho
257 
258 end subroutine calculate_stanley_density_scalar
259 
260 !> Calls the appropriate subroutine to calculate the density of sea water for 1-D array inputs.
261 !! If rho_ref is present, the anomaly with respect to rho_ref is returned.
262 subroutine calculate_density_array(T, S, pressure, rho, start, npts, EOS, rho_ref, scale)
263  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
264  real, dimension(:), intent(in) :: S !< Salinity [ppt]
265  real, dimension(:), intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
266  real, dimension(:), intent(inout) :: rho !< Density (in-situ if pressure is local) [kg m-3] or [R ~> kg m-3]
267  integer, intent(in) :: start !< Start index for computation
268  integer, intent(in) :: npts !< Number of point to compute
269  type(EOS_type), pointer :: EOS !< Equation of state structure
270  real, optional, intent(in) :: rho_ref !< A reference density [kg m-3]
271  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
272  !! in combination with scaling given by US [various]
273  integer :: j
274 
275  if (.not.associated(eos)) call mom_error(fatal, &
276  "calculate_density_array called with an unassociated EOS_type EOS.")
277 
278  select case (eos%form_of_EOS)
279  case (eos_linear)
280  call calculate_density_linear(t, s, pressure, rho, start, npts, &
281  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS, rho_ref)
282  case (eos_unesco)
283  call calculate_density_unesco(t, s, pressure, rho, start, npts, rho_ref)
284  case (eos_wright)
285  call calculate_density_wright(t, s, pressure, rho, start, npts, rho_ref)
286  case (eos_teos10)
287  call calculate_density_teos10(t, s, pressure, rho, start, npts, rho_ref)
288  case (eos_nemo)
289  call calculate_density_nemo(t, s, pressure, rho, start, npts, rho_ref)
290  case default
291  call mom_error(fatal, "calculate_density_array: EOS%form_of_EOS is not valid.")
292  end select
293 
294  if (present(scale)) then ; if (scale /= 1.0) then ; do j=start,start+npts-1
295  rho(j) = scale * rho(j)
296  enddo ; endif ; endif
297 
298 end subroutine calculate_density_array
299 
300 !> Calls the appropriate subroutine to calculate the density of sea water for 1-D array inputs
301 !! including the variance of T, S and covariance of T-S.
302 !! The calculation uses only the second order correction in a series as discussed
303 !! in Stanley et al., 2020.
304 !! If rho_ref is present, the anomaly with respect to rho_ref is returned.
305 subroutine calculate_stanley_density_array(T, S, pressure, Tvar, TScov, Svar, rho, start, npts, EOS, rho_ref, scale)
306  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
307  real, dimension(:), intent(in) :: S !< Salinity [ppt]
308  real, dimension(:), intent(in) :: pressure !< Pressure [Pa]
309  real, dimension(:), intent(in) :: Tvar !< Variance of potential temperature referenced to the surface [degC2]
310  real, dimension(:), intent(in) :: TScov !< Covariance of potential temperature and salinity [degC ppt]
311  real, dimension(:), intent(in) :: Svar !< Variance of salinity [ppt2]
312  real, dimension(:), intent(inout) :: rho !< Density (in-situ if pressure is local) [kg m-3]
313  integer, intent(in) :: start !< Start index for computation
314  integer, intent(in) :: npts !< Number of point to compute
315  type(EOS_type), pointer :: EOS !< Equation of state structure
316  real, optional, intent(in) :: rho_ref !< A reference density [kg m-3].
317  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
318  !! from kg m-3 to the desired units [R m3 kg-1]
319  ! Local variables
320  real, dimension(size(T)) :: d2RdTT, d2RdST, d2RdSS, d2RdSp, d2RdTp ! Second derivatives of density wrt T,S,p
321  integer :: j
322 
323  if (.not.associated(eos)) call mom_error(fatal, &
324  "calculate_density_array called with an unassociated EOS_type EOS.")
325 
326  select case (eos%form_of_EOS)
327  case (eos_linear)
328  call calculate_density_linear(t, s, pressure, rho, start, npts, &
329  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS, rho_ref)
330  call calculate_density_second_derivs_linear(t, s, pressure, d2rdss, d2rdst, &
331  d2rdtt, d2rdsp, d2rdtp, start, npts)
332  case (eos_wright)
333  call calculate_density_wright(t, s, pressure, rho, start, npts, rho_ref)
334  call calculate_density_second_derivs_wright(t, s, pressure, d2rdss, d2rdst, &
335  d2rdtt, d2rdsp, d2rdtp, start, npts)
336  case (eos_teos10)
337  call calculate_density_teos10(t, s, pressure, rho, start, npts, rho_ref)
338  call calculate_density_second_derivs_teos10(t, s, pressure, d2rdss, d2rdst, &
339  d2rdtt, d2rdsp, d2rdtp, start, npts)
340  case default
341  call mom_error(fatal, "calculate_stanley_density_array: EOS%form_of_EOS is not valid.")
342  end select
343 
344  ! Equation 25 of Stanley et al., 2020.
345  do j=start,start+npts-1
346  rho(j) = rho(j) &
347  + ( 0.5 * d2rdtt(j) * tvar(j) + ( d2rdst(j) * tscov(j) + 0.5 * d2rdss(j) * svar(j) ) )
348  enddo
349 
350  if (present(scale)) then ; if (scale /= 1.0) then ; do j=start,start+npts-1
351  rho(j) = scale * rho(j)
352  enddo ; endif ; endif
353 
354 end subroutine calculate_stanley_density_array
355 
356 !> Calls the appropriate subroutine to calculate the density of sea water for 1-D array inputs,
357 !! potentially limiting the domain of indices that are worked on.
358 !! If rho_ref is present, the anomaly with respect to rho_ref is returned.
359 subroutine calculate_density_1d(T, S, pressure, rho, EOS, dom, rho_ref, scale)
360  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
361  real, dimension(:), intent(in) :: S !< Salinity [ppt]
362  real, dimension(:), intent(in) :: pressure !< Pressure [R L2 T-2 ~> Pa]
363  real, dimension(:), intent(inout) :: rho !< Density (in-situ if pressure is local) [R ~> kg m-3]
364  type(EOS_type), pointer :: EOS !< Equation of state structure
365  integer, dimension(2), optional, intent(in) :: dom !< The domain of indices to work on, taking
366  !! into account that arrays start at 1.
367  real, optional, intent(in) :: rho_ref !< A reference density [kg m-3]
368  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
369  !! in combination with scaling given by US [various]
370  ! Local variables
371  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
372  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
373  real :: rho_unscale ! A factor to convert density from R to kg m-3 [kg m-3 R-1 ~> 1]
374  real :: rho_reference ! rho_ref converted to [kg m-3]
375  real, dimension(size(rho)) :: pres ! Pressure converted to [Pa]
376  integer :: i, is, ie, npts
377 
378  if (.not.associated(eos)) call mom_error(fatal, &
379  "calculate_density_1d called with an unassociated EOS_type EOS.")
380 
381  if (present(dom)) then
382  is = dom(1) ; ie = dom(2) ; npts = 1 + ie - is
383  else
384  is = 1 ; ie = size(rho) ; npts = 1 + ie - is
385  endif
386 
387  p_scale = eos%RL2_T2_to_Pa
388  rho_unscale = eos%R_to_kg_m3
389 
390  if ((p_scale == 1.0) .and. (rho_unscale == 1.0)) then
391  call calculate_density_array(t, s, pressure, rho, is, npts, eos, rho_ref=rho_ref)
392  elseif (present(rho_ref)) then ! This is the same as above, but with some extra work to rescale variables.
393  do i=is,ie ; pres(i) = p_scale * pressure(i) ; enddo
394  rho_reference = rho_unscale*rho_ref
395  call calculate_density_array(t, s, pres, rho, is, npts, eos, rho_ref=rho_reference)
396  else ! There is rescaling of variables, but rho_ref is not present. Passing a 0 value of rho_ref
397  ! changes answers at roundoff for some equations of state, like Wright and UNESCO.
398  do i=is,ie ; pres(i) = p_scale * pressure(i) ; enddo
399  call calculate_density_array(t, s, pres, rho, is, npts, eos)
400  endif
401 
402  rho_scale = eos%kg_m3_to_R
403  if (present(scale)) rho_scale = rho_scale * scale
404  if (rho_scale /= 1.0) then ; do i=is,ie
405  rho(i) = rho_scale * rho(i)
406  enddo ; endif
407 
408 end subroutine calculate_density_1d
409 
410 !> Calls the appropriate subroutine to calculate the density of sea water for 1-D array inputs
411 !! including the variance of T, S and covariance of T-S,
412 !! potentially limiting the domain of indices that are worked on.
413 !! The calculation uses only the second order correction in a series as discussed
414 !! in Stanley et al., 2020.
415 !! If rho_ref is present, the anomaly with respect to rho_ref is returned.
416 subroutine calculate_stanley_density_1d(T, S, pressure, Tvar, TScov, Svar, rho, EOS, dom, rho_ref, scale)
417  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
418  real, dimension(:), intent(in) :: S !< Salinity [ppt]
419  real, dimension(:), intent(in) :: pressure !< Pressure [R L2 T-2 ~> Pa]
420  real, dimension(:), intent(in) :: Tvar !< Variance of potential temperature [degC2]
421  real, dimension(:), intent(in) :: TScov !< Covariance of potential temperature and salinity [degC ppt]
422  real, dimension(:), intent(in) :: Svar !< Variance of salinity [ppt2]
423  real, dimension(:), intent(inout) :: rho !< Density (in-situ if pressure is local) [R ~> kg m-3]
424  type(EOS_type), pointer :: EOS !< Equation of state structure
425  integer, dimension(2), optional, intent(in) :: dom !< The domain of indices to work on, taking
426  !! into account that arrays start at 1.
427  real, optional, intent(in) :: rho_ref !< A reference density [kg m-3]
428  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
429  !! in combination with scaling given by US [various]
430  ! Local variables
431  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
432  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
433  real, dimension(size(rho)) :: pres ! Pressure converted to [Pa]
434  real, dimension(size(T)) :: d2RdTT, d2RdST, d2RdSS, d2RdSp, d2RdTp ! Second derivatives of density wrt T,S,p
435  integer :: i, is, ie, npts
436 
437  if (.not.associated(eos)) call mom_error(fatal, &
438  "calculate_density_1d called with an unassociated EOS_type EOS.")
439 
440  if (present(dom)) then
441  is = dom(1) ; ie = dom(2) ; npts = 1 + ie - is
442  else
443  is = 1 ; ie = size(rho) ; npts = 1 + ie - is
444  endif
445 
446  p_scale = eos%RL2_T2_to_Pa
447  do i=is,ie
448  pres(i) = p_scale * pressure(i)
449  enddo
450 
451  select case (eos%form_of_EOS)
452  case (eos_linear)
453  call calculate_density_linear(t, s, pres, rho, 1, npts, &
454  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS, rho_ref)
455  call calculate_density_second_derivs_linear(t, s, pres, d2rdss, d2rdst, &
456  d2rdtt, d2rdsp, d2rdtp, 1, npts)
457  case (eos_wright)
458  call calculate_density_wright(t, s, pres, rho, 1, npts, rho_ref)
459  call calculate_density_second_derivs_wright(t, s, pres, d2rdss, d2rdst, &
460  d2rdtt, d2rdsp, d2rdtp, 1, npts)
461  case (eos_teos10)
462  call calculate_density_teos10(t, s, pres, rho, 1, npts, rho_ref)
463  call calculate_density_second_derivs_teos10(t, s, pres, d2rdss, d2rdst, &
464  d2rdtt, d2rdsp, d2rdtp, 1, npts)
465  case default
466  call mom_error(fatal, "calculate_stanley_density_scalar: EOS is not valid.")
467  end select
468 
469  ! Equation 25 of Stanley et al., 2020.
470  do i=is,ie
471  rho(i) = rho(i) &
472  + ( 0.5 * d2rdtt(i) * tvar(i) + ( d2rdst(i) * tscov(i) + 0.5 * d2rdss(i) * svar(i) ) )
473  enddo
474 
475  rho_scale = eos%kg_m3_to_R
476  if (present(scale)) rho_scale = rho_scale * scale
477  if (rho_scale /= 1.0) then ; do i=is,ie
478  rho(i) = rho_scale * rho(i)
479  enddo ; endif
480 
481 end subroutine calculate_stanley_density_1d
482 
483 !> Calls the appropriate subroutine to calculate the specific volume of sea water
484 !! for 1-D array inputs.
485 subroutine calculate_spec_vol_array(T, S, pressure, specvol, start, npts, EOS, spv_ref, scale)
486  real, dimension(:), intent(in) :: T !< potential temperature relative to the surface [degC]
487  real, dimension(:), intent(in) :: S !< salinity [ppt]
488  real, dimension(:), intent(in) :: pressure !< pressure [Pa]
489  real, dimension(:), intent(inout) :: specvol !< in situ specific volume [kg m-3]
490  integer, intent(in) :: start !< the starting point in the arrays.
491  integer, intent(in) :: npts !< the number of values to calculate.
492  type(EOS_type), pointer :: EOS !< Equation of state structure
493  real, optional, intent(in) :: spv_ref !< A reference specific volume [m3 kg-1]
494  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale specific
495  !! volume in combination with scaling given by US [various]
496 
497  real, dimension(size(specvol)) :: rho ! Density [kg m-3]
498  integer :: j
499 
500  if (.not.associated(eos)) call mom_error(fatal, &
501  "calculate_spec_vol_array called with an unassociated EOS_type EOS.")
502 
503  select case (eos%form_of_EOS)
504  case (eos_linear)
505  call calculate_spec_vol_linear(t, s, pressure, specvol, start, npts, &
506  eos%rho_T0_S0, eos%drho_dT, eos%drho_dS, spv_ref)
507  case (eos_unesco)
508  call calculate_spec_vol_unesco(t, s, pressure, specvol, start, npts, spv_ref)
509  case (eos_wright)
510  call calculate_spec_vol_wright(t, s, pressure, specvol, start, npts, spv_ref)
511  case (eos_teos10)
512  call calculate_spec_vol_teos10(t, s, pressure, specvol, start, npts, spv_ref)
513  case (eos_nemo)
514  call calculate_density_nemo(t, s, pressure, rho, start, npts)
515  if (present(spv_ref)) then
516  specvol(:) = 1.0 / rho(:) - spv_ref
517  else
518  specvol(:) = 1.0 / rho(:)
519  endif
520  case default
521  call mom_error(fatal, "calculate_spec_vol_array: EOS%form_of_EOS is not valid.")
522  end select
523 
524  if (present(scale)) then ; if (scale /= 1.0) then ; do j=start,start+npts-1
525  specvol(j) = scale * specvol(j)
526  enddo ; endif ; endif
527 
528 end subroutine calculate_spec_vol_array
529 
530 !> Calls the appropriate subroutine to calculate specific volume of sea water
531 !! for scalar inputs.
532 subroutine calc_spec_vol_scalar(T, S, pressure, specvol, EOS, spv_ref, scale)
533  real, intent(in) :: T !< Potential temperature referenced to the surface [degC]
534  real, intent(in) :: S !< Salinity [ppt]
535  real, intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
536  real, intent(out) :: specvol !< In situ? specific volume [m3 kg-1] or [R-1 ~> m3 kg-1]
537  type(EOS_type), pointer :: EOS !< Equation of state structure
538  real, optional, intent(in) :: spv_ref !< A reference specific volume [m3 kg-1] or [R-1 m3 kg-1]
539  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale specific
540  !! volume in combination with scaling given by US [various]
541 
542  real, dimension(1) :: Ta, Sa, pres, spv ! Rescaled single element array versions of the arguments.
543  real :: spv_reference ! spv_ref converted to [m3 kg-1]
544  real :: spv_scale ! A factor to convert specific volume from m3 kg-1 to the desired units [kg R-1 m-3 ~> 1]
545 
546  if (.not.associated(eos)) call mom_error(fatal, &
547  "calc_spec_vol_scalar called with an unassociated EOS_type EOS.")
548 
549  pres(1) = eos%RL2_T2_to_Pa*pressure
550  ta(1) = t ; sa(1) = s
551 
552  if (present(spv_ref)) then
553  spv_reference = eos%kg_m3_to_R*spv_ref
554  call calculate_spec_vol_array(ta, sa, pres, spv, 1, 1, eos, spv_reference)
555  else
556  call calculate_spec_vol_array(ta, sa, pres, spv, 1, 1, eos)
557  endif
558  specvol = spv(1)
559 
560  spv_scale = eos%R_to_kg_m3
561  if (present(scale)) spv_scale = spv_scale * scale
562  if (spv_scale /= 1.0) then
563  specvol = spv_scale * specvol
564  endif
565 
566 end subroutine calc_spec_vol_scalar
567 
568 !> Calls the appropriate subroutine to calculate the specific volume of sea water for 1-D array
569 !! inputs, potentially limiting the domain of indices that are worked on.
570 subroutine calc_spec_vol_1d(T, S, pressure, specvol, EOS, dom, spv_ref, scale)
571  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
572  real, dimension(:), intent(in) :: S !< Salinity [ppt]
573  real, dimension(:), intent(in) :: pressure !< Pressure [R L2 T-2 ~> Pa]
574  real, dimension(:), intent(inout) :: specvol !< In situ specific volume [R-1 ~> m3 kg-1]
575  type(EOS_type), pointer :: EOS !< Equation of state structure
576  integer, dimension(2), optional, intent(in) :: dom !< The domain of indices to work on, taking
577  !! into account that arrays start at 1.
578  real, optional, intent(in) :: spv_ref !< A reference specific volume [R-1 ~> m3 kg-1]
579  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale
580  !! output specific volume in combination with
581  !! scaling given by US [various]
582  ! Local variables
583  real, dimension(size(specvol)) :: pres ! Pressure converted to [Pa]
584  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
585  real :: spv_unscale ! A factor to convert specific volume from R-1 to m3 kg-1 [m3 kg-1 R ~> 1]
586  real :: spv_scale ! A factor to convert specific volume from m3 kg-1 to the desired units [kg m-3 R-1 ~> 1]
587  real :: spv_reference ! spv_ref converted to [m3 kg-1]
588  integer :: i, is, ie, npts
589 
590  if (.not.associated(eos)) call mom_error(fatal, &
591  "calc_spec_vol_1d called with an unassociated EOS_type EOS.")
592 
593  if (present(dom)) then
594  is = dom(1) ; ie = dom(2) ; npts = 1 + ie - is
595  else
596  is = 1 ; ie = size(specvol) ; npts = 1 + ie - is
597  endif
598 
599  p_scale = eos%RL2_T2_to_Pa
600  spv_unscale = eos%kg_m3_to_R
601 
602  if ((p_scale == 1.0) .and. (spv_unscale == 1.0)) then
603  call calculate_spec_vol_array(t, s, pressure, specvol, is, npts, eos, spv_ref)
604  elseif (present(spv_ref)) then ! This is the same as above, but with some extra work to rescale variables.
605  do i=is,ie ; pres(i) = p_scale * pressure(i) ; enddo
606  spv_reference = spv_unscale*spv_ref
607  call calculate_spec_vol_array(t, s, pres, specvol, is, npts, eos, spv_reference)
608  else ! There is rescaling of variables, but spv_ref is not present. Passing a 0 value of spv_ref
609  ! changes answers at roundoff for some equations of state, like Wright and UNESCO.
610  do i=is,ie ; pres(i) = p_scale * pressure(i) ; enddo
611  call calculate_spec_vol_array(t, s, pres, specvol, is, npts, eos)
612  endif
613 
614  spv_scale = eos%R_to_kg_m3
615  if (present(scale)) spv_scale = spv_scale * scale
616  if (spv_scale /= 1.0) then ; do i=is,ie
617  specvol(i) = spv_scale * specvol(i)
618  enddo ; endif
619 
620 end subroutine calc_spec_vol_1d
621 
622 
623 !> Calls the appropriate subroutine to calculate the freezing point for scalar inputs.
624 subroutine calculate_tfreeze_scalar(S, pressure, T_fr, EOS, pres_scale)
625  real, intent(in) :: S !< Salinity [ppt]
626  real, intent(in) :: pressure !< Pressure [Pa] or [other]
627  real, intent(out) :: T_fr !< Freezing point potential temperature referenced
628  !! to the surface [degC]
629  type(EOS_type), pointer :: EOS !< Equation of state structure
630  real, optional, intent(in) :: pres_scale !< A multiplicative factor to convert pressure into Pa
631 
632  ! Local variables
633  real :: p_scale ! A factor to convert pressure to units of Pa.
634 
635  if (.not.associated(eos)) call mom_error(fatal, &
636  "calculate_TFreeze_scalar called with an unassociated EOS_type EOS.")
637 
638  p_scale = 1.0 ; if (present(pres_scale)) p_scale = pres_scale
639 
640  select case (eos%form_of_TFreeze)
641  case (tfreeze_linear)
642  call calculate_tfreeze_linear(s, p_scale*pressure, t_fr, eos%TFr_S0_P0, &
643  eos%dTFr_dS, eos%dTFr_dp)
644  case (tfreeze_millero)
645  call calculate_tfreeze_millero(s, p_scale*pressure, t_fr)
646  case (tfreeze_teos10)
647  call calculate_tfreeze_teos10(s, p_scale*pressure, t_fr)
648  case default
649  call mom_error(fatal, "calculate_TFreeze_scalar: form_of_TFreeze is not valid.")
650  end select
651 
652 end subroutine calculate_tfreeze_scalar
653 
654 !> Calls the appropriate subroutine to calculate the freezing point for a 1-D array.
655 subroutine calculate_tfreeze_array(S, pressure, T_fr, start, npts, EOS, pres_scale)
656  real, dimension(:), intent(in) :: S !< Salinity [ppt]
657  real, dimension(:), intent(in) :: pressure !< Pressure [Pa] or [other]
658  real, dimension(:), intent(inout) :: T_fr !< Freezing point potential temperature referenced
659  !! to the surface [degC]
660  integer, intent(in) :: start !< Starting index within the array
661  integer, intent(in) :: npts !< The number of values to calculate
662  type(EOS_type), pointer :: EOS !< Equation of state structure
663  real, optional, intent(in) :: pres_scale !< A multiplicative factor to convert pressure into Pa.
664 
665  ! Local variables
666  real, dimension(size(pressure)) :: pres ! Pressure converted to [Pa]
667  real :: p_scale ! A factor to convert pressure to units of Pa.
668  integer :: j
669 
670  if (.not.associated(eos)) call mom_error(fatal, &
671  "calculate_TFreeze_scalar called with an unassociated EOS_type EOS.")
672 
673  p_scale = 1.0 ; if (present(pres_scale)) p_scale = pres_scale
674 
675  if (p_scale == 1.0) then
676  select case (eos%form_of_TFreeze)
677  case (tfreeze_linear)
678  call calculate_tfreeze_linear(s, pressure, t_fr, start, npts, &
679  eos%TFr_S0_P0, eos%dTFr_dS, eos%dTFr_dp)
680  case (tfreeze_millero)
681  call calculate_tfreeze_millero(s, pressure, t_fr, start, npts)
682  case (tfreeze_teos10)
683  call calculate_tfreeze_teos10(s, pressure, t_fr, start, npts)
684  case default
685  call mom_error(fatal, "calculate_TFreeze_scalar: form_of_TFreeze is not valid.")
686  end select
687  else
688  do j=start,start+npts-1 ; pres(j) = p_scale * pressure(j) ; enddo
689  select case (eos%form_of_TFreeze)
690  case (tfreeze_linear)
691  call calculate_tfreeze_linear(s, pres, t_fr, start, npts, &
692  eos%TFr_S0_P0, eos%dTFr_dS, eos%dTFr_dp)
693  case (tfreeze_millero)
694  call calculate_tfreeze_millero(s, pres, t_fr, start, npts)
695  case (tfreeze_teos10)
696  call calculate_tfreeze_teos10(s, pres, t_fr, start, npts)
697  case default
698  call mom_error(fatal, "calculate_TFreeze_scalar: form_of_TFreeze is not valid.")
699  end select
700  endif
701 
702 end subroutine calculate_tfreeze_array
703 
704 !> Calls the appropriate subroutine to calculate density derivatives for 1-D array inputs.
705 subroutine calculate_density_derivs_array(T, S, pressure, drho_dT, drho_dS, start, npts, EOS, scale)
706  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
707  real, dimension(:), intent(in) :: S !< Salinity [ppt]
708  real, dimension(:), intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
709  real, dimension(:), intent(inout) :: drho_dT !< The partial derivative of density with potential
710  !! temperature [kg m-3 degC-1] or [R degC-1 ~> kg m-3 degC-1]
711  real, dimension(:), intent(inout) :: drho_dS !< The partial derivative of density with salinity,
712  !! in [kg m-3 ppt-1] or [R degC-1 ~> kg m-3 ppt-1]
713  integer, intent(in) :: start !< Starting index within the array
714  integer, intent(in) :: npts !< The number of values to calculate
715  type(EOS_type), pointer :: EOS !< Equation of state structure
716  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
717  !! in combination with scaling given by US [various]
718 
719  ! Local variables
720  integer :: j
721 
722  if (.not.associated(eos)) call mom_error(fatal, &
723  "calculate_density_derivs called with an unassociated EOS_type EOS.")
724 
725  select case (eos%form_of_EOS)
726  case (eos_linear)
727  call calculate_density_derivs_linear(t, s, pressure, drho_dt, drho_ds, eos%Rho_T0_S0, &
728  eos%dRho_dT, eos%dRho_dS, start, npts)
729  case (eos_unesco)
730  call calculate_density_derivs_unesco(t, s, pressure, drho_dt, drho_ds, start, npts)
731  case (eos_wright)
732  call calculate_density_derivs_wright(t, s, pressure, drho_dt, drho_ds, start, npts)
733  case (eos_teos10)
734  call calculate_density_derivs_teos10(t, s, pressure, drho_dt, drho_ds, start, npts)
735  case (eos_nemo)
736  call calculate_density_derivs_nemo(t, s, pressure, drho_dt, drho_ds, start, npts)
737  case default
738  call mom_error(fatal, "calculate_density_derivs_array: EOS%form_of_EOS is not valid.")
739  end select
740 
741  if (present(scale)) then ; if (scale /= 1.0) then ; do j=start,start+npts-1
742  drho_dt(j) = scale * drho_dt(j)
743  drho_ds(j) = scale * drho_ds(j)
744  enddo ; endif ; endif
745 
746 end subroutine calculate_density_derivs_array
747 
748 
749 !> Calls the appropriate subroutine to calculate density derivatives for 1-D array inputs.
750 subroutine calculate_density_derivs_1d(T, S, pressure, drho_dT, drho_dS, EOS, dom, scale)
751  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
752  real, dimension(:), intent(in) :: S !< Salinity [ppt]
753  real, dimension(:), intent(in) :: pressure !< Pressure [R L2 T-2 ~> Pa]
754  real, dimension(:), intent(inout) :: drho_dT !< The partial derivative of density with potential
755  !! temperature [R degC-1 ~> kg m-3 degC-1]
756  real, dimension(:), intent(inout) :: drho_dS !< The partial derivative of density with salinity
757  !! [R degC-1 ~> kg m-3 ppt-1]
758  type(EOS_type), pointer :: EOS !< Equation of state structure
759  integer, dimension(2), optional, intent(in) :: dom !< The domain of indices to work on, taking
760  !! into account that arrays start at 1.
761  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
762  !! in combination with scaling given by US [various]
763  ! Local variables
764  real, dimension(size(drho_dT)) :: pres ! Pressure converted to [Pa]
765  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
766  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
767  integer :: i, is, ie, npts
768 
769  if (.not.associated(eos)) call mom_error(fatal, &
770  "calculate_density_derivs called with an unassociated EOS_type EOS.")
771 
772  if (present(dom)) then
773  is = dom(1) ; ie = dom(2) ; npts = 1 + ie - is
774  else
775  is = 1 ; ie = size(drho_dt) ; npts = 1 + ie - is
776  endif
777 
778  p_scale = eos%RL2_T2_to_Pa
779 
780  if (p_scale == 1.0) then
781  call calculate_density_derivs_array(t, s, pressure, drho_dt, drho_ds, is, npts, eos)
782  else
783  do i=is,ie ; pres(i) = p_scale * pressure(i) ; enddo
784  call calculate_density_derivs_array(t, s, pres, drho_dt, drho_ds, is, npts, eos)
785  endif
786 
787  rho_scale = eos%kg_m3_to_R
788  if (present(scale)) rho_scale = rho_scale * scale
789  if (rho_scale /= 1.0) then ; do i=is,ie
790  drho_dt(i) = rho_scale * drho_dt(i)
791  drho_ds(i) = rho_scale * drho_ds(i)
792  enddo ; endif
793 
794 end subroutine calculate_density_derivs_1d
795 
796 
797 !> Calls the appropriate subroutines to calculate density derivatives by promoting a scalar
798 !! to a one-element array
799 subroutine calculate_density_derivs_scalar(T, S, pressure, drho_dT, drho_dS, EOS, scale)
800  real, intent(in) :: T !< Potential temperature referenced to the surface [degC]
801  real, intent(in) :: S !< Salinity [ppt]
802  real, intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
803  real, intent(out) :: drho_dT !< The partial derivative of density with potential
804  !! temperature [kg m-3 degC-1] or [R degC-1 ~> kg m-3 degC-1]
805  real, intent(out) :: drho_dS !< The partial derivative of density with salinity,
806  !! in [kg m-3 ppt-1] or [R ppt-1 ~> kg m-3 ppt-1]
807  type(EOS_type), pointer :: EOS !< Equation of state structure
808  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
809  !! in combination with scaling given by US [various]
810  ! Local variables
811  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
812  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
813  integer :: j
814 
815  if (.not.associated(eos)) call mom_error(fatal, &
816  "calculate_density_derivs called with an unassociated EOS_type EOS.")
817 
818  p_scale = eos%RL2_T2_to_Pa
819 
820  select case (eos%form_of_EOS)
821  case (eos_linear)
822  call calculate_density_derivs_linear(t, s, p_scale*pressure, drho_dt, drho_ds, &
823  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS)
824  case (eos_wright)
825  call calculate_density_derivs_wright(t, s, p_scale*pressure, drho_dt, drho_ds)
826  case (eos_teos10)
827  call calculate_density_derivs_teos10(t, s, p_scale*pressure, drho_dt, drho_ds)
828  case default
829  call mom_error(fatal, "calculate_density_derivs_scalar: EOS%form_of_EOS is not valid.")
830  end select
831 
832  rho_scale = eos%kg_m3_to_R
833  if (present(scale)) rho_scale = rho_scale * scale
834  if (rho_scale /= 1.0) then
835  drho_dt = rho_scale * drho_dt
836  drho_ds = rho_scale * drho_ds
837  endif
838 
839 end subroutine calculate_density_derivs_scalar
840 
841 !> Calls the appropriate subroutine to calculate density second derivatives for 1-D array inputs.
842 subroutine calculate_density_second_derivs_array(T, S, pressure, drho_dS_dS, drho_dS_dT, drho_dT_dT, &
843  drho_dS_dP, drho_dT_dP, start, npts, EOS, scale)
844  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
845  real, dimension(:), intent(in) :: S !< Salinity [ppt]
846  real, dimension(:), intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
847  real, dimension(:), intent(inout) :: drho_dS_dS !< Partial derivative of beta with respect to S
848  !! [kg m-3 ppt-2] or [R ppt-2 ~> kg m-3 ppt-2]
849  real, dimension(:), intent(inout) :: drho_dS_dT !< Partial derivative of beta with respect to T
850  !! [kg m-3 ppt-1 degC-1] or [R ppt-1 degC-1 ~> kg m-3 ppt-1 degC-1]
851  real, dimension(:), intent(inout) :: drho_dT_dT !< Partial derivative of alpha with respect to T
852  !! [kg m-3 degC-2] or [R degC-2 ~> kg m-3 degC-2]
853  real, dimension(:), intent(inout) :: drho_dS_dP !< Partial derivative of beta with respect to pressure
854  !! [kg m-3 ppt-1 Pa-1] or [R ppt-1 Pa-1 ~> kg m-3 ppt-1 Pa-1]
855  real, dimension(:), intent(inout) :: drho_dT_dP !< Partial derivative of alpha with respect to pressure
856  !! [kg m-3 degC-1 Pa-1] or [R degC-1 Pa-1 ~> kg m-3 degC-1 Pa-1]
857  integer, intent(in) :: start !< Starting index within the array
858  integer, intent(in) :: npts !< The number of values to calculate
859  type(EOS_type), pointer :: EOS !< Equation of state structure
860  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
861  !! in combination with scaling given by US [various]
862  ! Local variables
863  real, dimension(size(pressure)) :: pres ! Pressure converted to [Pa]
864  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
865  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
866  real :: I_p_scale ! The inverse of the factor to convert pressure to units of Pa [R L2 T-2 Pa-1 ~> 1]
867  integer :: j
868 
869  if (.not.associated(eos)) call mom_error(fatal, &
870  "calculate_density_derivs called with an unassociated EOS_type EOS.")
871 
872  p_scale = eos%RL2_T2_to_Pa
873 
874  if (p_scale == 1.0) then
875  select case (eos%form_of_EOS)
876  case (eos_linear)
877  call calculate_density_second_derivs_linear(t, s, pressure, drho_ds_ds, drho_ds_dt, &
878  drho_dt_dt, drho_ds_dp, drho_dt_dp, start, npts)
879  case (eos_wright)
880  call calculate_density_second_derivs_wright(t, s, pressure, drho_ds_ds, drho_ds_dt, &
881  drho_dt_dt, drho_ds_dp, drho_dt_dp, start, npts)
882  case (eos_teos10)
883  call calculate_density_second_derivs_teos10(t, s, pressure, drho_ds_ds, drho_ds_dt, &
884  drho_dt_dt, drho_ds_dp, drho_dt_dp, start, npts)
885  case default
886  call mom_error(fatal, "calculate_density_derivs: EOS%form_of_EOS is not valid.")
887  end select
888  else
889  do j=start,start+npts-1 ; pres(j) = p_scale * pressure(j) ; enddo
890  select case (eos%form_of_EOS)
891  case (eos_linear)
892  call calculate_density_second_derivs_linear(t, s, pres, drho_ds_ds, drho_ds_dt, &
893  drho_dt_dt, drho_ds_dp, drho_dt_dp, start, npts)
894  case (eos_wright)
895  call calculate_density_second_derivs_wright(t, s, pres, drho_ds_ds, drho_ds_dt, &
896  drho_dt_dt, drho_ds_dp, drho_dt_dp, start, npts)
897  case (eos_teos10)
898  call calculate_density_second_derivs_teos10(t, s, pres, drho_ds_ds, drho_ds_dt, &
899  drho_dt_dt, drho_ds_dp, drho_dt_dp, start, npts)
900  case default
901  call mom_error(fatal, "calculate_density_derivs: EOS%form_of_EOS is not valid.")
902  end select
903  endif
904 
905  rho_scale = eos%kg_m3_to_R
906  if (present(scale)) rho_scale = rho_scale * scale
907  if (rho_scale /= 1.0) then ; do j=start,start+npts-1
908  drho_ds_ds(j) = rho_scale * drho_ds_ds(j)
909  drho_ds_dt(j) = rho_scale * drho_ds_dt(j)
910  drho_dt_dt(j) = rho_scale * drho_dt_dt(j)
911  drho_ds_dp(j) = rho_scale * drho_ds_dp(j)
912  drho_dt_dp(j) = rho_scale * drho_dt_dp(j)
913  enddo ; endif
914 
915  if (p_scale /= 1.0) then
916  i_p_scale = 1.0 / p_scale
917  do j=start,start+npts-1
918  drho_ds_dp(j) = i_p_scale * drho_ds_dp(j)
919  drho_dt_dp(j) = i_p_scale * drho_dt_dp(j)
920  enddo
921  endif
922 
923 end subroutine calculate_density_second_derivs_array
924 
925 !> Calls the appropriate subroutine to calculate density second derivatives for scalar nputs.
926 subroutine calculate_density_second_derivs_scalar(T, S, pressure, drho_dS_dS, drho_dS_dT, drho_dT_dT, &
927  drho_dS_dP, drho_dT_dP, EOS, scale)
928  real, intent(in) :: T !< Potential temperature referenced to the surface [degC]
929  real, intent(in) :: S !< Salinity [ppt]
930  real, intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
931  real, intent(out) :: drho_dS_dS !< Partial derivative of beta with respect to S
932  !! [kg m-3 ppt-2] or [R ppt-2 ~> kg m-3 ppt-2]
933  real, intent(out) :: drho_dS_dT !< Partial derivative of beta with respect to T
934  !! [kg m-3 ppt-1 degC-1] or [R ppt-1 degC-1 ~> kg m-3 ppt-1 degC-1]
935  real, intent(out) :: drho_dT_dT !< Partial derivative of alpha with respect to T
936  !! [kg m-3 degC-2] or [R degC-2 ~> kg m-3 degC-2]
937  real, intent(out) :: drho_dS_dP !< Partial derivative of beta with respect to pressure
938  !! [kg m-3 ppt-1 Pa-1] or [R ppt-1 Pa-1 ~> kg m-3 ppt-1 Pa-1]
939  real, intent(out) :: drho_dT_dP !< Partial derivative of alpha with respect to pressure
940  !! [kg m-3 degC-1 Pa-1] or [R degC-1 Pa-1 ~> kg m-3 degC-1 Pa-1]
941  type(EOS_type), pointer :: EOS !< Equation of state structure
942  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale density
943  !! in combination with scaling given by US [various]
944  ! Local variables
945  real :: rho_scale ! A factor to convert density from kg m-3 to the desired units [R m3 kg-1 ~> 1]
946  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
947  real :: I_p_scale ! The inverse of the factor to convert pressure to units of Pa [R L2 T-2 Pa-1 ~> 1]
948 
949  if (.not.associated(eos)) call mom_error(fatal, &
950  "calculate_density_derivs called with an unassociated EOS_type EOS.")
951 
952  p_scale = eos%RL2_T2_to_Pa
953 
954  select case (eos%form_of_EOS)
955  case (eos_linear)
956  call calculate_density_second_derivs_linear(t, s, p_scale*pressure, drho_ds_ds, drho_ds_dt, &
957  drho_dt_dt, drho_ds_dp, drho_dt_dp)
958  case (eos_wright)
959  call calculate_density_second_derivs_wright(t, s, p_scale*pressure, drho_ds_ds, drho_ds_dt, &
960  drho_dt_dt, drho_ds_dp, drho_dt_dp)
961  case (eos_teos10)
962  call calculate_density_second_derivs_teos10(t, s, p_scale*pressure, drho_ds_ds, drho_ds_dt, &
963  drho_dt_dt, drho_ds_dp, drho_dt_dp)
964  case default
965  call mom_error(fatal, "calculate_density_derivs: EOS%form_of_EOS is not valid.")
966  end select
967 
968  rho_scale = eos%kg_m3_to_R
969  if (present(scale)) rho_scale = rho_scale * scale
970  if (rho_scale /= 1.0) then
971  drho_ds_ds = rho_scale * drho_ds_ds
972  drho_ds_dt = rho_scale * drho_ds_dt
973  drho_dt_dt = rho_scale * drho_dt_dt
974  drho_ds_dp = rho_scale * drho_ds_dp
975  drho_dt_dp = rho_scale * drho_dt_dp
976  endif
977 
978  if (p_scale /= 1.0) then
979  i_p_scale = 1.0 / p_scale
980  drho_ds_dp = i_p_scale * drho_ds_dp
981  drho_dt_dp = i_p_scale * drho_dt_dp
982  endif
983 
984 end subroutine calculate_density_second_derivs_scalar
985 
986 !> Calls the appropriate subroutine to calculate specific volume derivatives for an array.
987 subroutine calculate_spec_vol_derivs_array(T, S, pressure, dSV_dT, dSV_dS, start, npts, EOS)
988  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
989  real, dimension(:), intent(in) :: S !< Salinity [ppt]
990  real, dimension(:), intent(in) :: pressure !< Pressure [Pa]
991  real, dimension(:), intent(inout) :: dSV_dT !< The partial derivative of specific volume with potential
992  !! temperature [m3 kg-1 degC-1]
993  real, dimension(:), intent(inout) :: dSV_dS !< The partial derivative of specific volume with salinity
994  !! [m3 kg-1 ppt-1]
995  integer, intent(in) :: start !< Starting index within the array
996  integer, intent(in) :: npts !< The number of values to calculate
997  type(EOS_type), pointer :: EOS !< Equation of state structure
998 
999  ! Local variables
1000  real, dimension(size(T)) :: press ! Pressure converted to [Pa]
1001  real, dimension(size(T)) :: rho ! In situ density [kg m-3]
1002  real, dimension(size(T)) :: dRho_dT ! Derivative of density with temperature [kg m-3 degC-1]
1003  real, dimension(size(T)) :: dRho_dS ! Derivative of density with salinity [kg m-3 ppt-1]
1004  integer :: j
1005 
1006  if (.not.associated(eos)) call mom_error(fatal, &
1007  "calculate_spec_vol_derivs_array called with an unassociated EOS_type EOS.")
1008 
1009  select case (eos%form_of_EOS)
1010  case (eos_linear)
1011  call calculate_specvol_derivs_linear(t, s, pressure, dsv_dt, dsv_ds, start, &
1012  npts, eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS)
1013  case (eos_unesco)
1014  call calculate_density_unesco(t, s, pressure, rho, start, npts)
1015  call calculate_density_derivs_unesco(t, s, pressure, drho_dt, drho_ds, start, npts)
1016  do j=start,start+npts-1
1017  dsv_dt(j) = -drho_dt(j)/(rho(j)**2)
1018  dsv_ds(j) = -drho_ds(j)/(rho(j)**2)
1019  enddo
1020  case (eos_wright)
1021  call calculate_specvol_derivs_wright(t, s, pressure, dsv_dt, dsv_ds, start, npts)
1022  case (eos_teos10)
1023  call calculate_specvol_derivs_teos10(t, s, pressure, dsv_dt, dsv_ds, start, npts)
1024  case (eos_nemo)
1025  call calculate_density_nemo(t, s, pressure, rho, start, npts)
1026  call calculate_density_derivs_nemo(t, s, pressure, drho_dt, drho_ds, start, npts)
1027  do j=start,start+npts-1
1028  dsv_dt(j) = -drho_dt(j)/(rho(j)**2)
1029  dsv_ds(j) = -drho_ds(j)/(rho(j)**2)
1030  enddo
1031  case default
1032  call mom_error(fatal, "calculate_spec_vol_derivs_array: EOS%form_of_EOS is not valid.")
1033  end select
1034 
1035 end subroutine calculate_spec_vol_derivs_array
1036 
1037 !> Calls the appropriate subroutine to calculate specific volume derivatives for 1-d array inputs,
1038 !! potentially limiting the domain of indices that are worked on.
1039 subroutine calc_spec_vol_derivs_1d(T, S, pressure, dSV_dT, dSV_dS, EOS, dom, scale)
1040  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
1041  real, dimension(:), intent(in) :: S !< Salinity [ppt]
1042  real, dimension(:), intent(in) :: pressure !< Pressure [R L2 T-2 ~> Pa]
1043  real, dimension(:), intent(inout) :: dSV_dT !< The partial derivative of specific volume with potential
1044  !! temperature [R-1 degC-1 ~> m3 kg-1 degC-1]
1045  real, dimension(:), intent(inout) :: dSV_dS !< The partial derivative of specific volume with salinity
1046  !! [R-1 ppt-1 ~> m3 kg-1 ppt-1]
1047  type(EOS_type), pointer :: EOS !< Equation of state structure
1048  integer, dimension(2), optional, intent(in) :: dom !< The domain of indices to work on, taking
1049  !! into account that arrays start at 1.
1050  real, optional, intent(in) :: scale !< A multiplicative factor by which to scale specific
1051  !! volume in combination with scaling given by US [various]
1052 
1053  ! Local variables
1054  real, dimension(size(dSV_dT)) :: press ! Pressure converted to [Pa]
1055  real :: spv_scale ! A factor to convert specific volume from m3 kg-1 to the desired units [kg R-1 m-3 ~> 1]
1056  real :: p_scale ! A factor to convert pressure to units of Pa [Pa T2 R-1 L-2 ~> 1]
1057  integer :: i, is, ie, npts
1058 
1059  if (.not.associated(eos)) call mom_error(fatal, &
1060  "calculate_spec_vol_derivs_1d called with an unassociated EOS_type EOS.")
1061 
1062  if (present(dom)) then
1063  is = dom(1) ; ie = dom(2) ; npts = 1 + ie - is
1064  else
1065  is = 1 ; ie = size(dsv_dt) ; npts = 1 + ie - is
1066  endif
1067  p_scale = eos%RL2_T2_to_Pa
1068 
1069  if (p_scale == 1.0) then
1070  call calculate_spec_vol_derivs_array(t, s, pressure, dsv_dt, dsv_ds, is, npts, eos)
1071  else
1072  do i=is,ie ; press(i) = p_scale * pressure(i) ; enddo
1073  call calculate_spec_vol_derivs_array(t, s, press, dsv_dt, dsv_ds, is, npts, eos)
1074  endif
1075 
1076  spv_scale = eos%R_to_kg_m3
1077  if (present(scale)) spv_scale = spv_scale * scale
1078  if (spv_scale /= 1.0) then ; do i=is,ie
1079  dsv_dt(i) = spv_scale * dsv_dt(i)
1080  dsv_ds(i) = spv_scale * dsv_ds(i)
1081  enddo ; endif
1082 
1083 end subroutine calc_spec_vol_derivs_1d
1084 
1085 
1086 !> Calls the appropriate subroutine to calculate the density and compressibility for 1-D array
1087 !! inputs. If US is present, the units of the inputs and outputs are rescaled.
1088 subroutine calculate_compress_array(T, S, press, rho, drho_dp, start, npts, EOS)
1089  real, dimension(:), intent(in) :: T !< Potential temperature referenced to the surface [degC]
1090  real, dimension(:), intent(in) :: S !< Salinity [PSU]
1091  real, dimension(:), intent(in) :: press !< Pressure [Pa] or [R L2 T-2 ~> Pa]
1092  real, dimension(:), intent(inout) :: rho !< In situ density [kg m-3] or [R ~> kg m-3]
1093  real, dimension(:), intent(inout) :: drho_dp !< The partial derivative of density with pressure
1094  !! (also the inverse of the square of sound speed)
1095  !! [s2 m-2] or [T2 L-2]
1096  integer, intent(in) :: start !< Starting index within the array
1097  integer, intent(in) :: npts !< The number of values to calculate
1098  type(EOS_type), pointer :: EOS !< Equation of state structure
1099 
1100  ! Local variables
1101  real, dimension(size(press)) :: pressure ! Pressure converted to [Pa]
1102  integer :: i, is, ie
1103 
1104  if (.not.associated(eos)) call mom_error(fatal, &
1105  "calculate_compress called with an unassociated EOS_type EOS.")
1106 
1107  is = start ; ie = is + npts - 1
1108  do i=is,ie ; pressure(i) = eos%RL2_T2_to_Pa * press(i) ; enddo
1109 
1110  select case (eos%form_of_EOS)
1111  case (eos_linear)
1112  call calculate_compress_linear(t, s, pressure, rho, drho_dp, start, npts, &
1113  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS)
1114  case (eos_unesco)
1115  call calculate_compress_unesco(t, s, pressure, rho, drho_dp, start, npts)
1116  case (eos_wright)
1117  call calculate_compress_wright(t, s, pressure, rho, drho_dp, start, npts)
1118  case (eos_teos10)
1119  call calculate_compress_teos10(t, s, pressure, rho, drho_dp, start, npts)
1120  case (eos_nemo)
1121  call calculate_compress_nemo(t, s, pressure, rho, drho_dp, start, npts)
1122  case default
1123  call mom_error(fatal, "calculate_compress: EOS%form_of_EOS is not valid.")
1124  end select
1125 
1126  if (eos%kg_m3_to_R /= 1.0) then ; do i=is,ie
1127  rho(i) = eos%kg_m3_to_R * rho(i)
1128  enddo ; endif
1129  if (eos%L_T_to_m_s /= 1.0) then ; do i=is,ie
1130  drho_dp(i) = eos%L_T_to_m_s**2 * drho_dp(i)
1131  enddo ; endif
1132 
1133 end subroutine calculate_compress_array
1134 
1135 !> Calculate density and compressibility for a scalar. This just promotes the scalar to an array
1136 !! with a singleton dimension and calls calculate_compress_array. If US is present, the units of
1137 !! the inputs and outputs are rescaled.
1138 subroutine calculate_compress_scalar(T, S, pressure, rho, drho_dp, EOS)
1139  real, intent(in) :: T !< Potential temperature referenced to the surface [degC]
1140  real, intent(in) :: S !< Salinity [ppt]
1141  real, intent(in) :: pressure !< Pressure [Pa] or [R L2 T-2 ~> Pa]
1142  real, intent(out) :: rho !< In situ density [kg m-3] or [R ~> kg m-3]
1143  real, intent(out) :: drho_dp !< The partial derivative of density with pressure (also the
1144  !! inverse of the square of sound speed) [s2 m-2] or [T2 L-2]
1145  type(EOS_type), pointer :: EOS !< Equation of state structure
1146 
1147  ! Local variables
1148  real, dimension(1) :: Ta, Sa, pa, rhoa, drho_dpa
1149 
1150  if (.not.associated(eos)) call mom_error(fatal, &
1151  "calculate_compress called with an unassociated EOS_type EOS.")
1152  ta(1) = t ; sa(1) = s; pa(1) = pressure
1153 
1154  call calculate_compress_array(ta, sa, pa, rhoa, drho_dpa, 1, 1, eos)
1155  rho = rhoa(1) ; drho_dp = drho_dpa(1)
1156 
1157 end subroutine calculate_compress_scalar
1158 
1159 
1160 !> This subroutine returns a two point integer array indicating the domain of i-indices
1161 !! to work on in EOS calls based on information from a hor_index type
1162 function eos_domain(HI, halo) result(EOSdom)
1163  type(hor_index_type), intent(in) :: hi !< The horizontal index structure
1164  integer, optional, intent(in) :: halo !< The halo size to work on; missing is equivalent to 0.
1165  integer, dimension(2) :: eosdom !< The index domain that the EOS will work on, taking into account
1166  !! that the arrays inside the EOS routines will start at 1.
1167 
1168  ! Local variables
1169  integer :: halo_sz
1170 
1171  halo_sz = 0 ; if (present(halo)) halo_sz = halo
1172 
1173  eosdom(1) = hi%isc - (hi%isd-1) - halo_sz
1174  eosdom(2) = hi%iec - (hi%isd-1) + halo_sz
1175 
1176 end function eos_domain
1177 
1178 
1179 !> Calls the appropriate subroutine to calculate analytical and nearly-analytical
1180 !! integrals in pressure across layers of geopotential anomalies, which are
1181 !! required for calculating the finite-volume form pressure accelerations in a
1182 !! non-Boussinesq model. There are essentially no free assumptions, apart from the
1183 !! use of Boole's rule to do the horizontal integrals, and from a truncation in the
1184 !! series for log(1-eps/1+eps) that assumes that |eps| < 0.34.
1185 subroutine analytic_int_specific_vol_dp(T, S, p_t, p_b, alpha_ref, HI, EOS, &
1186  dza, intp_dza, intx_dza, inty_dza, halo_size, &
1187  bathyP, dP_tiny, useMassWghtInterp)
1188  type(hor_index_type), intent(in) :: hi !< The horizontal index structure
1189  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1190  intent(in) :: t !< Potential temperature referenced to the surface [degC]
1191  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1192  intent(in) :: s !< Salinity [ppt]
1193  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1194  intent(in) :: p_t !< Pressure at the top of the layer [R L2 T-2 ~> Pa] or [Pa]
1195  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1196  intent(in) :: p_b !< Pressure at the bottom of the layer [R L2 T-2 ~> Pa] or [Pa]
1197  real, intent(in) :: alpha_ref !< A mean specific volume that is subtracted out
1198  !! to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1]
1199  !! The calculation is mathematically identical with different values of
1200  !! alpha_ref, but this reduces the effects of roundoff.
1201  type(eos_type), pointer :: eos !< Equation of state structure
1202  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1203  intent(inout) :: dza !< The change in the geopotential anomaly across
1204  !! the layer [L2 T-2 ~> m2 s-2] or [m2 s-2]
1205  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1206  optional, intent(inout) :: intp_dza !< The integral in pressure through the layer of the
1207  !! geopotential anomaly relative to the anomaly at the bottom of the
1208  !! layer [R L4 T-4 ~> Pa m2 s-2] or [Pa m2 s-2]
1209  real, dimension(HI%IsdB:HI%IedB,HI%jsd:HI%jed), &
1210  optional, intent(inout) :: intx_dza !< The integral in x of the difference between the
1211  !! geopotential anomaly at the top and bottom of the layer divided by
1212  !! the x grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]
1213  real, dimension(HI%isd:HI%ied,HI%JsdB:HI%JedB), &
1214  optional, intent(inout) :: inty_dza !< The integral in y of the difference between the
1215  !! geopotential anomaly at the top and bottom of the layer divided by
1216  !! the y grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]
1217  integer, optional, intent(in) :: halo_size !< The width of halo points on which to calculate dza.
1218  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1219  optional, intent(in) :: bathyp !< The pressure at the bathymetry [R L2 T-2 ~> Pa] or [Pa]
1220  real, optional, intent(in) :: dp_tiny !< A miniscule pressure change with
1221  !! the same units as p_t [R L2 T-2 ~> Pa] or [Pa]
1222  logical, optional, intent(in) :: usemasswghtinterp !< If true, uses mass weighting
1223  !! to interpolate T/S for top and bottom integrals.
1224  ! Local variables
1225  real :: pres_scale ! A unit conversion factor from the rescaled units of pressure to Pa [Pa T2 R-1 L-2 ~> 1]
1226  real :: sv_scale ! A multiplicative factor by which to scale specific
1227  ! volume from m3 kg-1 to the desired units [kg m-3 R-1 ~> 1]
1228 
1229  if (.not.associated(eos)) call mom_error(fatal, &
1230  "int_specific_vol_dp called with an unassociated EOS_type EOS.")
1231 
1232  ! We should never reach this point with quadrature. EOS_quadrature indicates that numerical
1233  ! integration be used instead of analytic. This is a safety check.
1234  if (eos%EOS_quadrature) call mom_error(fatal, "EOS_quadrature is set!")
1235 
1236  select case (eos%form_of_EOS)
1237  case (eos_linear)
1238  call int_spec_vol_dp_linear(t, s, p_t, p_b, alpha_ref, hi, eos%kg_m3_to_R*eos%Rho_T0_S0, &
1239  eos%kg_m3_to_R*eos%dRho_dT, eos%kg_m3_to_R*eos%dRho_dS, dza, &
1240  intp_dza, intx_dza, inty_dza, halo_size, &
1241  bathyp, dp_tiny, usemasswghtinterp)
1242  case (eos_wright)
1243  call int_spec_vol_dp_wright(t, s, p_t, p_b, alpha_ref, hi, dza, intp_dza, intx_dza, &
1244  inty_dza, halo_size, bathyp, dp_tiny, usemasswghtinterp, &
1245  sv_scale=eos%R_to_kg_m3, pres_scale=eos%RL2_T2_to_Pa)
1246  case default
1247  call mom_error(fatal, "No analytic integration option is available with this EOS!")
1248  end select
1249 
1250 end subroutine analytic_int_specific_vol_dp
1251 
1252 !> This subroutine calculates analytical and nearly-analytical integrals of
1253 !! pressure anomalies across layers, which are required for calculating the
1254 !! finite-volume form pressure accelerations in a Boussinesq model.
1255 subroutine analytic_int_density_dz(T, S, z_t, z_b, rho_ref, rho_0, G_e, HI, EOS, dpa, &
1256  intz_dpa, intx_dpa, inty_dpa, bathyT, dz_neglect, useMassWghtInterp)
1257  type(hor_index_type), intent(in) :: hi !< Ocean horizontal index structure
1258  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1259  intent(in) :: t !< Potential temperature referenced to the surface [degC]
1260  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1261  intent(in) :: s !< Salinity [ppt]
1262  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1263  intent(in) :: z_t !< Height at the top of the layer in depth units [Z ~> m]
1264  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1265  intent(in) :: z_b !< Height at the bottom of the layer [Z ~> m]
1266  real, intent(in) :: rho_ref !< A mean density [R ~> kg m-3] or [kg m-3], that is
1267  !! subtracted out to reduce the magnitude of each of the
1268  !! integrals.
1269  real, intent(in) :: rho_0 !< A density [R ~> kg m-3] or [kg m-3], that is used
1270  !! to calculate the pressure (as p~=-z*rho_0*G_e)
1271  !! used in the equation of state.
1272  real, intent(in) :: g_e !< The Earth's gravitational acceleration
1273  !! [L2 Z-1 T-2 ~> m s-2] or [m2 Z-1 s-2 ~> m s-2]
1274  type(eos_type), pointer :: eos !< Equation of state structure
1275  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1276  intent(inout) :: dpa !< The change in the pressure anomaly
1277  !! across the layer [R L2 T-2 ~> Pa] or [Pa]
1278  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1279  optional, intent(inout) :: intz_dpa !< The integral through the thickness of the
1280  !! layer of the pressure anomaly relative to the
1281  !! anomaly at the top of the layer [R L2 Z T-2 ~> Pa m]
1282  real, dimension(HI%IsdB:HI%IedB,HI%jsd:HI%jed), &
1283  optional, intent(inout) :: intx_dpa !< The integral in x of the difference between
1284  !! the pressure anomaly at the top and bottom of the
1285  !! layer divided by the x grid spacing [R L2 T-2 ~> Pa]
1286  real, dimension(HI%isd:HI%ied,HI%JsdB:HI%JedB), &
1287  optional, intent(inout) :: inty_dpa !< The integral in y of the difference between
1288  !! the pressure anomaly at the top and bottom of the
1289  !! layer divided by the y grid spacing [R L2 T-2 ~> Pa]
1290  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed), &
1291  optional, intent(in) :: bathyt !< The depth of the bathymetry [Z ~> m]
1292  real, optional, intent(in) :: dz_neglect !< A miniscule thickness change [Z ~> m]
1293  logical, optional, intent(in) :: usemasswghtinterp !< If true, uses mass weighting to
1294  !! interpolate T/S for top and bottom integrals.
1295  ! Local variables
1296  real :: rho_scale ! A multiplicative factor by which to scale density from kg m-3 to the
1297  ! desired units [R m3 kg-1 ~> 1]
1298  real :: pres_scale ! A multiplicative factor to convert pressure into Pa [Pa T2 R-1 L-2 ~> 1]
1299 
1300  if (.not.associated(eos)) call mom_error(fatal, &
1301  "int_density_dz called with an unassociated EOS_type EOS.")
1302 
1303  ! We should never reach this point with quadrature. EOS_quadrature indicates that numerical
1304  ! integration be used instead of analytic. This is a safety check.
1305  if (eos%EOS_quadrature) call mom_error(fatal, "EOS_quadrature is set!")
1306 
1307  select case (eos%form_of_EOS)
1308  case (eos_linear)
1309  rho_scale = eos%kg_m3_to_R
1310  if (rho_scale /= 1.0) then
1311  call int_density_dz_linear(t, s, z_t, z_b, rho_ref, rho_0, g_e, hi, &
1312  rho_scale*eos%Rho_T0_S0, rho_scale*eos%dRho_dT, rho_scale*eos%dRho_dS, &
1313  dpa, intz_dpa, intx_dpa, inty_dpa, bathyt, dz_neglect, usemasswghtinterp)
1314  else
1315  call int_density_dz_linear(t, s, z_t, z_b, rho_ref, rho_0, g_e, hi, &
1316  eos%Rho_T0_S0, eos%dRho_dT, eos%dRho_dS, &
1317  dpa, intz_dpa, intx_dpa, inty_dpa, bathyt, dz_neglect, usemasswghtinterp)
1318  endif
1319  case (eos_wright)
1320  rho_scale = eos%kg_m3_to_R
1321  pres_scale = eos%RL2_T2_to_Pa
1322  if ((rho_scale /= 1.0) .or. (pres_scale /= 1.0)) then
1323  call int_density_dz_wright(t, s, z_t, z_b, rho_ref, rho_0, g_e, hi, &
1324  dpa, intz_dpa, intx_dpa, inty_dpa, bathyt, &
1325  dz_neglect, usemasswghtinterp, rho_scale, pres_scale)
1326  else
1327  call int_density_dz_wright(t, s, z_t, z_b, rho_ref, rho_0, g_e, hi, &
1328  dpa, intz_dpa, intx_dpa, inty_dpa, bathyt, &
1329  dz_neglect, usemasswghtinterp)
1330  endif
1331  case default
1332  call mom_error(fatal, "No analytic integration option is available with this EOS!")
1333  end select
1334 
1335 end subroutine analytic_int_density_dz
1336 
1337 !> Returns true if the equation of state is compressible (i.e. has pressure dependence)
1338 logical function query_compressible(EOS)
1339  type(eos_type), pointer :: eos !< Equation of state structure
1340 
1341  if (.not.associated(eos)) call mom_error(fatal, &
1342  "query_compressible called with an unassociated EOS_type EOS.")
1343 
1344  query_compressible = eos%compressible
1345 end function query_compressible
1346 
1347 !> Initializes EOS_type by allocating and reading parameters
1348 subroutine eos_init(param_file, EOS, US)
1349  type(param_file_type), intent(in) :: param_file !< Parameter file structure
1350  type(eos_type), pointer :: eos !< Equation of state structure
1351  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
1352  optional :: us
1353  ! Local variables
1354 #include "version_variable.h"
1355  character(len=40) :: mdl = "MOM_EOS" ! This module's name.
1356  character(len=40) :: tmpstr
1357 
1358  if (.not.associated(eos)) call eos_allocate(eos)
1359 
1360  ! Read all relevant parameters and write them to the model log.
1361  call log_version(param_file, mdl, version, "")
1362 
1363  call get_param(param_file, mdl, "EQN_OF_STATE", tmpstr, &
1364  "EQN_OF_STATE determines which ocean equation of state "//&
1365  "should be used. Currently, the valid choices are "//&
1366  '"LINEAR", "UNESCO", "WRIGHT", "NEMO" and "TEOS10". '//&
1367  "This is only used if USE_EOS is true.", default=eos_default)
1368  select case (uppercase(tmpstr))
1369  case (eos_linear_string)
1370  eos%form_of_EOS = eos_linear
1371  case (eos_unesco_string)
1372  eos%form_of_EOS = eos_unesco
1373  case (eos_wright_string)
1374  eos%form_of_EOS = eos_wright
1375  case (eos_teos10_string)
1376  eos%form_of_EOS = eos_teos10
1377  case (eos_nemo_string)
1378  eos%form_of_EOS = eos_nemo
1379  case default
1380  call mom_error(fatal, "interpret_eos_selection: EQN_OF_STATE "//&
1381  trim(tmpstr) // "in input file is invalid.")
1382  end select
1383  call mom_mesg('interpret_eos_selection: equation of state set to "' // &
1384  trim(tmpstr)//'"', 5)
1385 
1386  if (eos%form_of_EOS == eos_linear) then
1387  eos%Compressible = .false.
1388  call get_param(param_file, mdl, "RHO_T0_S0", eos%Rho_T0_S0, &
1389  "When EQN_OF_STATE="//trim(eos_linear_string)//", "//&
1390  "this is the density at T=0, S=0.", units="kg m-3", &
1391  default=1000.0)
1392  call get_param(param_file, mdl, "DRHO_DT", eos%dRho_dT, &
1393  "When EQN_OF_STATE="//trim(eos_linear_string)//", "//&
1394  "this is the partial derivative of density with "//&
1395  "temperature.", units="kg m-3 K-1", default=-0.2)
1396  call get_param(param_file, mdl, "DRHO_DS", eos%dRho_dS, &
1397  "When EQN_OF_STATE="//trim(eos_linear_string)//", "//&
1398  "this is the partial derivative of density with "//&
1399  "salinity.", units="kg m-3 PSU-1", default=0.8)
1400  endif
1401 
1402  call get_param(param_file, mdl, "EOS_QUADRATURE", eos%EOS_quadrature, &
1403  "If true, always use the generic (quadrature) code "//&
1404  "code for the integrals of density.", default=.false.)
1405 
1406  call get_param(param_file, mdl, "TFREEZE_FORM", tmpstr, &
1407  "TFREEZE_FORM determines which expression should be "//&
1408  "used for the freezing point. Currently, the valid "//&
1409  'choices are "LINEAR", "MILLERO_78", "TEOS10"', &
1410  default=tfreeze_default)
1411  select case (uppercase(tmpstr))
1412  case (tfreeze_linear_string)
1413  eos%form_of_TFreeze = tfreeze_linear
1414  case (tfreeze_millero_string)
1415  eos%form_of_TFreeze = tfreeze_millero
1416  case (tfreeze_teos10_string)
1417  eos%form_of_TFreeze = tfreeze_teos10
1418  case default
1419  call mom_error(fatal, "interpret_eos_selection: TFREEZE_FORM "//&
1420  trim(tmpstr) // "in input file is invalid.")
1421  end select
1422 
1423  if (eos%form_of_TFreeze == tfreeze_linear) then
1424  call get_param(param_file, mdl, "TFREEZE_S0_P0",eos%TFr_S0_P0, &
1425  "When TFREEZE_FORM="//trim(tfreeze_linear_string)//", "//&
1426  "this is the freezing potential temperature at "//&
1427  "S=0, P=0.", units="deg C", default=0.0)
1428  call get_param(param_file, mdl, "DTFREEZE_DS",eos%dTFr_dS, &
1429  "When TFREEZE_FORM="//trim(tfreeze_linear_string)//", "//&
1430  "this is the derivative of the freezing potential "//&
1431  "temperature with salinity.", &
1432  units="deg C PSU-1", default=-0.054)
1433  call get_param(param_file, mdl, "DTFREEZE_DP",eos%dTFr_dP, &
1434  "When TFREEZE_FORM="//trim(tfreeze_linear_string)//", "//&
1435  "this is the derivative of the freezing potential "//&
1436  "temperature with pressure.", &
1437  units="deg C Pa-1", default=0.0)
1438  endif
1439 
1440  if ((eos%form_of_EOS == eos_teos10 .OR. eos%form_of_EOS == eos_nemo) .AND. &
1441  eos%form_of_TFreeze /= tfreeze_teos10) then
1442  call mom_error(fatal, "interpret_eos_selection: EOS_TEOS10 or EOS_NEMO \n" //&
1443  "should only be used along with TFREEZE_FORM = TFREEZE_TEOS10 .")
1444  endif
1445 
1446  ! Unit conversions
1447  eos%m_to_Z = 1. ; if (present(us)) eos%m_to_Z = us%m_to_Z
1448  eos%kg_m3_to_R = 1. ; if (present(us)) eos%kg_m3_to_R = us%kg_m3_to_R
1449  eos%R_to_kg_m3 = 1. ; if (present(us)) eos%R_to_kg_m3 = us%R_to_kg_m3
1450  eos%RL2_T2_to_Pa = 1. ; if (present(us)) eos%RL2_T2_to_Pa = us%RL2_T2_to_Pa
1451  eos%L_T_to_m_s = 1. ; if (present(us)) eos%L_T_to_m_s = us%L_T_to_m_s
1452 
1453 end subroutine eos_init
1454 
1455 !> Manually initialized an EOS type (intended for unit testing of routines which need a specific EOS)
1456 subroutine eos_manual_init(EOS, form_of_EOS, form_of_TFreeze, EOS_quadrature, Compressible, &
1457  Rho_T0_S0, drho_dT, dRho_dS, TFr_S0_P0, dTFr_dS, dTFr_dp)
1458  type(eos_type), pointer :: eos !< Equation of state structure
1459  integer, optional, intent(in) :: form_of_eos !< A coded integer indicating the equation of state to use.
1460  integer, optional, intent(in) :: form_of_tfreeze !< A coded integer indicating the expression for
1461  !! the potential temperature of the freezing point.
1462  logical, optional, intent(in) :: eos_quadrature !< If true, always use the generic (quadrature)
1463  !! code for the integrals of density.
1464  logical, optional, intent(in) :: compressible !< If true, in situ density is a function of pressure.
1465  real , optional, intent(in) :: rho_t0_s0 !< Density at T=0 degC and S=0 ppt [kg m-3]
1466  real , optional, intent(in) :: drho_dt !< Partial derivative of density with temperature
1467  !! in [kg m-3 degC-1]
1468  real , optional, intent(in) :: drho_ds !< Partial derivative of density with salinity
1469  !! in [kg m-3 ppt-1]
1470  real , optional, intent(in) :: tfr_s0_p0 !< The freezing potential temperature at S=0, P=0 [degC]
1471  real , optional, intent(in) :: dtfr_ds !< The derivative of freezing point with salinity
1472  !! in [degC ppt-1]
1473  real , optional, intent(in) :: dtfr_dp !< The derivative of freezing point with pressure
1474  !! in [degC Pa-1]
1475 
1476  if (present(form_of_eos )) eos%form_of_EOS = form_of_eos
1477  if (present(form_of_tfreeze)) eos%form_of_TFreeze = form_of_tfreeze
1478  if (present(eos_quadrature )) eos%EOS_quadrature = eos_quadrature
1479  if (present(compressible )) eos%Compressible = compressible
1480  if (present(rho_t0_s0 )) eos%Rho_T0_S0 = rho_t0_s0
1481  if (present(drho_dt )) eos%drho_dT = drho_dt
1482  if (present(drho_ds )) eos%dRho_dS = drho_ds
1483  if (present(tfr_s0_p0 )) eos%TFr_S0_P0 = tfr_s0_p0
1484  if (present(dtfr_ds )) eos%dTFr_dS = dtfr_ds
1485  if (present(dtfr_dp )) eos%dTFr_dp = dtfr_dp
1486 
1487 end subroutine eos_manual_init
1488 
1489 !> Allocates EOS_type
1490 subroutine eos_allocate(EOS)
1491  type(eos_type), pointer :: eos !< Equation of state structure
1492 
1493  if (.not.associated(eos)) allocate(eos)
1494 end subroutine eos_allocate
1495 
1496 !> Deallocates EOS_type
1497 subroutine eos_end(EOS)
1498  type(eos_type), pointer :: eos !< Equation of state structure
1499 
1500  if (associated(eos)) deallocate(eos)
1501 end subroutine eos_end
1502 
1503 !> Set equation of state structure (EOS) to linear with given coefficients
1504 !!
1505 !! \note This routine is primarily for testing and allows a local copy of the
1506 !! EOS_type (EOS argument) to be set to use the linear equation of state
1507 !! independent from the rest of the model.
1508 subroutine eos_use_linear(Rho_T0_S0, dRho_dT, dRho_dS, EOS, use_quadrature)
1509  real, intent(in) :: rho_t0_s0 !< Density at T=0 degC and S=0 ppt [kg m-3]
1510  real, intent(in) :: drho_dt !< Partial derivative of density with temperature [kg m-3 degC-1]
1511  real, intent(in) :: drho_ds !< Partial derivative of density with salinity [kg m-3 ppt-1]
1512  logical, optional, intent(in) :: use_quadrature !< If true, always use the generic (quadrature)
1513  !! code for the integrals of density.
1514  type(eos_type), pointer :: eos !< Equation of state structure
1515 
1516  if (.not.associated(eos)) call mom_error(fatal, &
1517  "MOM_EOS.F90: EOS_use_linear() called with an unassociated EOS_type EOS.")
1518 
1519  eos%form_of_EOS = eos_linear
1520  eos%Compressible = .false.
1521  eos%Rho_T0_S0 = rho_t0_s0
1522  eos%dRho_dT = drho_dt
1523  eos%dRho_dS = drho_ds
1524  eos%EOS_quadrature = .false.
1525  if (present(use_quadrature)) eos%EOS_quadrature = use_quadrature
1526 
1527 end subroutine eos_use_linear
1528 
1529 
1530 !> Convert T&S to Absolute Salinity and Conservative Temperature if using TEOS10
1531 subroutine convert_temp_salt_for_teos10(T, S, HI, kd, mask_z, EOS)
1532  integer, intent(in) :: kd !< The number of layers to work on
1533  type(hor_index_type), intent(in) :: hi !< The horizontal index structure
1534  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed,kd), &
1535  intent(inout) :: t !< Potential temperature referenced to the surface [degC]
1536  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed,kd), &
1537  intent(inout) :: s !< Salinity [ppt]
1538  real, dimension(HI%isd:HI%ied,HI%jsd:HI%jed,kd), &
1539  intent(in) :: mask_z !< 3d mask regulating which points to convert.
1540  type(eos_type), pointer :: eos !< Equation of state structure
1541 
1542  integer :: i, j, k
1543  real :: gsw_sr_from_sp, gsw_ct_from_pt, gsw_sa_from_sp
1544  real :: p
1545 
1546  if (.not.associated(eos)) call mom_error(fatal, &
1547  "convert_temp_salt_to_TEOS10 called with an unassociated EOS_type EOS.")
1548 
1549  if ((eos%form_of_EOS /= eos_teos10) .and. (eos%form_of_EOS /= eos_nemo)) return
1550 
1551  do k=1,kd ; do j=hi%jsc,hi%jec ; do i=hi%isc,hi%iec
1552  if (mask_z(i,j,k) >= 1.0) then
1553  s(i,j,k) = gsw_sr_from_sp(s(i,j,k))
1554 ! Get absolute salinity from practical salinity, converting pressures from Pascal to dbar.
1555 ! If this option is activated, pressure will need to be added as an argument, and it should be
1556 ! moved out into module that is not shared between components, where the ocean_grid can be used.
1557 ! S(i,j,k) = gsw_sa_from_sp(S(i,j,k),pres(i,j,k)*1.0e-4,G%geoLonT(i,j),G%geoLatT(i,j))
1558  t(i,j,k) = gsw_ct_from_pt(s(i,j,k), t(i,j,k))
1559  endif
1560  enddo ; enddo ; enddo
1561 end subroutine convert_temp_salt_for_teos10
1562 
1563 !> Return value of EOS_quadrature
1564 logical function eos_quadrature(EOS)
1565  type(eos_type), pointer :: eos !< Equation of state structure
1566 
1567  eos_quadrature = eos%EOS_quadrature
1568 
1569 end function eos_quadrature
1570 
1571 !> Extractor routine for the EOS type if the members need to be accessed outside this module
1572 subroutine extract_member_eos(EOS, form_of_EOS, form_of_TFreeze, EOS_quadrature, Compressible, &
1573  Rho_T0_S0, drho_dT, dRho_dS, TFr_S0_P0, dTFr_dS, dTFr_dp)
1574  type(eos_type), pointer :: eos !< Equation of state structure
1575  integer, optional, intent(out) :: form_of_eos !< A coded integer indicating the equation of state to use.
1576  integer, optional, intent(out) :: form_of_tfreeze !< A coded integer indicating the expression for
1577  !! the potential temperature of the freezing point.
1578  logical, optional, intent(out) :: eos_quadrature !< If true, always use the generic (quadrature)
1579  !! code for the integrals of density.
1580  logical, optional, intent(out) :: compressible !< If true, in situ density is a function of pressure.
1581  real , optional, intent(out) :: rho_t0_s0 !< Density at T=0 degC and S=0 ppt [kg m-3]
1582  real , optional, intent(out) :: drho_dt !< Partial derivative of density with temperature
1583  !! in [kg m-3 degC-1]
1584  real , optional, intent(out) :: drho_ds !< Partial derivative of density with salinity
1585  !! in [kg m-3 ppt-1]
1586  real , optional, intent(out) :: tfr_s0_p0 !< The freezing potential temperature at S=0, P=0 [degC]
1587  real , optional, intent(out) :: dtfr_ds !< The derivative of freezing point with salinity
1588  !! [degC PSU-1]
1589  real , optional, intent(out) :: dtfr_dp !< The derivative of freezing point with pressure
1590  !! [degC Pa-1]
1591 
1592  if (present(form_of_eos )) form_of_eos = eos%form_of_EOS
1593  if (present(form_of_tfreeze)) form_of_tfreeze = eos%form_of_TFreeze
1594  if (present(eos_quadrature )) eos_quadrature = eos%EOS_quadrature
1595  if (present(compressible )) compressible = eos%Compressible
1596  if (present(rho_t0_s0 )) rho_t0_s0 = eos%Rho_T0_S0
1597  if (present(drho_dt )) drho_dt = eos%drho_dT
1598  if (present(drho_ds )) drho_ds = eos%dRho_dS
1599  if (present(tfr_s0_p0 )) tfr_s0_p0 = eos%TFr_S0_P0
1600  if (present(dtfr_ds )) dtfr_ds = eos%dTFr_dS
1601  if (present(dtfr_dp )) dtfr_dp = eos%dTFr_dp
1602 
1603 end subroutine extract_member_eos
1604 
1605 end module mom_eos
1606 
1607 !> \namespace mom_eos
1608 !!
1609 !! The MOM_EOS module is a wrapper for various equations of state (e.g. Linear,
1610 !! Wright, UNESCO) and provides a uniform interface to the rest of the model
1611 !! independent of which equation of state is being used.
mom_eos_unesco::calculate_density_unesco
Compute the in situ density of sea water (in [kg m-3]), or its anomaly with respect to a reference de...
Definition: MOM_EOS_UNESCO.F90:22
mom_tfreeze::calculate_tfreeze_teos10
Compute the freezing point conservative temperature [degC] from absolute salinity [g/kg] and pressure...
Definition: MOM_TFreeze.F90:33
mom_eos_linear
A simple linear equation of state for sea water with constant coefficients.
Definition: MOM_EOS_linear.F90:2
mom_eos::calculate_specific_vol_derivs
Calculate the derivatives of specific volume with temperature and salinity from T, S, and P.
Definition: MOM_EOS.F90:87
mom_eos_wright::calculate_density_second_derivs_wright
For a given thermodynamic state, return the second derivatives of density with various combinations o...
Definition: MOM_EOS_Wright.F90:50
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:68
mom_eos::calculate_tfreeze
Calculates the freezing point of sea water from T, S and P.
Definition: MOM_EOS.F90:98
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_eos_wright::calculate_density_derivs_wright
For a given thermodynamic state, return the derivatives of density with temperature and salinity...
Definition: MOM_EOS_Wright.F90:44
mom_eos_wright
The equation of state using the Wright 1997 expressions.
Definition: MOM_EOS_Wright.F90:2
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_eos_teos10::calculate_density_second_derivs_teos10
For a given thermodynamic state, return the second derivatives of density with various combinations o...
Definition: MOM_EOS_TEOS10.F90:47
mom_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_eos_unesco::calculate_spec_vol_unesco
Compute the in situ specific volume of sea water (in [m3 kg-1]), or an anomaly with respect to a refe...
Definition: MOM_EOS_UNESCO.F90:29
mom_eos_wright::calculate_density_wright
Compute the in situ density of sea water (in [kg m-3]), or its anomaly with respect to a reference de...
Definition: MOM_EOS_Wright.F90:32
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains. ...
Definition: MOM_hor_index.F90:16
mom_tfreeze
Freezing point expressions.
Definition: MOM_TFreeze.F90:2
mom_eos_wright::calculate_spec_vol_wright
Compute the in situ specific volume of sea water (in [m3 kg-1]), or an anomaly with respect to a refe...
Definition: MOM_EOS_Wright.F90:39
mom_eos_teos10::calculate_density_teos10
Compute the in situ density of sea water ([kg m-3]), or its anomaly with respect to a reference densi...
Definition: MOM_EOS_TEOS10.F90:28
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_eos::calculate_spec_vol
Calculates specific volume of sea water from T, S and P.
Definition: MOM_EOS.F90:75
mom_eos::calculate_compress
Calculates the compressibility of water from T, S, and P.
Definition: MOM_EOS.F90:103
mom_eos_linear::calculate_density_derivs_linear
For a given thermodynamic state, return the derivatives of density with temperature and salinity usin...
Definition: MOM_EOS_linear.F90:40
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_eos_linear::calculate_spec_vol_linear
Compute the specific volume of sea water (in m^3/kg), or its anomaly from a reference value...
Definition: MOM_EOS_linear.F90:34
mom_eos_teos10
The equation of state using the TEOS10 expressions.
Definition: MOM_EOS_TEOS10.F90:2
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_eos_nemo
The equation of state using the expressions of Roquet et al. that are used in NEMO.
Definition: MOM_EOS_NEMO.F90:2
mom_eos_teos10::calculate_density_derivs_teos10
For a given thermodynamic state, return the derivatives of density with conservative temperature and ...
Definition: MOM_EOS_TEOS10.F90:41
mom_eos::calculate_density_derivs
Calculate the derivatives of density with temperature and salinity from T, S, and P...
Definition: MOM_EOS.F90:81
mom_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_eos_linear::calculate_density_linear
Compute the density of sea water (in kg/m^3), or its anomaly from a reference density, using a simple linear equation of state from salinity (in psu), potential temperature (in deg C) and pressure [Pa].
Definition: MOM_EOS_linear.F90:27
mom_eos_nemo::calculate_density_derivs_nemo
For a given thermodynamic state, return the derivatives of density with conservative temperature and ...
Definition: MOM_EOS_NEMO.F90:34
mom_tfreeze::calculate_tfreeze_millero
Compute the freezing point potential temperature [degC] from salinity [PSU] and pressure [Pa] using t...
Definition: MOM_TFreeze.F90:27
mom_eos_linear::calculate_density_second_derivs_linear
For a given thermodynamic state, return the second derivatives of density with various combinations o...
Definition: MOM_EOS_linear.F90:47
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_eos_unesco
The equation of state using the Jackett and McDougall fits to the UNESCO EOS.
Definition: MOM_EOS_UNESCO.F90:2
mom_string_functions
Handy functions for manipulating strings.
Definition: MOM_string_functions.F90:2
mom_tfreeze::calculate_tfreeze_linear
Compute the freezing point potential temperature [degC] from salinity [ppt] and pressure [Pa] using a...
Definition: MOM_TFreeze.F90:18
mom_eos_nemo::calculate_density_nemo
Compute the in situ density of sea water ([kg m-3]), or its anomaly with respect to a reference densi...
Definition: MOM_EOS_NEMO.F90:28
mom_eos_teos10::calculate_spec_vol_teos10
Compute the in situ specific volume of sea water (in [m3 kg-1]), or an anomaly with respect to a refe...
Definition: MOM_EOS_TEOS10.F90:35
mom_eos::calculate_density_second_derivs
Calculates the second derivatives of density with various combinations of temperature, salinity, and pressure from T, S and P.
Definition: MOM_EOS.F90:93
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_eos::eos_type
A control structure for the equation of state.
Definition: MOM_EOS.F90:108