:py:mod:`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, 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` :returns: array Heating rate in each atmospheric layer (W/kg). .. py:function:: turbulent_diffusion_numba(timestep, Nlay, p_lay, p_lev, dmass, t_lay_ov_mu, g, Kzz, tracer_array, 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). :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. :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:: 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, 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) !======================================================================!