MOM6
mom_eos_wright::calculate_density_derivs_wright Interface Reference

Detailed Description

For a given thermodynamic state, return the derivatives of density with temperature and salinity.

Definition at line 44 of file MOM_EOS_Wright.F90.

Private functions

subroutine calculate_density_derivs_scalar_wright (T, S, pressure, drho_dT, drho_dS)
 The scalar version of calculate_density_derivs which promotes scalar inputs to a 1-element array and then demotes the output back to a scalar. More...
 
subroutine calculate_density_derivs_array_wright (T, S, pressure, drho_dT, drho_dS, start, npts)
 For a given thermodynamic state, return the thermal/haline expansion coefficients. More...
 

Functions and subroutines

◆ calculate_density_derivs_array_wright()

subroutine mom_eos_wright::calculate_density_derivs_wright::calculate_density_derivs_array_wright ( real, dimension(:), intent(in)  T,
real, dimension(:), intent(in)  S,
real, dimension(:), intent(in)  pressure,
real, dimension(:), intent(out)  drho_dT,
real, dimension(:), intent(out)  drho_dS,
integer, intent(in)  start,
integer, intent(in)  npts 
)
private

For a given thermodynamic state, return the thermal/haline expansion coefficients.

Parameters
[in]tPotential temperature relative to the surface [degC].
[in]sSalinity [PSU].
[in]pressurepressure [Pa].
[out]drho_dtThe partial derivative of density with potential temperature [kg m-3 degC-1].
[out]drho_dsThe partial derivative of density with salinity, in [kg m-3 PSU-1].
[in]startThe starting point in the arrays.
[in]nptsThe number of values to calculate.

Definition at line 200 of file MOM_EOS_Wright.F90.

200  real, intent(in), dimension(:) :: T !< Potential temperature relative to the
201  !! surface [degC].
202  real, intent(in), dimension(:) :: S !< Salinity [PSU].
203  real, intent(in), dimension(:) :: pressure !< pressure [Pa].
204  real, intent(out), dimension(:) :: drho_dT !< The partial derivative of density with potential
205  !! temperature [kg m-3 degC-1].
206  real, intent(out), dimension(:) :: drho_dS !< The partial derivative of density with salinity,
207  !! in [kg m-3 PSU-1].
208  integer, intent(in) :: start !< The starting point in the arrays.
209  integer, intent(in) :: npts !< The number of values to calculate.
210 
211  ! Local variables
212  real :: al0, p0, lambda, I_denom2
213  integer :: j
214 
215  do j=start,start+npts-1
216  al0 = (a0 + a1*t(j)) + a2*s(j)
217  p0 = (b0 + b4*s(j)) + t(j) * (b1 + t(j)*((b2 + b3*t(j))) + b5*s(j))
218  lambda = (c0 +c4*s(j)) + t(j) * (c1 + t(j)*((c2 + c3*t(j))) + c5*s(j))
219 
220  i_denom2 = 1.0 / (lambda + al0*(pressure(j) + p0))
221  i_denom2 = i_denom2 *i_denom2
222  drho_dt(j) = i_denom2 * &
223  (lambda* (b1 + t(j)*(2.0*b2 + 3.0*b3*t(j)) + b5*s(j)) - &
224  (pressure(j)+p0) * ( (pressure(j)+p0)*a1 + &
225  (c1 + t(j)*(c2*2.0 + c3*3.0*t(j)) + c5*s(j)) ))
226  drho_ds(j) = i_denom2 * (lambda* (b4 + b5*t(j)) - &
227  (pressure(j)+p0) * ( (pressure(j)+p0)*a2 + (c4 + c5*t(j)) ))
228  enddo
229 

◆ calculate_density_derivs_scalar_wright()

subroutine mom_eos_wright::calculate_density_derivs_wright::calculate_density_derivs_scalar_wright ( real, intent(in)  T,
real, intent(in)  S,
real, intent(in)  pressure,
real, intent(out)  drho_dT,
real, intent(out)  drho_dS 
)
private

The scalar version of calculate_density_derivs which promotes scalar inputs to a 1-element array and then demotes the output back to a scalar.

Parameters
[in]tPotential temperature relative to the surface [degC].
[in]sSalinity [PSU].
[in]pressurepressure [Pa].
[out]drho_dtThe partial derivative of density with potential temperature [kg m-3 degC-1].
[out]drho_dsThe partial derivative of density with salinity, in [kg m-3 PSU-1].

Definition at line 235 of file MOM_EOS_Wright.F90.

235  real, intent(in) :: T !< Potential temperature relative to the surface [degC].
236  real, intent(in) :: S !< Salinity [PSU].
237  real, intent(in) :: pressure !< pressure [Pa].
238  real, intent(out) :: drho_dT !< The partial derivative of density with potential
239  !! temperature [kg m-3 degC-1].
240  real, intent(out) :: drho_dS !< The partial derivative of density with salinity,
241  !! in [kg m-3 PSU-1].
242 
243  ! Local variables needed to promote the input/output scalars to 1-element arrays
244  real, dimension(1) :: T0, S0, P0
245  real, dimension(1) :: drdt0, drds0
246 
247  t0(1) = t
248  s0(1) = s
249  p0(1) = pressure
250  call calculate_density_derivs_array_wright(t0, s0, p0, drdt0, drds0, 1, 1)
251  drho_dt = drdt0(1)
252  drho_ds = drds0(1)
253 

The documentation for this interface was generated from the following file: