!        Generated by TAPENADE     (INRIA, Ecuador team)
!  tapenade 3.x
!
!  Differentiation of rhow in forward (tangent) mode:
!   variations   of useful results: rhow
!   with respect to varying inputs: p t
! MIT License
!
! Copyright (c) 2020 SHEMAT-Suite
!
! Permission is hereby granted, free of charge, to any person obtaining a copy
! of this software and associated documentation files (the "Software"), to deal
! in the Software without restriction, including without limitation the rights
! to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
! copies of the Software, and to permit persons to whom the Software is
! furnished to do so, subject to the following conditions:
!
! The above copyright notice and this permission notice shall be included in all
! copies or substantial portions of the Software.
!
! THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
! IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
! FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
! AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
! LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
! OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
! SOFTWARE.
!> @brief rhow calculates the density of pure water [kg/m3]
!> @param[in] p pressure [MPa]
!> @param[in] t temperature [degC]
!> @return rhow density of pure water [kg/m3]
!> @details
!> rhow(i,j,k,ismpl) calculates the density in (in kg/m^3) of pure
!> water, given temperature (t, in degC), and pressure (p,in Pa)\n \n
!>
!> Source: \n
!>
!> Magri, F. (2005), Mechanismus und Fluiddynamik der
!> Salzwasserzirkulation im Norddeutschen Becken: Ergebnisse
!> thermohaliner numerischer Simulationen, Doktorarbeit,
!> Geoforschungszentrum Potsdam \n \n
!>
!> Specifically: Equation (2.12) and Table 2-1
DOUBLE PRECISION FUNCTION g_RHOW(p, g_p, t, g_t, rhow)
  IMPLICIT NONE
! Local temperature (degC)
  DOUBLE PRECISION, INTENT(IN) :: t
  DOUBLE PRECISION, INTENT(IN) :: g_t
! Local pressure (MPa)
  DOUBLE PRECISION, INTENT(IN) :: p
  DOUBLE PRECISION, INTENT(IN) :: g_p
! Pressure in kPa
  DOUBLE PRECISION :: pk
  DOUBLE PRECISION :: g_pk
! Monomials in pressure and temperature
  DOUBLE PRECISION :: pk2, t2, t4, t6
  DOUBLE PRECISION :: g_pk2, g_t2, g_t4, g_t6
! Coefficients from Magri2005 Equation (2.12)
  DOUBLE PRECISION :: a, b, c, d, e, f, g
  DOUBLE PRECISION :: g_a, g_b, g_c, g_d, g_e, g_f, g_g
! Coefficients from Magri2005 Table 2-1
  DOUBLE PRECISION, PARAMETER :: a0=9.99792877961606d02
  DOUBLE PRECISION, PARAMETER :: a1=5.07605113140940d-04
  DOUBLE PRECISION, PARAMETER :: a2=-5.28425478164183d-10
  DOUBLE PRECISION, PARAMETER :: b0=5.13864847162196d-02
  DOUBLE PRECISION, PARAMETER :: b1=-3.61991396354483d-06
  DOUBLE PRECISION, PARAMETER :: b2=7.97204102509724d-12
  DOUBLE PRECISION, PARAMETER :: c0=-7.53557031774437d-03
  DOUBLE PRECISION, PARAMETER :: c1=6.32712093275576d-08
  DOUBLE PRECISION, PARAMETER :: c2=-1.66203631393248d-13
  DOUBLE PRECISION, PARAMETER :: d0=4.60380647957350d-05
  DOUBLE PRECISION, PARAMETER :: d1=-5.61299059722121d-10
  DOUBLE PRECISION, PARAMETER :: d2=1.80924436489400d-15
  DOUBLE PRECISION, PARAMETER :: e0=-2.26651454175013d-07
  DOUBLE PRECISION, PARAMETER :: e1=3.36874416675978d-12
  DOUBLE PRECISION, PARAMETER :: e2=-1.30352149261326d-17
  DOUBLE PRECISION, PARAMETER :: f0=6.14889851856743d-10
  DOUBLE PRECISION, PARAMETER :: f1=-1.06165223196756d-14
  DOUBLE PRECISION, PARAMETER :: f2=4.75014903737416d-20
  DOUBLE PRECISION, PARAMETER :: g0=-7.39221950969522d-13
  DOUBLE PRECISION, PARAMETER :: g1=1.42790422913922d-17
  DOUBLE PRECISION, PARAMETER :: g2=-7.13130230531541d-23
  DOUBLE PRECISION :: rhow
! Unit [kPa] is needed after Magri(2005); the forwarded p is in [MPa]
  g_pk = 1000*g_p
  pk = p*1000
! Monomials in pressure and temperature
  g_pk2 = 2*pk*g_pk
  pk2 = pk*pk
  g_t2 = 2*t*g_t
  t2 = t*t
  g_t4 = 2*t2*g_t2
  t4 = t2*t2
  g_t6 = t2*g_t4 + t4*g_t2
  t6 = t4*t2
! Coefficients from Magri2005, Equation (2.12)
  g_a = a1*g_pk + a2*g_pk2
  a = a0 + a1*pk + a2*pk2
  g_b = b1*g_pk + b2*g_pk2
  b = b0 + b1*pk + b2*pk2
  g_c = c1*g_pk + c2*g_pk2
  c = c0 + c1*pk + c2*pk2
  g_d = d1*g_pk + d2*g_pk2
  d = d0 + d1*pk + d2*pk2
  g_e = e1*g_pk + e2*g_pk2
  e = e0 + e1*pk + e2*pk2
  g_f = f1*g_pk + f2*g_pk2
  f = f0 + f1*pk + f2*pk2
  g_g = g1*g_pk + g2*g_pk2
  g = g0 + g1*pk + g2*pk2
! Water density [kg/m3]
  g_rhow = g_a + t*g_b + (b+d*t2+f*t4)*g_t + t2*g_c + c*g_t2&
&   + t*(t2*g_d+d*g_t2) + t4*g_e + e*g_t4 + t*(t4*g_f+f*g_t4&
&   ) + t6*g_g + g*g_t6
  rhow = a + b*t + c*t2 + d*t2*t + e*t4 + f*t4*t + g*t6
  RETURN
END FUNCTION g_RHOW

