<commit_before>""" Collection of activation parameterizations.

"""

import numpy as np
from scipy.special import erfc

from thermo import ( es, rho_air, ka, ka_cont, dv, dv_cont,
                     sigma_w, kohler_crit )
import constants as c

def _unpack_aerosols(aerosols):
    """ Convert a list of :class:`AerosolSpecies` into lists of aerosol properties.

    Parameters
    ----------
    aerosols : list of :class:`AerosolSpecies`

    Returns
    -------
    dictionary of lists of aerosol properties

    """

    species, mus, sigmas, kappas, Ns = [], [], [], [], []
    for a in aerosols:
        species.append(a.species)
        mus.append(a.distribution.mu)
        sigmas.append(a.distribution.sigma)
        Ns.append(a.distribution.N)
        kappas.append(a.kappa)

    species = np.asarray(species)
    mus     = np.asarray(mus)
    sigmas  = np.asarray(sigmas)
    kappas  = np.asarray(kappas)
    Ns      = np.asarray(Ns)

    return dict(species=species, mus=mus, sigmas=sigmas, Ns=Ns, kappas=kappas)

def lognormal_activation(smax, mu, sigma, N, kappa, sgi=None, T=None, approx=True):
    """ Compute the activated number/fraction from a lognormal mode

    Parameters
    ----------
    smax : float
        Maximum parcel supersaturation
    mu, sigma, N : floats
        Lognormal mode parameters; ``mu`` should be in meters
    kappa : float
        Hygroscopicity of material in aerosol mode
    sgi :float, optional
        Modal critical supersaturation; if not provided, this method will
        go ahead and compute them, but a temperature ``T`` must also be passed
    T : float, optional
        Parcel temperature; only necessary if no ``sgi`` was passed
    approx : boolean, optional (default=False)
        If computing modal critical supersaturations, use the approximated
        Kohler theory

    Returns
    -------
    N_act, act_frac : floats
        Activated number concentration and fraction for the given mode

    """

    if not sgi:
        assert T
        _, sgi = kohler_crit(T, mu, kappa, approx)

    ui = 2.*np.log(sgi/smax)/(3.*np.sqrt(2.)*np.log(sigma))
    N_act = 0.5*N*erfc(ui)
    act_frac = N_act/N

    return N_act, act_frac

def binned_activation(Smax, T, rs, aerosol):
    """ Compute the activation statistics of a given aerosol, its transient
    size distribution, and updraft characteristics. Following Nenes et al, 2001
    also compute the kinetic limitation statistics for the aerosol.

    Parameters
    ----------
    Smax : float
        Environmental maximum supersaturation.
    T : float
        Environmental temperature.
    rs : array of floats
        Wet radii of aerosol/droplet population.
    aerosol : :class:`AerosolSpecies`
        The characterization of the dry aerosol.

    Returns
    -------
    eq, kn: floats
        Activated fractions
    alpha : float
        N_kn / N_eq
    phi : float
        N_unact / N_kn

    """
    kappa = aerosol.kappa
    r_drys = aerosol.r_drys
    Nis = aerosol.Nis
    N_tot = np.sum(Nis)

    # Coerce the reference wet droplet sizes (rs) to an array if they were passed as
    # a Series
    if hasattr(rs, 'values'):
        rs = rs.values

    r_crits, s_crits = zip(*[kohler_crit(T, r_dry, kappa, True) for r_dry in r_drys])
    s_crits = np.array(s_crits)
    r_crits = np.array(r_crits)

    ## Equilibrium calculation - all aerosol whose critical supersaturation is
    ##                           less than the environmental supersaturation
    activated_eq = (Smax >= s_crits)
    N_eq = np.sum(Nis[activated_eq])
    eq_frac = N_eq/N_tot

    ## Kinetic calculation - find the aerosol with the smallest dry radius which has
    ##                       grown past its critical radius, and count this aerosol and
    ##                       all larger ones as activated. This will include some large
    ##                       particles which haven't grown to critical size yet.
    is_r_large = rs >= r_crits
    if not np.any(is_r_large):
        N_kn, kn_frac = 0., 0.
        phi = 1.0
    else:
        smallest_ind = np.where(is_r_large)[0][0]
        N_kn = np.sum(Nis[smallest_ind:])
        kn_frac = N_kn/N_tot

        ## Unactivated - all droplets smaller than their critical size
        droplets = range(smallest_ind, len(Nis))
        Nis_drops = Nis[droplets]
        r_crits_drops = r_crits[droplets]
        rs_drops = rs[droplets]
        too_small = (rs_drops < r_crits_drops)
<commit_msg>Adding user-control of 'approx' flag for kohler theory in binned_activation<commit_after>""" Collection of activation parameterizations.

"""

import numpy as np
from scipy.special import erfc

