Add non-constant surface tension

[claresinger]

I'm interested in modifying this line with calculation of the Kelvin term to allow for variable surface tension (const.sgm not constant anymore). https://github.com/atmos-cloud-sim-uj/PySDM/blob/81243955ae257038c3a427d2618d32972fd1de02/PySDM/backends/numba/numba_helpers.py#L101

I want instead to replace it with an expression for the surface tension that allows for bulk-surface partitioning of surface-active organic species. Something along the lines of this compressed film model (https://doi.org/10.1126/science.aad4889) where the surface tension is a function of the wet radius, dry radius, organic fraction, and temperature. It seems like rw, rd, and T are available in this scope already (need to just pass r and rd to A(T)). I would then just need to add another attribute for f_org to describe the fraction of the aerosol particle that is organic. Do you foresee any issues doing this?

[slayoo]

We have an ongoing effort to enable either proper injection or singleton-based global setting for the constants and formulae catalogue. Stay tuned... :)

[slayoo]

OK, the relevant "machinery" is part of the https://github.com/atmos-cloud-sim-uj/PySDM/pull/489 A draft of compressed_film surface tension function is here: https://github.com/slayoo/PySDM/blob/koehler_refactor/PySDM/physics/surface_tension/compressed_film.py

@claresinger please suggest a minimal test case we could start with?

[claresinger]

@slayoo This is very exciting. Thank you!

I'm looping in fellow Caltech grad student @rxward (Ryan Ward). He has been doing relevant laboratory experiments and has some Julia code written here: https://github.com/rxward/ccn_study/blob/master/ccn_tools.jl. It's based on Forestieri et al (2018) and Ruehl et al (2016). Some of the examples in these papers are probably good test cases, but I will defer to Ryan on this.

[slayoo]

Hello @rxward! Getting closer to merging all this new stuff. To keep you updated, all the "physics" options such as the type of surface tension representation are now parameters to the PySDM.physics.Formulae class __init__ method. An instance of Formulae is then injected to the backend. This is then used wherever needed, and automagically propagates to GUIs, netCDF output, etc, so e.g. in the 2d-kinematic GUI you get: and in the netCDF output:

$ ncdump -c ~/Downloads/tmpewr6ur_s.nc | grep -2 "global attr"
        Z:units = "metres" ;

// global attributes:
        :dt = 5.f ;
        :formulae = "condensation_coordinate: VolumeLogarithm, diffusion_kinetics: FuchsSutugin, drop_growth: MaxwellMason, fastmath: bool, hydrostatics: Default, hygroscopicity: KappaKoehlerLeadingTerms, latent_heat: Kirchhoff, particle_advection: ExplicitInSpace, saturation_vapour_pressure: FlatauWalkoCotton, state_variable_triplet: RhodThdQv, surface_tension: Constant, trivia: Trivia, ventilation: Neglect" ;

[claresinger]

Hi @slayoo! The GUI looks great. After talking with @rxward, we think the parcel example from Lowe et al. (2019) Fig 2 / Table 1 would be a great test case and example to implement in the code.

In terms of changes to the code. I think we’ll need to add an attribute of the superdroplets f_org, the fraction of organics in the aerosol. We will also need to specify what type of organic molecule(s) by specifying the delta_min and sig_org chemical parameters for the compressed film model. We propose carrying these parameters in a structure called chem, but you can decide how is best to do this. Finally, we propose renaming sgm in constants.py to sgm_water.

A rough sketch of the sigma function is here:

# Function to implement Lowe et al. (2019) simplified compressed film framework.
def sigma(T, v_wet, v_dry, f_org):

    # chemical parameters
    delta_min = chem.delta_min  # minimum organic film thickness [nm]
    sgm_org = chem.sgm_org      # organic surface tension [N m-1]

    # convert volumes to diameters
    d_coat = 2 * ((3*v_dry) / (4*np.pi))**(1/3)
    d_seed = (d_coat**3 - f_org * d_coat**3)**(1/3)
    d_wet = 2 * ((3*v_wet) / (4*np.pi))**(1/3)

    # calculate the total volume of organic
    v_beta = v_dry - np.pi/6 * d_seed**3

    # calculate the min shell volume, v_delta by Eq. (S3.4.3)
    v_delta = np.pi/6 * (d_wet**3 - (d_wet**3 - 2*delta_min)**3)

    # calculate the coverage parameter using Eq. (S3.4.2)
    c_beta = np.min(v_beta/v_delta, 1)

    # calculate sigma
    sgm = (1-c_beta)*const.sgm_water + c_beta*sgm_org

    return sgm

[rxward]

Hi @slayoo! Just wanted to provide some updates here. Specifically, we want to add some information into the chem struct that we proposed for carrying chemical parameters. It's likely that we want to have multiple species, and weren't sure if it would be best to hold all of this info in a single class or if multiple (one per species) would be better. An example could look something like:

# define the attrs for each species
class Species:
    def __init__(self, delta_min, sigma_org):
        self.delta_min = delta_min
        self.sigma_org = sigma_org

# add each species to the chem dict as a class      
chem = {
'species1':Species(val1,val2),
'species2':Species(val3,val4)
}

# call any species/attribute combo
chem['species1'].delta_min

Of course, there are many ways to do this, so this I just an example. Also, we were wondering where f_org would be stored as an additional attribute. This is a physical property of the aerosol rather than a chemical property, however (i.e, what fraction of the aerosol is made of organics vs inorganic), so perhaps it would need to be carried somewhere else.

Finally, it bares mentioning that this compressed film framework is only an approximation, presented by the Lowe et al. work. We will adapt it to be a bit more computationally robust once we've tested against the Lowe paper. This should be trivial though.

[slayoo]

@rxward @claresinger Apologies for taking so long to reply, and big big thanks for providing all the above information! I'll work on it within the next few days.

[slayoo]

Re f_org, elaborating on what I've just wrote in a comment to #570, I think we will need to have two attributes: inorganic mass and organic mass within each droplet - to make it work well with the collision handling code. These would be both based on the ExtensiveAttribute class. The f_org would then inherit from DerivedAttribute. I'm very sorry for not replying promptly here, and big thanks for working on it!

[slayoo]

One more thought: to make condensation "feel" the f_org, it will likely be best to make a similar dry radius handling as for aqueous chemistry where the dry radius is computed on the fly based on other attributes defining the composition of the particle, see: https://github.com/atmos-cloud-sim-uj/PySDM/blob/83649d5af8a949f355debbb8743765c520810cc7/PySDM/attributes/impl/mapper.py#L23

[slayoo]

One another issue is how to initialise the organic/inorganic fraction. In the aqueous chemistry example, the initialisation is done e.g. here: https://github.com/atmos-cloud-sim-uj/PySDM-examples/blob/18c23c4d2a3c1394f81259067e57914004879a96/PySDM_examples/Kreidenweis_et_al_2003/simulation.py#L24

[slayoo]

concerning, the proposed Species class, IIUC the needs here, all the information should best be carried by particle attributes extensive or derived ones

[slayoo]

Finally, I guess the best way to start would be to draft a new example using the Parcel environment and try adding new attributes step by step and plotting their values using new or existing products. How to best iterate @rxward? Perhaps we could set up a screen-sharing call and draft the new example notebook together?

[slayoo]

Well, given #586 and https://github.com/atmos-cloud-sim-uj/PySDM-examples/pull/42 are almost done, let me close this one. Thanks a lot @claresinger and @rxward !