For a given thermodynamic state, return the derivatives of density with temperature and salinity.
Definition at line 44 of file MOM_EOS_Wright.F90.
|
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...
|
|
◆ 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] | t | Potential temperature relative to the surface [degC]. |
[in] | s | Salinity [PSU]. |
[in] | pressure | pressure [Pa]. |
[out] | drho_dt | The partial derivative of density with potential temperature [kg m-3 degC-1]. |
[out] | drho_ds | The partial derivative of density with salinity, in [kg m-3 PSU-1]. |
[in] | start | The starting point in the arrays. |
[in] | npts | The number of values to calculate. |
Definition at line 200 of file MOM_EOS_Wright.F90.
200 real,
intent(in),
dimension(:) :: T
202 real,
intent(in),
dimension(:) :: S
203 real,
intent(in),
dimension(:) :: pressure
204 real,
intent(out),
dimension(:) :: drho_dT
206 real,
intent(out),
dimension(:) :: drho_dS
208 integer,
intent(in) :: start
209 integer,
intent(in) :: npts
212 real :: al0, p0, lambda, I_denom2
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))
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)) ))
◆ 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] | t | Potential temperature relative to the surface [degC]. |
[in] | s | Salinity [PSU]. |
[in] | pressure | pressure [Pa]. |
[out] | drho_dt | The partial derivative of density with potential temperature [kg m-3 degC-1]. |
[out] | drho_ds | The 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
236 real,
intent(in) :: S
237 real,
intent(in) :: pressure
238 real,
intent(out) :: drho_dT
240 real,
intent(out) :: drho_dS
244 real,
dimension(1) :: T0, S0, P0
245 real,
dimension(1) :: drdt0, drds0
250 call calculate_density_derivs_array_wright(t0, s0, p0, drdt0, drds0, 1, 1)
The documentation for this interface was generated from the following file:
- /glade/work/altuntas/cesm.sandboxes/cesm2_2_alpha_X_mom/components/mom/MOM6/src/equation_of_state/MOM_EOS_Wright.F90