:py:mod:`exo_k.atm_evolution.atm_evol` ====================================== .. py:module:: exo_k.atm_evolution.atm_evol .. autoapi-nested-parse:: @author: jeremy leconte Module Contents --------------- .. py:class:: Atm_evolution(bg_vmr=None, verbose=False, **kwargs) Bases: :py:obj:`object` Model of atmospheric evolution. Uses exo_k.Atm class to compute radiative transfer Initializes atmospheric profiles. Most arguments are passed directly to exo_k.Atm class through **kwargs .. warning:: Layers are counted from the top down (increasing pressure order). All methods follow the same convention. .. py:property:: time Yields current time in seconds .. py:property:: time_hist Yields the array of the times for the last call to evolve (in seconds) .. py:method:: set_options(reset_rad_model=False, check_keys=True, Kzz=None, cp=None, verbose=False, **kwargs) This method is used to store the global options in the `Settings` object. Arguments are all passed through **kwargs. Sometimes, one needs to reset the radiative model to take into account some modifications (like the databases). Normally, this should be automatic, but you can force it with `reset_rad_model=True` .. py:method:: initialize_condensation(condensing_species=None, **kwargs) This method initializes the condensation module by listing all the condensing vapors and linking them to their condensed form. For each vapor-condensate pair, a :class:`Condensible_species` object is created with the thermodynamical data provided. Here is an example of dictionary to provide as input to include CH4 condensation ``` condensing_species={'ch4':{'Latent_heat_vaporization': 5.25e5, 'cp_vap': 2.232e3, 'Mvap': 16.e-3, 'T_ref': 90., 'Psat_ref': 0.11696e5}} ``` .. py:method:: setup_radiative_model(k_database=None, k_database_stellar=None, cia_database=None, cia_database_stellar=None, gas_vmr=None, **kwargs) This method initializes the exo_k.Atm object that will be used to carry out radiative transfer calculations. This is where the radiative data used are chosen and transmitted to the radiative transfer module, along with many other parameters including the incoming stellar flux (`flux_top_dw`), the blackbody temperature of the star (`Tstar`), the If a `k_database_stellar` is provided, then this is this database that will be used to treat the scattering and absorption of incoming radiation. In that case, `k_database` will be used to treat the emission of the atmosphere. The effective cos(zenith angle) for the incoming stellar radiation can then be specified independently with the `mu0_stellar` keyword. If no `k_database_stellar` is provided, `k_database` will be used to treat both the atmospheric emission and the stellar radiation. Running a model with `k_database_stellar=k_database` yields the same results at twice the cost. :param k_database: radiative database for the molecules in the atmospheres. :type k_database: `exo_k.Kdatabase` objects :param k_database_stellar: radiative database for the molecules in the atmospheres. :type k_database_stellar: `exo_k.Kdatabase` objects :param cia_database: radiative database for the CIA ofmolecules in the atmospheres. :type cia_database: `exo_k.CIA_database` object :param cia_database_stellar: radiative database for the CIA ofmolecules in the atmospheres. :type cia_database_stellar: `exo_k.CIA_database` object .. py:method:: compute_average_fluxes() Use the averaged heating rates to compute the various fluxes (W/m^2) at the level interfaces. These fluxes are positive when the energy flows upward. To be consistent with radiative fluxes, the first array value corresponds to the top of atmosphere and should be 0 in most cases. The last value corresponds to the flux between the deepest layer (considered to be the surface) and the layer just above. .. py:method:: evolve(N_timestep=1, N_kernel=10000, timestep_factor=1.0, dT_max=100.0, verbose=False, check_cons=False, **kwargs) The time variable used in the model is tau=t/cp. The equation we are solving in each layer is thus .. math:: c_p \frac{d T}{d t}= \frac{d T}{d tau} = \sum_i H_i For a given timestep `dtau`, the physical time elapsed in second can be computed using `dt=dtau*cp` To work, the heating rates (H) must be computed in W/kg. This also means that if one needs the physical rate of change of another quantity (like dq/dt) from the delta q over a time step, one needs to do `dq/dt = delta q / (timestep * cp)` :param N_timestep: Number of timesteps to perform. :type N_timestep: :class:`int` :param N_kernel: Maximal number of timesteps between two computations of the radiative kernel. :type N_kernel: :class:`int` :param timestep_factor: Multiplicative factor applied to timestep computed automatically by the radiative module. timestep_factor > 1 can lead to unstabilities. :type timestep_factor: :class:`float` :param dT_max: Maximum temperature increment in a single timestep. :type dT_max: :class:`float` .. py:method:: equilibrate(Fnet_tolerance=None, N_iter_max=10, N_timestep_ini=100, N_timestep_max=1000000, verbose=False, **kwargs) Evolves an atmosphere until it is at equilibrium. Equilibrium is assumed to be reached when the net top of atmosphere flux remains within +-Fnet_tolerance of the internal flux for a whole evolution step. The number of timesteps per evolution step in multiplied by two at each iteration, starting from N_timestep_ini, until the limit of N_timestep_max is reached. :param Fnet_tolerance: Tolerance on net flux in W/m^2 to identify convergence. :type Fnet_tolerance: :class:`float` :param N_iter_max: Max number of successive calls to evolve :type N_iter_max: :class:`int` :param N_timestep_ini: Initial number of timesteps in a single evolution step :type N_timestep_ini: :class:`int` :param N_timestep_max: Max number of timesteps in a single evolution step :type N_timestep_max: :class:`int` .. py:method:: moist_convective_adjustment(timestep, Htot, moist_inhibition=True, verbose=False) This method computes the vapor and temperature tendencies do to moist convectoin in saturated layers. The tracer array in modified in place. :param timestep: physical timestep of the current step (in s/cp). :type timestep: :class:`float` :param Htot: Total heating rate (in W/kg) of all physical processes already computed :type Htot: :class:`array`, :class:`np.ndarray` :returns: H_madj: array, np.ndarray Heating rate due to large scale condensation (W/kg) .. py:method:: molecular_diffusion(timestep, Htot, atm, cp) Mixes energy following a diffusion equation with a constant Dmol parameter (self.Dmol in m^2/s). :param timestep: physical timestep of the current step (in s/cp). (needs to be converted before it is sent to `turbulent diffusion) :type timestep: :class:`float` :param Htot: Total heating rate (in W/kg) of all physical processes already computed :type Htot: :class:`array`, :class:`np.ndarray` :param atm: The Atm object used in the radiative transfer which contains many state variables. :type atm: :class:`Atm` object .. py:method:: condensation(timestep, Htot, verbose=False) This method computes the vapor and temperature tendencies do to large scale condensation in saturated layers. The tracer array in modified in place. :param timestep: physical timestep of the current step (in s/cp). :type timestep: :class:`float` :param Htot: Total heating rate (in W/kg) of all physical processes already computed :type Htot: :class:`array`, :class:`np.ndarray` :returns: H_cond: array, np.ndarray Heating rate due to large scale condensation (W/kg) .. py:method:: rainout(timestep, Htot, verbose=False) This method computes rainout. Condensates are carried down and reevaporated whenever there is "room" in an unsaturated layer. The option `evap_coeff` acts has an efficiency factor. `evap_coeff=1` is the efficient evaporation limit. When `evap_coeff<1` the maximum amount of condensates that can be reevaporated in a single layer is multiplied by `evap_coeff` All condensates are finaly evaporated in the last layer or when T > Tboil. The tracer array is modified in place. :param timestep: physical timestep of the current step (in s/cp). :type timestep: :class:`float` :param Htot: Total heating rate (in W/kg) of all physical processes already computed :type Htot: :class:`array`, :class:`np.ndarray` :returns: H_rain: array, np.ndarray Heating rate due to re evaporation (W/kg) .. py:method:: compute_hybrid_coordinates() Compute hybrid coordinates as in GCM. This will be used when surface pressure changes. Convention : sigma = (p-ptop)/(psurf-ptop) For each layer/level, the pressure is p = sigma * psurf + gamma .. py:method:: compute_mass_flux(dvapor_mass, sum_dvapor_mass) Computes the mass flux through the hybrid coordinate interfaces (kg/s/m^2; positive upward). (see Methods in Leconte et al. (2013)) W[k] is the mass flux between layer k-1 et k. :param sum_dvapor_mass: Total mass of vapor added to the atmosphere in the last timestep. :type sum_dvapor_mass: :class:`float` :param dvapor_mass: mass of vapor added to each layer. :type dvapor_mass: :class:`array`, :class:`np.ndarray` For the moment, W[0] = W[Nlay] .. py:method:: mass_redistribution(qarray_before_condensation) Update new mass and new pressure of a layer due to the evaporation and condensation of a given species, for more details see Methods, Leconte et al., 2013 (Nature) .. py:method:: radiative_acceleration(timestep=0.0, acceleration_mode=0, verbose=False, **kwargs) "Computes an acceleration factor and a new heating rate to speed up convergence :param acceleration_mode: 0: no acceleration 1 or 3: acceleration limited to radiative zones. 2 or 4: acceleration in convective zones as well. :type acceleration_mode: :class:`int` * 1: acceleration limited to radiative zones. The largest radiative timescale in non-radiative zones is used as reference radiative timsescale to compute acceleration. * 2: Same as mode 1 + acceleration in convective zones * 3 acceleration limited to radiative zones. The smallest radiative timescale in radiative zones is used as reference radiative timsescale to compute acceleration. * 4: Same as mode 3 + acceleration in convective zones .. py:method:: heating_rate(physical_process) Returns heating rates in W/kg per layer averaged over last call to evolve. Possible physical_processes are rad, cond, conv, rain, madj, tot .. py:method:: net_flux(physical_process) Returns net_flux in W/m^2 averaged over last call to evolve. Possible physical_processes are rad, cond, conv, rain, madj, tot .. py:method:: qsat(mol) Returns the saturation specific concentration of molecule mol (kg/kg) .. py:method:: write_pickle(filename, data_reduction_level=1) Saves the instance in a pickle file :param filename: Path to pickle file :type filename: :class:`str` :param data_reduction_level: Level of data to delete. 0: keep everything (results in big files). 1: removes some arrays, should not affect subsequent evolution. 2: removes the k and cia databases. The radiative model will need to be reset. This can be done with `set_options(k_database=, cia_database=, reset_rad_model=True)` after re-loading the `Atm_evolution` instance. :type data_reduction_level: :class:`int`