# -*- coding: utf-8 -*-
"""
@author: jeremy leconte
"""
import numpy as np
import numba
import exo_k.util.cst as cst
[docs]
class Condensing_species(object):
def __init__(self, Latent_heat_vaporization = 0., cp_vap = 0., Mvap = 0.,
T_ref = 0., Psat_ref = 0., cp_cond = None):
"""
Parameters
----------
Latent_heat_vaporization: float
specific Latent heat vaporization (J/kg)
cp_vap: float
Specific heat capacity of the vapor (J/kg/K)
Mvap: float
Molar mass of vapor (kg/mol)
T_ref: float
Reference temperature
Psat_ref: float
Saturation vapor pressure at the reference temperature (Pa)
cp_cond: float (optional)
Specific heat capacity of the condensate (J/kg/K).
Assumed equal to cp_vap if not provided.
"""
self.Latent_heat_vaporization = Latent_heat_vaporization
self.cp_vap = cp_vap
self.T_ref = T_ref
self.Psat_ref = Psat_ref
self.Mvap = Mvap
self.Rvap = cst.RGP / self.Mvap
if cp_cond is None:
self.cp_cond = cp_vap
else:
self.cp_cond = cp_cond
self.delta_cp = self.cp_vap - self.cp_cond
self.delta_cp_R = self.delta_cp / self.Rvap
self.LovR = self.Latent_heat_vaporization / self.Rvap
self.c1 = self.LovR / self.T_ref - self.delta_cp_R
self.c2 = self.delta_cp_R * self.T_ref - self.LovR
[docs]
def Lvap(self, T):
"""Latent heat at temperature T
Parameters
----------
T: array
Temperature in layers (K)
"""
return Lvap_T(T, self.Latent_heat_vaporization, self.T_ref, self.delta_cp)
[docs]
def Psat(self, T):
"""Saturation vapor pressure for the condensing species
Parameters
----------
T: array
Temperature in layers (K)
"""
return Psat_T(T, self.T_ref, self.Psat_ref, self.c1, self.c2, self.delta_cp_R)
[docs]
def Tsat(self, P):
"""Boiling temperature for the condensing species
NOT EXACT IF DELTA CP != 0.
Parameters
----------
P: array
Pressure in layers (Pa)
"""
return Tsat_P(P, self.Psat_ref, self.c1, self.c2)
[docs]
def qsat(self, psat, p, epsilon):
"""Saturation vapor mass mixing ratio for the condensing species
Parameters
----------
psat: array
saturation vapor pressure
p: array
pressure at the layer center
epsilon: float or array
Ratio of the molar mass of the vapor over the background molar mass
"""
return Qsat(psat, p, epsilon)
[docs]
def dPsat_dT(self, T):
"""Saturation vapor pressure derivative for the condensing species
Parameters
----------
T: array
Temperature in layers (K)
"""
return dPsat_dT(T, self.Latent_heat_vaporization, self.T_ref,
self.Psat_ref, self.Rvap, self.delta_cp, self.delta_cp_R, self.c1, self.c2)
[docs]
def dlnPsat_dlnT(self, T):
"""Saturation vapor pressure for the condensing species and its derivative
Parameters
----------
T: array
Temperature in layers (K)
"""
return dlnPsat_dlnT(T, self.Latent_heat_vaporization, self.T_ref,
self.delta_cp, self.Rvap)
[docs]
def moist_adiabat(self, T, P, cp, Mgas):
"""Computes the threshold thermal gradient (d ln T / d ln P) for a moist atmosphere.
Parameters
----------
T: array
Temperature in layers (K)
P: array
pressure at the layer center
cp: float
specific heat capacity at constant pressure of the background gas
Mgas: float or array
Molar mass of the background atmosphere
Returns
-------
array
Moist adiabat lapse rate
Lvap: array
Latent heat at temperature T
psat: array
Saturation vapor pressure for the condensing species (Pa)
qsat: array
Saturation vapor mass mixing ratio for the condensing species
dqsat_dt: array
Derivative of qsat with respect to temperature at fixed pressure
q_crit: array
Critical mass mixing ratio for the inhibition of moist convection
(Eq. 17 of Leconte et al. 2017)
"""
return moist_adiabat(T, P, cp, Mgas, self.cp_vap, self.Mvap, self.Rvap,
self.Latent_heat_vaporization, self.T_ref, self.Psat_ref,
self.delta_cp, self.delta_cp_R, self.c1, self.c2)
[docs]
def compute_condensation_parameters(self, T, P, Mgas):
"""Computes necessary quantities to compute
large scale condensation.
Parameters
----------
T: array
Temperature in layers (K)
P: array
pressure at the layer center
Mgas: float or array
Molar mass of the background atmosphere
"""
return compute_condensation_parameters(T, P, Mgas, self.Mvap,
self.Rvap, self.Latent_heat_vaporization, self.T_ref,
self.Psat_ref, self.delta_cp, self.delta_cp_R, self.c1, self.c2)
def __repr__(self):
"""Method to output parameters
"""
output="""
Lvap : {lvap}
cp (vap) : {cp}
T0 (K): {t0}
Psat(T0) (Pa): {p0}""".format(lvap=self.Latent_heat_vaporization,
cp=self.cp_vap, t0=self.T_ref, p0=self.Psat_ref)
return output
[docs]
class Condensation_Thermodynamical_Parameters(object):
def __init__(self, Latent_heat_vaporization = 0., cp_vap = 0., Mvap = 0.,
T_ref = 0., Psat_ref = 0., cp_cond = None):
"""
Parameters
----------
Latent_heat_vaporization: float
specific Latent heat vaporization (J/kg)
cp_vap: float
Specific heat capacity of the vapor (J/kg/K)
Mvap: float
Molar mass of vapor (kg/mol)
T_ref: float
Reference temperature
Psat_ref: float
Saturation vapor pressure at the reference temperature (Pa)
cp_cond: float (optional)
Specific heat capacity of the condensate (J/kg/K).
Assumed equal to cp_vap if not provided.
"""
Rvap = cst.RGP / Mvap
if cp_cond is None:
cp_cond = cp_vap
else:
cp_cond = cp_cond
delta_cp = cp_vap - cp_cond
delta_cp_R = delta_cp / Rvap
LovR = Latent_heat_vaporization / Rvap
c1 = LovR / T_ref - delta_cp_R
c2 = delta_cp_R * T_ref - LovR
self.th_params = np.array([cp_vap, Mvap, Rvap,
Latent_heat_vaporization, T_ref, Psat_ref,
delta_cp, delta_cp_R, c1, c2])
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def Lvap_T(T, Latent_heat_vaporization, T_ref, delta_cp):
"""Latent heat at temperature T
Parameters
----------
T: array
Temperature in layers (K)
"""
return Latent_heat_vaporization * np.ones_like(T)
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def Psat_T(T, T_ref, Psat_ref, c1, c2, delta_cp_R):
"""GCM version for test
"""
Pref_solid_liquid=611.14
Trefvaporization=35.86
Trefsublimation=7.66
Tmin=8.
r3vaporization=17.269
r3sublimation=21.875
T_h2O_ice_liq=273.16
psat = np.where(T>T_h2O_ice_liq,
Pref_solid_liquid*np.exp(r3vaporization*(T-T_h2O_ice_liq)/(T-Trefvaporization)),# liquid / vapour
Pref_solid_liquid*np.exp(r3sublimation*(T-T_h2O_ice_liq)/(T-Trefsublimation)) # solid / vapour
)
return psat
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def Tsat_P(P, Psat_ref, c1, c2):
"""GCM version for test
"""
Pref_solid_liquid=611.14
Trefvaporization=35.86
Trefsublimation=7.66
Tmin=8.
r3vaporization=17.269
r3sublimation=21.875
T_h2O_ice_liq=273.16
Tsat = np.where(P<Pref_solid_liquid,
(T_h2O_ice_liq*r3sublimation- Trefsublimation*np.log(P/Pref_solid_liquid))/(r3sublimation-np.log(P/Pref_solid_liquid)),
(T_h2O_ice_liq*r3vaporization- Trefvaporization*np.log(P/Pref_solid_liquid))/(r3vaporization-np.log(P/Pref_solid_liquid))
)
return Tsat
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def Qsat(psat, p, epsilon):
"""Saturation vapor mass mixing ratio for the condensing species
Parameters
----------
psat: array
saturation vapor pressure
p: array
pressure at the layer center
epsilon: float or array
Ratio of the molar mass of the vapor over the background molar mass
"""
psat_tmp = np.core.umath.minimum(psat,p)
fac = p + (epsilon -1.) * psat_tmp
qsat = epsilon * psat_tmp / fac
return qsat
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def dPsat_dT(T, Latent_heat_vaporization, T_ref, Psat_ref, Rvap, delta_cp, delta_cp_R, c1, c2):
"""Saturation vapor pressure derivative for the condensing species
Parameters
----------
T: array
Temperature in layers (K)
"""
Lvap = Lvap_T(T, Latent_heat_vaporization, T_ref, delta_cp)
psat = Psat_T(T, T_ref, Psat_ref, c1, c2, delta_cp_R)
Trefvaporization=35.86
Trefsublimation=7.66
r3vaporization=17.269
r3sublimation=21.875
T_h2O_ice_liq=273.16
dpsat = np.where(T>T_h2O_ice_liq,
psat*r3vaporization*(T_h2O_ice_liq-Trefvaporization)/(T-Trefvaporization)**2,# liquid / vapour
psat*r3sublimation*(T_h2O_ice_liq-Trefsublimation)/(T-Trefsublimation)**2 # solid / vapour
)
return dpsat, psat, Lvap
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def dlnPsat_dlnT(T, Latent_heat_vaporization, T_ref, delta_cp, Rvap):
"""Saturation vapor pressure for the condensing species and its derivative
Parameters
----------
T: array
Temperature in layers (K)
"""
Lvap = Lvap_T(T, Latent_heat_vaporization, T_ref, delta_cp)
Trefvaporization=35.86
Trefsublimation=7.66
r3vaporization=17.269
r3sublimation=21.875
T_h2O_ice_liq=273.16
dlnpsat_dt = np.where(T>T_h2O_ice_liq,
r3vaporization*(T_h2O_ice_liq-Trefvaporization)/(T-Trefvaporization)**2,# liquid / vapour
r3sublimation*(T_h2O_ice_liq-Trefsublimation)/(T-Trefsublimation)**2 # solid / vapour
)
return dlnpsat_dt*T, Lvap
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def moist_adiabat(T, P, cp, Mgas, cp_vap, Mvap, Rvap, Latent_heat_vaporization, T_ref, Psat_ref, delta_cp, delta_cp_R, c1, c2):
"""Computes the threshold thermal gradient (d ln T / d ln P) for a moist atmosphere.
Parameters
----------
T: array
Temperature in layers (K)
P: array
pressure at the layer center
cp: float
specific heat capacity at constant pressure of the background gas
Mgas: float or array
Molar mass of the background atmosphere
Returns
-------
array
Moist adiabat lapse rate
Lvap: array
Latent heat at temperature T
psat: array
Saturation vapor pressure for the condensing species (Pa)
qsat: array
Saturation vapor mass mixing ratio for the condensing species
dqsat_dt: array
Derivative of qsat with respect to temperature at fixed pressure
q_crit: array
Critical mass mixing ratio for the inhibition of moist convection
(Eq. 17 of Leconte et al. 2017)
"""
epsilon = Mvap / Mgas
psat = Psat_T(T, T_ref, Psat_ref, c1, c2, delta_cp_R)
qsat = Qsat(psat, P, epsilon)
dlpsat_dlT, Lvap = dlnPsat_dlnT(T, Latent_heat_vaporization, T_ref, delta_cp, Rvap)
dqsat_dt = qsat**2 * P * dlpsat_dlT / (epsilon*psat*T)
qa=1.-qsat
qLvt = qsat * Lvap / T
# fac = qsat * cp_vap + qa * cp + qLvt * dlpsat_dlT
fac = cp + qLvt * dlpsat_dlT
q_crit = epsilon / ((epsilon - 1.) * dlpsat_dlT)
return ( cst.RGP / Mgas + qLvt ) / fac, Lvap, psat, qsat, dqsat_dt, q_crit # is missing p/p_a terms
# return ( cst.RGP / Mgas + qLvt ) / fac*0.8, Lvap, psat, qsat, dqsat_dt, q_crit # is missing p/p_a terms
# factor 0.8 to mimic 3D
[docs]
@numba.jit(nopython=True, fastmath=True, cache=True)
def compute_condensation_parameters(T, P, Mgas, Mvap, Rvap, Latent_heat_vaporization, T_ref, Psat_ref, delta_cp, delta_cp_R, c1, c2):
"""Computes necessary quantities to compute
large scale condensation.
Parameters
----------
T: array
Temperature in layers (K)
P: array
pressure at the layer center
Mgas: float or array
Molar mass of the background atmosphere
"""
epsilon = Mvap / Mgas
psat = Psat_T(T, T_ref, Psat_ref, c1, c2, delta_cp_R)
qsat = Qsat(psat, P, epsilon)
dlpsat_dlT, Lvap = dlnPsat_dlnT(T, Latent_heat_vaporization, T_ref, delta_cp, Rvap)
dqsat_dt = qsat**2 * P * dlpsat_dlT / (epsilon*psat*T)
return Lvap, qsat, dqsat_dt