exo_k.atm_evolution.convection ============================== .. py:module:: exo_k.atm_evolution.convection .. autoapi-nested-parse:: @author: jeremy leconte Module Contents --------------- .. py:function:: dry_convective_adjustment_numba(timestep, Nlay, t_lay, exner, dmass, tracer_array, Mgas, verbose=False) Computes the heating rates needed to adjust unstable regions of a given atmosphere to a convectively neutral T profile on a given timestep. :param timestep: Duration of the adjustment. If given in seconds, H=dT/timestep is in K/s. In the current model, timestep is in second over cp. This ensures that H=dT/timestep is in W/kg. :type timestep: :class:`float` :param Nlay: Number of atmospheric layers :type Nlay: :class:`int` :param t_lay: Temperatures of the atmospheric layers (K) :type t_lay: :class:`array`, :class:`np.ndarray` :param exner: Exner function computed at the layer centers ((p/psurf)**rcp) .. math:: \Pi=(p / p_{s})^{R/c_p} :type exner: :class:`array`, :class:`np.ndarray` :param dmass: Mass of gas in each layer (kg/m^2) :type dmass: :class:`array`, :class:`np.ndarray` :returns: array Heating rate in each atmospheric layer (W/kg). .. py:function:: convective_acceleration_numba(timestep, Nlay, H_rad, rad_layers, tau_rad, tau_rads, dmass, convective_acceleration_mode=0, convective_acceleration_factor=1.0, verbose=False) Computes a heating rate for the whole convective region to accelerate convergence :param timestep: Duration of the adjustment. If given in seconds, H=dT/timestep is in K/s. In the current model, timestep is in second over cp. This ensures that H=dT/timestep is in W/kg. :type timestep: :class:`float` :param Nlay: Number of atmospheric layers :type Nlay: :class:`int` :param H_rad: Radiative heating rate :type H_rad: :class:`array`, :class:`np.ndarray` :param rad_layers: Elements in the array are true if layer is purely radiative :type rad_layers: :class:`array`, :class:`np.ndarray` of :class:`bool` :param tau_rad: Baseline radiative timescale for the atmosphere. (e.g. the min of tau_rads) Should use the same units as timestep. :type tau_rad: :class:`array`, :class:`np.ndarray` :param tau_rads: Radiative timescale for each layer. Should use the same units as timestep. :type tau_rads: :class:`array`, :class:`np.ndarray` :param dmass: Mass of gas in each layer (kg/m^2) :type dmass: :class:`array`, :class:`np.ndarray` :param convective_acceleration_mode: convective_acceleration_mode=0 is the default mode =1 emulates an old (and bad) behavior. :type convective_acceleration_mode: :class:`int` :param convective_acceleration_factor: Multiplying factor for convection acceleration :type convective_acceleration_factor: :class:`int` :returns: array Heating rate in each atmospheric layer (W/kg). .. py:function:: turbulent_diffusion_numba(timestep, Nlay, p_lay, p_lev, dmass, t_lay, exner, t_lay_ov_mu, g, Kzz, tracer_array, cp, mix_potential_temp=False, verbose=False) Solves turbulent diffusion equation: .. math:: \rho frac{\partial q}{\partial t} = \frac{\partial F_{diff}}{\partial z} with a diffusive flux given by .. math:: F_{diff} = - \rho K_{zz} \frac{\partial q}{\partial z} The equation is solved with an implicit scheme assuming that there is no flux at the top and bottom boundaries (evaporation must be treated separately for now). For the potential temperature mixing, the equation in not completely nrj conserving. One should mix enthalpy as in GCM. :param timestep: Time step in seconds. :type timestep: :class:`float` :param Nlay: Number of atmospheric layers :type Nlay: :class:`int` :param t_lay_ov_mu: Temperatures of the atmospheric layers divided by the molar_mass in kg/mol :type t_lay_ov_mu: :class:`array`, :class:`np.ndarray` :param p_lay: Pressure at the layer centers (Pa) :type p_lay: :class:`array`, :class:`np.ndarray` :param p_lev: Presure at the Nlay+1 level boundaries (Pa) :type p_lev: :class:`array`, :class:`np.ndarray` :param dmass: Mass of gas in each layer (kg/m^2) :type dmass: :class:`array`, :class:`np.ndarray` :param g: Gravity (m/s^2) :type g: :class:`float` :param Kzz: Eddy mixing coefficient (m^2/s) :type Kzz: :class:`float` :param tracer_array: Array containing the mass mixing ratio of all tracers at each layer before the mixing :type tracer_array: :class:`array`, :class:`np.ndarray (Ntrac`, :class:`Nlay)` :returns: new_tracers: array, np.ndarray (Ntrac, Nlay) Array containing the mass mixing ratio of all tracers at each layer after the mixing .. py:function:: compute_condensation_numba(timestep, Nlay, tlay, play, cp, Mgas, qarray, idx_vap, idx_cond, thermo_parameters, latent_heating=True, condensation_timestep_reducer=None, verbose=False) Computes the heating rates needed to bring sursaturated regions back to saturation :param timestep: :param Nlay: Number of layers :type Nlay: :class:`float` :param tlay: Layer temperatures :type tlay: :class:`array`, :class:`np.ndarray` :param cp: Specific heat capacity :type cp: :class:`float` :param Mgas: Layer molar mass :type Mgas: :class:`array`, :class:`np.ndarray` :param qarray: Array of tracer specific concentrations :type qarray: :class:`array`, :class:`np.ndarray` :param idx_vap: Indices for condensing species in tracer array :type idx_vap: :class:`int` :param idx_cond: Indices for condensing species in tracer array :type idx_cond: :class:`int` :param thermo_parameters: Array of thermodynamical paramaters for condensing species. See `:class:`~exo_k.atm_evolution.condensation.Condensation_Thermodynamical_Parameters` :type thermo_parameters: :class:`array`, :class:`np.ndarray` :param latent_heating: Whether latent heating should be taken into account :type latent_heating: :class:`bool` :returns: H_cond: array, np.ndarray Heating rate in W/kg .. py:function:: moist_convective_adjustment_numba(timestep, Nlay, tlay, play, dmass, cp, Mgas, q_array, i_vap, i_cond, thermo_parameters, moist_inhibition=True, verbose=False) Computes the heating rates needed to adjust unstable regions of a given atmosphere to a moist adiabat. Based on *Description algorithmique d'un ensemble de paramétrisation physique - phylmd* Zhaoxin Laurent LI :param timestep: :param Nlay: Number of layers :type Nlay: :class:`float` :param tlay: Layer temperatures :type tlay: :class:`array`, :class:`np.ndarray` :param play: Pressure at layer centers :type play: :class:`array` :param dmass: mass of layers in kg/m^2 :type dmass: :class:`array`, :class:`np.ndarray` :param cp: specific heat capacity at constant pressure :type cp: :class:`float` :param q_array: mass mixing ratio of tracers :type q_array: :class:`array`, :class:`np.ndarray` :param i_vap: index of vapor tracer in qarray :type i_vap: :class:`int` :param i_cond: index of condensate tracer in qarray :type i_cond: :class:`int` :param qsat: Saturation mmr for each layer :type qsat: :class:`array`, :class:`np.ndarray` :param dqsat_dt: d qsat / dT in each layer :type dqsat_dt: :class:`array`, :class:`np.ndarray` :param Lvap: Latent heat of vaporization (can have different values in each layer if Lvap=f(T)) :type Lvap: :class:`array`, :class:`np.ndarray` :param dlnt_dlnp_moist: threshold thermal gradient (d ln T / d ln P) for a moist atmosphere computed at the layer centers. :type dlnt_dlnp_moist: :class:`array`, :class:`np.ndarray` :param q_crit: Critical mass mixing ratio for the inhibition of moist convection (Eq. 17 of Leconte et al. 2017) :type q_crit: :class:`array`, :class:`np.ndarray` :returns: H_madj: array, np.ndarray Heating rate in each atmospheric layer (W/kg). new_q: array, np.ndarray tracer mmr array after adjustment. new_t: array, np.ndarray Temperature of layers after adjustment. .. py:function:: moist_convective_adjustment_cloud_fraction_numba(timestep, Nlay, tlay, play, dmass, cp, Mgas, q_array, i_vap, i_cond, thermo_parameters, moist_inhibition=True, verbose=False, humidity_distribution_width=0.2) Computes the heating rates needed to adjust unstable regions of a given atmosphere to a moist adiabat. This approach with a cloud fraction is experimental and should be used with extreme caution. It does not seem to give good results when re-evaporation is very efficient because there is constant re-evaporation and condensation in the convective region. :param timestep: :param Nlay: Number of layers :type Nlay: :class:`float` :param tlay: Layer temperatures :type tlay: :class:`array` :param play: Pressure at layer centers :type play: :class:`array` :param dmass: mass of layers in kg/m^2 :type dmass: :class:`array` :param cp: specific heat capacity at constant pressure :type cp: :class:`float` :param q_array: mass mixing ratio of tracers :type q_array: :class:`array` :param i_vap: index of vapor tracer in qarray :type i_vap: :class:`int` :param i_cond: index of condensate tracer in qarray :type i_cond: :class:`int` :param qsat: Saturation mmr for each layer :type qsat: :class:`array` :param dqsat_dt: d qsat / dT in each layer :type dqsat_dt: :class:`array` :param Lvap: Latent heat of vaporization (can have different values in each layer if Lvap=f(T)) :type Lvap: :class:`array` :param dlnt_dlnp_moist: threshold thermal gradient (d ln T / d ln P) for a moist atmosphere computed at the layer centers. :type dlnt_dlnp_moist: :class:`array` :param q_crit: Critical mass mixing ratio for the inhibition of moist convection (Eq. 17 of Leconte et al. 2017) :type q_crit: :class:`array` :returns: H_madj: array Heating rate in each atmospheric layer (W/kg). new_q: array tracer mmr array after adjustment. new_t: array Temperature of layers after adjustment. .. py:function:: compute_rainout_numba(timestep, Nlay, tlay, play, dmass, cp, Mgas, qarray, idx_vap, idx_cond, thermo_parameters, evap_coeff, qvap_deep, q_cloud=0.0, latent_heating=False, total_cloud_fraction=1.0, humidity_distribution_width=0.2, verbose=False) Computes the heating rates needed to adjust unstable regions of a given atmosphere to a moist adiabat. :param timestep: :param Nlay: Number of layers :type Nlay: :class:`float` :param tlay: Layer temperatures :type tlay: :class:`array`, :class:`np.ndarray` .. py:function:: molecular_diffusion_numba(timestep, Nlay, p_lay, p_lev, dmass, tlay, mu, g, Dmol, verbose=False) Solves turbulent diffusion equation: .. math:: \rho frac{\partialT}{\partial t} = \frac{\partial F_{diff}}{\partial z} with a diffusive flux given by .. math:: F_{diff} = - \rho D_{mol} \frac{\partial T}{\partial z} The equation is solved with an implicit scheme assuming that there is no flux at the top and bottom boundaries (evaporation must be treated separately for now). :param timestep: Time step in seconds. :type timestep: :class:`float` :param Nlay: Number of atmospheric layers :type Nlay: :class:`int` :param t_lay_ov_mu: Temperatures of the atmospheric layers divided by the molar_mass in kg/mol :type t_lay_ov_mu: :class:`array`, :class:`np.ndarray` :param p_lay: Pressure at the layer centers (Pa) :type p_lay: :class:`array`, :class:`np.ndarray` :param p_lev: Presure at the Nlay+1 level boundaries (Pa) :type p_lev: :class:`array`, :class:`np.ndarray` :param dmass: Mass of gas in each layer (kg/m^2) :type dmass: :class:`array`, :class:`np.ndarray` :param g: Gravity (m/s^2) :type g: :class:`float` :param Dmol: molecular diffusion coefficient (m^2/s) :type Dmol: :class:`float` :returns: new_tlay: array, np.ndarray (Nlay) Array containing the temperature at each layer after the mixing .. py:function:: DTRIDGL(L, AF, BF, CF, DF) ! GCM2.0 Feb 2003 ! DOUBLE PRECISION VERSION OF TRIDGL DIMENSION AF(L),BF(L),CF(L),DF(L),XK(L) DIMENSION AS(2*L),DS(2*L) !* THIS SUBROUTINE SOLVES A SYSTEM OF TRIDIAGIONAL MATRIX !* EQUATIONS. THE FORM OF THE EQUATIONS ARE: !* A(I)*X(I-1) + B(I)*X(I) + C(I)*X(I+1) = D(I) !======================================================================!