MOM6
mom_eos_nemo::calculate_density_derivs_nemo Interface Reference

Detailed Description

For a given thermodynamic state, return the derivatives of density with conservative temperature and absolute salinity, the expressions derived for use with NEMO.

Definition at line 34 of file MOM_EOS_NEMO.F90.

Private functions

subroutine calculate_density_derivs_scalar_nemo (T, S, pressure, drho_dt, drho_ds)
 Wrapper to calculate_density_derivs_array for scalar inputs. More...
 
subroutine calculate_density_derivs_array_nemo (T, S, pressure, drho_dT, drho_dS, start, npts)
 For a given thermodynamic state, calculate the derivatives of density with conservative temperature and absolute salinity, using the expressions derived for use with NEMO. More...
 

Functions and subroutines

◆ calculate_density_derivs_array_nemo()

subroutine mom_eos_nemo::calculate_density_derivs_nemo::calculate_density_derivs_array_nemo ( 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, calculate the derivatives of density with conservative temperature and absolute salinity, using the expressions derived for use with NEMO.

Parameters
[in]tConservative temperature [degC].
[in]sAbsolute salinity [g kg-1].
[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 ppt-1].
[in]startThe starting point in the arrays.
[in]nptsThe number of values to calculate.

Definition at line 268 of file MOM_EOS_NEMO.F90.

268  real, intent(in), dimension(:) :: T !< Conservative temperature [degC].
269  real, intent(in), dimension(:) :: S !< Absolute salinity [g kg-1].
270  real, intent(in), dimension(:) :: pressure !< pressure [Pa].
271  real, intent(out), dimension(:) :: drho_dT !< The partial derivative of density with potential
272  !! temperature [kg m-3 degC-1].
273  real, intent(out), dimension(:) :: drho_dS !< The partial derivative of density with salinity,
274  !! in [kg m-3 ppt-1].
275  integer, intent(in) :: start !< The starting point in the arrays.
276  integer, intent(in) :: npts !< The number of values to calculate.
277 
278  ! Local variables
279  real :: zp,zt , zh , zs , zr0, zn , zn0, zn1, zn2, zn3
280  integer :: j
281 
282  do j=start,start+npts-1
283  !Conversions
284  zs = s(j) !gsw_sr_from_sp(S(j)) !Convert practical salinity to absolute salinity
285  zt = t(j) !gsw_ct_from_pt(S(j),T(j)) !Convert potantial temp to conservative temp
286  zp = pressure(j)* pa2db !Convert pressure from Pascal to decibar
287 
288  !The following algorithm was provided by Roquet in a private communication.
289  !It is not necessarily the algorithm used in NEMO ocean!
290  zp = zp * r1_p0 ! pressure (first converted to decibar)
291  zt = zt * r1_t0 ! temperature
292  zs = sqrt( abs( zs + rdeltas ) * r1_s0 ) ! square root salinity
293  !
294  ! alpha
295  zn3 = alp003
296  !
297  zn2 = alp012*zt + alp102*zs+alp002
298  !
299  zn1 = ((alp031*zt &
300  & + alp121*zs+alp021)*zt &
301  & + (alp211*zs+alp111)*zs+alp011)*zt &
302  & + ((alp301*zs+alp201)*zs+alp101)*zs+alp001
303  !
304  zn0 = ((((alp050*zt &
305  & + alp140*zs+alp040)*zt &
306  & + (alp230*zs+alp130)*zs+alp030)*zt &
307  & + ((alp320*zs+alp220)*zs+alp120)*zs+alp020)*zt &
308  & + (((alp410*zs+alp310)*zs+alp210)*zs+alp110)*zs+alp010)*zt &
309  & + ((((alp500*zs+alp400)*zs+alp300)*zs+alp200)*zs+alp100)*zs+alp000
310  !
311  zn = ( ( zn3 * zp + zn2 ) * zp + zn1 ) * zp + zn0
312  !
313  drho_dt(j) = -zn
314  !
315  ! beta
316  !
317  zn3 = bet003
318  !
319  zn2 = bet012*zt + bet102*zs+bet002
320  !
321  zn1 = ((bet031*zt &
322  & + bet121*zs+bet021)*zt &
323  & + (bet211*zs+bet111)*zs+bet011)*zt &
324  & + ((bet301*zs+bet201)*zs+bet101)*zs+bet001
325  !
326  zn0 = ((((bet050*zt &
327  & + bet140*zs+bet040)*zt &
328  & + (bet230*zs+bet130)*zs+bet030)*zt &
329  & + ((bet320*zs+bet220)*zs+bet120)*zs+bet020)*zt &
330  & + (((bet410*zs+bet310)*zs+bet210)*zs+bet110)*zs+bet010)*zt &
331  & + ((((bet500*zs+bet400)*zs+bet300)*zs+bet200)*zs+bet100)*zs+bet000
332  !
333  zn = ( ( zn3 * zp + zn2 ) * zp + zn1 ) * zp + zn0
334  !
335  drho_ds(j) = zn / zs
336  enddo
337 

◆ calculate_density_derivs_scalar_nemo()

subroutine mom_eos_nemo::calculate_density_derivs_nemo::calculate_density_derivs_scalar_nemo ( real, intent(in)  T,
real, intent(in)  S,
real, intent(in)  pressure,
real, intent(out)  drho_dt,
real, intent(out)  drho_ds 
)
private

Wrapper to calculate_density_derivs_array for scalar inputs.

Parameters
[in]tPotential temperature relative to the surface [degC].
[in]sSalinity [g kg-1].
[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 ppt-1].

Definition at line 342 of file MOM_EOS_NEMO.F90.

342  real, intent(in) :: T !< Potential temperature relative to the surface [degC].
343  real, intent(in) :: S !< Salinity [g kg-1].
344  real, intent(in) :: pressure !< Pressure [Pa].
345  real, intent(out) :: drho_dT !< The partial derivative of density with potential
346  !! temperature [kg m-3 degC-1].
347  real, intent(out) :: drho_dS !< The partial derivative of density with salinity,
348  !! in [kg m-3 ppt-1].
349  ! Local variables
350  real :: al0, p0, lambda
351  integer :: j
352  real, dimension(1) :: T0, S0, pressure0
353  real, dimension(1) :: drdt0, drds0
354 
355  t0(1) = t
356  s0(1) = s
357  pressure0(1) = pressure
358 
359  call calculate_density_derivs_array_nemo(t0, s0, pressure0, drdt0, drds0, 1, 1)
360  drho_dt = drdt0(1)
361  drho_ds = drds0(1)

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