diff --git a/docs/src/microphysics/saturation_adjustment.md b/docs/src/microphysics/saturation_adjustment.md index a221f520..1aef4e7a 100644 --- a/docs/src/microphysics/saturation_adjustment.md +++ b/docs/src/microphysics/saturation_adjustment.md @@ -7,19 +7,19 @@ Mixed-phase saturation adjustment is described by [Pressel2015](@citet). ## Moist static energy and total moisture mass fraction The saturation adjustment solver (specific to our anelastic formulation) takes four inputs: - * moist static energy ``e`` - * total moisture mass fraction ``qᵗ`` - * height ``z`` - * reference pressure ``pᵣ`` +* moist static energy ``e``, +* total moisture mass fraction ``qᵗ``, +* height ``z``, and +* reference pressure ``pᵣ``. Note that moist static energy density ``ρᵣ e`` and moisture density ``ρᵣ qᵗ`` -are prognostic variables for `Breeze.AtmosphereModel` when using `AnelasticFormulation`, +are prognostic variables for [`AtmosphereModel`](@ref) when using [`AnelasticFormulation`](@ref), where ``ρᵣ`` is the reference density. With warm-phase microphysics, the moist static energy ``e`` is related to temperature ``T``, height ``z``, and liquid mass fraction ``qˡ`` by ```math -e \equiv cᵖᵐ \, T + g z - ℒˡᵣ qˡ , +e ≡ cᵖᵐ \, T + g z - ℒˡᵣ qˡ , ``` where ``cᵖᵐ`` is the mixture heat capacity, ``g`` is gravitational acceleration, @@ -56,7 +56,7 @@ where ``qᵈ = 1 - qᵗ`` is the dry air mass fraction, ``qᵛ`` is the specific ``Rᵈ`` is the dry air gas constant, and ``Rᵛ`` is the vapor gas constant. The density ``ρ`` is expressed in terms of ``pᵣ`` under the anelastic approximation. -In saturated conditions, we have ``qᵛ ≡ qᵛ⁺`` by definition, which leads to the expression +In saturated conditions, we have ``qᵛ ≡ qᵛ⁺`` by definition, which leads to the expression ```math qᵛ⁺ = \frac{ρᵛ⁺}{ρ} = \frac{Rᵐ}{Rᵛ} \frac{pᵛ⁺}{pᵣ} = \frac{Rᵈ}{Rᵛ} \left ( 1 - qᵗ \right ) \frac{pᵛ⁺}{pᵣ} + qᵛ⁺ \frac{pᵛ⁺}{pᵣ} . @@ -69,14 +69,14 @@ _valid only in saturated conditions and under the assumptions of saturation adju qᵛ⁺ = \frac{Rᵈ}{Rᵛ} \left ( 1 - qᵗ \right ) \frac{pᵛ⁺}{pᵣ - pᵛ⁺} . ``` -This expression can also be found in [Pressel2015](@citet), equation 37. +This expression can also be found in paper by [Pressel2015](@citet), equation (37). ## The saturation adjustment algorithm We compute the saturation adjustment temperature by solving the nonlinear algebraic equation ```math -0 = r(T) \equiv T - \frac{1}{cᵖᵐ} \left [ e - g z + ℒˡᵣ \max(0, qᵗ - qᵛ⁺) \right ] \, +0 = r(T) ≡ T - \frac{1}{cᵖᵐ} \left [ e - g z + ℒˡᵣ \max(0, qᵗ - qᵛ⁺) \right ] \, ``` where ``r`` is the "residual", using a secant method. @@ -117,7 +117,7 @@ In equilibrium (and thus under the assumptions of saturation adjustment), the sp ``qᵛ = qᵛ⁺``, while the liquid mass fraction is ```@example microphysics -qˡ = qᵗ - qᵛ⁺ +qˡ = qᵗ - qᵛ⁺ ``` It is small but greater than zero → the typical situation in clouds on Earth. @@ -134,9 +134,9 @@ z = 0.0 e = cᵖᵐ * T + g * z - ℒˡᵣ * qˡ ``` -Moist static energy has units ``\mathrm{m^2 / s^2}``, or ``\mathrm{J}{kg}``. +Moist static energy has units ``\mathrm{m^2 / s^2}``, or ``\mathrm{J} / \mathrm{kg}``. Next we show that the saturation adjustment solver recovers the input temperature -by passing it an "unadjusted" moisture mass fraction into `compute_temperature`, +by passing it an "unadjusted" moisture mass fraction into [`Breeze.AtmosphereModels.compute_temperature`](@ref), ```@example microphysics using Breeze.AtmosphereModels: compute_temperature diff --git a/docs/src/thermodynamics.md b/docs/src/thermodynamics.md index 05e9e9d6..5bcd1d2f 100644 --- a/docs/src/thermodynamics.md +++ b/docs/src/thermodynamics.md @@ -83,7 +83,7 @@ By definition, all of the mass fractions sum up to unity, so that, using ``qᵗ = qᵛ + qˡ + qⁱ``, the dry air mass fraction can be diagnosed with ``qᵈ = 1 - qᵗ``. The sometimes tedious bookkeeping required to correctly diagnose the effective mixture properties -of moist air are facilitated by Breeze's handy `MoistureMassFractions` abstraction. +of moist air are facilitated by Breeze's handy [`MoistureMassFractions`](@ref Breeze.Thermodynamics.MoistureMassFractions) abstraction. For example, ```@example thermo @@ -328,7 +328,7 @@ thermo = ThermodynamicConstants() thermo.molar_gas_constant ``` -`ThermodynamicConstants`, which is central to Breeze's implementation of moist thermodynamics. +[`ThermodynamicConstants`](@ref), which is central to Breeze's implementation of moist thermodynamics. holds constants like the molar gas constant and molar masses, latent heats, gravitational acceleration, and more, ```@example thermo @@ -480,7 +480,7 @@ cˡ = thermo.liquid.heat_capacity Δcˡ = cᵖᵛ - cˡ ``` -This difference ``\Delta c^l ≈`` $(round(1885 - 4181, digits=1)) J/(kg⋅K) is negative because +This difference ``\Delta c^l ≈`` -2296 J/(kg⋅K) is negative because water vapor has a lower heat capacity than liquid water. ### Mixed-phase saturation vapor pressure @@ -535,7 +535,7 @@ pᵛᵐ⁺ = [saturation_vapor_pressure(Tⁱ, thermo, mixed_surface) for Tⁱ in using CairoMakie fig = Figure() -ax = Axis(fig[1, 1], xlabel="Temperature (ᵒK)", ylabel="Saturation vapor pressure pᵛ⁺ (Pa)", +ax = Axis(fig[1, 1], xlabel="Temperature (ᵒK)", ylabel="Saturation vapor pressure pᵛ⁺ (Pa)", yscale = log10, xticks=200:20:320) lines!(ax, T, pᵛˡ⁺, label="liquid", linewidth=2) lines!(ax, T, pᵛⁱ⁺, label="ice", linestyle=:dash, linewidth=2) @@ -588,7 +588,7 @@ using CairoMakie fig = Figure(size=(1000, 400)) # Panel 1: Saturation vapor pressure -ax1 = Axis(fig[1, 1], xlabel="Temperature (K)", ylabel="Saturation vapor pressure (Pa)", +ax1 = Axis(fig[1, 1], xlabel="Temperature (K)", ylabel="Saturation vapor pressure (Pa)", yscale=log10, title="Saturation vapor pressure") for (i, λ) in enumerate(λ_values) @@ -600,7 +600,7 @@ end axislegend(ax1, position=:lt) # Panel 2: Saturation specific humidity -ax2 = Axis(fig[1, 2], xlabel="Temperature (K)", ylabel="Saturation specific humidity (kg/kg)", +ax2 = Axis(fig[1, 2], xlabel="Temperature (K)", ylabel="Saturation specific humidity (kg/kg)", title="Saturation specific humidity") for (i, λ) in enumerate(λ_values) diff --git a/src/AtmosphereModels/anelastic_formulation.jl b/src/AtmosphereModels/anelastic_formulation.jl index 7d9cd2bc..81235bd2 100644 --- a/src/AtmosphereModels/anelastic_formulation.jl +++ b/src/AtmosphereModels/anelastic_formulation.jl @@ -23,6 +23,14 @@ import Oceananigans.TimeSteppers: compute_pressure_correction!, make_pressure_co ##### Formulation definition ##### +""" +$(TYPEDSIGNATURES) + +AnelasticFormulation is a dynamical formulation wherein the density and pressure are +small perturbations from a dry, hydrostatic, adiabatic `reference_state`. +The prognostic energy variable is the moist static energy density. +The energy density equation includes a buoyancy flux term, following [Pauluis2008](@citet). +""" struct AnelasticFormulation{R} reference_state :: R end @@ -41,16 +49,6 @@ end Base.show(io::IO, formulation::AnelasticFormulation) = print(io, "AnelasticFormulation") -function AnelasticFormulation(grid, state_constants, thermo) - pᵣ = Field{Nothing, Nothing, Center}(grid) - ρᵣ = Field{Nothing, Nothing, Center}(grid) - set!(pᵣ, z -> reference_pressure(z, state_constants, thermo)) - set!(ρᵣ, z -> reference_density(z, state_constants, thermo)) - fill_halo_regions!(pᵣ) - fill_halo_regions!(ρᵣ) - return AnelasticFormulation(state_constants, pᵣ, ρᵣ) -end - ##### ##### Thermodynamic state ##### diff --git a/src/AtmosphereModels/microphysics_interface.jl b/src/AtmosphereModels/microphysics_interface.jl index d38c1e72..6bcdd218 100644 --- a/src/AtmosphereModels/microphysics_interface.jl +++ b/src/AtmosphereModels/microphysics_interface.jl @@ -35,7 +35,7 @@ end """ $(TYPEDSIGNATURES) -Build and return `MoistureMassFractions` at `(i, j, k)` for the given `grid`, +Build and return [`MoistureMassFractions`](@ref) at `(i, j, k)` for the given `grid`, `microphysics`, `microphysical_fields`, and total moisture mass fraction `qᵗ`. Dispatch is provided for `::Nothing` microphysics here. Specific microphysics diff --git a/src/MoistAirBuoyancies.jl b/src/MoistAirBuoyancies.jl index 6d5712da..83857405 100644 --- a/src/MoistAirBuoyancies.jl +++ b/src/MoistAirBuoyancies.jl @@ -51,7 +51,7 @@ Adapt.adapt_structure(to, mb::MoistAirBuoyancy) = """ $(TYPEDSIGNATURES) -Return a MoistAirBuoyancy formulation that can be provided as input to an +Return a `MoistAirBuoyancy` formulation that can be provided as input to an `Oceananigans.NonhydrostaticModel`. !!! note "Required tracers" @@ -72,12 +72,12 @@ MoistAirBuoyancy: └── thermodynamics: ThermodynamicConstants{Float64} ``` -To build a model with MoistAirBuoyancy, we include potential temperature and total specific humidity +To build a model with `MoistAirBuoyancy`, we include potential temperature and total specific humidity tracers `θ` and `qᵗ` to the model. ```jldoctest mab model = NonhydrostaticModel(; grid, buoyancy, tracers = (:θ, :qᵗ)) - + # output NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0) ├── grid: 1×1×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 1×1×3 halo @@ -97,7 +97,7 @@ function MoistAirBuoyancy(grid; reference_state = ReferenceState(grid, thermodynamics; base_pressure, potential_temperature = reference_potential_temperature) - + return MoistAirBuoyancy(reference_state, thermodynamics) end @@ -168,9 +168,9 @@ The saturation equilibrium temperature satisfies the nonlinear relation with ``ℒˡᵣ`` the latent heat at the reference temperature ``Tᵣ``, ``cᵖᵐ`` the mixture specific heat, ``Π`` the Exner function, ``qˡ = \\max(0, qᵗ - qᵛ⁺)`` the condensate specific humidity, ``qᵗ`` is the -total specific humidity, ``qᵛ⁺`` is the saturation specific humidity. +total specific humidity, and ``qᵛ⁺`` is the saturation specific humidity. -The saturation equilibrium temperature is thus obtained by solving ``r(T)``, where +The saturation equilibrium temperature is thus obtained by solving ``r(T) = 0``, where ```math r(T) ≡ T - θ Π - ℒˡᵣ qˡ / cᵖᵐ . ``` @@ -211,7 +211,7 @@ Solution of ``r(T) = 0`` is found via the [secant method](https://en.wikipedia.o # since ``qᵛ⁺₁(T₁)`` underestimates the saturation specific humidity, # and therefore qˡ₁ is overestimated. This is similar to an approach # used in Pressel et al 2015. However, it doesn't work for large liquid fractions. - T₂ = T₁ + 1 + T₂ = T₁ + 1 #= ℒˡᵣ = thermo.liquid.reference_latent_heat @@ -248,7 +248,7 @@ end # This estimate assumes that the specific humidity is itself the saturation # specific humidity, eg ``qᵛ = qᵛ⁺``. Knowledge of the specific humidity -# is needed to compute the mixture gas constant, and thus density, +# is needed to compute the mixture gas constant, and thus density, # which in turn is needed to compute the _saturation_ specific humidity. # This consideration culminates in a new expression for saturation specific humidity # used below, and also written in Pressel et al 2015, equation 37. @@ -279,7 +279,7 @@ end cᵖᵐ = mixture_heat_capacity(q, thermo) qˡ = q.liquid θ = 𝒰.potential_temperature - return T - Π * θ - ℒˡᵣ * qˡ / cᵖᵐ + return T - Π * θ - ℒˡᵣ * qˡ / cᵖᵐ end ##### diff --git a/src/Thermodynamics/vapor_saturation.jl b/src/Thermodynamics/vapor_saturation.jl index ec60f018..96de6e35 100644 --- a/src/Thermodynamics/vapor_saturation.jl +++ b/src/Thermodynamics/vapor_saturation.jl @@ -6,7 +6,7 @@ Compute the [saturation vapor pressure](https://en.wikipedia.org/wiki/Vapor_pres using the Clausius-Clapeyron relation, ```math -dpᵛ⁺ / dT = pᵛ⁺ ℒᵝ(T) / (Rᵛ T^2) , +𝖽pᵛ⁺ / 𝖽T = pᵛ⁺ ℒᵝ(T) / (Rᵛ T^2) , ``` where the temperature-dependent latent heat of the surfaceis ``ℒᵝ(T)``. @@ -24,7 +24,7 @@ and the specific heat of phase ``β``. Note that we typically parameterize the latent heat interms of a reference temperature ``T = Tᵣ`` that is well above absolute zero. In that case, the latent heat is written - + ```math ℒᵝ = ℒᵝᵣ + Δcᵝ (T - Tᵣ), \\qquad \\text{and} \\qquad