Thermodynamics

ClimateMachine.ThermodynamicsModule
Thermodynamics

Moist thermodynamic functions, e.g., for air pressure (atmosphere equation of state), latent heats of phase transitions, saturation vapor pressures, and saturation specific humidities.

AbstractParameterSet's

Many functions defined in this module rely on CLIMAParameters.jl. CLIMAParameters.jl defines several functions (e.g., many planet parameters). For example, to compute the mole-mass ratio:

using CLIMAParameters.Planet: molmass_ratio
using CLIMAParameters: AbstractEarthParameterSet
struct EarthParameterSet <: AbstractEarthParameterSet end
param_set = EarthParameterSet()
_molmass_ratio = molmass_ratio(param_set)

Because these parameters are widely used throughout this module, param_set is an argument for many Thermodynamics functions.

source

Thermodynamic State Constructors

ClimateMachine.Thermodynamics.PhasePartitionType
PhasePartition

Represents the mass fractions of the moist air mixture.

Constructors

PhasePartition(q_tot::Real[, q_liq::Real[, q_ice::Real]])
PhasePartition(ts::ThermodynamicState)

See also PhasePartition_equil

Fields

  • tot

    total specific humidity

  • liq

    liquid water specific humidity (default: 0)

  • ice

    ice specific humidity (default: 0)

source
ClimateMachine.Thermodynamics.ThermodynamicStateType
ThermodynamicState{FT}

A thermodynamic state, which can be initialized for various thermodynamic formulations (via its sub-types). All ThermodynamicState's have access to functions to compute all other thermodynamic properties.

source
ClimateMachine.Thermodynamics.PhaseDryType
PhaseDry{FT} <: ThermodynamicState

A dry thermodynamic state (q_tot = 0).

Constructors

PhaseDry(param_set, e_int, ρ)

Fields

  • param_set

    parameter set, used to dispatch planet parameter function calls

  • e_int

    internal energy

  • ρ

    density of dry air

source
ClimateMachine.Thermodynamics.PhaseEquilType
PhaseEquil{FT} <: ThermodynamicState

A thermodynamic state assuming thermodynamic equilibrium (therefore, saturation adjustment may be needed).

Constructors

PhaseEquil(param_set, e_int, ρ, q_tot)

Fields

  • param_set

    parameter set, used to dispatch planet parameter function calls

  • e_int

    internal energy

  • ρ

    density of air (potentially moist)

  • q_tot

    total specific humidity

  • T

    temperature: computed via saturation_adjustment

source
ClimateMachine.Thermodynamics.PhaseNonEquilType

PhaseNonEquil{FT} <: ThermodynamicState

A thermodynamic state assuming thermodynamic non-equilibrium (therefore, temperature can be computed directly).

Constructors

PhaseNonEquil(param_set, e_int, q::PhasePartition, ρ)

Fields

  • param_set

    parameter set, used to dispatch planet parameter function calls

  • e_int

    internal energy

  • ρ

    density of air (potentially moist)

  • q

    phase partition

source
ClimateMachine.Thermodynamics.LiquidIcePotTempSHumEquilFunction
LiquidIcePotTempSHumEquil(param_set, θ_liq_ice, ρ, q_tot)

Constructs a PhaseEquil thermodynamic state from:

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • θ_liq_ice liquid-ice potential temperature
  • ρ (moist-)air density
  • q_tot total specific humidity
  • temperature_tol temperature tolerance for saturation adjustment
  • maxiter maximum iterations for saturation adjustment
source
ClimateMachine.Thermodynamics.LiquidIcePotTempSHumNonEquilFunction
LiquidIcePotTempSHumNonEquil(param_set, θ_liq_ice, ρ, q_pt)

Constructs a PhaseNonEquil thermodynamic state from:

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • θ_liq_ice liquid-ice potential temperature
  • ρ (moist-)air density
  • q_pt phase partition

and, optionally

  • potential_temperature_tol potential temperature for non-linear equation solve
  • maxiter maximum iterations for non-linear equation solve
source

Thermodynamic state methods

ClimateMachine.Thermodynamics.air_densityFunction
air_density(param_set, T, p[, q::PhasePartition])

