def dynamic_enthalpy_CT_exact(SA, CT, p): r"""Calculates the dynamic enthalpy of seawater from Absolute Salinity and Conservative Temperature and pressure. Dynamic enthalpy is defined as enthalpy minus potential enthalpy (Young, 2010). Parameters ---------- SA : array_like Absolute Salinity [g/kg] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p : array_like sea pressure [dbar] Returns ------- dynamic_enthalpy_CT_exact : array_like dynamic enthalpy [J/kg] See Also -------- TODO Notes ----- This function uses the full Gibbs function. There is an alternative to calling this function, namely dynamic_enthalpy(SA, CT, p), which uses the computationally efficient 48-term expression for density in terms of SA, CT and p (McDougall et al., 2011). Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See apendix A.30. .. [2] McDougall T.J., P.M. Barker, R. Feistel and D.R. Jackett, 2011: A computationally efficient 48-term expression for the density of seawater in terms of Conservative Temperature, and related properties of seawater. .. [3] Young, W.R., 2010: Dynamic enthalpy, Conservative Temperature, and the seawater Boussinesq approximation. Journal of Physical Oceanography, 40, 394-400. Modifications: 2011-04-05. Trevor McDougall and Paul Barker. """ t = t_from_CT(SA, CT, p) return enthalpy_t_exact(SA, t, p) - cp0 * CT
def rho_CT_exact(SA, CT, p): r"""Calculates in-situ density from Absolute Salinity and Conservative Temperature. Parameters ---------- SA : array_like Absolute Salinity [g/kg] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p : array_like sea pressure [dbar] Returns ------- rho_CT_exact : array_like in-situ density [kg/m**3] See Also -------- TODO Notes ----- The potential density with respect to reference pressure, p_ref, is obtained by calling this function with the pressure argument being p_ref (i.e. "rho_CT_exact(SA, CT, p_ref)"). This function uses the full Gibbs function. There is an alternative to calling this function, namely rho_CT(SA, CT, p), which uses the computationally efficient 48-term expression for density in terms of SA, CT and p (McDougall et al., 2011). Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See Eqn. (2.8.2). .. [2] McDougall T.J., P.M. Barker, R. Feistel and D.R. Jackett, 2011: A computationally efficient 48-term expression for the density of seawater in terms of Conservative Temperature, and related properties of seawater. Modifications: 2011-04-03. Trevor McDougall and Paul Barker. """ t = t_from_CT(SA, CT, p) return rho_t_exact(SA, t, p)
def alpha_CT_exact(SA, CT, p): r"""Calculates the thermal expansion coefficient of seawater with respect to Conservative Temperature from Absolute Salinity and Conservative Temperature. Parameters ---------- SA : array_like Absolute salinity [g kg :sup:`-1`] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p : array_like pressure [dbar] Returns ------- alpha_CT_exact : array_like thermal expansion coefficient [K :sup:`-1`] with respect to Conservative Temperature See Also -------- TODO Notes ----- This function uses the full Gibbs function. There is an alternative to calling this function, namely alpha_wrt_CT(SA, CT, p) which uses the computationally efficient 48-term expression for density in terms of SA, CT and p (McDougall et al., (2011)). Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See Eqn. (2.18.3). .. [2] McDougall T.J., P.M. Barker, R. Feistel and D.R. Jackett, 2011: A computationally efficient 48-term expression for the density of seawater in terms of Conservative Temperature, and related properties of seawater. Modifications: 2011-03-23. David Jackett, Trevor McDougall and Paul Barker. """ t = t_from_CT(SA, CT, p) return alpha_wrt_CT_t_exact(SA, t, p)
def specvol_anom_CT_exact(SA, CT, p): r"""Calculates specific volume anomaly from Absolute Salinity, Conservative Temperature and pressure. The reference value of Absolute Salinity is SSO and the reference value of Conservative Temperature is equal to 0 deg C. Parameters ---------- SA : array_like Absolute Salinity [g/kg] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p : array_like sea pressure [dbar] Returns ------- specvol_anom_CT_exact : array_like specific volume anomaly [m**3/kg] See Also -------- TODO Notes ----- This function uses the full Gibbs function. There is an alternative to calling this function, namely specvol_anom_CT(SA, CT, p), which uses the computationally efficient 48-term expression for density in terms of SA, CT and p (McDougall et al., 2011). Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See Eqn. (3.7.3). .. [2] McDougall T.J., P.M. Barker, R. Feistel and D.R. Jackett, 2011: A computationally efficient 48-term expression for the density of seawater in terms of Conservative Temperature, and related properties of seawater. Modifications: 2011-04-06. Trevor McDougall and Paul Barker. """ t = t_from_CT(SA, CT, p) return specvol_anom_t_exact(SA, t, p)
def rho_alpha_beta_CT_exact(SA, CT, p): r"""Calculates in-situ density, the appropriate thermal expansion coefficient and the appropriate saline contraction coefficient of seawater from Absolute Salinity and Conservative Temperature. See the individual functions rho_CT_exact, alpha_CT_exact, and beta_CT_exact. Retained for compatibility with the Matlab GSW toolbox. """ t = t_from_CT(SA, CT, p) rho_CT_exact = rho_t_exact(SA, t, p) alpha_CT_exact = alpha_wrt_CT_t_exact(SA, t, p) beta_CT_exact = beta_const_CT_t_exact(SA, t, p) return rho_CT_exact, alpha_CT_exact, beta_CT_exact
def sound_speed_CT_exact(SA, CT, p): r"""Calculates the speed of sound in seawater from Absolute Salinity and Conservative Temperature and pressure. Parameters ---------- SA : array_like Absolute salinity [g kg :sup:`-1`] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p : array_like pressure [dbar] Returns ------- sound_speed_CT_exact : array_like Speed of sound in seawater [m s :sup:`-1`] See Also -------- TODO Notes ----- TODO Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See Eqn. (2.17.1). Modifications: 2011-04-05. David Jackett, Paul Barker and Trevor McDougall. """ t = t_from_CT(SA, CT, p) return sound_speed_t_exact(SA, t, p)
def internal_energy_CT_exact(SA, CT, p): r"""Calculates the specific internal energy of seawater from Absolute Salinity, Conservative Temperature and pressure. Parameters ---------- SA : array_like Absolute Salinity [g/kg] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p : array_like sea pressure [dbar] Returns ------- internal_energy_CT_exact: array_like specific internal energy (u) [J/kg] See Also -------- TODO Notes ----- TODO Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See Eqn. (2.11.1). Modifications: 2011-04-05. Trevor McDougall. """ t = t_from_CT(SA, CT, p) return internal_energy_t_exact(SA, t, p)
def t_freezing(SA, p, saturation_fraction=1): r"""Calculates the in-situ temperature at which seawater freezes. Parameters ---------- SA : array_like Absolute Salinity [g/kg] p : array_like sea pressure [dbar] saturation_fraction : fraction between 0, 1. The saturation fraction of dissolved air in seawater. Default is 0 or completely saturated. Returns ------- t_freezing : array_like in-situ temperature at which seawater freezes [:math:`^\circ` C (ITS-90)] See Also -------- TODO Notes ----- TODO Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See sections 3.33 and 3.34. Modifications: 2011-11-03. Trevor McDougall, Paul Barker and Rainer Feistal. """ """This function, t_freezing, calculates the in-situ freezing temperature, t_freezing, of seawater by first evaluating a polynomial of the Conservative Temperature at which seawater freezes, CT_freezing, using the GSW function CT_freezing. The in-situ freezing temperature is then calculated using the GSW function t_from_CT. However, if one wanted to compute the in-situ freezing temperature directly from a single polynomial expression without first calculating the Conservative Temperature at the freezing point, the following lines of code achieve this. The error of the following fit is similar to that of the present function, t_freezing, and ranges between -8e-4 K and 3e-4 K when compared with the in-situ freezing temperature evaluated by Newton-Raphson iteration of the equality of the chemical potentials of water in seawater and in ice. (Note that the in-situ freezing temperature can be found by this exact method using the function sea_ice_freezingtemperature_si in the SIA library). SA_r = SA * 1e-2 x = np.sqrt(SA_r) p_r = p * 1e-4 t_freeze = T[0] + SA_r * (T[1] + x * (T[2] + x * (T[3] + x * (T[4] + x * (T[5] + T[6] * x))))) + p_r * (T[7] + p_r * (T[8] + T[9] * p_r)) + SA_r * p_r * (T[10] + p_r * (T[12] + p_r * (T[15] + T[21] * SA_r)) + SA_r * (T[13] + T[17] * p_r + T[19] * SA_r) + x * (T[11] + p_r * (T[14] + T[18] * p_r) + SA_r * (T[16] + T[20] * p_r + T[22] * SA_r))) Adjust for the effects of dissolved air t_freezing -= saturation_fraction * (1e-3) * (2.4 - SA / 70.33008) """ SA, p, saturation_fraction = np.broadcast_arrays(SA, p, saturation_fraction) if (SA < 0).any(): raise ValueError('SA must be non-negative!') if np.logical_or(saturation_fraction < 0, saturation_fraction > 1).any(): raise ValueError('Saturation_fraction MUST be between zero and one.') CT_freeze = CT_freezing(SA, p, saturation_fraction) t_freeze = t_from_CT(SA, CT_freeze, p) tmp = np.logical_or(p > 10000, SA > 120) out = np.logical_or(tmp, p + SA * 71.428571428571402 > 13571.42857142857) t_freeze[out] = np.ma.masked return t_freeze
def enthalpy_diff_CT_exact(SA, CT, p_shallow, p_deep): r"""Calculates the difference of the specific enthalpy of seawater between two different pressures, p_deep (the deeper pressure) and p_shallow (the shallower pressure), at the same values of SA and CT. The output (enthalpy_diff_CT_exact) is the specific enthalpy evaluated at (SA, CT, p_deep) minus the specific enthalpy at (SA,CT,p_shallow). parameters ---------- Parameters ---------- SA : array_like Absolute Salinity [g/kg] CT : array_like Conservative Temperature [:math:`^\circ` C (ITS-90)] p_shallow : array_like lower sea pressure [dbar] p_deep : array-like upper sea pressure [dbar] returns ------- enthalpy_diff_CT_exact : array_like difference of specific enthalpy [J/kg] (deep minus shallow) See Also -------- TODO Notes ----- This function uses the full Gibbs function. There is an alternative to calling this function, namely enthalpy_diff_CT(SA, CT, p), which uses the computationally efficient 48-term expression for density in terms of SA, CT and p (McDougall et al., 2011). Examples -------- TODO References ---------- .. [1] IOC, SCOR and IAPSO, 2010: The international thermodynamic equation of seawater - 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides No. 56, UNESCO (English), 196 pp. See Eqns (3.32.2). .. [2] McDougall T.J., P.M. Barker, R. Feistel and D.R. Jackett, 2011: A computationally efficient 48-term expression for the density of seawater in terms of Conservative Temperature, and related properties of seawater. Modifications: 2011-04-06. Trevor McDougall and Paul Barker. """ t_shallow = t_from_CT(SA, CT, p_shallow) t_deep = t_from_CT(SA, CT, p_deep) return (enthalpy_t_exact(SA, t_deep, p_deep) - enthalpy_t_exact(SA, t_shallow, p_shallow))
def sigma4_CT_exact(SA, CT): r"""Calculates potential density anomaly with reference pressure of 4000 dbar.""" t = t_from_CT(SA, CT, 4000.) return rho_t_exact(SA, t, 4000.) - 1000