from thermo import ( es, rho_air, ka, ka_cont, dv, dv_cont,
                     sigma_w, kohler_crit )
import constants as c

def _unpack_aerosols(aerosols):
    """ Convert a list of :class:`AerosolSpecies` into lists of aerosol properties.

    Parameters
    ----------
    aerosols : list of :class:`AerosolSpecies`

    Returns
    -------
    dictionary of lists of aerosol properties

    """

    species, mus, sigmas, kappas, Ns = [], [], [], [], []
    for a in aerosols:
        species.append(a.species)
        mus.append(a.distribution.mu)
        sigmas.append(a.distribution.sigma)
        Ns.append(a.distribution.N)
        kappas.append(a.kappa)

    species = np.asarray(species)
    mus     = np.asarray(mus)
    sigmas  = np.asarray(sigmas)
    kappas  = np.asarray(kappas)
    Ns      = np.asarray(Ns)

    return dict(species=species, mus=mus, sigmas=sigmas, Ns=Ns, kappas=kappas)

def lognormal_activation(smax, mu, sigma, N, kappa, sgi=None, T=None, approx=True):
    """ Compute the activated number/fraction from a lognormal mode

    Parameters
    ----------
    smax : float
        Maximum parcel supersaturation
    mu, sigma, N : floats
        Lognormal mode parameters; ``mu`` should be in meters
    kappa : float
        Hygroscopicity of material in aerosol mode
    sgi :float, optional
        Modal critical supersaturation; if not provided, this method will
        go ahead and compute them, but a temperature ``T`` must also be passed
    T : float, optional
        Parcel temperature; only necessary if no ``sgi`` was passed
    approx : boolean, optional (default=False)
        If computing modal critical supersaturations, use the approximated
        Kohler theory

    Returns
    -------
    N_act, act_frac : floats
        Activated number concentration and fraction for the given mode

    """

    if not sgi:
        assert T
        _, sgi = kohler_crit(T, mu, kappa, approx)

    ui = 2.*np.log(sgi/smax)/(3.*np.sqrt(2.)*np.log(sigma))
    N_act = 0.5*N*erfc(ui)
    act_frac = N_act/N

    return N_act, act_frac

def binned_activation(Smax, T, rs, aerosol, approx=False):
    """ Compute the activation statistics of a given aerosol, its transient
    size distribution, and updraft characteristics. Following Nenes et al, 2001
    also compute the kinetic limitation statistics for the aerosol.

    Parameters
    ----------
    Smax : float
        Environmental maximum supersaturation.
    T : float
        Environmental temperature.
    rs : array of floats
        Wet radii of aerosol/droplet population.
    aerosol : :class:`AerosolSpecies`
        The characterization of the dry aerosol.
    approx : boolean
        Approximate Kohler theory rather than include detailed calculation (default False)

    Returns
    -------
    eq, kn: floats
        Activated fractions
    alpha : float
        N_kn / N_eq
    phi : float
        N_unact / N_kn

    """
    kappa = aerosol.kappa
    r_drys = aerosol.r_drys
    Nis = aerosol.Nis
    N_tot = np.sum(Nis)

    # Coerce the reference wet droplet sizes (rs) to an array if they were passed as
    # a Series
    if hasattr(rs, 'values'):
        rs = rs.values

    r_crits, s_crits = zip(*[kohler_crit(T, r_dry, kappa, approx) for r_dry in r_drys])
    s_crits = np.array(s_crits)
    r_crits = np.array(r_crits)

    ## Equilibrium calculation - all aerosol whose critical supersaturation is
    ##                           less than the environmental supersaturation
    activated_eq = (Smax >= s_crits)
    N_eq = np.sum(Nis[activated_eq])
    eq_frac = N_eq/N_tot

    ## Kinetic calculation - find the aerosol with the smallest dry radius which has
    ##                       grown past its critical radius, and count this aerosol and
    ##                       all larger ones as activated. This will include some large
    ##                       particles which haven't grown to critical size yet.
    is_r_large = rs >= r_crits
    if not np.any(is_r_large):
        N_kn, kn_frac = 0., 0.
        phi = 1.0
    else:
        smallest_ind = np.where(is_r_large)[0][0]
        N_kn = np.sum(Nis[smallest_ind:])
        kn_frac = N_kn/N_tot

        ## Unactivated - all droplets smaller than their critical size
        droplets = range(smallest_ind, len(Nis))
        Nis_drops = Nis[droplets]
        r_crits_drops = r_crits[droplets]
        rs_drops = rs[droplets]
        too_small = (rs_drops < r_crits_drops)

{'fraction_non_alphanumeric': 0.09260478384902349, 'fraction_numerical': 0.0032916392363396972, 'mean_word_length': 2.3647102251753416, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4196965', 'n_tokens_mistral': 2974, 'n_tokens_neox': 2751, 'n_words': 874}