.
Non-equilibrium cloud formation
CloudMicrophysics.MicrophysicsNonEq
— ModuleNon-equilibrium bulk microphysics scheme, which includes:
- condensation and evaporation of cloud liquid water and deposition and sublimation of cloud ice (relaxation to equilibrium)
CloudMicrophysics.MicrophysicsNonEq.τ_relax
— Functionτ_relax(liquid)
τ_relax(ice)
liquid
orice
- a type for cloud liquid water or ice
Returns the relaxation timescale for condensation and evaporation of cloud liquid water or the relaxation timescale for sublimation and deposition of cloud ice.
CloudMicrophysics.MicrophysicsNonEq.conv_q_vap_to_q_liq_ice
— Functionconv_q_vap_to_q_liq_ice(params, q_sat, q_)
params
- a struct with cloud water or ice free parametersq_sat
- liquid or ice equilibrium specific contentsq_
` - current liquid or ice specific contants
Returns the cloud water tendency due to condensation and evaporation or cloud ice tendency due to sublimation and vapor deposition. The tendency is obtained assuming a relaxation to equilibrium with a constant timescale and is based on the difference between specific contents in equilibrium at the current temperature and the current cloud condensate.
CloudMicrophysics.MicrophysicsNonEq.conv_q_vap_to_q_liq_ice_MM2015
— Functionconv_q_vap_to_q_liq_ice_MM2015(params, tps, q_tot, q_liq, q_ice, q_rai, q_sno, ρ, T)
conv_q_vap_to_q_liq_ice_MM2015(params, tps, q::PhasePartition, q_rai, q_sno, ρ, T)
params
- a struct with cloud water or ice free parameterstps
- thermodynamics parameters structq_tot
,q_liq
,q_ice
,q_rai
,q_sno
- specific contents of total water, cloud liquid and ice, rain and snow,q
- current PhasePartitionρ
- air density [kg/m3]T
- air temperature [K]
Returns the cloud water tendency due to condensation and evaporation or cloud ice tendency due to sublimation and vapor deposition. The formulation is based on Morrison and Grabowski 2008 and Morrison and Milbrandt 2015. It does NOT screen for small or negative values for humidities, so we suggest applying a limiter of choice to this function when running it in a model.
CloudMicrophysics.MicrophysicsNonEq.terminal_velocity
— Functionterminal_velocity(sediment, vel, ρ, q)
sediment
- a struct with sedimentation type (cloud liquid or ice)vel
- a struct with terminal velocity parametersρₐ
- air densityq
- cloud liquid or ice specific content
Returns the mass weighted average terminal velocity assuming a monodisperse size distribution with prescribed number concentration. The fall velocity of individual particles is parameterized following Chen et. al 2022, DOI: 10.1016/j.atmosres.2022.106171
0-moment precipitation microphysics
CloudMicrophysics.Microphysics0M
— ModuleMicrophysics0M
Zero-moment bulk microphysics scheme that instantly removes moisture above certain threshold. This is equivalent to instanteneous conversion of cloud condensate into precipitation and precipitation fallout with infinite terminal velocity.
CloudMicrophysics.Microphysics0M.remove_precipitation
— Functionremove_precipitation(params_0M, q_liq, q_ice; q_vap_sat)
remove_precipitation(params_0M, q; q_vap_sat)
params_0M
- a struct with 0-moment parametersq
- current Thermodynamics.PhasePartition orq_liq
andq_ice
specific contentsq_vap_sat
- specific humidity at saturation
Returns the q_tot
tendency due to the removal of precipitation. The tendency is obtained assuming a relaxation with a constant timescale to a state with precipitable water removed. The threshold for when to remove q_tot
is defined either by the condensate specific content or supersaturation. The thresholds and the relaxation timescale are defined Parameters0M
struct.
1-moment precipitation microphysics
CloudMicrophysics.Microphysics1M
— ModuleOne-moment bulk microphysics scheme, which includes:
- terminal velocity of precipitation
- autoconversion of cloud liquid water into rain and of cloud ice into snow
- accretion due to collisions between categories of condensed species
- evaporation and sublimation of hydrometeors
- melting of snow into rain
CloudMicrophysics.Microphysics1M.get_v0
— Functionget_v0(v, ρ)
v
- a struct with bulk 1-moment terminal velocity parametersρ
- air density (only forRain
)
Returns the proportionality coefficient in terminal velocity(r/r0).
CloudMicrophysics.Microphysics1M.get_n0
— Functionget_n0(pdf, q_sno, ρ)
pdf
- a struct with parameters for snow, ice, and rain size distributionq_sno
- snow specific content (only forSnow
)ρ
- air density (only forSnow
)
Returns the intercept parameter of the assumed Marshall-Palmer distribution
CloudMicrophysics.Microphysics1M.lambda_inverse
— Functionlambda_inverse(pdf, mass, q, ρ)
pdf
,mass
- structs with particle size distribution and mass parametersq
- specific content of rain, ice or snowρ
- air density
Returns the inverse of rate parameter of the assumed size distribution of particles (rain drops, ice crystals, snow crystals). The value is clipped at 1e-8.
CloudMicrophysics.Microphysics1M.terminal_velocity
— Functionterminal_velocity(precip, vel, ρ, q)
precip
- a struct with precipitation type (rain or snow)vel
- a struct with terminal velocity parametersρ
- air densityq
- rain or snow specific content
Returns the mass weighted average terminal velocity assuming a Marshall-Palmer (1948) distribution of particles. Fall velocity of individual rain drops is parameterized:
- assuming an empirical power-law relations for
velocity == Blk1MVelType
- following Chen et. al 2022, DOI: 10.1016/j.atmosres.2022.106171, for
velocity == Chen2022VelType
CloudMicrophysics.Microphysics1M.conv_q_liq_to_q_rai
— Functionconv_q_liq_to_q_rai(acnv, q_liq, smooth_transition)
acnv
- 1M autoconversion parametersq_liq
- liquid water specific contentsmooth_transition
- a flag to switch on smoothing
Returns the q_rai tendency due to collisions between cloud droplets (autoconversion), parametrized following Kessler (1995).
CloudMicrophysics.Microphysics1M.conv_q_ice_to_q_sno_no_supersat
— Functionconv_q_ice_to_q_sno_no_supersat(acnv, q_ice, smooth_transition)
acnv
- 1M autoconversion parametersq_ice
- cloud ice specific contentsmooth_transition
- a flag to switch on smoothing
Returns the q_sno tendency due to autoconversion from ice. This is a simplified version of a snow autoconversion rate that can be used in simulations where there is no supersaturation (for example in TC.jl when using saturation adjustment).
CloudMicrophysics.Microphysics1M.conv_q_ice_to_q_sno
— Functionconv_q_ice_to_q_sno(ice, aps, tps, q_tot, q_liq, q_ice, q_rai, q_sno, ρ, T)
conv_q_ice_to_q_sno(ice, aps, tps, q::PhasePartition, q_rai, q_sno, ρ, T)
ice
- a struct with ice parametersaps
- a struct with air propertiestps
- a struct with thermodynamics parametersq_tot
- total water specific contentq_liq
- cloud water specific contentq_ice
- ice specific contentq_rai
- rain specific contentq_sno
- snow specific contentq
- Thermodynamics.PhasePartitionρ
- air densityT
- air temperature
Returns the q_sno tendency due to autoconversion from ice. Parameterized following Harrington et al. (1996) and Kaul et al. (2015).
CloudMicrophysics.Microphysics1M.accretion
— Functionaccretion(cloud, precip, vel, ce, q_clo, q_pre, ρ)
cloud
- type for cloud water or cloud iceprecip
- type for rain or snowvel
- a struct with terminal velocity parametersce
- collision efficiency parametersq_clo
- cloud water or cloud ice specific contentq_pre
- rain water or snow specific contentρ
- air density
Returns the source of precipitating water (rain or snow) due to collisions with cloud water (liquid or ice).
CloudMicrophysics.Microphysics1M.accretion_rain_sink
— Functionaccretion_rain_sink(rain, ice, vel, ce, q_ice, q_rai, ρ)
rain
- rain type parametersice
- ice type parametersvel
- terminal velocity parameters for raince
- collision efficiency parametersq_ice
- cloud ice specific contentq_rai
- rain water specific contentρ
- air density
Returns the sink of rain water (partial source of snow) due to collisions with cloud ice.
CloudMicrophysics.Microphysics1M.accretion_snow_rain
— Functionaccretion_snow_rain(ce, type_i, type_j, blk1m_type_i, blk1m_type_j, q_i, q_j, ρ)
ce
- collision efficiency parametersi
- snow for temperatures below freezing or rain for temperatures above freezingj
- rain for temperatures below freezing or snow for temperatures above freezingtype_i
,type_j
- a type for snow or rainblk1mveltype_ti
,blk1mveltype_tj
- 1M terminal velocity parametersq_
- specific content of snow or rainρ
- air density
Returns the accretion rate between rain and snow. Collisions between rain and snow result in snow at temperatures below freezing and in rain at temperatures above freezing.
CloudMicrophysics.Microphysics1M.evaporation_sublimation
— Functionevaporation_sublimation(params, vel, aps, tps, q_tot, q_liq, q_ice, q_rai, q_snow, ρ, T)
evaporation_sublimation(params, vel, aps, tps, q::PhasePartition, q_rai, q_sno, ρ, T)
params
- a struct with rain or snow parametersvel
- a struct with terminal velocity parametersaps
- a struct with air parameterstps
- a struct with thermodynamics parametersq_tot
- total water specific contentq_liq
- liquid water specific contentq_ice
- ice specific contentq_rai
- rain specific contentq_sno
- snow specific contentq
- Thermodynamics.PhasePartitionρ
- air densityT
- air temperature
Returns the tendency due to rain evaporation or snow sublimation.
CloudMicrophysics.Microphysics1M.snow_melt
— Functionsnow_melt(snow, vel, aps, tps, q_sno, ρ, T)
snow
- snow parametersvel
- terminal velocity parametersaps
- air propertiestps
- thermodynamics parametersq_sno
- snow water specific contentρ
- air densityT
- air temperature
Returns the tendency due to snow melt.
2-moment precipitation microphysics
CloudMicrophysics.Microphysics2M
— ModuleDouble-moment bulk microphysics parametrizations including:
- autoconversion, accretion, self-collection, breakup, mean terminal velocity of raindrops, and rain evaporation rates from Seifert and Beheng 2006.
- number concentration adjustment from Horn 2012.
- additional double-moment bulk microphysics autoconversion and accretion rates from: Khairoutdinov and Kogan 2000, Beheng 1994, Tripoli and Cotton 1980, and Liu and Daum 2004.
Size distributions
Parameters
CloudMicrophysics.Microphysics2M.pdf_cloud_parameters
— Functionpdf_cloud_parameters(pdf_c, q, ρₐ, N)
Return the parameters of the size distribution of cloud particles in terms of diameter.
The size distribution is given by:
n(D) = N₀c * D^νcD * exp(-λc * D^μcD)
where
νcD = 3νc + 2
μcD = 3μc
Arguments
pdf_c
: Size distribution parameters for cloud droplets,CMP.CloudParticlePDF_SB2006
q
: Liquid mass content [kg/kg]ρₐ
: Air density [kg/m³]N
: Number concentration of the particle [1/m³]
Returns
(logN₀c, λc, νcD, μcD)
: Parameters of the generalized gamma distribution in terms of diameter
CloudMicrophysics.Microphysics2M.pdf_rain_parameters
— Functionpdf_rain_parameters(pdf_r, qᵣ, ρₐ, Nᵣ)
Return the parameters of the rain drop diameter distribution
n_r(D) = N_0 * exp(- D / Dr_mean)
where
D
is the diameter of the raindrop,N_0
[1/m³] is the number concentration of raindrops,Dr_mean
[m] is the mean diameter of the raindrops.
Note: in SB2006, Eq. (83) the distribution is given as:
f(D) = N_0 * exp(- λ_r D)
where λ_r ≡ 1 / Dr_mean
[1/m] is the inverse of the mean diameter of the raindrops.
Arguments
pdf_r
: struct containing size distribution parameters for rain. Can either beCMP.RainParticlePDF_SB2006_notlimited
orCMP.RainParticlePDF_SB2006_limited
. For the latter, the values forN_0
,Dr_mean
, andxr_mean
are limited to be within provided ranges.qᵣ
: mass of rain water [kg]ρₐ
: air density [kg/m³]Nᵣ
: number of rain drops [1/m³]
Returns
- A
NamedTuple
with the fields(; N₀r, Dr_mean, xr_mean)
CloudMicrophysics.Microphysics2M.pdf_cloud_parameters_mass
— Functionpdf_cloud_parameters_mass(pdf_c, q, ρₐ, N)
Return the parameters of the size distribution of cloud particles in terms of mass.
See log_pdf_cloud_parameters_mass
for more details.
CloudMicrophysics.Microphysics2M.pdf_rain_parameters_mass
— Functionpdf_rain_parameters_mass(pdf_r::CMP.RainParticlePDF_SB2006, qᵣ, ρₐ, Nᵣ)
Return the parameters of the rain drop diameter distribution in terms of mass.
As a function of diameter, the size distribution is given by:
n(D) = N₀r * exp(-D / Dr_mean)
In terms of mass (x
), the size distribution is given by:
f(x) = n(D(x)) * ∂D∂x(x)
= N₀ * exp(-D(x) / Dr_mean) * 2 / (π * ρw) * x^(-2/3)
= N₀ * 2 / (π * ρw) * x^(-2/3) * exp(- (6 / (π * ρw))^(1/3) / Dr_mean * x^(1/3))
where
D(x) = (6x / (π * ρw))^(1/3)
is the diameter of a raindrop of massx
.∂D∂x(x) = (6 / (π * ρw))^(1/3) * x^(-2/3)
is the derivative of the diameter with respect to the mass.
If we write the general form of the size distribution as:
f(x) = A * x^ν * exp(-B * x^μ)
then we have that:
A = N₀ * 2 / (π * ρw)
B = (6 / (π * ρw))^(1/3) / Dr_mean
ν = -2/3
μ = 1/3
CloudMicrophysics.Microphysics2M.log_pdf_cloud_parameters_mass
— Functionlog_pdf_cloud_parameters_mass(pdf_c, q, ρₐ, N)
Return the log of the parameters of the generalized gamma distribution of the form
f(x) = A * x^ν * exp(-B * x^μ), [Eq. (79) in Seifert and Beheng 2006, but using the symbol `B` instead of `λ`]
where
B = [ x̄ Γ(z₁) / Γ(z₂) ]^(-μ)
A = μ N B^(z₁) / Γ(z₁)
z₁ = (ν + 1) / μ
z₂ = (ν + 2) / μ
That is,
log(B) = - μ [ log(x̄) + logΓ(z₁) - logΓ(z₂) ]
log(A) = log(μ) + log(N) + z₁ * log(B) - logΓ(z₁)
Arguments
pdf_c
: Size distribution parameters for cloud droplets,CMP.CloudParticlePDF_SB2006
q
: Liquid mass content [kg/kg]ρₐ
: Air density [kg/m³]N
: Number concentration of the particle [1/m³]
Returns
(logA, logB)
: Log of the parameters of the generalized gamma distribution
Size distributions
CloudMicrophysics.DistributionTools.size_distribution
— Methodsize_distribution(pdf::CMP.RainParticlePDF_SB2006, q, ρₐ, N)
Return n(D)
, a function that computes the size distribution for rain particles at diameter D
Arguments
pdf
: Rain size distribution parameters,CMP.RainParticlePDF_SB2006
q
: Rain water specific content [kg/kg]ρₐ
: Density of air [kg/m³]N
: Rain water number concentration [1/m³]
CloudMicrophysics.DistributionTools.size_distribution
— Methodsize_distribution(pdf::CMP.CloudParticlePDF_SB2006, q, ρₐ, N)
Return n(D)
, a function that computes the size distribution for cloud particles at diameter D
The size distribution is given by:
n(D) = N₀c * D^(3νc + 2) * exp(-λc * D^(3μc))
Arguments
pdf
: Cloud size distribution parameters,CMP.CloudParticlePDF_SB2006
q
: Cloud water specific content [kg/kg]ρₐ
: Density of air [kg/m³]N
: Cloud water number concentration [1/m³]
CloudMicrophysics.Microphysics2M.size_distribution_value
— Functionsize_distribution_value(pdf, q, ρₐ, N, D)
Return the size distribution value for a cloud or rain particle of diameter D
.
See size_distribution
for more details.
CloudMicrophysics.Microphysics2M.get_size_distribution_bounds
— Functionget_size_distribution_bounds(pdf, q, ρₐ, N, p)
Return the minimum and maximum diameters of a cloud or rain particle such that the size distribution is within (1 - p) to (p) probability of the true size distribution.
Arguments
pdf
: Size distribution parameters for cloud or rain,CMP.RainParticlePDF_SB2006
orCMP.CloudParticlePDF_SB2006
q
: mass mixing ratio of cloud or rain waterρₐ
: density of airN
: number mixing ratio of cloud or rainp
: probability level (0 ≤ p ≤ 1)
Returns
D_min, D_max
: minimum and maximum diameters of a cloud or rain particle such that the size distribution is within (1 - p) to (p) probability of the true size distribution. All inputs and output diameters are in base SI units. The bounds are calculated through quantile functions of the size distribution.
Rates
CloudMicrophysics.Microphysics2M.LiqRaiRates
— TypeA structure containing the rates of change of the specific contents and number densities of liquid and rain water.
CloudMicrophysics.Microphysics2M.autoconversion
— Functionautoconversion(acnv, pdf_c, q_liq, q_rai, ρ, N_liq)
Compute autoconversion rates
Arguments
acnv
: Autoconversion parameters,CMP.AcnvSB2006
pdf_c
: Cloud size distribution parameters,CMP.CloudParticlePDF_SB2006
q_liq
: Cloud water specific content [kg/kg]q_rai
: Rain water specific content [kg/kg]ρ
: Air density [kg/m³]N_liq
: Cloud droplet number density [1/m³]
Returns
LiqRaiRates
withq_liq
,N_liq
,q_rai
,N_rai
tendencies due to collisions between cloud droplets (autoconversion)
CloudMicrophysics.Microphysics2M.accretion
— Functionaccretion(accr, q_liq, q_rai, ρ, N_liq)
Compute accretion rate
Arguments
accr
: Accretion parameters,CMP.AccrSB2006
q_liq
: Cloud water specific content [kg/kg]q_rai
: Rain water specific content [kg/kg]ρ
: Air density [kg/m³]N_liq
: Cloud droplet number density [1/m³]
Returns
LiqRaiRates
withq_liq
,N_liq
,q_rai
,N_rai
tendencies due to collisions between raindrops and cloud droplets (accretion)
accretion(accretion_scheme, q_liq, q_rai, ρ)
accretion_scheme
- type for 2-moment rain accretion parameterizationq_liq
- cloud water specific contentq_rai
- rain water specific contentρ
- air density (forKK2000Type
andBeheng1994Type
)
Returns the accretion rate of rain, parametrized following
- Khairoutdinov and Kogan (2000) for
scheme == KK2000Type
- Beheng (1994) for
scheme == B1994Type
- Tripoli and Cotton (1980) for
scheme == TC1980Type
CloudMicrophysics.Microphysics2M.liquid_self_collection
— Functionliquid_self_collection(acnv, pdf_c, q_liq, ρ, dN_liq_dt_au)
Compute liquid self-collection rate
Arguments
acnv
: 2-moment autoconversion parameterization,CMP.AcnvSB2006
pdf_c
: Cloud size distribution parameters,CMP.CloudParticlePDF_SB2006
q_liq
: Cloud water specific content [kg/kg]ρ
: Air density [kg/m³]dN_liq_dt_au
: Rate of change of cloud droplets number density due to autoconversion [1/m³/s]
Returns
- The cloud droplets number density tendency due to collisions of cloud droplets that produce larger cloud droplets (self-collection)
CloudMicrophysics.Microphysics2M.autoconversion_and_liquid_self_collection
— Functionautoconversion_and_liquid_self_collection(scheme, q_liq, q_rai, ρ, N_liq)
Compute autoconversion and liquid self-collection rates
Arguments
scheme
: 2-moment rain autoconversion parameterization,CMP.SB2006
q_liq
: Cloud water specific content [kg/kg]q_rai
: Rain water specific content [kg/kg]ρ
: Air density [kg/m³]N_liq
: Cloud droplet number density [1/m³]
Returns
(au, sc)
: ANamedTuple
containing the autoconversion rate and the liquid self-collection rate.
CloudMicrophysics.Microphysics2M.rain_self_collection
— Functionrain_self_collection(pdf, self, q_rai, ρ, N_rai)
Compute the rain self-collection rate
Arguments
pdf
: Rain size distribution parameters,CMP.RainParticlePDF_SB2006
self
: Rain self-collection parameters,CMP.SelfColSB2006
q_rai
: Rain water specific content [kg/kg]ρ
: Air density [kg/m³]N_rai
: Raindrops number density [1/m³]
Returns
- The raindrops number density tendency due to collisions of raindrops that produce larger raindrops (self-collection).
CloudMicrophysics.Microphysics2M.rain_breakup
— Functionrain_breakup(pdf, brek, q_rai, ρ, N_rai, dN_rai_dt_sc)
Compute the raindrops number density tendency due to breakup of raindrops
Arguments
pdf
: Rain size distribution parameters,CMP.RainParticlePDF_SB2006
brek
: Rain breakup parameters,CMP.BreakupSB2006
q_rai
: Rain water specific contentρ
: Air densityN_rai
: Raindrops number densitydN_rai_dt_sc
: Rate of change of raindrops number density due to self-collection
Returns
- The raindrops number density tendency due to breakup of raindrops that produce smaller raindrops
CloudMicrophysics.Microphysics2M.rain_self_collection_and_breakup
— Functionrain_self_collection_and_breakup(params, q_rai, ρ, N_rai)
Compute the raindrops self-collection and breakup rates.
Arguments
params
: 2-moment rain size distribution parameters,CMP.SB2006
including raindrop size distribution, self collection, and breakup parametersq_rai
: Rain water specific contentρ
: Air densityN_rai
: Raindrops number density
Returns
(sc, br)
: ANamedTuple
containing the raindrops self-collection and breakup rates, respectively.
CloudMicrophysics.Microphysics2M.rain_terminal_velocity
— Functionrain_terminal_velocity(SB2006, vel, q_rai, ρ, N_rai)
Compute the raindrops terminal velocity.
Arguments
pdf_r
: Rain size distribution parameters,CMP.RainParticlePDF_SB2006
vel
: Terminal velocity parameters,CMP.Chen2022VelTypeRain
q_rai
: Rain water specific contentρ
: Air densityN_rai
: Raindrops number density
Returns
A tuple containing the number and mass weigthed mean fall velocities of raindrops in [m/s]. Assuming an exponential size distribution from Seifert and Beheng 2006 for scheme == SB2006Type
Fall velocity of individual rain drops is parameterized:
- assuming an empirical relation similar to Rogers (1993) for
velo_scheme == SB2006VelType
- following Chen et. al 2022, DOI: 10.1016/j.atmosres.2022.106171 for
velo_scheme == Chen2022Type
CloudMicrophysics.Microphysics2M.rain_evaporation
— Functionrain_evaporation(evap, aps, tps, q, q_rai, ρ, N_rai, T)
evap
- evaporation parameterization schemeaps
- air propertiestps
- thermodynamics parametersq
- phase partitionq_rai
- rain specific contentρ
- air densityN_rai
- raindrops number densityT
- air temperature
Returns a named tuple containing the tendency of raindrops number density and rain water specific content due to rain rain_evaporation, assuming a power law velocity relation for fall velocity of individual drops and an exponential size distribution, for scheme == SB2006Type
CloudMicrophysics.Microphysics2M.conv_q_liq_to_q_rai
— Functionconv_q_liq_to_q_rai(acnv, q_liq, ρ, N_d; smooth_transition)
acnv
- 2-moment rain autoconversion parameterizationq_liq
- cloud water specific contentρ
- air densityN_d
- prescribed cloud droplet number concentration
Returns the q_rai tendency due to collisions between cloud droplets (autoconversion), parametrized following:
- Khairoutdinov and Kogan (2000) for
scheme == KK2000Type
- Beheng (1994) for
scheme == B1994Type
- Tripoli and Cotton (1980) for
scheme == TC1980Type
- Liu and Daum (2004) for
scheme ==LD2004Type
The Beheng1994Type
, TC1980Type
and LD2004Type
of schemes additionally accept smooth_transition
flag that smoothes their thershold behaviour if set to true
. The default value is false
.
CloudMicrophysics.Microphysics2M.number_increase_for_mass_limit
— Functionnumber_increase_for_mass_limit(numadj, x_max, q, ρ, N)
Compute the tendency (rate of change) of number concentration N
required to ensure that the mean particle mass x = ρq / N
does not exceed the upper limit x_max
. Returns a positive tendency when the mean mass is too high (x > x_max
), and zero otherwise. The method is based on Horn (2012, DOI: 10.5194/gmd-5-345-2012).
Arguments
numadj
: Number concentration adjustment parameters (CMP.NumberAdjustmentHorn2012)q
: Mass mixing ratio [kg/kg]ρ
: Air density [kg/m³]N
: Number concentration [1/m³]x_max
: Maximum allowed mean particle mass [kg]
Returns
- The rate of change of number concentration [1/(m³·s)] needed to bring the mean mass within the upper bound.
CloudMicrophysics.Microphysics2M.number_decrease_for_mass_limit
— Functionnumber_decrease_for_mass_limit(numadj, x_min, q, ρ, N)
Compute the tendency (rate of change) of number concentration N
required to ensure that the mean particle mass x = ρq / N
does not fall below the lower limit x_min
. Returns a negative tendency when the mean mass is too low (x < x_min
), and zero otherwise. The method is based on Horn (2012, DOI: 10.5194/gmd-5-345-2012).
Arguments
numadj
: Number concentration adjustment parameters (CMP.NumberAdjustmentHorn2012)q
: Mass mixing ratio [kg/kg]ρ
: Air density [kg/m³]N
: Number concentration [1/m³]x_min
: Minimum allowed mean particle mass [kg]
Returns
- The rate of change of number concentration [1/(m³·s)] needed to bring the mean mass within the lower bound.
Distribution tools for 2-moment microphysics
CloudMicrophysics.DistributionTools
— ModuleDistributionTools
A module containing tools for working with size distributions.
Currently, it contains tools for working with the generalized gamma distribution and the exponential distribution.
Mostly used for the 2-moment microphysics.
CloudMicrophysics.DistributionTools.generalized_gamma_quantile
— Functiongeneralized_gamma_quantile(ν, μ, B, Y)
Calculate the quantile (inverse cumulative distribution function) for a generalized gamma distribution parameterized in the form:
g(x) = A ⋅ x^ν ⋅ exp(-B ⋅ x^μ)
Arguments
ν, μ, B
: The PDF parametersY
: The probability level (0 ≤ Y ≤ 1) for which to compute the quantile
Returns
x
: The value x such that P(X ≤ x) = Y
CloudMicrophysics.DistributionTools.generalized_gamma_cdf
— Functiongeneralized_gamma_cdf(ν, μ, B, x)
Calculate the cumulative distribution function (CDF) for a generalized gamma distribution parameterized in the form:
g(x) = A * x^ν * exp(-B * x^μ)
The CDF gives the probability P(X ≤ x).
Arguments
ν, μ, B
: The PDF parametersx
: The value at which to evaluate the CDF
Returns
p
: The probability P(X ≤ x)
CloudMicrophysics.DistributionTools.exponential_cdf
— Functionexponentialcdf(Dmean, D)
Calculate the cumulative distribution function (CDF) for an exponential distribution parameterized in the form:
n(D) = N₀ * exp(-D / D_mean)
where N₀ is a normalizing constant such that the total probability is 1.
Arguments
D_mean
: The mean value of the distributionD
: The point at which to evaluate the CDF (must be ≥ 0)
Returns
p
: The probability P(X ≤ D)
CloudMicrophysics.DistributionTools.exponential_quantile
— Functionexponential_quantile(D_mean, Y)
Calculate the quantile (inverse cumulative distribution function) for an exponential distribution parameterized in the form:
n(D) = N₀ * exp(-D / D_mean)
where N₀ is a normalizing constant such that the total probability is 1.
Arguments
Y
: The probability level (0 ≤ Y ≤ 1) for which to compute the quantileD_mean
: The mean value of the distribution
Returns
D
: The value D such that P(X ≤ D) = Y
P3 scheme
CloudMicrophysics.P3Scheme
— ModulePredicted particle properties scheme P3 for ice [22].
See online docs for more information.
Construct parameterization set
CloudMicrophysics.Parameters.ParametersP3
— TypeParametersP3{FT}
Parameters for P3 bulk microphysics scheme.
From Morrison and Milbrandt (2015) [22]
Fields
mass
: Mass-size relation, e.g.MassPowerLaw
area
: Area-size relation, e.g.AreaPowerLaw
slope
: Slope relation, e.g.SlopePowerLaw
orSlopeConstant
vent
: Ventilation relation, e.g.VentilationFactor
ρ_i
: Cloud ice density [kg m⁻³
]ρ_l
: Cloud liquid water density [kg m⁻³
]T_freeze
: Water freeze temperature [K
]
CloudMicrophysics.Parameters.ParametersP3
— MethodParametersP3(FT)
Create a ParametersP3
object from a floating point type FT
.
Examples
julia> import CloudMicrophysics.Parameters as CMP
julia> CMP.ParametersP3(Float64)
ParametersP3{Float64}
├── mass: MassPowerLaw
│ ├── α_va = 0.018537721864540644 [kg μm^(-β_va)]
│ └── β_va = 1.9 [-]
├── area: AreaPowerLaw
│ ├── γ = 0.2285 [μm^(2-σ)]
│ └── σ = 1.88 [-]
├── slope: SlopePowerLaw
│ ├── a = 0.00191 [m^b]
│ ├── b = 0.8 [-]
│ ├── c = 2.0 [-]
│ └── μ_max = 6.0 [-]
├── vent: VentilationFactor
│ ├── vent_a = 0.78 [-]
│ └── vent_b = 0.308 [-]
├── ρ_i = 916.7 [kg m⁻³]
├── ρ_l = 1000.0 [kg m⁻³]
└── T_freeze = 273.15 [K]
CloudMicrophysics.Parameters.ParametersP3
— MethodParametersP3(toml_dict::CP.AbstractTOMLDict; [slope_law = :powerlaw])
Create a ParametersP3
object from a ClimaParams
TOML dictionary.
Arguments
toml_dict::CP.AbstractTOMLDict
: AClimaParams
TOML dictionaryslope_law
: Slope law to use (:constant
or, by default,:powerlaw
)
Sub-parameterizations
CloudMicrophysics.Parameters.MassPowerLaw
— TypeMassPowerLaw{FT}
Parameters for mass(size) relation.
From measurements of mass grown by vapor diffusion and aggregation in midlatitude cirrus by Brown and Francis (1995) [23]
A part of the ParametersP3
parameter set.
The BF1995_mass_coeff_alpha
parameter is provided in units of [g μm^(-β_va)
] but the α_va
field is stored in SI-like units of [kg m^(-β_va)
] for consistency with the rest of the code.
Fields
α_va
: Coefficient in mass(size) relation [kg m^(-β_va)
]β_va
: Coefficient in mass(size) relation [-
]
CloudMicrophysics.Parameters.AreaPowerLaw
— TypeAreaPowerLaw{FT}
Parameters for area(size) relation.
\[A(D) = γ D^σ\]
where γ
and σ
are coefficients in area(size) for ice side plane, column, bullet, and planar polycrystal aggregates. Values are from Mitchell (1996) [24]
A part of the ParametersP3
parameter set.
Fields
γ
: Scale [μm^(2-σ)
]σ
: Power [-
]
CloudMicrophysics.Parameters.SlopeLaw
— TypeSlopeLaw{FT}
The top-level super-type for slope parameterizations.
See SlopePowerLaw
and SlopeConstant
for concrete implementations.
CloudMicrophysics.Parameters.SlopePowerLaw
— TypeSlopePowerLaw{FT}
Slope parameter μ as a power law in shape parameter λ:
\[μ(λ) = a λ^b - c\]
and is limited to:
\[0 ≤ μ ≤ μ_{max}\]
See also Eq. 3 in Morrison and Milbrandt (2015) [22]
A part of the ParametersP3
parameter set.
Fields
a
: Scale [m^b
]b
: Power [-
]c
: Offset [-
]μ_max
: Upper limiter [-
]
CloudMicrophysics.Parameters.SlopeConstant
— TypeSlopeConstant{FT}
Slope parameter μ as a constant:
\[μ(λ) = μ_{const}\]
A part of the ParametersP3
parameter set.
Fields
μ
: Slope parameter μ [-
]
CloudMicrophysics.Parameters.VentilationFactor
— TypeVentilationFactor{FT}
Parameters for ventilation factor:
\[F(D) = a_{vent} + b_{vent} N_{Sc}^{1/3} N_{Re}(D)^{1/2}\]
where N_{Sc}
is the Schmidt number and N_{Re}(D)
is the Reynolds number for a particle with diameter D
.
From Seifert and Beheng (2006) [10], see also Eq. (13-61) in Pruppacher and Klett (2010) [62]
A part of the ParametersP3
parameter set.
Fields
vent_a
: Ventilation factor a [-
]vent_b
: Ventilation factor b [-
]
Obtain particle state
CloudMicrophysics.P3Scheme.P3State
— TypeP3State{FT}
State of the P3 scheme.
This struct bundles the P3 parameterizations params
, the provided rime state (F_rim
, ρ_rim
), and the derived threshold variables (D_th
, D_gr
, D_cr
, ρ_g
).
To obtain a P3State
object, use the get_state
function.
Fields
params
:CMP.ParametersP3
objectL_ice
: The ice mass concentration [kg/m³]N_ice
: The ice number concentration [1/m³]F_rim
: Rime mass fractionρ_rim
: Rime density
CloudMicrophysics.P3Scheme.get_state
— Functionget_state(params; L_ice, N_ice, F_rim, ρ_rim)
Create a P3State
from CMP.ParametersP3
and rime state parameters.
Arguments
params
:CMP.ParametersP3
object
Keyword arguments
L_ice
: ice mass concentration [kg/m³]N_ice
: ice number concentration [1/m³]F_rim
: rime mass fraction [-],F_rim = L_rim / L_ice
ρ_rim
: rime density [kg/m³],ρ_rim = L_rim / B_rim
Examples
```jldoctest julia> import CloudMicrophysics.Parameters as CMP, CloudMicrophysics.P3Scheme as P3
julia> FT = Float32;
julia> params = CMP.ParametersP3(FT);
julia> state = P3.getstate(params; Frim = FT(0.5), ρrim = FT(916.7)) P3State{Float32} ├── params = {MassPowerLaw, AreaPowerLaw, SlopePowerLaw, VentilationFactor} ├── Frim = 0.5 [-] └── ρ_rim = 916.7 [kg/m^3] ```
State relationships
CloudMicrophysics.P3Scheme.get_thresholds_ρ_g
— Functionget_thresholds_ρ_g(state::P3State)
get_thresholds_ρ_g(params::CMP.ParametersP3, F_rim, ρ_rim)
Returns
(; D_th, D_gr, D_cr, ρ_g)
: The thresholds for the size distribution, and the density of total (deposition + rime) ice mass for graupel [kg/m³]
See get_D_th
, get_D_gr
, get_D_cr
, and get_ρ_g
for more details.
CloudMicrophysics.P3Scheme.get_ρ_d
— Functionget_ρ_d(mass::MassPowerLaw, F_rim, ρ_rim)
get_ρ_d(state::P3State)
Exact solution for the density of the unrimed portion of the particle as function of the rime mass fraction F_rim
, mass power law parameters mass
, and rime density ρ_rim
.
Arguments
mass
:CMP.MassPowerLaw
parametersF_rim
: rime mass fractionρ_rim
: rime density
Returns
ρ_d
: density of the unrimed portion of the particle [kg/m³]
Examples
julia> import CloudMicrophysics.Parameters as CMP,
ClimaParams as CP,
CloudMicrophysics.P3Scheme as P3
julia> FT = Float64;
julia> mass = CMP.MassPowerLaw(CP.create_toml_dict(FT));
julia> F_rim, ρ_rim = FT(0.5), FT(916.7);
julia> ρ_d = P3.get_ρ_d(mass, F_rim, ρ_rim)
488.9120789986412
CloudMicrophysics.P3Scheme.get_ρ_g
— Functionget_ρ_g(F_rim, ρ_rim, ρ_d)
get_ρ_g(mass::MassPowerLaw, F_rim, ρ_rim)
Return the density of total (deposition + rime) ice mass for graupel [kg/m³]
Arguments
F_rim
: rime mass fraction (L_rim / L_ice
) [-]ρ_rim
: rime density (L_rim / B_rim
) [kg/m³]ρ_d
: density of the unrimed portion of the particle [kg/m³], seeget_ρ_d
Returns
ρ_g
: density of total (deposition + rime) ice mass for graupel [kg/m³]
Notes:
See Eq. 16 in [22].
CloudMicrophysics.P3Scheme.get_D_th
— Functionget_D_th(mass::MassPowerLaw, ρ_i)
Return the critical size separating spherical and nonspherical ice [meters]
See Eq. 8 in [22].
CloudMicrophysics.P3Scheme.get_D_gr
— Functionget_D_gr(mass::MassPowerLaw, ρ_g)
Return the size of equal mass for graupel and unrimed ice [meters]
See Eq. 15 in [22].
CloudMicrophysics.P3Scheme.get_D_cr
— Functionget_D_cr(mass::MassPowerLaw, F_rim, ρ_g)
Return the size of equal mass for graupel and partially rimed ice [meters]
See Eq. 14 in [22].
Derived quantities
Main particle methods
These methods are a function of the particle diameter, D
.
CloudMicrophysics.P3Scheme.ice_mass
— Functionice_mass(state, D)
ice_mass(params, F_rim, ρ_rim, D)
Return the mass of a particle based on where it falls in the particle-size-based properties regime.
Arguments
state
: TheP3State
, orparams, F_rim, ρ_rim
: TheCMP.ParametersP3
, rime mass fraction, and rime density,D
: maximum particle dimension [m]
CloudMicrophysics.P3Scheme.ice_density
— Functionice_density(state::P3State, D)
ice_density(params::CMP.ParametersP3, F_rim, ρ_rim, D)
Return the density of a particle at diameter D
Arguments
state
: TheP3State
, orparams, F_rim, ρ_rim
: TheCMP.ParametersP3
, rime mass fraction, and rime density,D
: maximum particle dimension [m]
Notes:
The density of nonspherical particles is assumed to be the particle mass divided by the volume of a sphere with the same D [22]. Needed for aspect ratio calculation, so we assume zero liquid fraction.
CloudMicrophysics.P3Scheme.∂ice_mass_∂D
— Function∂ice_mass_∂D(state::P3State, D)
∂ice_mass_∂D(params::CMP.ParametersP3, F_rim, ρ_rim, D)
Return the derivative of the ice mass with respect to the particle diameter.
CloudMicrophysics.P3Scheme.ice_area
— Functionice_area(state::P3State, D)
ice_area(params::CMP.ParametersP3, F_rim, ρ_rim, D)
Return the cross-sectional area of a particle based on where it falls in the particle-size-based properties regime.
Arguments
state
:P3State
object, orparams, F_rim, ρ_rim
: TheCMP.ParametersP3
, rime mass fraction, and rime density,D
: maximum particle dimension [m]
CloudMicrophysics.P3Scheme.ϕᵢ
— Functionϕᵢ(state::P3State, D)
ϕᵢ(params::CMP.ParametersP3, F_rim, ρ_rim, D)
Returns the aspect ratio (ϕ) for an ice particle
Arguments
state
: TheP3State
, orparams, F_rim, ρ_rim
: TheCMP.ParametersP3
, rime mass fraction, and rime density,D
: maximum dimension of ice particle [m]
Notes
The density of nonspherical particles is assumed to be equal to the particle mass divided by the volume of a spherical particle with the same D_max [22]. Assuming zero liquid fraction and oblate shape.
Distribution parameters
CloudMicrophysics.P3Scheme.get_μ
— Functionget_μ(slope::CMP.SlopeLaw, logλ)
get_μ(state::P3State, logλ)
Compute the slope parameter μ
Arguments
slope
:CMP.SlopeLaw
object, orstate
:P3State
object, orparams
:CMP.ParametersP3
objectlogλ
: The log of the slope parameter [log(1/m)]
CloudMicrophysics.P3Scheme.get_logN₀
— Functionget_logN₀(N_ice, μ, logλ)
Compute log(N₀)
given the state
, N
, and logλ
,
N = N₀ ∫ G(D) dD
log N₀ = log N - log(∫G(D) dD)
= log(N) - log( ∫D^μ e^{-λD} dD )
= log(N) - M⁰
Arguments
N_ice
: The number concentration [1/m³]μ
: The shape parameter [-
]logλ
: The log of the slope parameter [log(1/m)]
CloudMicrophysics.P3Scheme.get_distribution_logλ
— Functionget_distribution_logλ(params, L_ice, N_ice, F_rim, ρ_rim; [logλ_min, logλ_max])
get_distribution_logλ(state; [logλ_min, logλ_max])
Solve for the distribution parameters given the state, and the mass (L
) and number (N
) concentrations.
The assumed distribution is of the form
\[N′(D) = N₀ D^μ e^{-λD}\]
where N′(D)
is the number concentration at diameter D
and μ
is the slope parameter. The slope parameter is parameterized, e.g. CMP.SlopePowerLaw
or CMP.SlopeConstant
.
This algorithm solves for logλ = log(λ)
and log_N₀ = log(N₀)
given L_ice
and N_ice
by solving the equations:
\[\begin{align*} \log(L) &= \log ∫_0^∞ m(D) N′(D)\ \mathrm{d}D, \\ \log(N) &= \log ∫_0^∞ N′(D)\ \mathrm{d}D, \\ \end{align*}\]
where m(D)
is the mass of a particle at diameter D
(see ice_mass
). The procedure is decribed in detail in the P3 docs.
Arguments
state
:P3State
objectparams
:CMP.ParametersP3
objectL_ice
: The mass concentration [kg/m³]N_ice
: The number concentration [1/m³]F_rim
: The rime mass fractionρ_rim
: The rime density
Keyword arguments
logλ_min
: The minimum value of the search bounds [log(1/m)], default islog(1e1)
logλ_max
: The maximum value of the search bounds [log(1/m)], default islog(1e7)
Distribution relationships
CloudMicrophysics.P3Scheme.logN′ice
— FunctionlogN′ice(state, logλ)
Compute the log of the ice particle number concentration at diameter D
given the distribution dist
CloudMicrophysics.DistributionTools.size_distribution
— Methodsize_distribution(state::P3State, logλ)
Return n(D)
, a function that computes the size distribution for ice particles at diameter D
Arguments
state
: TheP3State
logλ
: The log of the slope parameter [log(1/m)]
CloudMicrophysics.P3Scheme.loggamma_inc_moment
— Functionloggamma_inc_moment(D₁, D₂, μ, logλ, [k = 0], [scale = 1])
Compute log(Iᵏ)
where Iᵏ
is the following integral:
``I^k = ∫_{D₁}^{D₂} G(D) D^k dD``
$G(D) ≡ D^μ e^{-λD}$ is the (unnormalized) gamma kernel, and k
is an arbitrary exponent.
If scale
is provided, log(scale ⋅ Iᵏ)
is returned.
With appropriate scaling, we can compute useful quantities like:
- the
k
-th moment of the ice PSD, $M^k = N₀ I^k$ - combined power law and moment weighted integrals, $∫_{D₁}^{D₂} (aD^b) D^n K(D) dD ≡ a I^(b + n)$
Arguments
D₁
: The minimum diameter [m
]D₂
: The maximum diameter [m
]μ
: The PSD shape parameter [-
]logλ
: The log of the PSD slope parameter [log(1/m)
]k
: An arbitrary exponent [-
], default is0
scale
: The scale factor [-
], default is1
Extended help
Implementation details
We can write ∫_D₁^D₂ G(D) D^k dD
, where G(D) = D^μ e^{-λD}
as: ∫_D₁^∞ G(D) D^k dD - ∫_D₂^∞ G(D) D^k dD
with the transformation x = λD
, and z = μ+k+1
, each term can be written as: ∫_{Dᵢ}^∞ G(D) D^k dD = ∫_{λDᵢ}^∞ x^z e^{-x} dx / λ^z = Γ(z, λDᵢ) / λ^z
where Γ(z, λDᵢ) = q ⋅ Γ(z)
and q
is the incomplete gamma function ratio given by (_, q) = SF.gamma_inc(z, x)
. This means that the integral ∫_{Dᵢ}^∞ G(D) D^k dD
is computed as: Γ(z) ⋅ q / λ^z
The full integral from D₁
to D₂
is then: Γ(z) ⋅ (q_D₁ - q_D₂) / λ^z
In log-space, this is: - z log(λ) + logΓ(z) + log(q_D₁ - q_D₂)
CloudMicrophysics.P3Scheme.loggamma_moment
— Functionloggamma_moment(μ, logλ; [k = 0], [scale = 1])
Compute log(scale ⋅ ∫_0^∞ G(D) D^k dD)
, where G(D) ≡ D^μ e^{-λD}
is the (unnormalized) gamma kernel, k
is an arbitrary exponent, and scale
is a scale factor.
Arguments
μ
: The PSD shape parameter [-
]logλ
: The log of the PSD slope parameter [log(1/m)
]
Keyword arguments
k
: An arbitrary exponent [-
], default is0
scale
: The scale factor [-
], default is1
.
The implementation follows the same logic as loggamma_inc_moment
, but with D₁ = 0
and D₂ = ∞
, which implies q_D₁ = 1
and q_D₂ = 0
.
CloudMicrophysics.P3Scheme.logmass_gamma_moment
— Functionlogmass_gamma_moment(state, logλ; [n=0])
Compute log(∫_0^∞ Dⁿ m(D) N′(D) dD)
given the state
and logλ
. This is the log of the n
-th moment of the mass-weighted PSD.
Arguments
state
:P3State
objectμ
: The shape parameter [-
]logλ
: The log of the slope parameter [log(1/m)]
Keyword arguments
n
: The order of the moment, default is0
Note:
- For
n = 0
, this evaluates tolog(L/N₀)
- For
n = 1
, this evaluates to the (unnormalized) mass-weighted mean particle size, seeD_m
CloudMicrophysics.P3Scheme.logLdivN
— FunctionlogLdivN(state, logλ)
Compute log(L/N)
given the state
and logλ
Arguments
state
:P3State
objectlogλ
: The log of the slope parameter [log(1/m)]
Derived integral quantities
These methods integrate over the particle size distribution.
CloudMicrophysics.P3Scheme.D_m
— FunctionD_m(state::P3State, logλ)
Compute the mass weighted mean particle size [m]
Parameters
state
:P3State
objectlogλ
: The log of the slope parameter [log(1/m)]
CloudMicrophysics.P3Scheme.ice_particle_terminal_velocity
— Functionice_particle_terminal_velocity(state, velocity_params, ρₐ; [use_aspect_ratio])
Returns the terminal velocity of a single ice particle as a function of its size (maximum dimension, D
) using the Chen 2022 parametrization.
Arguments
state
: AP3State
velocity_params
: ACMP.Chen2022VelType
with terminal velocity parametersρₐ
: Air density [kg/m³]
Keyword arguments
use_aspect_ratio
: Bool flag set totrue
if we want to consider the effects of particle aspect ratio on its terminal velocity (default:true
)
CloudMicrophysics.P3Scheme.ice_terminal_velocity_number_weighted
— Functionice_terminal_velocity_number_weighted(state, logλ velocity_params, ρₐ; [use_aspect_ratio], [∫kwargs...])
Return the terminal velocity of the number-weighted mean ice particle size.
Arguments
state
: AP3State
logλ
: The log of the slope parameter [log(1/m)]velocity_params
: ACMP.Chen2022VelType
with terminal velocity parametersρₐ
: Air density [kg/m³]
Keyword arguments
use_aspect_ratio
: Bool flag set totrue
if we want to consider the effects of particle aspect ratio on its terminal velocity (default:true
)∫kwargs...
: Additional optional keyword arguments to pass to∫fdD
See also ice_terminal_velocity_mass_weighted
CloudMicrophysics.P3Scheme.ice_terminal_velocity_mass_weighted
— Functionice_terminal_velocity_mass_weighted(state, logλ, velocity_params, ρₐ; [use_aspect_ratio], [∫kwargs...])
Return the terminal velocity of the mass-weighted mean ice particle size.
Arguments
state
: AP3State
logλ
: The log of the slope parameter [log(1/m)]velocity_params
: ACMP.Chen2022VelType
with terminal velocity parametersρₐ
: Air density [kg/m³]
Keyword arguments
use_aspect_ratio
: Bool flag set totrue
if we want to consider the effects of particle aspect ratio on its terminal velocity (default:true
)∫kwargs...
: Keyword arguments passed to∫fdD
CloudMicrophysics.P3Scheme.het_ice_nucleation
— Functionhet_ice_nucleation(pdf_c, p3, tps, q, N, T, ρₐ, p, aerosol)
- aerosol - aerosol parameters (supported types: desert dust, illite, kaolinite)
- tps - thermodynamics parameters
- qₚ - phase partition
- N_liq - cloud water number concentration
- RH - relative humidity
- T - temperature
- ρₐ - air density
- dt - model time step
Returns a named tuple with ice number concentration and ice content hetergoeneous freezing rates from cloud droplets.
CloudMicrophysics.P3Scheme.ice_melt
— Functionice_melt(state, logλ, Chen2022, aps, tps, Tₐ, ρₐ, dt; ∫kwargs...)
Arguments
state
: aP3State
objectlogλ
: the log of the slope parameter [log(1/m)]Chen2022
: struct containing Chen 2022 velocity parametersaps
: air propertiestps
: thermodynamics parametersTₐ
: temperature (K)ρₐ
: air densitydt
: model time step (for limiting the tendnecy)
Keyword arguments
∫kwargs
: Named tuple of keyword arguments passed to∫fdD
Returns the melting rate of ice (QIMLT in Morrison and Mildbrandt (2015)).
Supporting integral methods
CloudMicrophysics.P3Scheme.∫fdD
— Function∫fdD(f, dist; [p = 1e-6], kwargs...)
Integrate the function f
over the size distribution dist
This function is useful for integrating functions over the size distribution dist
. It is a light wrapper around QGK.quadgk
that automatically inserts appropriate size distribution thresholds as integration limits.
This method calls ∫fdD_error
, which returns both the value of the integral and the estimated error. Since the error is typically not of interest, this method only returns the value of the integral.
Arguments
f
: The function to integratedist
: The distribution object, passed tointegral_bounds
to set the integration limits.p
: The integration bounds are set to thep
-th and1-p
-th quantiles of the size distribution. Default:p = 1e-6
(i.e. 99.9998% of the size distribution is integrated).kwargs
: Additional optional keyword arguments to pass toQGK.quadgk
rtol
: The relative tolerance for the integration, default:rtol = sqrt(eps(FT))
atol
: The absolute tolerance for the integration, default:atol = 0
maxevals
: The maximum number of function evaluations, default:maxevals = 10^7
order
: The order of the quadrature rule, default:order = 7
Returns
value
: The value of the integral
To achieve highest accuracy, which can be challenging when integrating over the size distribution, it is recommended to increase the order
of the quadrature rule and set rtol = 0
. Experimentally, order = 44
has been found to be sufficient for most cases.
For convenience, passing accurate = true
will set rtol = 0
and order = 44
.
Examples
julia> import CloudMicrophysics.Parameters as CMP,
CloudMicrophysics.P3Scheme as P3
julia> params = CMP.ParametersP3(Float64);
julia> state = P3.get_state(params; F_rim = 0.0, ρ_rim = 400.0, L_ice = 0.002, N_ice = 1000.0);
julia> logλ = P3.get_distribution_logλ(state);
julia> f(D) = D^3 * P3.size_distribution(state, logλ)(D); # Define a function to integrate
julia> P3.∫fdD(f, state, logλ; p = 0.01) # Integrate the function
0.0008519464332296608
julia> P3.∫fdD(state, logλ; p = 0.01) do D # Integrate with a `do`-block
P3.ice_mass(state, D) * P3.size_distribution(state, logλ)(D)
end
0.0017027833723511712
CloudMicrophysics.P3Scheme.∫fdD_error
— Function∫fdD_error(f, dist; p, [moment_order = 0], [accurate = false], kwargs...)
Integrate the function f
over the size distribution dist
Returns
value
: The value of the integralerror
: The estimated error of the integral
Notes
See ∫fdD
, which only returns the value of the integral and not the error, for details.
CloudMicrophysics.P3Scheme.integral_bounds
— Functionintegral_bounds(state::P3State, logλ; p, moment_order = 0)
Compute the integration bounds for the P3 size distribution,
Mⁿ = ∫_a^b Dⁿ * N′(D) dD = N₀ * ∫_a^b Dⁿ D^μ * exp(-λ * D) dD
where Mⁿ
is the n
-th moment of the size distribution. Here n ≡ moment_order
and a
and b
are the integration bounds. For a proper moment, a=0
and b=∞
. For the numerical integration, a
and b
are determined by this function.
Arguments
state
:P3State
objectlogλ
: The log of the slope parameter [log(1/m)]p
: The integration bounds are set to thep
-th and1-p
-th quantiles of the size distribution.moment_order
: For integrands proportional to moments of the size distribution,moment_order
can be used to indicate the order of the moment. May provide more accurate bounds; thus more accurate integration. Default:moment_order = 0
.
Returns
bnds
: The integration bounds (aTuple
), for use in [QGK.quadgk
].
Aerosol model
CloudMicrophysics.AerosolModel
— ModuleA container for information on aerosol size distribution
and chemical properties.
The size distribution is a sum of lognormal internally mixed modes.
The chemical composition can be expressed using kappa parameter
or hygroscopicity parameter B.
CloudMicrophysics.AerosolModel.Mode_B
— TypeMode_B
Represents the sizes and chemical composition of aerosol particles in one size distribution mode. The mode is assumed to be made up of internally mixed components and follow a lognormal size distribution. The chemical composition of aerosol particles in this mode is described using the parameters from Abdul-Razzak and Ghan 2000.
CloudMicrophysics.AerosolModel.Mode_κ
— TypeMode_κ
Represents the sizes and chemical composition of aerosol particles in one size distribution mode. The mode is assumed to be made up of internally mixed components and follow a lognormal size distribution. The chemical composition of aerosol particles in this mode is described using the parameters from Petters and Kreidenweis 2007.
CloudMicrophysics.AerosolModel.AerosolDistribution
— TypeAerosolDistribution
Represents the aerosol size distribution as a tuple with different modes. All modes have to either be of type ModeB (Abdul-Razzak and Ghan 2000) or of type Modeκ (Petters and Kreidenweis 2007).
Constructors
AerosolDistribution(modes::T)
Aerosol activation
CloudMicrophysics.AerosolActivation
— ModuleAerosol activation scheme, which includes:
- mean hygroscopicity for each mode of the aerosol size distribution
- critical supersaturation for each mode of the aerosol size distribution
- maximum supersaturation
- total number of particles actived
- total mass of particles actived
CloudMicrophysics.AerosolActivation.mean_hygroscopicity_parameter
— Functionmean_hygroscopicity_parameter(ap, ad)
ap
- a struct with aerosol activation parametersad
- a struct with aerosol distribution (B or κ based)
Returns a tuple of hygroscopicity parameters (one tuple element for each aerosol size distribution mode). The tuple is computed either as mass-weighted B parameters (Abdul-Razzak and Ghan 2000) or volume weighted kappa parameters (Petters and Kreidenweis 2007). Implemented via a dispatch based on aerosol distribution mode type.
CloudMicrophysics.AerosolActivation.max_supersaturation
— Functionmax_supersaturation(ap, ad, aip, tps, T, p, w, q)
ap
- a struct with aerosol activation parametersad
- a struct with aerosol distributionaip
- a struct with air parameterstps
- a struct with thermodynamics parametersT
- air temperaturep
- air pressurew
- vertical velocityq
- phase partition
Returns the maximum supersaturation.
CloudMicrophysics.AerosolActivation.N_activated_per_mode
— FunctionN_activated_per_mode(ap, ad, aip, tps, T, p, w, q)
ap
- a struct with aerosol activation parametersad
- aerosol distribution structaip
- a struct with air parameterstps
- a struct with thermodynamics parametersT
- air temperaturep
- air pressurew
- vertical velocityq
- phase partition
Returns the number of activated aerosol particles in each aerosol size distribution mode.
CloudMicrophysics.AerosolActivation.M_activated_per_mode
— FunctionM_activated_per_mode(ap, ad, aip, tps, T, p, w, q)
ap
- a struct with aerosol activation parametersad
- a struct with aerosol distribution parametersaip
- a struct with air parameterstps
- a struct with thermodynamics parametersT
- air temperaturep
- air pressurew
- vertical velocityq
- phase partition
Returns the mass of activated aerosol particles per mode of the aerosol size distribution.
CloudMicrophysics.AerosolActivation.total_N_activated
— Functiontotal_N_activated(ap, ad, aip, tps, T, p, w, q)
ap
- a struct with aerosol activation parametersad
- aerosol distribution structaip
- a struct with air propertiestps
- a struct with thermodynamics parametersT
- air temperaturep
- air pressurew
- vertical velocityq
- phase partition
Returns the total number of activated aerosol particles.
CloudMicrophysics.AerosolActivation.total_M_activated
— Functiontotal_M_activated(ap, ad, aip, tps, T, p, w, q)
ap
- a struct with aerosol activation parametersad
- aerosol distribution structaip
- a struct with air propertiestps
- a struct with thermodynamics parametersT
- air temperaturep
- air pressurew
- vertical velocityq
- phase partition
Returns the total mass of activated aerosol particles.
Artifact calling
CloudMicrophysics.ArtifactCalling
— ModuleCall artifacts from Artifacts.toml
CloudMicrophysics.ArtifactCalling.AIDA_ice_nucleation
— FunctionAIDA_ice_nucleation(data_file_name)
data_file_name
- name of the data file on Caltech box.
Returns the filepath of the data file in Caltech box.
Heterogeneous ice nucleation
CloudMicrophysics.HetIceNucleation
— ModuleParameterization for heterogenous ice nucleation.
CloudMicrophysics.HetIceNucleation.dust_activated_number_fraction
— Functiondust_activated_number_fraction(dust, ip, Si, T)
dust
- a struct with dust parametersip
- a struct with ice nucleation parametersSi
- ice saturation ratioT
- air temperature [K]
Returns the number fraction of mineral dust particles acting as deposition nuclei (n ice nuclei / n dust particles). From Mohler et al 2006 Table 2 (averages from different measurements excluding those where a was not measured)
CloudMicrophysics.HetIceNucleation.MohlerDepositionRate
— FunctionMohlerDepositionRate(dust, ip, S_i, T, dSi_dt, N_aer)
dust
- a struct with dust parametersip
- a struct with ice nucleation parametersSi
- ice saturationT
- ambient temperaturedSi_dt
- change in ice saturation over timeN_aer
- number of unactivated aerosols
Returns the ice nucleation rate from deposition. From Mohler et al 2006 equation 5.
CloudMicrophysics.HetIceNucleation.deposition_J
— Functiondeposition_J(dust, Δa_w)
dust
- a struct with dust parametersΔa_w
- change in water activity [unitless].
Returns the deposition nucleation rate coefficient, J
, in m^-2 s^-1 for different minerals in liquid droplets. The free parameters m
and c
are derived from China et al (2017) see DOI: 10.1002/2016JD025817 Returns zero for unsupported aerosol types.
CloudMicrophysics.HetIceNucleation.ABIFM_J
— FunctionABIFM_J(dust, Δa_w)
dust
- a struct with dust parametersΔa_w
- change in water activity [unitless].
Returns the immersion freezing nucleation rate coefficient, J
, in m^-2 s^-1 for different minerals in liquid droplets. The free parameters m
and c
are taken from Knopf & Alpert 2013 see DOI: 10.1039/C3FD00035D Returns zero for unsupported aerosol types.
CloudMicrophysics.HetIceNucleation.P3_deposition_N_i
— FunctionP3_deposition_N_i(ip, T)
ip
- a struct with ice nucleation parameters,T
- air temperature [K].
Returns the number of ice nucleated via deposition nucleation with units of m^-3. From Thompson et al 2004 eqn 2 as used in Morrison & Milbrandt 2015.
CloudMicrophysics.HetIceNucleation.P3_het_N_i
— FunctionP3_het_N_i(ip, T, N_l, B, V_l, a, Δt)
ip
- a struct with ice nucleation parameters,T
- air temperature [K],N_l
- number of droplets [m^-3],B
- water-type dependent parameter [cm^-3 s^-1],V_l
- volume of droplets to be frozen [m^3],a
- empirical parameter [C^-1],Δt
- timestep.
Returns the number of ice nucleated within Δt via heterogeneous freezing with units of m^-3. From Pruppacher & Klett 1997 eqn (9-51) as used in Morrison & Milbrandt 2015.
CloudMicrophysics.HetIceNucleation.INP_concentration_frequency
— FunctionINP_concentration_frequency(params,INPC,T)
params
- a struct with INPC(T) distribution parametersINPC
- concentration of ice nucleating particles [m^-3]T
- air temperature [K]
Returns the relative frequency of a given INP concentration, depending on the temperature. Based on Frostenberg et al., 2023. See DOI: 10.5194/acp-23-10883-2023
CloudMicrophysics.HetIceNucleation.INP_concentration_mean
— FunctionINP_concentration_mean(T)
T
- air temperature [K]
Returns the logarithm of mean INP concentration (in m^-3), depending on the temperature. Based on the function μ(T) in Frostenberg et al., 2023. See DOI: 10.5194/acp-23-10883-2023
Homogeneous ice nucleation
CloudMicrophysics.HomIceNucleation
— ModuleParameterization for homogeneous ice nucleation
CloudMicrophysics.HomIceNucleation.homogeneous_J_cubic
— Functionhomogeneous_J_cubic(ip, Δa_w)
ip
- a struct with ice nucleation parametersΔa_w
- change in water activity [-].
Returns the homogeneous freezing nucleation rate coefficient, J
, in m^-3 s^-1 for sulphuric acid solutions. Parameterization based on Koop 2000, DOI: 10.1038/35020537.
CloudMicrophysics.HomIceNucleation.homogeneous_J_linear
— Functionhomogeneous_J_linear(ip, Δa_w)
ip
- a struct with ice nucleation parametersΔa_w
- change in water activity [-].
Returns the homogeneous freezing nucleation rate coefficient, J
, in m^-3 s^-1 for sulphuric acid solutions. Parameterization derived from a linear fit of the Koop 2000 parameterization, DOI: 10.1038/35020537.
Cloud diagnostics
CloudMicrophysics.CloudDiagnostics
— ModuleCloud diagnostics
- radar reflectivity (1-moment and 2-moment)
- effective radius (1-moment and 2-moment)
CloudMicrophysics.CloudDiagnostics.radar_reflectivity_1M
— Functionradar_reflectivity_1M(precip, q, ρ)
precip
- struct with 1-moment microphysics rain free parametersq
- specific content of rainρ
- air density
Returns logarithmic radar reflectivity for the 1-moment microphysics based on the assumed rain particle size distribution. Normalized by the reflectivty of 1 millimiter drop in a volume of 1m3. The values are clipped at -150 dBZ.
CloudMicrophysics.CloudDiagnostics.radar_reflectivity_2M
— Functionradar_reflectivity_2M(structs, q_liq, q_rai, N_liq, N_rai, ρ_air)
- `structs` - structs microphysics 2-moment with SB2006 cloud droplets
and raindrops size distributions parameters
- `q_liq` - cloud water specific content
- `q_rai` - rain water specific content
- `N_liq` - cloud droplet number density
- `N_rai` - rain droplet number density
- `ρ_air` - air density
Returns logarithmic radar reflectivity for the 2-moment microphysics SB2006 based on the assumed cloud and rain particle size distribuions. Normalized by the reflectivty of 1 millimiter drop in a volume of 1m3. The values are clipped at -150 dBZ.
CloudMicrophysics.CloudDiagnostics.effective_radius_const
— Functioneffective_radius_const(cloud_params)
- cloud_params - a struct with cloud liquid or cloud ice parameters
Returns a constant assumed effective radius for clouds
CloudMicrophysics.CloudDiagnostics.effective_radius_Liu_Hallet_97
— Functioneffective_radius_Liu_Hallet_97(q_liq, q_rai, N_liq, N_rai, ρ_air, ρ_w)
- `q_liq` - cloud water specific content
- `q_rai` - rain water specific content
- `N_liq` - cloud droplet number density
- `N_rai` - rain droplet number density
- `ρ_air` - air density
- `ρ_w` - water density
Returns effective radius using the "1/3" power law from Liu and Hallett (1997). If not provided by the user, it is assumed that there is no rain present and that the cloud droplet number concentration is 100 1/cm3.
CloudMicrophysics.CloudDiagnostics.effective_radius_2M
— Functioneffective_radius_2M(structs, q_liq, q_rai, N_liq, N_rai, ρ_air)
- `structs` - structs with SB2006 cloud droplets and raindrops
size distribution parameters
- `q_liq` - cloud water specific content
- `q_rai` - rain water specific content
- `N_liq` - cloud droplet number density
- `N_rai` - rain droplet number density
- `ρ_air` - air density
Returns effective radius for the 2-moment microphysics scheme. Computed based on the assumed cloud and rain particle size distributions.
Common utility functions
CloudMicrophysics.Common
— ModuleModule for functions shared by different parameterizations.
CloudMicrophysics.Common.G_func
— FunctionG_func(air_props, tps, T, Liquid())
G_func(air_props, tps, T, Ice())
air_props
- struct with air parameterstps
- struct with thermodynamics parametersT
- air temperatureLiquid()
,Ice()
- liquid or ice phase to dispatch over.
Utility function combining thermal conductivity and vapor diffusivity effects.
CloudMicrophysics.Common.logistic_function
— Functionlogistic_function(x, x_0, k)
x
- independent variablex_0
- threshold value for xk
- growth rate of the curve, characterizing steepness of the transition
Returns the value of the logistic function for smooth transitioning at thresholds. This is a normalized curve changing from 0 to 1 while x varies from 0 to Inf (for positive k). For x < 0 the value at x = 0 (zero) is returned. For x_0 = 0 H(x) is returned.
CloudMicrophysics.Common.logistic_function_integral
— Functionlogistic_function_integral(x, x_0, k)
x
- independent variablex_0
- threshold value for xk
- growth rate of the logistic function, characterizing steepness of the transition
Returns the value of the indefinite integral of the logistic function, for smooth transitioning of piecewise linear profiles at thresholds. This curve smoothly transition from y = 0 for 0 < x < x0 to y = x - x0 for x_0 < x.
CloudMicrophysics.Common.H2SO4_soln_saturation_vapor_pressure
— FunctionH2SO4_soln_saturation_vapor_pressure(prs, x, T)
prs
- a struct with H2SO4 solution free parametersx
- wt percent sulphuric acid [unitless]T
- air temperature [K].
Returns the saturation vapor pressure above a sulphuric acid solution droplet in Pa. x
is, for example, 0.1 if droplets are 10 percent sulphuric acid by weight
CloudMicrophysics.Common.a_w_xT
— Functiona_w_xT(H2SO4_prs, tps, x, T)
H2SO4_prs
- a struct with H2SO4 solution free parameterstps
- a struct with thermodynamics parametersx
- wt percent sulphuric acid [unitless]T
- air temperature [K].
Returns water activity of H2SO4 containing droplet. x
is, for example, 0.1 if droplets are 10 percent sulphuric acid by weight.
CloudMicrophysics.Common.a_w_eT
— Functiona_w_eT(tps, e, T)
tps
- struct with thermodynamics parameterse
- partial pressure of water [Pa]T
- air temperature [K].
Returns water activity of pure water droplet. Valid when droplet is in equilibrium with surroundings.
CloudMicrophysics.Common.a_w_ice
— Functiona_w_ice(tps, T)
tps
- struct with thermodynamics parametersT
- air temperature [K].
Returns water activity of ice.
CloudMicrophysics.Common.Chen2022_vel_coeffs
— FunctionChen2022_vel_coeffs(coeffs, ρₐ)
Chen2022_vel_coeffs(coeffs, ρₐ, ρᵢ)
Compute the coefficients for the Chen 2022 terminal velocity parametrization.
Arguments
coeffs
: a struct with terminal velocity free parametersCMP.Chen2022VelTypeRain
: Fetch from Table B1CMP.Chen2022VelTypeSmallIce
: Fetch from Table B2CMP.Chen2022VelTypeLargeIce
: Fetch from Table B4
ρₐ
: air density [kg/m³]ρᵢ
: apparent density of ice particles [kg/m³], only used forCMP.Chen2022VelTypeSmallIce
andCMP.Chen2022VelTypeLargeIce
See [12] for more details.
CloudMicrophysics.Common.Chen2022_monodisperse_pdf
— FunctionChen2022_monodisperse_pdf(a, b, c)
Arguments
a
,b
,c
: free parameters defined in [12]
Returns
pdf(D)
: The monodisperse particle distribution function as a function of diameter,D
, in [m/s].
CloudMicrophysics.Common.Chen2022_exponential_pdf
— FunctionChen2022_exponential_pdf(a, b, c, λ_inv, k)
- a, b, c, - free parameters defined in Chen etl al 2022
- λ_inv - inverse of the size distribution parameter
- k - size distribution moment for which we compute the bulk fall speed
Returns the addends of the bulk fall speed of rain or ice particles following Chen et al 2022 DOI: 10.1016/j.atmosres.2022.106171 in [m/s]. Assuming exponential size distribution and hence μ=0.
CloudMicrophysics.Common.volume_sphere_D
— Functionvolume_sphere_D(D)
Calculate the volume of a sphere with diameter D.
\[V = D^3 * π / 6\]
See also volume_sphere_R
.
CloudMicrophysics.Common.volume_sphere_R
— Functionvolume_sphere_R(R)
Calculate the volume of a sphere with radius R.
\[V = (2R)^3 * π / 6\]
See also volume_sphere_D
.
CloudMicrophysics.Common.ventilation_factor
— Functionventilation_factor(vent, aps, v_term)
Returns a function that computes the ventilation factor for a particle as a function of its diameter, D
.
The ventilation factor parameterizes the increase in the mass and heat exchange for falling particles.
Arguments
vent
: Ventilation parameterization constants,CMP.VentilationFactor
aps
: Parameters with air properties,CMP.AirProperties
v_term
: A functionv_term(D)
that returns the terminal velocity of a particle with diameterD
Returns
F_v(D)
: The ventilation factor as a function of diameter,D
See e.g. [10] Eq. (24) for the definition of the ventilation factor.
Parameters
CloudMicrophysics.Parameters
— ModuleParameters
A module for CloudMicrophysics.jl free parameters.
CloudMicrophysics.Parameters.ParametersType
— TypeParametersType{FT}
The top-level super-type for all cloud microphysics parameters
CloudMicrophysics.Parameters.AerosolType
— TypeAerosolType{FT}
The top-level super-type for all aerosol properties
CloudMicrophysics.Parameters.AerosolDistributionType
— TypeAerosolDistributionType
The top-level super-type for all aerosol distribution types
CloudMicrophysics.Parameters.CloudCondensateType
— TypeCloudCondensateType{FT}
The top-level super-type for cloud condensate types (liquid and ice)
CloudMicrophysics.Parameters.PrecipitationType
— TypePrecipitationType{FT}
The top-level super-type for precipitation types (rain and snow)
CloudMicrophysics.Parameters.TerminalVelocityType
— TypeTerminalVelocityType{FT}
The top-level super-type for terminal velocity parameterizations
CloudMicrophysics.Parameters.Precipitation2MType
— TypePrecipitation2MType
The top-level super-type for 2-moment precipitation parameterizations
CloudMicrophysics.Parameters.AirProperties
— TypeAirProperties{FT}
Parameters with air properties.
Fields
K_therm
: thermal conductivity of air [w/m/K]D_vapor
: diffusivity of water vapor [m2/s]ν_air
: kinematic viscosity of air [m2/s]
CloudMicrophysics.Parameters.WaterProperties
— TypeWaterProperties{FT}
Parameters with water properties.
Fields
ρw
: density of liquid water [kg/m3]ρi
: density of ice [kg/m3]
CloudMicrophysics.Parameters.ArizonaTestDust
— TypeArizonaTestDust{FT}
Parameters for Arizona Test Dust from Mohler et al, 2006. DOI: 10.5194/acp-6-3007-2006
Fields
S₀_warm
: S₀ for T > T_thr [-]S₀_cold
: S₀ for T < T_thr [-]a_warm
: a for T > T_thr [-]a_cold
: a for T < T_thr [-]deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.DesertDust
— TypeDesertDust{FT}
Parameters for desert dust from Knopf and Alpert 2013 DOI: 10.1039/C3FD00035D and from Mohler et al, 2006 DOI: 10.5194/acp-6-3007-2006
Fields
S₀_warm
: S₀ for T > T_thr [-]S₀_cold
: S₀ for T < T_thr [-]a_warm
: a for T > T_thr [-]a_cold
: a for T < T_thr [-]ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.AsianDust
— TypeAsianDust{FT}
Parameters for Asian Dust
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.MiddleEasternDust
— TypeMiddleEasternDust{FT}
Parameters for Middle Eastern Dust
Fields
ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.SaharanDust
— TypeSaharanDust{FT}
Parameters for Saharan Dust
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]
CloudMicrophysics.Parameters.Illite
— TypeIllite{FT}
Parameters for illite from Knopf and Alpert 2013 DOI: 10.1039/C3FD00035D
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.Kaolinite
— TypeKaolinite{FT}
Parameters for kaolinite from Knopf and Alpert 2013 DOI: 10.1039/C3FD00035D and China et al 2017 DOI: 10.1002/2016JD025817
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.Feldspar
— TypeFeldspar{FT}
Parameters for Feldspar from Alpert et al 2022 DOI: 10.1039/D1EA00077B
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]
CloudMicrophysics.Parameters.Ferrihydrite
— TypeFerrihydrite{FT}
Parameters for Ferrihydrite from Alpert et al 2022 DOI: 10.1039/D1EA00077B
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]
CloudMicrophysics.Parameters.Dust
— TypeDust{FT}
Parameters for generic dust.
Fields
deposition_m
: m coefficient for deposition nucleation J [-]deposition_c
: c coefficient for deposition nucleation J [-]ABIFM_m
: m coefficient for immersion freezing J [-]ABIFM_c
: c coefficient for immersion freezing J [-]
CloudMicrophysics.Parameters.Seasalt
— TypeSeasalt{FT}
Parameters for seasalt
Fields
M
: molar mass [kg/mol]ρ
: density [kg/m3]ϕ
: osmotic coefficient [-]ν
: ion number [-]ϵ
: water soluble mass fraction [-]κ
: hygroscopicity parameter [-]
CloudMicrophysics.Parameters.Sulfate
— TypeSulfate{FT}
Parameters for sulfate aerosol
Fields
M
: molar mass [kg/mol]ρ
: density [kg/m3]ϕ
: osmotic coefficient [-]ν
: ion number [-]ϵ
: water soluble mass fraction [-]κ
: hygroscopicity parameter [-]
CloudMicrophysics.Parameters.AerosolActivationParameters
— TypeAerosolActivationParameters{FT}
Parameters for Abdul-Razzak and Ghan 2000 aerosol activation scheme DOI: 10.1029/1999JD901161
Fields
M_w
: molar mass of water [kg/mol]R
: gas constant [J/mol/K]ρ_w
: cloud water density [kg/m3]σ
: surface tension of water [N/m]g
: gravitational_acceleration [m/s2]f1
: scaling coefficient in Abdul-Razzak and Ghan 2000 [-]f2
: scaling coefficient in Abdul-Razzak and Ghan 2000 [-]g1
: scaling coefficient in Abdul-Razzak and Ghan 2000 [-]g2
: scaling coefficient in Abdul-Razzak and Ghan 2000 [-]p1
: power of (zeta / eta) in Abdul-Razzak and Ghan 2000 [-]p2
: power of (S_m^2 / (zeta + 3 * eta)) in Abdul-Razzak and Ghan 2000 [-]
CloudMicrophysics.Parameters.IceNucleationParameters
— TypeIceNucleationParameters{FT, DEP, HOM, P3_type}
Parameters for ice nucleation
Fields
deposition
homogeneous
p3
CloudMicrophysics.Parameters.Frostenberg2023
— TypeFrostenberg2023{FT}
Parameters for frequency distribution of INP concentration DOI: 10.5194/acp-23-10883-2023
Fields
σ
: standard deviationa
: coefficientb
: coefficient
CloudMicrophysics.Parameters.H2SO4SolutionParameters
— TypeH2SO4SolutionParameters{FT}
Parameters for water activity of H2SO4 solutions from Luo et al 1995. DOI: 10.1029/94GL02988
Fields
T_max
: max temperature for which the parameterization is valid [K]T_min
: min temperature for which the parameterization is valid [K]w_2
: coefficient [-]c1
: coefficient [-]c2
: coefficient [-]c3
: coefficient [-]c4
: coefficient [-]c5
: coefficient [-]c6
: coefficient [-]c7
: coefficient [-]
CloudMicrophysics.Parameters.Mohler2006
— TypeMohler2006{FT}
Parameters for ice nucleation from Mohler et al 2006 DOI: 10.5194/acp-6-3007-2006
Fields
Sᵢ_max
: max allowed supersaturation [-]T_thr
: threshold temperature [K]
CloudMicrophysics.Parameters.Koop2000
— TypeKoop2000{FT}
Parameters for ice nucleation from Koop et al 2000 DOI: 10.1038/35020537
Fields
Δa_w_min
: min Δaw [-]Δa_w_max
: max Δaw [-]c₁
: coefficient [-]c₂
: coefficient [-]c₃
: coefficient [-]c₄
: coefficient [-]linear_c₁
: coefficient [-]linear_c₂
: coefficient [-]
CloudMicrophysics.Parameters.H2S04NucleationParameters
— TypeH2S04NucleationParameters{FT}
Parameters for pure sulfuric acid nucleation from Dunne et al 1016 DOI:10.1126/science.aaf2649
Fields
p_b_n
p_b_i
u_b_n
u_b_i
v_b_n
v_b_i
w_b_n
w_b_i
p_t_n
p_t_i
u_t_n
u_t_i
v_t_n
v_t_i
w_t_n
w_t_i
p_A_n
p_A_i
a_n
a_i
CloudMicrophysics.Parameters.OrganicNucleationParameters
— TypeOrganicNucleationParameters{FT}
Parameters for pure organic nucleation from Kirkby 2016 DOI: 10.1038/nature17953
Fields
a_1
a_2
a_3
a_4
a_5
Y_MTO3
Y_MTOH
k_MTO3
k_MTOH
exp_MTO3
exp_MTOH
CloudMicrophysics.Parameters.MixedNucleationParameters
— TypeMixedNucleationParameters{FT}
Parameters for mixed organic and sulfuric acid nucleation from Riccobono et al 2014 DOI:10.1126/science.1243527
Fields
k_H2SO4org
k_MTOH
exp_MTOH
CloudMicrophysics.Parameters.Parameters0M
— TypeParameters0M{FT}
Parameters for zero-moment bulk microphysics scheme
Fields
τ_precip
: precipitation timescale [s]qc_0
: condensate specific content precipitation threshold [-]S_0
: supersaturation precipitation threshold [-]
CloudMicrophysics.Parameters.ParticlePDFSnow
— TypeParticlePDFSnow{FT}
A struct with snow size distribution parameters
Fields
μ
: snow size distribution coefficient [1/m4]ν
: snow size distribution coefficient [-]
CloudMicrophysics.Parameters.ParticlePDFIceRain
— TypeParticlePDFIceRain{FT}
A struct with snow size distribution parameters
Fields
n0
: Size distribution coefficient [1/m4]
CloudMicrophysics.Parameters.ParticleMass
— TypeParticleMass{FT}
A struct with coefficients of the assumed mass(size) relationship for particles
m(r) = m0 χm (r/r0)^(me + Δm)
Fields
r0
: particle length scale [m]m0
: mass size relation coefficient [kg]me
: mass size relation coefficient [-]Δm
: mass size relation coefficient [-]χm
: mass size relation coefficient [-]
CloudMicrophysics.Parameters.ParticleArea
— TypeParticleArea{FT}
A struct with coefficients of the assumed crosssectionarea(size) relationship for particles
a(r) = a0 χa (r/r0)^(ae + Δa)
Fields
a0
: cross section size relation coefficient [m2]ae
: cross section size relation coefficient [-]Δa
: cross section size relation coefficient [-]χa
: cross section size relation coefficient [-]
CloudMicrophysics.Parameters.SnowAspectRatio
— TypeSnowAspectRatio{FT}
A struct with aspect ratio coefficients
Fields
ϕ
: aspect ratio [-]κ
: power law coefficient in terminal velocity parameterization from Chen et al 2022 [-]
CloudMicrophysics.Parameters.Ventilation
— TypeVentilation{FT}
A struct with ventilation coefficients
Fields
a
: ventilation coefficienta
[-]b
: ventilation coefficientb
[-]
CloudMicrophysics.Parameters.Acnv1M
— TypeAcnv1M{FT}
A struct with autoconversion parameters
Fields
τ
: autoconversion timescale [s]q_threshold
: condensate specific content autoconversion threshold [-]k
: threshold smooth transition steepness [-]
CloudMicrophysics.Parameters.CloudLiquid
— TypeCloudLiquid{FT}
The parameters and type for cloud liquid water condensate
Fields
τ_relax
: condensation evaporation non_equil microphysics relaxation timescale [s]ρw
: water density [kg/m3]r_eff
: effective radius [m]
CloudMicrophysics.Parameters.CloudIce
— TypeCloudIce{FT, MS}
The parameters and type for cloud ice condensate
Fields
pdf
: a struct with size distribution parametersmass
: a struct with mass size relation parametersr0
: particle length scale [m]r_ice_snow
: ice snow threshold radius [m]τ_relax
: deposition sublimation non_equil microphysics relaxation timescale [s]ρᵢ
: cloud ice apparent density [kg/m3]r_eff
: effective radius [m]
CloudMicrophysics.Parameters.Rain
— TypeRain{FT, MS, AR, VT, AC}
The parameters and type for rain
Fields
pdf
: a struct with size distribution parametersmass
: a struct with mass size relation parametersarea
: a struct with cross section size relation parameersvent
: a struct with ventilation coefficientsacnv1M
: a struct with cloud water to rain autoconversion parametersr0
: particle length scale [m]
CloudMicrophysics.Parameters.Snow
— TypeSnow{FT, PD, MS, AR, VT, AP, AC}
The parameters and type for snow
Fields
pdf
: a struct with size distribution parametersmass
: a struct with mass size relation parametersarea
: a struct with cross section size relation parameersvent
: a struct with ventilation coefficientsaspr
: a struct with aspect ratio parametersacnv1M
: a struct with ice to snow autoconversion parametersr0
: particle length scale [m]T_freeze
: freezing temperature of water [K]ρᵢ
: snow apparent density [kg/m3]
CloudMicrophysics.Parameters.CollisionEff
— TypeCollisionEff{FT}
Collision efficiency parameters for the 1-moment scheme
Fields
e_liq_rai
: cloud liquid-rain collision efficiency [-]e_liq_sno
: cloud liquid-snow collision efficiency [-]e_ice_rai
: cloud ice-rain collision efficiency [-]e_ice_sno
: cloud ice-snow collision efficiency [-]e_rai_sno
: rain-snow collision efficiency [-]
CloudMicrophysics.Parameters.KK2000
— TypeKK2000
The type and parameters for 2-moment precipitation formation by Khairoutdinov and Kogan (2000)
DOI:10.1175/1520-0493(2000)128<0229:ANCPPI>2.0.CO;2
Fields
acnv
: Autoconversion parametersaccr
: Accretion parameters
CloudMicrophysics.Parameters.AcnvKK2000
— TypeAcnvKK2000
Khairoutdinov and Kogan (2000) autoconversion parameters
Fields
A
: Autoconversion coefficient Aa
: Autoconversion coefficient ab
: Autoconversion coefficient bc
: Autoconversion coefficient c
CloudMicrophysics.Parameters.AccrKK2000
— TypeAccrKK2000
Khairoutdinov and Kogan (2000) accretion parameters
Fields
A
: Accretion coefficient Aa
: Accretion coefficient ab
: Accretion coefficient b
CloudMicrophysics.Parameters.B1994
— TypeB1994
The type and parameter for 2-moment precipitation formation by Beheng (1994) DOI: 10.1016/0169-8095(94)90020-5
Fields
acnv
: Autoconversion coeff Caccr
: Autoconversion coeff a
CloudMicrophysics.Parameters.AcnvB1994
— TypeAcnvB1994
Beheng (1994) autoconversion parameters
Fields
C
: Autoconversion coeff Ca
: Autoconversion coeff ab
: Autoconversion coeff bc
: Autoconversion coeff cN_0
: Autoconversion coeff N_0d_low
: Autoconversion coeff d_lowd_high
: Autoconversion coeff d_highk
: Threshold for smooth tranistion steepness
CloudMicrophysics.Parameters.AccrB1994
— TypeAccrB1994
Beheng (1994) accretion parameters
Fields
A
: Accretion coefficient A
CloudMicrophysics.Parameters.TC1980
— TypeTC1980
The type and parameters for 2-moment precipitation formation by Tripoli and Cotton (1980)
DOI: 10.1175/1520-0450(1980)019<1037:ANIOSF>2.0.CO;2
Fields
acnv
: Autoconversion parametersaccr
: Accretion parameters
CloudMicrophysics.Parameters.AcnvTC1980
— TypeAcnvTC1980
Tripoli and Cotton (1980) autoconversion parameters
Fields
a
: Autoconversion coefficient ab
: Autoconversion coefficient bD
: Autoconversion coefficient Dr_0
: Autoconversion coefficient r_0me_liq
: Autoconversion coefficient me_liqm0_liq_coeff
: Autoconversion coefficient m0liqcoeffk
: Threshold for smooth tranistion steepness
CloudMicrophysics.Parameters.AccrTC1980
— TypeAccrTC1980
Tripoli and Cotton (1980) accretion parameters
Fields
A
: Accretion coefficient A
CloudMicrophysics.Parameters.LD2004
— TypeLD2004
The type and parameters for 2-moment precipitation formation by Liu and Daum (2004)
DOI: 10.1175/1520-0469(2004)061<1539:POTAPI>2.0.CO;2
Fields
R_6C_0
: Autoconversion coefficient R6C0E_0
: Autoconversion coefficient E_0ρ_w
: liquid water density [kg/m3]k
: Threshold for smooth tranistion steepness
CloudMicrophysics.Parameters.VarTimescaleAcnv
— TypeVarTimescaleAcnv
The type for 2-moment precipitation formation based on the 1-moment parameterization with variable time scale Azimi et al (2023)
Fields
τ
: Timescale [s]α
: Powerlaw coefficient [-]
CloudMicrophysics.Parameters.SB2006
— TypeSB2006
The type and parameters for 2-moment precipitation formation by Seifert and Beheng (2006). The pdf_r type choses between running with or without limiters on raindrop size distribution parameters
DOI: 10.1007/s00703-005-0112-4
Fields
pdf_c
: Cloud particle size distribution parameterspdf_r
: Rain particle size distribution parametersacnv
: Autoconversion parametersaccr
: Accretion parametersself
: Rain selfcollection parametersbrek
: Rain breakup parametersevap
: Rain evaporation parametersnumadj
: Number concentration adjustment parameter
CloudMicrophysics.Parameters.RainParticlePDF_SB2006
— TypeRainParticlePDF_SB2006
Abstract type for the size distribution parameters of rain particles
See RainParticlePDF_SB2006_limited
and RainParticlePDF_SB2006_notlimited
for the concrete types.
CloudMicrophysics.Parameters.RainParticlePDF_SB2006_limited
— TypeRainParticlePDF_SB2006_limited
Rain size distribution parameters from SB2006 including the limiters on drop maximum mass and the size distribution coefficinets N0 and lambda
Fields
νr
: Raindrop size distribution coefficient νrμr
: Raindrop size distribution coefficient μrxr_min
: Raindrop minimal massxr_max
: Raindrop maximum massN0_min
: Raindrop size distribution coefficient N0 minN0_max
: Raindrop size distribution coefficient N0 maxλ_min
: Raindrop size distribution coefficient lambda minλ_max
: Raindrop size distribution coefficient lambda maxρw
: Cloud liquid water density [kg/m3]ρ0
: Reference air density [kg/m3]
CloudMicrophysics.Parameters.RainParticlePDF_SB2006_notlimited
— TypeRainParticlePDF_SB2006
Rain size distribution parameters from SB2006 but without the limiters
Fields
νr
: Raindrop size distribution coefficient νrμr
: Raindrop size distribution coefficient μrxr_min
: Raindrop minimum massxr_max
: Raindrop maximum massρw
: Cloud liquid water density [kg/m3]ρ0
: Reference air density [kg/m3]
CloudMicrophysics.Parameters.CloudParticlePDF_SB2006
— TypeCloudParticlePDF_SB2006
Cloud droplets size distribution parameters from SB2006
Fields
νc
: Cloud droplet size distribution coefficient νcμc
: Cloud droplet size distribution coefficient μcxc_min
: Cloud droplets minimum massxc_max
: Cloud droplets maximum massρw
: Cloud liquid water density [kg/m3]
CloudMicrophysics.Parameters.AcnvSB2006
— TypeAcnvSB2006
Autoconversion parameters from SB2006
Fields
kcc
: Collection kernel coefficientx_star
: Minimum mass of rain dropletsρ0
: Reference air density [kg/m3]A
: Autoconversion correcting function coeff Aa
: Autoconversion correcting function coeff ab
: Autoconversion correcting function coeff b
CloudMicrophysics.Parameters.AccrSB2006
— TypeAccrSB2006
Accretion parameters from SB2006
Fields
kcr
: Collection kernel coefficient Kcrτ0
: Accretion correcting function coefficient τ_0ρ0
: Reference air density [kg/m3]c
: Accretion correcting function coefficient c
CloudMicrophysics.Parameters.SelfColSB2006
— TypeSelfColSB2006
Rain selfcollection parameters from SB2006
Fields
krr
: Collection kernel coefficient krrκrr
: Collection kernel coefficient kappa rrd
: Raindrop self collection coefficient d
CloudMicrophysics.Parameters.BreakupSB2006
— TypeBreakupSB2006
Rain breakup parameters from SB2006
Fields
Deq
: Raindrop equilibrium mean diamaterDr_th
: Raindrop breakup mean diamater thresholdkbr
: Raindrops breakup coefficient kbrκbr
: Raindrops breakup coefficient kappa br
CloudMicrophysics.Parameters.EvaporationSB2006
— TypeEvaporationSB2006
Rain evaporation parameters from SB2006
Fields
av
: Ventillation coefficient abv
: Ventillation coefficient bα
: Rain evapoartion coefficient αβ
: Rain evapoartion coefficient βρ0
: Reference air density [kg/m3]
CloudMicrophysics.Parameters.NumberAdjustmentHorn2012
— TypeNumberAdjustmentHorn2012
Number concentration adjustment parameter from Horn (2012, DOI: 10.5194/gmd-5-345-2012)
Fields
τ
: Number concentration adjustment timescale [s]
CloudMicrophysics.Parameters.Blk1MVelType
— TypeBlk1MVelType
The type for precipitation terminal velocity from the simple 1-moment scheme (defined for rain and snow)
Fields
rain
snow
CloudMicrophysics.Parameters.Blk1MVelTypeRain
— TypeBlk1MVelTypeRain
The type for precipitation terminal velocity from the simple 1-moment scheme for rain
Fields
r0
: particle length scale [m]ve
: rain terminal velocity size relation coefficient [-]Δv
: rain terminal velocity size relation coefficient [-]χv
: rain terminal velocity size relation coefficient [-]ρw
: cloud water density [kg/m3]C_drag
: rain drop drag coefficient [-]grav
: gravitational acceleration [m/s2]
CloudMicrophysics.Parameters.Blk1MVelTypeSnow
— TypeBlk1MVelTypeSnow
The type for precipitation terminal velocity from the simple 1-moment scheme for snow
Fields
r0
: particle length scale [m]ve
: snow terminal velocity size relation coefficient [-]Δv
: snow terminal velocity size relation coefficient [-]χv
: snow terminal velocity size relation coefficient [-]v0
: snow terminal velocity size relation coefficient [m/s]
CloudMicrophysics.Parameters.SB2006VelType
— TypeSB2006VelType
The type for precipitation terminal velocity from Seifert and Beheng 2006 (Defined only for rain)
Fields
ρ0
aR
bR
cR
CloudMicrophysics.Parameters.Chen2022VelType
— TypeChen2022VelType
The type for precipitation terminal velocity from Chen et al 2022 DOI: 10.1016/j.atmosres.2022.106171 (defined for rain, snow and cloud ice)
Fields
rain
small_ice
large_ice
CloudMicrophysics.Parameters.Chen2022VelTypeSmallIce
— TypeChen2022VelTypeSmallIce
The type for precipitation terminal velocity from Chen et al 2022 for small ice. See Table B3 for parameter definitions. DOI: 10.1016/j.atmosres.2022.106171
Fields
A
B
C
E
F
G
cutoff
: cutoff for small vs large ice particle dimension [m]
CloudMicrophysics.Parameters.Chen2022VelTypeLargeIce
— TypeChen2022VelTypeLargeIce
The type for precipitation terminal velocity from Chen et al 2022 for large ice. See Table B4 for parameter definitions. DOI: 10.1016/j.atmosres.2022.106171
Fields
A
B
C
E
F
G
H
cutoff
: cutoff for small vs large ice particle dimension [m]
CloudMicrophysics.Parameters.Chen2022VelTypeRain
— TypeChen2022VelTypeRain
The type for precipitation terminal velocity from Chen et al 2022 for rain. See Table B1 for parameter definitions. DOI: 10.1016/j.atmosres.2022.106171
Fields
ρ0
a
a3_pow
b
b_ρ
c
Precipitation susceptibility
CloudMicrophysics.PrecipitationSusceptibility.precipitation_susceptibility_autoconversion
— Functionprecipitation_susceptibility_autoconversion(param_set, scheme, q_liq, q_rai, ρ, N_liq)
scheme
- type for 2-moment rain autoconversion parameterizationq_liq
- cloud water specific contentq_rai
- rain water specific contentρ
- air densityN_liq
- cloud droplet number density
Returns the precipitation susceptibility rates due to autoconversion as a precip_susceptibility_rates
object, using automatic differentiation. Works for any 2-moment scheme, as long as autoconversion is defined for it.
CloudMicrophysics.PrecipitationSusceptibility.precipitation_susceptibility_accretion
— Functionprecipitation_susceptibility_accretion(param_set, scheme, q_liq, q_rai, ρ, N_liq)
scheme
- type for 2-moment rain autoconversion parameterizationq_liq
- cloud water specific contentq_rai
- rain water specific contentρ
- air densityN_liq
- cloud droplet number density
Returns the precipitation susceptibility rates due to accretion as a precip_susceptibility_rates
object, using automatic differentiation. Works for any 2-moment scheme, as long as accretion is defined for it.