Balance Laws

The balance law

ClimateMachine.BalanceLaws.BalanceLawType
abstract type BalanceLaw end

An abstract type representing a PDE balance law of the form:

\[\frac{dq}{dt} = \nabla \cdot F_1(q, a, t) + \nabla \cdot F_2(q, \nabla g, h, a, t) + S(q, \nabla g, a, t)\]

where:

  • $q$ is the prognostic state,
  • $a$ is the auxiliary state,
  • $g = G(q, a, t)$ is the gradient state (variables of which we compute the gradient),
  • $h$ is the hyperdiffusive state.

Subtypes of BalanceLaw should define the following interfaces:

If vars(bl, ::GradientFlux, FT) is non-empty, then the following should be defined:

If vars(bl, ::Hyperdiffusive, FT) is non-empty, then the following should be defined:

Additional functions:

source

Tendency types and methods

ClimateMachine.BalanceLaws.eq_tendsFunction
eq_tends(::AbstractPrognosticVariable, ::BalanceLaw, ::AbstractTendencyType)

A tuple of TendencyDefs given

  • AbstractPrognosticVariable prognostic variable
  • AbstractTendencyType tendency type
  • BalanceLaw balance law

i.e., a tuple of TendencyDefs corresponding to F₁, F₂, or S for a single prognostic variable in:

`∂_t Yᵢ + (∇•F₁(Y))ᵢ + (∇•F₂(Y,G)))ᵢ = (S(Y,G))ᵢ`
source
ClimateMachine.BalanceLaws.prognostic_varsFunction
prognostic_vars(::BalanceLaw)

A tuple of AbstractPrognosticVariables given the BalanceLaw.

i.e., a tuple of AbstractPrognosticVariables corresponding to the column-vector Yᵢ in:

`∂_t Yᵢ + (∇•F₁(Y))ᵢ + (∇•F₂(Y,G)))ᵢ = (S(Y,G))ᵢ`
source
ClimateMachine.BalanceLaws.get_prog_stateFunction
var, name = get_prog_state(state::Union{Vars, Grad}, pv::AbstractPrognosticVariable)

Returns a tuple of two elements. var is a Vars or Grad object, and name is a Symbol. They should be linked such that getproperty(var, name) returns the corresponding prognostic variable type pv.

Example

get_prog_state(state, ::TotalMoisture) = (state.moisture, :ρq_tot)
var, name = get_prog_state(state, TotalMoisture())
@test getproperty(var, name) == state.moisture.ρq_tot
source
ClimateMachine.BalanceLaws.precomputeFunction
precompute(bl, args, ::AbstractTendencyType)

A nested NamedTuple of precomputed (cached) values and or objects. This is useful for "expensive" point-wise quantities that are used in multiple tendency terms. For example, computing a quantity that requires iteration.

source
ClimateMachine.BalanceLaws.prognostic_var_source_mapFunction
prognostic_var_source_map(driver_sources::Tuple)

A DispatchedTuple, given a Tuple of the driver/experiment sources.

Note

prognostic_vars, which returns a Tuple of prognostic variable types, must be defined for boundary condition types.

source
ClimateMachine.BalanceLaws.show_tendenciesFunction
show_tendencies(
    bl;
    include_module = false,
    table_complete = false,
)

Show a table of the tendencies for each prognostic variable for the given balance law.

Arguments

  • include_module[ = false] will print not remove the module where each prognostic variable is defined (e.g., Atmos.Mass).
  • table_complete[ = false] will print a warning (if false) that the tendency table is incomplete.

Requires definitions for

for the balance law.

source

Methods for fluxes and sources

ClimateMachine.BalanceLaws.ΣfluxesFunction
Σfluxes(fluxes::NTuple, bl, args)

Sum of the fluxes where

  • fluxes is an NTuple{N, TendencyDef{Flux{O}}} where {N, O}
  • bl is the balance law
  • args are the arguments passed to the individual flux functions
source
ClimateMachine.BalanceLaws.ΣsourcesFunction
Σsources(sources::NTuple, bl, args)

Sum of the sources where

  • sources is an NTuple{N, TendencyDef{Source}} where {N}
  • bl is the balance law
  • args are the arguments passed to the individual source functions
source

State variable types

ClimateMachine.BalanceLaws.AuxiliaryType
Auxiliary <: AbstractStateType

Auxiliary variables help serve several purposes:

  • Pre-compute and store "expensive" variables, for example, quantities computed in vertical integrals.
  • Diagnostic exports
source

Interface

Variable specification methods

ClimateMachine.BalanceLaws.vars_stateFunction
BalanceLaws.vars_state(::BL, ::AbstractStateType, FT)

Defines the state variables of a BalanceLaw subtype BL with floating point type FT.

For each AbstractStateType, this should return a NamedTuple type, with element type either FT, an SArray with element type FT or another NamedTuple satisfying the same property.

For convenience, we recommend using the VariableTemplates.@vars macro.

Example

struct MyBalanceLaw <: BalanceLaw end

