Example Microphysics Implementation

This page walks through the implementation of a simple microphysics scheme to illustrate the developer interface. We implement a scheme that represents droplet and ice particle nucleation with constant-rate relaxation of specific humidity to saturation.

Defining the Scheme

First, we define a struct to hold the scheme parameters:

using Breeze

struct ExplicitMicrophysics{FT}
    vapor_to_liquid :: FT
    vapor_to_ice :: FT
end

Prognostic Field Names

This scheme is fully prognostic, meaning we carry vapor, liquid, and ice density as prognostic variables:

import Breeze.AtmosphereModels: prognostic_field_names

prognostic_field_names(::ExplicitMicrophysics) = (:ρqᵛ, :ρqˡ, :ρqⁱ)

prognostic_field_names (generic function with 13 methods)

Materializing Fields

When we materialize the microphysics fields, we must include all prognostic fields plus any diagnostic fields needed:

using Oceananigans: CenterField
import Breeze.AtmosphereModels: materialize_microphysical_fields

function materialize_microphysical_fields(::ExplicitMicrophysics, grid, boundary_conditions)
    # Prognostic fields (density-weighted)
    ρqᵛ = CenterField(grid; boundary_conditions=boundary_conditions.ρqᵛ)
    ρqˡ = CenterField(grid; boundary_conditions=boundary_conditions.ρqˡ)
    ρqⁱ = CenterField(grid; boundary_conditions=boundary_conditions.ρqⁱ)

    # Diagnostic field (specific humidity)
    qᵛ = CenterField(grid)

    return (; ρqᵛ, ρqˡ, ρqⁱ, qᵛ)
end

materialize_microphysical_fields (generic function with 9 methods)

Note we include the diagnostic field qᵛ (vapor mass fraction) in addition to the three prognostic fields.

Building the Microphysical State

The microphysical state encapsulates local values for tendency computation:

using Breeze.AtmosphereModels: AbstractMicrophysicalState
import Breeze.AtmosphereModels: microphysical_state

struct ExplicitMicrophysicsState{FT} <: AbstractMicrophysicalState{FT}
    qᵛ :: FT
    qˡ :: FT
    qⁱ :: FT
end

function microphysical_state(::ExplicitMicrophysics, ρ, μ::NamedTuple, 𝒰, velocities)
    # Convert density-weighted prognostics to specific quantities
    # velocities is required for interface compatibility (used by some schemes for aerosol activation)
    qᵛ = μ.ρqᵛ / ρ
    qˡ = μ.ρqˡ / ρ
    qⁱ = μ.ρqⁱ / ρ
    return ExplicitMicrophysicsState(qᵛ, qˡ, qⁱ)
end

microphysical_state (generic function with 14 methods)

Computing Tendencies

Tendencies are computed from the microphysical state. Each prognostic variable needs a tendency method:

import Breeze.AtmosphereModels: microphysical_tendency
using Breeze.Thermodynamics: temperature, saturation_specific_humidity,
                              PlanarLiquidSurface, PlanarIceSurface

# Tendency for liquid water density
@inline function microphysical_tendency(em::ExplicitMicrophysics, ::Val{:ρqˡ}, ρ, ℳ, 𝒰, constants)
    T = temperature(𝒰, constants)
    q⁺ˡ = saturation_specific_humidity(T, ρ, constants, PlanarLiquidSurface())
    τᵛˡ = em.vapor_to_liquid
    # Relaxation toward liquid saturation
    return ρ * (ℳ.qᵛ - q⁺ˡ) / τᵛˡ
end

# Tendency for ice density
@inline function microphysical_tendency(em::ExplicitMicrophysics, ::Val{:ρqⁱ}, ρ, ℳ, 𝒰, constants)
    T = temperature(𝒰, constants)
    q⁺ⁱ = saturation_specific_humidity(T, ρ, constants, PlanarIceSurface())
    τᵛⁱ = em.vapor_to_ice
    # Relaxation toward ice saturation
    return ρ * (ℳ.qᵛ - q⁺ⁱ) / τᵛⁱ
end

# Tendency for vapor density (conservation: what's lost to liquid/ice)
@inline function microphysical_tendency(em::ExplicitMicrophysics, ::Val{:ρqᵛ}, ρ, ℳ, 𝒰, constants)
    Sˡ = microphysical_tendency(em, Val(:ρqˡ), ρ, ℳ, 𝒰, constants)
    Sⁱ = microphysical_tendency(em, Val(:ρqⁱ), ρ, ℳ, 𝒰, constants)
    return -Sˡ - Sⁱ
end

microphysical_tendency (generic function with 18 methods)

Updating Auxiliary Fields

The update_microphysical_auxiliaries! function writes diagnostic fields:

import Breeze.AtmosphereModels: update_microphysical_auxiliaries!

@inline function update_microphysical_auxiliaries!(μ, i, j, k, grid,
                                                    ::ExplicitMicrophysics,
                                                    ℳ::ExplicitMicrophysicsState,
                                                    ρ, 𝒰, constants)
    @inbounds μ.qᵛ[i, j, k] = 𝒰.moisture_mass_fractions.vapor
    return nothing
end

update_microphysical_auxiliaries! (generic function with 9 methods)

Computing Moisture Fractions

The moisture_fractions function partitions total moisture:

using Breeze.Thermodynamics: MoistureMassFractions
import Breeze.AtmosphereModels: moisture_fractions

@inline function moisture_fractions(::ExplicitMicrophysics, ℳ::ExplicitMicrophysicsState, qᵗ)
    return MoistureMassFractions(ℳ.qᵛ, ℳ.qˡ, ℳ.qⁱ)
end

moisture_fractions (generic function with 16 methods)

Thermodynamic Adjustment

This is a fully prognostic scheme with no saturation adjustment, so we simply return the state unchanged:

import Breeze.AtmosphereModels: maybe_adjust_thermodynamic_state

@inline maybe_adjust_thermodynamic_state(𝒰, ::ExplicitMicrophysics, qᵗ, constants) = 𝒰

maybe_adjust_thermodynamic_state (generic function with 8 methods)

Summary

To implement a new microphysics scheme, we need to:

1. Define a struct to hold scheme parameters
2. Implement prognostic_field_names to list density-weighted prognostic variables
3. Implement materialize_microphysical_fields to create prognostic and diagnostic fields
4. Define a state type inheriting from AbstractMicrophysicalState
5. Implement microphysical_state to build state from prognostic scalars
6. Implement microphysical_tendency for each prognostic variable
7. Implement update_microphysical_auxiliaries! to update diagnostic fields
8. Implement moisture_fractions to partition moisture
9. Implement maybe_adjust_thermodynamic_state (trivial for non-equilibrium schemes)

For schemes with sedimentation, you would also implement velocity fields and the microphysical_velocities function.