The Code

Submodules

virga.calc_mie module

virga.calc_mie.calc_new_mieff(wave_in, nn, kk, radius, rup, fort_calc_mie=False)
virga.calc_mie.fort_mie_calc(RO, RFR, RFI, THET, JX, R, RE2, TMAG2, WVNO)

Given the refractive indices at a certain wavelength this module calculates the Mie scattering by a stratified sphere.The basic code used was that described in the report: ” Subroutines for computing the parameters of the electromagnetic radiation scattered by a sphere ” J.V. Dave, I B M Scientific Center, Palo Alto , California. Report NO. 320 - 3236 .. MAY 1968 .

Parameters:
  • RO (float) – Outer Shell Radius (cm)

  • RFR (float) – Real refractive index of shell layer (in the form n= RFR-i*RFI)

  • RFI (float) – Imaginary refractive index of shell layer (in the form n= RFR-i*RFI)

  • THET (ndarray) – Angle in degrees between the directions of the incident and the scattered radiation.

  • JX (integer) – Total number of THET for which calculations are required

  • R (float) – Radius of core (cm)`

  • RE2 (float) – Real refractive index of core (in the form n= RE2-i*TMAG2)

  • TMAG2 (float) – Imaginary refractive index of core (in the form n= RE2-i*TMAG2)

  • WVNO (float) – Wave-number corresponding to the wavelength. (cm^-1)

Returns:

  • QEXT (float) – Efficiency factor for extinction,VAN DE HULST,P.14 ‘ 127

  • QSCAT (float) – Efficiency factor for scattering,VAN DE HULST,P.14 ‘ 127

  • CTBQRS (float) – Average(cos(theta))*QSCAT,VAN DE HULST,P.14 ‘ 127

  • ISTATUS (integer) – Convergence indicator, 0 if converged, -1 if otherwise.

virga.calc_mie.get_r_grid(r_min=1e-05, n_radii=40)

Get spacing of radii to run Mie code

r_minfloat

Minimum radius to compute (cm)

n_radiiint

Number of radii to compute

virga.calc_mie.get_refrind(igas, directory)

Reads reference files with wavelength, and refractory indecies. This function relies on input files being structured as a 4 column file with columns: index, wavelength (micron), nn, kk

Parameters:
  • igas (str) – Gas name

  • directory (str) – Directory were reference files are located.

virga.gas_properties module

virga.gas_properties.Al2O3(mw_atmos, mh=1, gas_mmr=None)

Defines properties for Al2O3 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 2.51e-6

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • Mean molecular weight of gas,

  • Gas mass mixing ratio

  • Density of gas cgs

Notes

virga.gas_properties.CH4(mw_atmos, mh=1, gas_mmr=None)

Defines properties for CH4 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio. None points to the default value of : 4.9e-4

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • mean molecular weight of gas,

  • gas mass mixing ratio

  • density of gas cgs

virga.gas_properties.CaAl12O19(mw_atmos, mh=1, gas_mmr=None)

Defines properties for CaAl12O19 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 8.80e-7

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • Mean molecular weight of gas,

  • Gas mass mixing ratio

  • Density of gas cgs

virga.gas_properties.CaTiO3(mw_atmos, mh=1, gas_mmr=None)

Defines properties for CaTiO3 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 2.51e-6

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • Mean molecular weight of gas,

  • Gas mass mixing ratio

  • Density of gas cgs

virga.gas_properties.Cr(mw_atmos, mh=1, gas_mmr=None)

Defines properties for Cr as condensible

mw_atmosfloat

Mean molecular weight of the atmosphere amu

gas_mmrfloat , optional

Gas mass mixing ratio None points to the default value of : 1xSolar = 8.87e-7 10xSolar = 8.6803E-06 50xSolar = 4.1308E-05

mhfloat , optional

Metallicity, Default is 1=1xSolar

Mean molecular weight of gas, Gas mass mixing ratio Density of gas cgs

Proc. 16, 379 (2010)

virga.gas_properties.Fe(mw_atmos, mh=1, gas_mmr=None)

Defines properties for Fe as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 5.78e-5

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • mean molecular weight of gas,

  • gas mass mixing ratio

  • density of gas cgs

virga.gas_properties.H2O(mw_atmos, mh=1, gas_mmr=None)

Defines properties for H2O as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 7.54e-4

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • mean molecular weight of gas,

  • gas mass mixing ratio

  • density of gas cgs

virga.gas_properties.KCl(mw_atmos, mh=1, gas_mmr=None)

Defines properties for KCl as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 1xSolar = 2.55E-07 10xSolar = 2.1829E-06 50xSolar = 8.1164E-06

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • Mean molecular weight of gas,

  • Gas mass mixing ratio

  • Density of gas cgs

virga.gas_properties.Mg2SiO4(mw_atmos, mh=1, gas_mmr=None)

Defines properties for Mg2SiO4 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 3.58125e-05

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • Mean molecular weight of gas,

  • Gas mass mixing ratio

  • Density of gas cgs

virga.gas_properties.MgSiO3(mw_atmos, mh=1, gas_mmr=None)

Defines properties for MgSiO3 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 2.75e-3

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • Mean molecular weight of gas,

  • Gas mass mixing ratio

  • Density of gas cgs

virga.gas_properties.MnS(mw_atmos, mh=1, gas_mmr=None)

Defines properties for MnS as condensible

mw_atmosfloat

Mean molecular weight of the atmosphere amu

gas_mmrfloat , optional

Gas mass mixing ratio None points to the default value of : 6.32e-7

mhfloat , optional

Metallicity, Default is 1=1xSolar

Mean molecular weight of gas, Gas mass mixing ratio Density of gas cgs

Proc. 16, 379 (2010)

virga.gas_properties.NH3(mw_atmos, mh=1, gas_mmr=None)

Defines properties for NH3 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio None points to the default value of : 1.34e-4

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • mean molecular weight of gas,

  • gas mass mixing ratio

  • density of gas cgs

virga.gas_properties.Na2S(mw_atmos, mh=1, gas_mmr=None)

Defines properties for Na2S as condensible

mw_atmosfloat

Mean molecular weight of the atmosphere amu

gas_mmrfloat , optional

Gas mass mixing ratio None points to the default value of : 3.97e-6

mhfloat , optional

Metallicity, Default is 1=1xSolar

Mean molecular weight of gas, Gas mass mixing ratio Density of gas cgs

Proc. 16, 379 (2010)

virga.gas_properties.SiO2(mw_atmos, mh=1, gas_mmr=None)

Defines properties for TiO2 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio. None points to the default value of : 1.69e-7

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • mean molecular weight of gas,

  • gas mass mixing ratio

  • density of gas cgs

virga.gas_properties.TiO2(mw_atmos, mh=1, gas_mmr=None)

Defines properties for TiO2 as condensible

Parameters:
  • mw_atmos (float) – Mean molecular weight of the atmosphere amu

  • gas_mmr (float , optional) – Gas mass mixing ratio. None points to the default value of : 1.69e-7

  • mh (float , optional) – Metallicity, Default is 1=1xSolar

Returns:

  • mean molecular weight of gas,

  • gas mass mixing ratio

  • density of gas cgs

virga.gas_properties.ZnS(mw_atmos, mh=1, gas_mmr=None)

Defines properties for ZnS as condensible

mw_atmosfloat

Mean molecular weight of the atmosphere amu

gas_mmrfloat , optional

Gas mass mixing ratio None points to the default value of : 8.40e-8

mhfloat , optional

Metallicity, Default is 1=1xSolar

Mean molecular weight of gas, Gas mass mixing ratio Density of gas cgs

Proc. 16, 379 (2010)

virga.justdoit module

class virga.justdoit.Atmosphere(condensibles, fsed=0.5, b=1, eps=0.01, mh=1, mmw=2.2, sig=2.0, param='const', verbose=False, supsat=0, gas_mmr=None)

Bases: object

compute(directory=None, as_dict=True)
Parameters:
  • atmo (class) – Atmosphere class

  • directory (str, optional) – Directory string that describes where refrind files are

  • as_dict (bool) – Default = True, option to view full output as dictionary

Returns:

  • dict – When as_dict=True. Dictionary output that contains full output. See tutorials for explanation of all output.

  • opd, w0, g0 – Extinction per layer, single scattering abledo, asymmetry parameter, All are ndarrays that are nlayer by nwave

constants()
get_atmo_parameters()

Defines all the atmospheric parameters needed for the calculation

Note: Some of this is repeated in the layer() function. This is done on purpose because layer is used to get a “virtual” layer off the user input’s grid. These parameters are also needed, though, to get the initial mixing parameters from the user. Therefore, we put them in both places.

get_kz_mixl(df, constant_kz, latent_heat, convective_overshoot, kz_min)

Computes kz profile and mixing length given user input. In brief the options are:

  1. Input Kz

  2. Input constant kz

3) Input convective heat flux (supply chf in df) 3a) Input convective heat flux, correct for latent heat (supply chf in df and set latent_heat=True) and/or 3b) Input convective heat flux, correct for convective overshoot (supply chf, convective_overshoot=1/3) 4) Set kz_min to prevent kz from going too low (any of the above and set kz_min~1e5)

Parameters:
  • df (dataframe or dict) – Dataframe from input with “pressure”(bars),”temperature”(K). MUST have at least two columns with names “pressure” and “temperature”. Optional columns include the eddy diffusion “kz” in cm^2/s CGS units, and the convective heat flux ‘chf’ also in cgs (e.g. sigma_csg T^4)

  • constant_kz (float) – Constant value for kz, if kz is supplied in df or filename, it will inheret that value and not use this constant_value Default = None

  • latent_heat (bool) – optional, Default = False. The latent heat factors into the mixing length. When False, the mixing length goes as the scale height When True, the mixing length is scaled by the latent heat

  • convective_overshoot (float) – Optional, Default is None. But the default value used in Ackerman & Marley 2001 is 1./3. If you are unsure of what to pick, start there. This is only used when the This is ONLY used when a chf (convective heat flux) is supplied

  • kz_min (float) – Minimum Kz value. This will reset everything below kz_min to kz_min. Default = 1e5 cm2/s

gravity(gravity=None, gravity_unit=None, radius=None, radius_unit=None, mass=None, mass_unit=None)

Get gravity based on mass and radius, or gravity inputs

Parameters:
  • gravity (float) – (Optional) Gravity of planet

  • gravity_unit (astropy.unit) – (Optional) Unit of Gravity

  • radius (float) – (Optional) radius of planet MUST be specified for thermal emission!

  • radius_unit (astropy.unit) – (Optional) Unit of radius

  • mass (float) – (Optional) mass of planet

  • mass_unit (astropy.unit) – (Optional) Unit of mass

kz(df=None, constant_kz=None, chf=None, kz_min=100000.0, latent_heat=False)

Define Kz in CGS. Should be on same grid as pressure. This overwrites whatever was defined in get_pt ! Users can define kz by:

  1. Defining a DataFrame with keys ‘pressure’ (in bars), and ‘kz’

  2. Defining constant kz

  3. Supplying a convective heat flux and prescription for latent_heat

Parameters:
  • df (pandas.DataFrame, dict) – Dataframe or dictionary with ‘kz’ as one of the fields.

  • constant_kz (float) – Constant value for kz in units of cm^2/s

  • chf (ndarray) – Convective heat flux in cgs units (e.g. sigma T^4). This will be used to compute the kzz using the methodology of Gierasch and Conrath (1985)

  • latent_heat (bool) – optional, Default = False. The latent heat factors into the mixing length. When False, the mixing length goes as the scale height When True, the mixing length is scaled by the latent heat This is ONLY used when a chf (convective heat flux) is supplied

ptk(df=None, filename=None, kz_min=100000.0, constant_kz=None, latent_heat=False, convective_overshoot=None, Teff=None, alpha_pressure=None, **pd_kwargs)

Read in file or define dataframe.

Parameters:
  • df (dataframe or dict) – Dataframe with “pressure”(bars),”temperature”(K). MUST have at least two columns with names “pressure” and “temperature”. Optional columns include the eddy diffusion “kz” in cm^2/s CGS units, and the convective heat flux ‘chf’ also in cgs (e.g. sigma_csg T^4)

  • filename (str) – Filename read in. Will be read in with pd.read_csv and should result in two named headers “pressure”(bars),”temperature”(K). Optional columns include the eddy diffusion “kz” in cm^2/s CGS units, and the convective heat flux ‘chf’ also in cgs (e.g. sigma_csg T^4) Use pd_kwargs to ensure file is read in properly.

  • kz_min (float, optional) – Minimum Kz value. This will reset everything below kz_min to kz_min. Default = 1e5 cm2/s

  • constant_kz (float, optional) – Constant value for kz, if kz is supplied in df or filename, it will inheret that value and not use this constant_value Default = None

  • latent_heat (bool) – optional, Default = False. The latent heat factors into the mixing length. When False, the mixing length goes as the scale height When True, the mixing length is scaled by the latent heat

  • convective_overshoot (float) – Optional, Default is None. But the default value used in Ackerman & Marley 2001 is 1./3. If you are unsure of what to pick, start there. This is only used when the This is ONLY used when a chf (convective heat flux) is supplied

  • Teff (float, optional) – Effective temperature. If None (default), Teff set to temperature at 1 bar

  • alpha_pressure (float) – Pressure at which we want fsed=alpha for variable fsed calculation

  • pd_kwargs (kwargs) – Pandas key words for file read in. If reading old style eddysed files, you would need: skiprows=3, sep=’s+’, header=None, names=[“ind”,”pressure”,”temperature”,”kz”]

virga.justdoit.available()

Print all available gas condensates

virga.justdoit.brown_dwarf()
virga.justdoit.calc_mie_db(gas_name, dir_refrind, dir_out, rmin=1e-08, rmax=0.054239131, nradii=60, fort_calc_mie=False)

Function that calculations new Mie database using PyMieScatt. :param gas_name: List of names of gasses. Or a single gas name.

See pyeddy.available() to see which ones are currently available.

Parameters:
  • dir_refrind (str) – Directory where you store optical refractive index files that will be created.

  • dir_out (str) – Directory where you want to store Mie parameter files. Will be stored as gas_name.Mieff. BEWARE FILE OVERWRITES.

  • rmin (float , optional) – (Default=1e-8) Units of cm. The minimum radius to compute Mie parameters for. Usually 0.1 microns is small enough. However, if you notice your mean particle radius is on the low end, you may compute your grid to even lower particle sizes.

  • nradii (int, optional) – (Default=60) number of radii points to compute grid on. 40 grid points for exoplanets/BDs is generally sufficient.

Returns:

  • Q extinction, Q scattering, asymmetry * Q scattering, radius grid (cm), wavelength grid (um)

  • The Q “efficiency factors” are = cross section / geometric cross section of particle

virga.justdoit.calc_optics(nwave, qc, qt, rg, reff, ndz, radius, dr, qext, qscat, cos_qscat, sig, rmin, nrad, verbose=False)

Calculate spectrally-resolved profiles of optical depth, single-scattering albedo, and asymmetry parameter.

Parameters:
  • nwave (int) – Number of wave points

  • qc (ndarray) – Condensate mixing ratio

  • qt (ndarray) – Gas + condensate mixing ratio

  • rg (ndarray) – Geometric mean radius of condensate

  • reff (ndarray) – Effective (area-weighted) radius of condensate (cm)

  • ndz (ndarray) – Column density of particle concentration in layer (#/cm^2)

  • radius (ndarray) – Radius bin centers (cm)

  • dr (ndarray) – Width of radius bins (cm)

  • qscat (ndarray) – Scattering efficiency

  • qext (ndarray) – Extinction efficiency

  • cos_qscat (ndarray) – qscat-weighted <cos (scattering angle)>

  • sig (float) – Width of the log normal particle distribution

  • verbose (bool) – print out warnings or not

Returns:

  • opd (ndarray) – extinction optical depth due to all condensates in layer

  • w0 (ndarray) – single scattering albedo

  • g0 (ndarray) – asymmetry parameter = Q_scat wtd avg of <cos theta>

  • opd_gas (ndarray) – cumulative (from top) opd by condensing vapor as geometric conservative scatterers

virga.justdoit.calc_optics_user_r_dist(wave_in, ndz, radius, radius_unit, r_distribution, qext, qscat, cos_qscat, verbose=False)

Calculate spectrally-resolved profiles of optical depth, single-scattering albedo, and asymmetry parameter for a user-input particle radius distribution

Parameters:
  • wave_in (ndarray) – your wavelength grid in microns

  • ndz (float) –

    Column density of total particle concentration (#/cm^2)

    Note: set to whatever, it’s your free knob —- this does not directly translate to something physical because it’s for all particles in your slab May have to use values of 1e8 or so

  • radius (ndarray) – Radius bin values - the range of particle sizes of interest. Maybe measured in the lab, Ensure radius_unit is specified

  • radius_unit (astropy.unit.Units) – Astropy compatible unit

  • qscat (ndarray) – Scattering efficiency

  • qext (ndarray) – Extinction efficiency

  • cos_qscat (ndarray) – qscat-weighted <cos (scattering angle)>

  • r_distribution (ndarray) – the radius distribution in each bin. Maybe measured from the lab, generated from microphysics, etc. Should integrate to 1.

  • verbose (bool) – print out warnings or not

Returns:

  • opd (ndarray) – extinction optical depth due to all condensates in layer

  • w0 (ndarray) – single scattering albedo

  • g0 (ndarray) – asymmetry parameter = Q_scat wtd avg of <cos theta>

virga.justdoit.calc_qc(gas_name, supsat, t_layer, p_layer, r_atmos, r_cloud, q_below, mixl, dz_layer, gravity, mw_atmos, mfp, visc, rho_p, w_convect, fsed, b, eps, param, z_bot, z_layer, z_alpha, z_min, sig, mh, rmin, nrad, og_vfall=True, z_cld=None)

Calculate condensate optical depth and effective radius for a layer, assuming geometric scatterers.

gas_namestr

Name of condensate

supsatfloat

Super saturation factor

t_layerfloat

Temperature of layer mid-pt (K)

p_layerfloat

Pressure of layer mid-pt (dyne/cm^2)

r_atmosfloat

specific gas constant for atmosphere (erg/K/g)

r_cloudfloat

specific gas constant for cloud species (erg/K/g)

q_belowfloat

total mixing ratio (vapor+condensate) below layer (g/g)

mxlfloat

convective mixing length scale (cm): no less than 1/10 scale height

dz_layerfloat

Thickness of layer cm

gravityfloat

Gravity of planet cgs

mw_atmosfloat

Molecular weight of the atmosphere

mfpfloat

atmospheric mean free path (cm)

viscfloat

atmospheric viscosity (dyne s/cm^2)

rho_pfloat

density of condensed vapor (g/cm^3)

w_convectfloat

convective velocity scale (cm/s)

fsedfloat

Sedimentation efficiency coefficient (unitless)

bfloat

Denominator of exponential in sedimentation efficiency (if param is ‘exp’)

eps: float

Minimum value of fsed function (if param=exp)

z_botfloat

Altitude at bottom of layer

z_layerfloat

Altitude of midpoint of layer (cm)

z_alphafloat

Altitude at which fsed=alpha for variable fsed calculation

paramstr

fsed parameterisation ‘const’ (constant), ‘exp’ (exponential density derivation)

sigfloat

Width of the log normal particle distrubtion

mhfloat

Metallicity NON log solar (1 = 1x solar)

rminfloat

Minium radius on grid (cm)

nradint

Number of radii on Mie grid

Returns:

  • qt_top (float) – gas + condensate mixing ratio at top of layer(g/g)

  • qc_layer (float) – condenstate mixing ratio (g/g)

  • qt_layer (float) – gas + condensate mixing ratio (g/g)

  • rg_layer (float) – geometric mean radius of condensate cm

  • reff_layer (float) – droplet effective radius (second moment of size distrib, cm)

  • ndz_layer (float) – number column density of condensate (cm^-3)

virga.justdoit.compute(atmo, directory=None, as_dict=True, og_solver=True, direct_tol=1e-15, refine_TP=True, og_vfall=True, analytical_rg=True, do_virtual=True)

Top level program to run eddysed. Requires running Atmosphere class before running this.

Parameters:
  • atmo (class) – Atmosphere class

  • directory (str, optional) – Directory string that describes where refrind files are

  • as_dict (bool, optional) – Default = False. Option to view full output as dictionary

  • og_solver (bool, optional) –

    Default=True. BETA. Contact developers before changing to False.

    Option to change mmr solver (True = original eddysed, False = new direct solver)

  • direct_tol (float , optional) – Only used if og_solver =False. Default = True. Tolerance for direct solver

  • refine_TP (bool, optional) – Only used if og_solver =False. Option to refine temperature-pressure profile for direct solver

  • og_vfall (bool, optional) – Option to use original A&M or new Khan-Richardson method for finding vfall

  • analytical_rg (bool, optional) – Only used if og_solver =False. Option to use analytical expression for rg, or alternatively deduce rg from calculation Calculation option will be most useful for future inclusions of alternative particle size distributions

  • do_virtual (bool) – If the user adds an upper bound pressure that is too low. There are cases where a cloud wants to form off the grid towards higher pressures. This enables that.

Returns:

  • opd, w0, g0 – Extinction per layer, single scattering abledo, asymmetry parameter, All are ndarrays that are nlayer by nwave

  • dict – Dictionary output that contains full output. See tutorials for explanation of all output.

virga.justdoit.condensation_t(gas_name, mh, mmw, pressure=array([1.00000000e-06, 2.63665090e-06, 6.95192796e-06, 1.83298071e-05, 4.83293024e-05, 1.27427499e-04, 3.35981829e-04, 8.85866790e-04, 2.33572147e-03, 6.15848211e-03, 1.62377674e-02, 4.28133240e-02, 1.12883789e-01, 2.97635144e-01, 7.84759970e-01, 2.06913808e+00, 5.45559478e+00, 1.43844989e+01, 3.79269019e+01, 1.00000000e+02]), gas_mmr=None)

Find condensation curve for any planet given a pressure. These are computed based on pressure vapor curves defined in pvaps.py.

Default is to compute condensation temperature on a pressure grid

Parameters:
  • gas_name (str) – Name of gas, which is case sensitive. See print_available to see which gases are available.

  • mh (float) – Metallicity in NOT log units. Solar =1

  • mmw (float) – Mean molecular weight of the atmosphere. Solar = 2.2

  • pressure (ndarray, list, float, optional) – Grid of pressures (bars) to compute condensation temperatures on. Default = np.logspace(-3,2,20)

Returns:

pressure (bars), condensation temperature (Kelvin)

Return type:

ndarray, ndarray

virga.justdoit.create_dict(qc, qt, rg, reff, ndz, opd, w0, g0, opd_gas, wave, pressure, temperature, gas_names, mh, mmw, fsed, sig, nrad, rmin, z, dz_layer, mixl, kz, scale_h, z_cld)
virga.justdoit.eddysed(t_top, p_top, t_mid, p_mid, condensibles, gas_mw, gas_mmr, rho_p, mw_atmos, gravity, kz, mixl, fsed, b, eps, scale_h, z_top, z_alpha, z_min, param, mh, sig, rmin, nrad, d_molecule, eps_k, c_p_factor, og_vfall=True, do_virtual=True, supsat=0, verbose=False)

Given an atmosphere and condensates, calculate size and concentration of condensates in balance between eddy diffusion and sedimentation.

Parameters:
  • t_top (ndarray) – Temperature at each layer (K)

  • p_top (ndarray) – Pressure at each layer (dyn/cm^2)

  • t_mid (ndarray) – Temperature at each midpoint (K)

  • p_mid (ndarray) – Pressure at each midpoint (dyn/cm^2)

  • condensibles (ndarray or list of str) – List or array of condensible gas names

  • gas_mw (ndarray) – Array of gas mean molecular weight from gas_properties

  • gas_mmr (ndarray) – Array of gas mixing ratio from gas_properties

  • rho_p (float) – density of condensed vapor (g/cm^3)

  • mw_atmos (float) – Mean molecular weight of the atmosphere

  • gravity (float) – Gravity of planet cgs

  • kz (float or ndarray) – Kzz in cgs, either float or ndarray depending of whether or not it is set as input

  • fsed (float) – Sedimentation efficiency coefficient, unitless

  • b (float) – Denominator of exponential in sedimentation efficiency (if param is ‘exp’)

  • eps (float) – Minimum value of fsed function (if param=exp)

  • scale_h (float) – Scale height of the atmosphere

  • z_top (float) – Altitude at each layer

  • z_alpha (float) – Altitude at which fsed=alpha for variable fsed calculation

  • param (str) – fsed parameterisation ‘const’ (constant), ‘exp’ (exponential density derivation)

  • mh (float) – Atmospheric metallicity in NON log units (e.g. 1 for 1x solar)

  • sig (float) – Width of the log normal particle distribution

  • d_molecule (float) – diameter of atmospheric molecule (cm) (Rosner, 2000) (3.711e-8 for air, 3.798e-8 for N2, 2.827e-8 for H2) Set in Atmosphere constants

  • eps_k (float) – Depth of the Lennard-Jones potential well for the atmosphere Used in the viscocity calculation (units are K) (Rosner, 2000)

  • c_p_factor (float) – specific heat of atmosphere (erg/K/g) . Usually 7/2 for ideal gas diatomic molecules (e.g. H2, N2). Technically does slowly rise with increasing temperature

  • og_vfall (bool , optional) – optional, default = True. True does the original fall velocity calculation. False does the updated one which runs a tad slower but is more consistent. The main effect of turning on False is particle sizes in the upper atmosphere that are slightly bigger.

  • do_virtual (bool,optional) – optional, Default = True which adds a virtual layer if the species condenses below the model domain.

  • supsat (float, optional) – Default = 0 , Saturation factor (after condensation)

Returns:

  • qc (ndarray) – condenstate mixing ratio (g/g)

  • qt (ndarray) – gas + condensate mixing ratio (g/g)

  • rg (ndarray) – geometric mean radius of condensate cm

  • reff (ndarray) – droplet effective radius (second moment of size distrib, cm)

  • ndz (ndarray) – number column density of condensate (cm^-3)

  • qc_path (ndarray) – vertical path of condensate

virga.justdoit.get_mie(gas, directory)

Get Mie parameters from old ass formatted files

virga.justdoit.get_r_grid(r_min=1e-08, n_radii=60)

Warning

Original code from A&M code. Discontinued function. See ‘get_r_grid’.

Get spacing of radii to run Mie code

r_minfloat

Minimum radius to compute (cm)

n_radiiint

Number of radii to compute

virga.justdoit.get_r_grid_w_max(r_min=1e-08, r_max=0.054239131, n_radii=60)

Get spacing of radii to run Mie code

r_minfloat

Minimum radius to compute (cm)

r_maxfloat

Maximum radius to compute (cm)

n_radiiint

Number of radii to compute

virga.justdoit.get_refrind(igas, directory)

Reads reference files with wavelength, and refractory indecies. This function relies on input files being structured as a 4 column file with columns: index, wavelength (micron), nn, kk

Parameters:
  • igas (str) – Gas name

  • directory (str) – Directory were reference files are located.

Return type:

wavelength, real part, imaginary part

virga.justdoit.hot_jupiter()
virga.justdoit.layer(gas_name, rho_p, t_layer, p_layer, t_top, t_bot, p_top, p_bot, kz, mixl, gravity, mw_atmos, gas_mw, q_below, supsat, fsed, b, eps, z_top, z_bot, z_alpha, z_min, param, sig, mh, rmin, nrad, d_molecule, eps_k, c_p_factor, og_vfall, z_cld)

Calculate layer condensate properties by iterating on optical depth in one model layer (convering on optical depth over sublayers)

gas_namestr

Name of condenstante

rho_pfloat

density of condensed vapor (g/cm^3)

t_layerfloat

Temperature of layer mid-pt (K)

p_layerfloat

Pressure of layer mid-pt (dyne/cm^2)

t_topfloat

Temperature at top of layer (K)

t_botfloat

Temperature at botton of layer (K)

p_topfloat

Pressure at top of layer (dyne/cm2)

p_botfloat

Pressure at botton of layer

kzfloat

eddy diffusion coefficient (cm^2/s)

mixlfloat

Mixing length (cm)

gravityfloat

Gravity of planet cgs

mw_atmosfloat

Molecular weight of the atmosphere

gas_mwfloat

Gas molecular weight

q_belowfloat

total mixing ratio (vapor+condensate) below layer (g/g)

supsatfloat

Super saturation factor

fsedfloat

Sedimentation efficiency coefficient (unitless)

bfloat

Denominator of exponential in sedimentation efficiency (if param is ‘exp’)

eps: float

Minimum value of fsed function (if param=exp)

z_topfloat

Altitude at top of layer

z_botfloat

Altitude at bottom of layer

z_alphafloat

Altitude at which fsed=alpha for variable fsed calculation

paramstr

fsed parameterisation ‘const’ (constant), ‘exp’ (exponential density derivation)

sigfloat

Width of the log normal particle distribution

mhfloat

Metallicity NON log soar (1=1xSolar)

rminfloat

Minium radius on grid (cm)

nradint

Number of radii on Mie grid

d_moleculefloat

diameter of atmospheric molecule (cm) (Rosner, 2000) (3.711e-8 for air, 3.798e-8 for N2, 2.827e-8 for H2) Set in Atmosphere constants

eps_kfloat

Depth of the Lennard-Jones potential well for the atmosphere Used in the viscocity calculation (units are K) (Rosner, 2000)

c_p_factorfloat

specific heat of atmosphere (erg/K/g) . Usually 7/2 for ideal gas diatomic molecules (e.g. H2, N2). Technically does slowly rise with increasing temperature

og_vfallbool

Use original or new vfall calculation

Returns:

  • qc_layer (ndarray) – condenstate mixing ratio (g/g)

  • qt_layer (ndarray) – gas + condensate mixing ratio (g/g)

  • rg_layer (ndarray) – geometric mean radius of condensate cm

  • reff_layer (ndarray) – droplet effective radius (second moment of size distrib, cm)

  • ndz_layer (ndarray) – number column density of condensate (cm^-3)

  • q_below (ndarray) – total mixing ratio (vapor+condensate) below layer (g/g)

virga.justdoit.picaso_format(opd, w0, g0, pressure=None, wavenumber=None)

Gets virga output to picaso format

Parameters:
  • opd (ndarray) – array from virga of extinction optical depths per layer

  • w0 (ndarray) – array from virga of single scattering albedo

  • g0 (ndarray) – array from virga of asymmetry

  • pressure (array) – pressure array in bars

  • wavenumber (array) – wavenumber arry in cm^(-1)

virga.justdoit.picaso_format_custom_wavenumber_grid(opd, w0, g0, wavenumber_grid)

This is currently redundant with picaso_format now that picaso_format reads in wavenumber grid. Keeping for now, but will discontinue soon.

virga.justdoit.picaso_format_slab(p_bottom, opd, w0, g0, wavenumber_grid, pressure_grid, p_top=None, p_decay=None)

Sets up a PICASO-readable dataframe that inserts a wavelength dependent aerosol layer at the user’s given pressure bounds, i.e., a wavelength-dependent slab of clouds or haze.

Parameters:
  • p_bottom (float) – the cloud/haze base pressure the upper bound of pressure (i.e., lower altitude bound) to set the aerosol layer. (Bars)

  • opd (ndarray) – wavelength-dependent optical depth of the aerosol

  • w0 (ndarray) – wavelength-dependent single scattering albedo of the aerosol

  • g0 (ndarray) – asymmetry parameter = Q_scat wtd avg of <cos theta>

  • wavenumber_grid (ndarray) – wavenumber grid in (cm^-1)

  • pressure_grid (ndarray) – bars, user-defined pressure grid for the model atmosphere

  • p_top (float) – bars, the cloud/haze-top pressure This cuts off the upper cloud region as a step function. You must specify either p_top or p_decay.

  • p_decay (ndarray) – noramlized to 1, unitless array the same size as pressure_grid which specifies a height dependent optical depth. The usual format of p_decay is a fsed like exponential decay ~np.exp(-fsed*z/H)

Return type:

Dataframe of aerosol layer with pressure (in levels - non-physical units!), wavenumber, opd, w0, and g0 to be read by PICASO

virga.justdoit.recommend_gas(pressure, temperature, mh, mmw, plot=False, returnplot=False, legend='inside', **plot_kwargs)

Recommends condensate species for a users calculation.

Parameters:
  • pressure (ndarray, list) – Pressure grid for user’s pressure-temperature profile. Unit=bars

  • temperature (ndarray, list) – Temperature grid (should be on same grid as pressure input). Unit=Kelvin

  • mh (float) – Metallicity in NOT log units. Solar =1

  • mmw (float) – Mean molecular weight of the atmosphere. Solar = 2.2

  • plot (bool, optional) – Default is False. Plots condensation curves against PT profile to demonstrate why it has chose the gases it did.

  • plot_kwargs (kwargs) – Plotting kwargs for bokeh figure

Returns:

pressure (bars), condensation temperature (Kelvin)

Return type:

ndarray, ndarray

virga.justdoit.temperate_neptune()
virga.justdoit.warm_neptune()

virga.justplotit module

virga.justplotit.all_optics(out)

Maps of the wavelength dependent single scattering albedo and cloud opacity and asymmetry parameter as a function of altitude.

Parameters:

out (dict) – Dictionary output from pyeddy run

Return type:

Three bokeh plots with the single scattering, optical depth, and assymetry maps

virga.justplotit.all_optics_1d(out, wave_range, return_output=False, legend=None, colors=('#0072B2', '#E69F00', '#F0E442', '#009E73', '#56B4E9', '#D55E00', '#CC79A7', '#000000'), **kwargs)

Plots 1d profiles of optical depth per layer, single scattering, and asymmetry averaged over the user input wave_range.

Parameters:
  • out (list or dict) – Either a list of output dictionaries or a single dictionary output from .compute(as_dict=True)

  • wave_range (list) – min and max wavelength in microns

  • return_output (bool) – Default is just to return a figure but you can also return all the 1d profiles

  • legend (bool) – Default is none. Legend for each component of out

  • **kwargs (keyword arguments) – Key word arguments will be supplied to each bokeh figure function

virga.justplotit.condensate_mmr(out, gas=None, compare=False, legend=None, **kwargs)

Condensate mean mass mixing ratio

Parameters:
  • out (dict) – Dictionary output from pyeddy run

  • color (method) – Method from bokeh.palletes. e.g. bokeh.palletes.virids, bokeh.palletes.magma

  • **kwargs (kwargs) – Kwargs for bokeh.figure()

virga.justplotit.find_nearest_1d(array, value)
virga.justplotit.fsed_from_output(out, labels, y_axis='pressure', color_indx=0, cld_bounds=False, gas_indx=None, **kwargs)
virga.justplotit.opd_by_gas(out, gas=None, color=<function magma>, compare=False, legend=None, **kwargs)

Optical depth for conservative geometric scatteres separated by gas. E.g. [Fig 7 in Morley+2012](https://arxiv.org/pdf/1206.4313.pdf)

Parameters:
  • out (dict) – Dictionary output from pyeddy run

  • color (method) – Method from bokeh.palletes. e.g. bokeh.palletes.virids, bokeh.palletes.magma

  • **kwargs (kwargs) – Kwargs for bokeh.figure()

virga.justplotit.plot_format(df)

Function to reformat plots

virga.justplotit.plot_fsed(pressure, z, scale_height, alpha, beta, epsilon=0.01, pres_alpha=None, **kwargs)
virga.justplotit.pressure_fig(**kwargs)
virga.justplotit.pt(out, with_condensation=True, return_condensation=False, **kwargs)

Plot PT profiles with condensation curves. Thick lines in this pt plot signify those gases that were turned on in the run, NOT those that were recommended.

Parameters:
  • out (dict) – Dictionary output from pyeddy run

  • with_condensation (bool) – Plots condensation curves of gases. Also plots those that were turned on in the calculation with line_width=5. All others are set to 1.

  • return_condensation (bool) – If true, it returns list of condenation temps for each gas, pressure grid, and a list of the gas names

  • **kwargs (kwargs) – Kwargs for bokeh.figure()

Returns:

  • if return_condensation (fig, cond temperature curves, ,cond pressure grid , name of all available gases)

  • else (fig)

virga.justplotit.radii(out, gas=None, at_pressure=0.001, compare=False, legend=None, p1w=300, p1h=300, p2w=300, p2h=300, color_indx=0)

Plots the particle radii profile along with the distribution, at a certain pressure.

Parameters:
  • out (dict) – Dictionary output from pyeddy run

  • pressure (float,optional) – Pressure level to plot full distribution (bars)

virga.optics module

virga.optics.get_refrind(condensible, directory='~/Documents/eddysed/input/optics')

Read old style refrind files

Parameters:

condensible (str) – Condensible name (e.g. Al2O3)

Return type:

micron_wave, nn, kk as ndarrays

virga.optics.init_optics(condensibles, nrad=40, rmin=1e-10, read_mie=True)

Setup up a particle size grid and calculate single-particle scattering and absorption efficiencies and other parameters to be used by calc_optics()

Parameters:
  • do_optics (bool) – (True/False) Calculate optics (T) or use pre computed files (F)

  • read_mie (bool) – (True/False) Read in Mie coefficients from pre compute files gas_name.mieff

  • condensibles (list of str) – Name in str of all condensible gases e.g. [‘H2O’,’CH4’]

  • nrad (int) – Number of radius grid points

  • rmin (float) – Minimum number of radius grid (cm)

Returns:

  • wave (array) – Wavelength bin centers (cm)

  • radius (array) – Radius bin centers (cm)

  • dr (array) – Widths of radius bins (cm)

  • qscat (array) – Scattering efficiency

  • qext (array) – Extinction efficiency

  • cos_qscat (array) – qscat * acerage <cos (scattering angle)>

virga.pvaps module

virga.pvaps.Al2O3(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.CH4(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.CaAl12O19(temp, p, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • p (float) – Pressure (dyne/cm^2)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.CaTiO3(temp, p, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • p (float) – Pressure (dyne/cm^2)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.Cr(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.Fe(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.H2O(temp, do_buck=True, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

  • do_buck (bool) – True means use Buck 1981 expresssion, False means use Wexler’s

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.KCl(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.Mg2SiO4(temp, p, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • p (float, ndarray) – Pressure dyne/cm^2

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.MgSiO3(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.MnS(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.NH3(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.Na2S(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.SiO2(temp, mh=1)

PLACEHODER Not suited yet for public use Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.pvaps.TiO2(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

virga.pvaps.ZnS(temp, mh=1)

Computes vapor pressure curve

Parameters:
  • temp (float, ndarray) – Temperature (K)

  • mh (float) – NON log metallicity relative to solar (1=1Xsolar)

Return type:

vapor pressure in dyne/cm^2

Notes

virga.roo_functions module

virga.root_functions.advdiff(qt, ad_qbelow=None, ad_qvs=None, ad_mixl=None, ad_dz=None, ad_rainf=None, zb=None, b=None, eps=None, param='const')

Calculate divergence from advective-diffusive balance for condensate in a model layer

All units are cgs

  1. Ackerman Feb-2000

Parameters:
  • qt (float) – total mixing ratio of condensate + vapor (g/g)

  • ad_qbelow (float) – total mixing ratio of vapor in underlying layer (g/g)

  • ad_qvs (float) – saturation mixing ratio (g/g)

  • ad_mixl (float) – convective mixing length (cm)

  • ad_dz (float) – layer thickness (cm)

  • ad_rainf (float) – rain efficiency factor

  • zb (float) – altitude at bottom of layer

  • b (float) – denominator of fsed exponential (if param is ‘exp’)

  • param (str) – fsed parameterisation ‘const’ (constant), ‘exp’ (exponential density derivation)

Returns:

ad_qc – mixing ratio of condensed condensate (g/g)

Return type:

float

virga.root_functions.find_cond_t(t_test, p_test=None, mh=None, mmw=None, gas_name=None, gas_mmr=None)

Root function used used to find condenstation temperature. E.g. the temperature when

log p_vap = log partial pressure of gas

Parameters:
  • t_test (float) – Temp (K)

  • p_test (float) – Pressure bars

  • mh (float) – NON log mh .. aka MH=1 for solar

  • mmw (float) – mean molecular weight (2.2 for solar)

  • gas_name (str) – gas name, case sensitive

virga.root_functions.find_rg(rg, fsed, rw, alpha, s, loc=0.0, dist='lognormal')

Root function used used to find the geometric mean radius of lognormal size distribution.

Useful if we consider more complicated size distributions for which analytical expressions for moments are not easily obtained.

Parameters:
  • rg (float) – Geometric mean radius

  • fsed (float) – Sedimentation efficiency

  • rw (float) – Fall velocity particle radius

  • alpha (float) – Exponent in power-law approximation for particle fall-speed

  • s (float) – s = log(sigma) where sigma is the geometric std dev of lognormal distn

  • loc (float) – Shift in lognormal distribution

  • dist (str) – Continuous random variable for particle size distribution

virga.root_functions.force_balance(vf, r, grav, mw_atmos, mfp, visc, t, p, rhop, gas_kinetics=True)

” Define force balance for spherical particles falling in atmosphere, namely equate gravitational and viscous drag forces.

Viscous drag assumed to be quadratic (Benchaita et al. 1983) and is a function of the Reynolds number dependent drag coefficient. Drag coefficient taken from Khan-Richardson model (Richardson et al. 2002) and is valid for 1e-2 < Re < 1e5.

Parameters:
  • vf (float) – particle sedimentation velocity (cm/s)

  • r (float) – particle radius (cm)

  • grav (float) – acceleration of gravity (cm/s^2)

  • mw_atmos (float) – atmospheric molecular weight (g/mol)

  • mfp (float) – atmospheric molecular mean free path (cm)

  • visc (float) – atmospheric dynamic viscosity (dyne s/cm^2) see Eqn. B2 in A&M

  • t (float) – atmospheric temperature (K)

  • p (float) – atmospheric pressure (dyne/cm^2)

  • rhop (float) – density of particle (g/cm^3)

virga.root_functions.moment(n, s, loc, scale, dist='lognormal')

Calculate moment of size distribution. Will be extended to include more than just lognormal

Parameters:
  • n (float) – nth moment to be calculated

  • s (float) – std dev

  • loc (float) – Shift in distribution

  • scale (float) – Scale distribution

  • dist (str) – Continuous random variable for particle size distribution

Returns:

moment_out – nth moment of distribution

Return type:

float

virga.root_functions.qvs_below_model(p_test, qv_dtdlnp=None, qv_p=None, qv_t=None, qv_factor=None, qv_gas_name=None, mh=None, q_below=None)

Calculates the pressure of saturation mixing ratio for a gas that is extrapolated below the model domain. This is specifically used if the gas has saturated below the model grid

virga.root_functions.solve_force_balance(solve_for, temp, grav, mw_atmos, mfp, visc, t, p, rhop, lo, hi)
This can be used to equate gravitational and viscous drag forces to find either
  1. fallspeed for a spherical particle of given radius, or

  2. the radius of a particle falling with a particular speed

at one layer in an atmosphere.

Parameters:
  • solve_for (str) – property we are solving force balance for either “vfall” or “rw”

  • temp (float) –

    the other property needed to solve the force balance for solve_for ie if solve_for=”vfall”, temp will be a radius (cm)

    if solve_for=”rw”, temp will be w_convect (cm/s)

  • grav (float) – acceleration of gravity (cm/s^2)

  • mw_atmos (float) – atmospheric molecular weight (g/mol)

  • mfp (float) – atmospheric molecular mean free path (cm)

  • visc (float) – atmospheric dynamic viscosity (dyne s/cm^2) see Eqn. B2 in A&M

  • t (float) – atmospheric temperature (K)

  • p (float) – atmospheric pressure (dyne/cm^2)

  • rhop (float) – density of particle (g/cm^3)

  • lo (float) – lower bound for root-finder

  • hi (float) – upper bound for root-finder

Returns:

soln.root – the value of the parameter solve_for

Return type:

float

virga.root_functions.vfall(r, grav, mw_atmos, mfp, visc, t, p, rhop)

Calculate fallspeed for a spherical particle at one layer in an atmosphere, depending on Reynolds number for Stokes flow.

For Re_Stokes < 1, use Stokes velocity with slip correction For Re_Stokes > 1, use fit to Re = exp( b1*x + b2*x^2 )

where x = log( Cd Re^2 / 24 ) where b2 = -0.1 (curvature term) and b1 from fit between Stokes at Re=1, Cd=24 and Re=1e3, Cd=0.45

and Precipitation, Reidel, Holland, 1978) and Carlson, Rossow, and Orton (J. Atmos. Sci. 45, p. 2066, 1988)

all units are cgs

  1. Ackerman Feb-2000

Parameters:
  • r (float) – particle radius (cm)

  • grav (float) – acceleration of gravity (cm/s^2)

  • mw_atmos (float) – atmospheric molecular weight (g/mol)

  • mfp (float) – atmospheric molecular mean free path (cm)

  • visc (float) – atmospheric dynamic viscosity (dyne s/cm^2) see Eqn. B2 in A&M

  • t (float) – atmospheric temperature (K)

  • p (float) – atmospheric pressure (dyne/cm^2)

  • rhop (float) – density of particle (g/cm^3)

virga.root_functions.vfall_aggregrates(r, grav, mw_atmos, t, p, rhop, D=2, Ragg=1)

Calculate fallspeed for a particle at one layer in an atmosphere, assuming low Reynolds number and in the molecular regime (Epstein drag), i.e., where Knudsen number (mfp/r) is high (>>1).

User chooses whether particle is spherical or an aggregate. If aggregrate, the monomer size is given by r and the outer effective radius of the aggregrate is Ragg.

all units are cgs

Parameters:
  • r (float) – monomer particle radius (cm)

  • mw_atmos (float) – atmospheric molecular weight (g/mol)

  • t (float) – atmospheric temperature (K)

  • p (float) – atmospheric pressure (dyne/cm^2)

  • rhop (float) – density of particle (g/cm^3)

  • grav (float) – acceleration of gravity (cm/s^2)

  • Ragg (float) – aggregate particle effective radius (cm). (Defaults to 1 for spherical monomers).

  • D (float) – fractal number (Default is 2 because function reduces to monomers at this value).

virga.root_functions.vfall_aggregrates_ohno(r, grav, mw_atmos, mfp, t, p, rhop, ad_qc, D=2)

Calculates fallspeed for a fractal aggegrate particle as performed by Ohno et al., 2020, with an outer effective radius of R_agg made up of monomers of size r, at one layer in an atmosphere. ***Requires characteristic radius ragg to have D > or equal 2, i.e., not very fluffy aggregrates.

Essentially a more explicit version of the regular “vfall_aggregates” function, and is not dependent on being in free molecular (Esptein) regime

all units are cgs

Parameters:
  • r (float) – monomer particle radius (cm)

  • grav (float) – acceleration of gravity (cm/s^2)

  • rhop (float) – density of aggregrate particle (g/cm^3)

  • t (float) – atmospheric temperature (K)

  • p (float) – atmospheric pressure (dyne/cm^2)

  • mw_atmos (float) – atmospheric molecular weight (g/mol)

  • mfp (float) – atmospheric molecular mean free path (cm)

  • ad_qc (float) – mixing ratio of condensed condensate (g/g)

  • D (float) – fractal number (Default is 2).

  • gas_mw (float) – condensate gas mean molecular weight (g/mol)

virga.root_functions.vfall_find_root(r, grav=None, mw_atmos=None, mfp=None, visc=None, t=None, p=None, rhop=None, w_convect=None)

This is used to find F(X)-y=0 where F(X) is the fall speed for a spherical particle and y is the convective velocity scale. When the two are balanced we arrive at our correct particle radius.

Therefore, it is the same as the vfall function but with the subtraction of the convective velocity scale (cm/s) w_convect.

Module contents