The (moist-)air density from the equation of state (ideal gas law) where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T air temperature
  • p pressure

and, optionally,

source
air_density(ts::ThermodynamicState)

The (moist-)air density, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.air_pressureFunction
air_pressure(param_set, T, ρ[, q::PhasePartition])

The air pressure from the equation of state (ideal gas law) where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T air temperature
  • ρ (moist-)air density

and, optionally,

source
air_pressure(ts::ThermodynamicState)

The air pressure from the equation of state (ideal gas law), given a thermodynamic state ts.

source
air_pressure(param_set, T::FT, T∞::FT, p∞::FT, ::DryAdiabaticProcess)

The air pressure for an isentropic process, where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • T∞ ambient temperature
  • p∞ ambient pressure
source
ClimateMachine.Thermodynamics.air_temperatureFunction
air_temperature(param_set, e_int, q::PhasePartition)

The air temperature, where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • e_int internal energy per unit mass

and, optionally,

source
air_temperature(ts::ThermodynamicState)

The air temperature, given a thermodynamic state ts.

source
air_temperature(param_set, p::FT, θ::FT, Φ::FT, ::DryAdiabaticProcess)

The air temperature for an isentropic process, where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • p pressure
  • θ potential temperature
source
ClimateMachine.Thermodynamics.air_temperature_from_liquid_ice_pottemp_non_linearFunction
air_temperature_from_liquid_ice_pottemp_non_linear(param_set, θ_liq_ice, ρ, q::PhasePartition)

Computes temperature T given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • θ_liq_ice liquid-ice potential temperature
  • ρ (moist-)air density
  • tol absolute tolerance for non-linear equation iterations. Can be one of:
    • SolutionTolerance() to stop when |x_2 - x_1| < tol
    • ResidualTolerance() to stop when |f(x)| < tol
  • maxiter maximum iterations for non-linear equation solve

and, optionally,

by finding the root of T - air_temperature_from_liquid_ice_pottemp_given_pressure(param_set, θ_liq_ice, air_pressure(param_set, T, ρ, q), q) = 0

source
ClimateMachine.Thermodynamics.cp_mFunction
cp_m(param_set, [q::PhasePartition])

The isobaric specific heat capacity of moist air given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details

and, optionally

source
cp_m(ts::ThermodynamicState)

The isobaric specific heat capacity of moist air, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.cv_mFunction
cv_m(param_set, [q::PhasePartition])

The isochoric specific heat capacity of moist air given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details

and, optionally

source
cv_m(ts::ThermodynamicState)

The isochoric specific heat capacity of moist air, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.dry_pottempFunction
dry_pottemp(param_set, T, ρ[, q::PhasePartition])

The dry potential temperature where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density

and, optionally,

source
dry_pottemp(ts::ThermodynamicState)

The dry potential temperature, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.exnerFunction
exner(param_set, T, ρ[, q::PhasePartition)])

The Exner function where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density

and, optionally,

source
exner(ts::ThermodynamicState)

The Exner function, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.gas_constantsFunction
(R_m, cp_m, cv_m, γ_m) = gas_constants(param_set, [q::PhasePartition])

Wrapper to compute all gas constants at once, where optionally,

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • q PhasePartition. Without this argument, the results are for dry air.

The function returns a tuple of

source
(R_m, cp_m, cv_m, γ_m) = gas_constants(ts::ThermodynamicState)

Wrapper to compute all gas constants at once, given a thermodynamic state ts.

The function returns a tuple of

source
ClimateMachine.Thermodynamics.internal_energyFunction
internal_energy(param_set, T[, q::PhasePartition])

The internal energy per unit mass, given a thermodynamic state ts or

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature

and, optionally,

source
internal_energy(ts::ThermodynamicState)

The internal energy per unit mass, given a thermodynamic state ts.

source
internal_energy(ρ::FT, ρe::FT, ρu::AbstractVector{FT}, e_pot::FT)

The internal energy per unit mass, given

  • ρ (moist-)air density
  • ρe total energy per unit volume
  • ρu momentum vector
  • e_pot potential energy (e.g., gravitational) per unit mass
source
ClimateMachine.Thermodynamics.internal_energy_satFunction
internal_energy_sat(param_set, T, ρ, q_tot, phase_type)