BalanceLaws.vars_state(::MyBalanceLaw, ::Prognostic, FT) =
    @vars(x::FT, y::SVector{3, FT})
BalanceLaws.vars_state(::MyBalanceLaw, ::Auxiliary, FT) =
    @vars(components::@vars(a::FT, b::FT))
source

Initial condition methods

ClimateMachine.BalanceLaws.init_state_auxiliary!Function
init_state_auxiliary!(
    ::BL,
    statearray_auxiliary,
    geom::LocalGeometry,
)

Sets the initial state of the auxiliary variables state_auxiliary at each node for a BalanceLaw subtype BL. By default this calls nodal_init_state_auxiliary!.

source
init_state_auxiliary!(::HBModel)

sets the initial value for auxiliary variables (those that aren't related to vertical integrals) dispatches to oceaninitaux! which is defined in a problem file such as SimpleBoxProblem.jl

source
No need to init, initialize by full model
source
init_state_auxiliary!(
    m::AtmosModel,
    aux::Vars,
    grid,
    direction
)

Initialise auxiliary variables for each AtmosModel subcomponent. Store Cartesian coordinate information in aux.coord.

source
init_state_auxiliary!(
    bl::BalanceLaw,
    f!,
    statearray_auxiliary,
    grid,
    direction;
    state_temporary = nothing
)

Apply f!(bl, state_auxiliary, tmp, geom) at each node, storing the result in statearray_auxiliary, where tmp are the values at the corresponding node in state_temporary and geom contains the geometry information.

source

Source term kernels

Integral kernels

Gradient/Laplacian kernels

ClimateMachine.BalanceLaws.compute_gradient_flux!Function
compute_gradient_flux!(
    ::BL,
    state_gradient_flux::Vars,
    ∇transformstate::Grad,
    state_prognostic::Vars,
    state_auxiliary::Vars,
    t::Real
)

Transformation of gradients to the diffusive variables for a BalanceLaw subtype BL. This should be a linear function of ∇transformstate

source

Boundary conditions

ClimateMachine.BalanceLaws.boundary_conditionsFunction
boundary_conditions(::BL)

Define the set of boundary conditions for the balance law BL. This should return a tuple, where a boundary tagged with the integer i will use the ith element of the tuple.

source
ClimateMachine.BalanceLaws.boundary_state!Function
boundary_state!(
    ::NumericalFluxGradient,
    ::BC
    ::BL,
    state_prognostic⁺::Vars,
    state_auxiliary⁺::Vars,
    normal⁻,
    state_prognostic⁻::Vars,
    state_auxiliary⁻::Vars,
    t
)
boundary_state!(
    ::NumericalFluxFirstOrder,
    ::BC
    ::BL,
    state_prognostic⁺::Vars,
    state_auxiliary⁺::Vars,
    normal⁻,
    state_prognostic⁻::Vars,
    state_auxiliary⁻::Vars,
    t
)
boundary_state!(
    ::NumericalFluxSecondOrder,
    ::BC
    ::BL,
    state_prognostic⁺::Vars,
    state_gradient_flux⁺::Vars,
    state_auxiliary⁺:
    Vars, normal⁻,
    state_prognostic⁻::Vars,
    state_gradient_flux⁻::Vars,
    state_auxiliary⁻::Vars,
    t
)

Specify the opposite (+ side) face for the boundary condition type BC with balance law BL.

  • NumericalFluxGradient numerical flux (internal method)
  • NumericalFluxFirstOrder first-order unknowns
  • NumericalFluxSecondOrder second-order unknowns
source

Auxiliary kernels

ClimateMachine.BalanceLaws.wavespeedFunction
wavespeed(
    ::BL,
    n⁻,
    state_prognostic::Vars,
    state_auxiliary::Vars,
    t::Real,
    direction
)

Wavespeed in the direction n⁻ for a BalanceLaw subtype BL. This is required to be defined if using a RusanovNumericalFlux numerical flux.

source
ClimateMachine.BalanceLaws.update_auxiliary_state!Function
update_auxiliary_state!(
    dg::DGModel,
    m::BalanceLaw,
    statearray_aux,
    t::Real,
    elems::UnitRange,
    [diffusive=false]
)

Hook to update the auxiliary state variables before calling any other functions.

By default, this calls nodal_update_auxiliary_state! at each node.

If diffusive=true, then state_gradflux is also passed to nodal_update_auxiliary_state!.

source
update_auxiliary_state!(::HBModel)

Applies the vertical filter to the zonal and meridional velocities to preserve numerical incompressibility Applies an exponential filter to θ to anti-alias the non-linear advective term

Doesn't actually touch the aux variables any more, but we need a better filter interface than this anyways

source
ClimateMachine.BalanceLaws.update_auxiliary_state_gradient!Function
update_auxiliary_state_gradient!(
    dg::DGModel,
    m::BalanceLaw,
    statearray_aux,
    t::Real,
    elems::UnitRange,
    [diffusive=false]
)

Hook to update the auxiliary state variables after the gradient computation.

By default, this calls nothing.

If diffusive=true, then state_gradflux is also passed to nodal_update_auxiliary_state!.

source