The internal energy per unit mass in thermodynamic equilibrium at saturation where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density
  • q_tot total specific humidity
  • phase_type a thermodynamic state type
source
internal_energy_sat(ts::ThermodynamicState)

The internal energy per unit mass in thermodynamic equilibrium at saturation, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.latent_heat_sublimFunction
latent_heat_sublim(param_set, T::FT) where {FT<:Real}

The specific latent heat of sublimation where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
source
latent_heat_sublim(ts::ThermodynamicState)

The specific latent heat of sublimation given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.latent_heat_vaporFunction
latent_heat_vapor(param_set, T::FT) where {FT<:Real}

The specific latent heat of vaporization where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
source
latent_heat_vapor(ts::ThermodynamicState)

The specific latent heat of vaporization given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.liquid_fractionFunction
liquid_fraction(param_set, T, phase_type[, q])

The fraction of condensate that is liquid where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • phase_type a thermodynamic state type

PhaseNonEquil behavior

If q.liq or q.ice are nonzero, the liquid fraction is computed from them.

ThermodynamicState

Otherwise, phase equilibrium is assumed so that the fraction of liquid is a function that is 1 above T_freeze and goes to zero below T_freeze.

source
liquid_fraction(ts::ThermodynamicState)

The fraction of condensate that is liquid given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.liquid_ice_pottempFunction
liquid_ice_pottemp(param_set, T, ρ, q::PhasePartition)

The liquid-ice potential temperature where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density

and, optionally,

source
liquid_ice_pottemp(ts::ThermodynamicState)

The liquid-ice potential temperature, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.liquid_ice_pottemp_satFunction
liquid_ice_pottemp_sat(param_set, T, ρ, phase_type[, q::PhasePartition])

The saturated liquid ice potential temperature where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density
  • phase_type a thermodynamic state type

and, optionally,

source
liquid_ice_pottemp_sat(param_set, T, ρ, phase_type, q_tot)

The saturated liquid ice potential temperature where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density
  • phase_type a thermodynamic state type
  • q_tot total specific humidity
source
liquid_ice_pottemp_sat(ts::ThermodynamicState)

The liquid potential temperature given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.q_vap_saturationFunction
q_vap_saturation(param_set, T, ρ, phase_type[, q::PhasePartition])

Compute the saturation specific humidity, given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density
  • phase_type a thermodynamic state type

and, optionally,

If the PhasePartition q is given, the saturation specific humidity is that of a mixture of liquid and ice, computed in a thermodynamically consistent way from the weighted sum of the latent heats of the respective phase transitions (Pressel et al., JAMES, 2015). That is, the saturation vapor pressure and from it the saturation specific humidity are computed from a weighted mean of the latent heats of vaporization and sublimation, with the weights given by the fractions of condensates q.liq/(q.liq + q.ice) and q.ice/(q.liq + q.ice) that are liquid and ice, respectively.

If the PhasePartition q is not given, or has zero liquid and ice specific humidities, the saturation specific humidity is that over a mixture of liquid and ice, with the fraction of liquid given by temperature dependent liquid_fraction(param_set, T, phase_type) and the fraction of ice by the complement 1 - liquid_fraction(param_set, T, phase_type).

source
q_vap_saturation(ts::ThermodynamicState)

Compute the saturation specific humidity, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.q_vap_saturation_genericFunction
q_vap_saturation_generic(param_set, T, ρ[, phase=Liquid()])

Compute the saturation specific humidity over a plane surface of condensate, given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density

and, optionally,

  • Liquid() indicating condensate is liquid
  • Ice() indicating condensate is ice
source
ClimateMachine.Thermodynamics.relative_humidityFunction
relative_humidity(param_set, T, p, phase_type, q::PhasePartition)

The relative humidity, given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • p pressure
  • phase_type a thermodynamic state type

and, optionally,

source
relative_humidity(ts::ThermodynamicState)

The relative humidity, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.saturation_adjustmentFunction
saturation_adjustment(
    param_set,
    e_int,
    ρ,
    q_tot,
    phase_type,
    maxiter,
    temperature_tol
)

Compute the temperature that is consistent with

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • e_int internal energy
  • ρ (moist-)air density
  • q_tot total specific humidity
  • phase_type a thermodynamic state type
  • maxiter maximum iterations for non-linear equation solve
  • temperature_tol temperature tolerance

by finding the root of

e_int - internal_energy_sat(param_set, T, ρ, q_tot, phase_type) = 0

using Newtons method with analytic gradients.

See also saturation_adjustment.

source
ClimateMachine.Thermodynamics.saturation_excessFunction
saturation_excess(param_set, T, ρ, phase_type, q::PhasePartition)

The saturation excess in equilibrium where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density
  • phase_type a thermodynamic state type
  • q PhasePartition

The saturation excess is the difference between the total specific humidity q.tot and the saturation specific humidity in equilibrium, and it is defined to be nonzero only if this difference is positive.

source
saturation_excess(ts::ThermodynamicState)

Compute the saturation excess in equilibrium, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.saturation_vapor_pressureFunction
saturation_vapor_pressure(param_set, T, Liquid())

Return the saturation vapor pressure over a plane liquid surface given

  • T temperature

  • param_set an AbstractParameterSet, see the Thermodynamics for more details

    saturation_vapor_pressure(param_set, T, Ice())

Return the saturation vapor pressure over a plane ice surface given

  • T temperature

  • param_set an AbstractParameterSet, see the Thermodynamics for more details

    saturation_vapor_pressure(param_set, T, LH_0, Δcp)

Compute the saturation vapor pressure over a plane surface by integration of the Clausius-Clapeyron relation.

The Clausius-Clapeyron relation

`dlog(p_v_sat)/dT = [LH_0 + Δcp * (T-T_0)]/(R_v*T^2)`

is integrated from the triple point temperature T_triple, using Kirchhoff's relation

`L = LH_0 + Δcp * (T - T_0)`

for the specific latent heat L with constant isobaric specific heats of the phases. The linear dependence of the specific latent heat on temperature T allows analytic integration of the Clausius-Clapeyron relation to obtain the saturation vapor pressure p_v_sat as a function of the triple point pressure press_triple.

source
ClimateMachine.Thermodynamics.soundspeed_airFunction
soundspeed_air(param_set, T[, q::PhasePartition])

The speed of sound in unstratified air, where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature

and, optionally,

source
soundspeed_air(ts::ThermodynamicState)

The speed of sound in unstratified air given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.supersaturationFunction
supersaturation(param_set, q, ρ, T, Liquid())
supersaturation(param_set, q, ρ, T, Ice())
  • param_set - abstract set with earth parameters
  • q - phase partition
  • ρ - air density,
  • T - air temperature
  • Liquid(), Ice() - liquid or ice phase to dispatch over.

Returns supersaturation (qv/qv_sat -1) over water or ice.

source
ClimateMachine.Thermodynamics.total_energyFunction
total_energy(param_set, e_kin, e_pot, T[, q::PhasePartition])

The total energy per unit mass, given

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • e_kin kinetic energy per unit mass
  • e_pot potential energy per unit mass
  • T temperature

and, optionally,

source
total_energy(e_kin, e_pot, ts::ThermodynamicState)

The total energy per unit mass given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.total_specific_humidityFunction
total_specific_humidity(ts::ThermodynamicState)
total_specific_humidity(param_set, T, p, relative_humidity)

Total specific humidity given

  • ts a thermodynamic state

or

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • p pressure
  • relative_humidity relative humidity (can exceed 1 when there is super saturation/condensate)
source
ClimateMachine.Thermodynamics.virtual_pottempFunction
virtual_pottemp(param_set, T, ρ[, q::PhasePartition])

The virtual potential temperature where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density

and, optionally,

source
virtual_pottemp(ts::ThermodynamicState)

The virtual potential temperature, given a thermodynamic state ts.

source
ClimateMachine.Thermodynamics.virtual_temperatureFunction
virtual_temperature(param_set, T, ρ[, q::PhasePartition])

The virtual temperature where

  • param_set an AbstractParameterSet, see the Thermodynamics for more details
  • T temperature
  • ρ (moist-)air density

and, optionally,

source
virtual_temperature(ts::ThermodynamicState)

The virtual temperature, given a thermodynamic state ts.

source

Dispatch types