Discontinuous Galerkin Methods

ClimateMachine.DGMethods.DGModelType
(dg::DGModel)(tendency, state_prognostic, _, t, α, β)

Uses spectral element discontinuous Galerkin in all direction to compute the tendency.

Computes the tendency terms compatible with IncrementODEProblem

tendency .= α .* dQdt(state_prognostic, p, t) .+ β .* tendency

The 4-argument form will just compute

tendency .= dQdt(state_prognostic, p, t)
source
DGModel(driver_config; kwargs...)

Initialize a DGModel given a DriverConfiguration and keyword arguments supported by DGModel.

source
ClimateMachine.DGMethods.ESDGModelType
ESDGModel

Contain type and functor that is used to evaluated the tendency for a entropy stable / conservative DGSEM discretization. Major fundamental difference between this and the more vanilla DGSEM is that the first order flux derivatives in the balance laws are evaluated using "flux-differencing". Namely, the following identities are used:

\[ ∂x f(q(x)) = 2∂x F(q(x), q(y))|_{x = y}, A(q(x)) ∂x q(x) = 2∂x D(q(x), q(y))|_{x = y},\]

where the numerical conservative flux F and numerical fluctuation flux D satisfy the following consistency and symmetry properties

\[ F(q, p) = F(p, q), F(q, q) = f(q), D(q, p) = B(q, p)(q - p), 2B(q, q) = A(p).\]

For the scheme to be entropy stable (and not just consistent) other properties of the numerical flux are also required. In particular, consider a balance laws of the form

\[ ∂t q + ∑_{j=1:d} (∂xj fj(q) + Aj(q) ∂xj q) = g(q, x, t),\]

where q is the state vector, fj is the conservative flux, and Aj nonconservative variable coefficient matrix, and g is the production field. Let there exists a scalar companion balance law of the form

\[ ∂t η(q) + ∑_{j=1:d} ∂xj ζj(q) = Π(q, x, t), Π(q, x, t) = β(q)^T g(q, x, t), β(q) = ∂q η(q).\]

Then for the scheme to be entropy stable it is requires that the numerical flux H(q, p) = F(q, p) + D(q, p) satisfy the following Tadmor-shuffle:

\[ β(q)^T Hj(q, p) - β(p)^T Hj(p, q) <= ψj(q) - ψ(p), ψj(q) = β(q)^T fj(q) - ζj(q);\]

when the equality is satisfied the scheme is called entropy conservative. For balance laws without a nonconservative term, ψj is the entropy potential.

source
ClimateMachine.DGMethods.DGFVModelType
(dgfvm::DGFVModel)(tendency, state_prognostic, _, t, α, β)

Uses spectral element discontinuous Galerkin in the horizontal and finite volume in the vertical to compute the tendency.

Computes the tendency terms compatible with IncrementODEProblem

tendency .= α .* dQdt(state_prognostic, p, t) .+ β .* tendency

The 4-argument form will just compute

tendency .= dQdt(state_prognostic, p, t)
source
ClimateMachine.DGMethods.remainder_DGModelFunction
remainder_DGModel(
    maindg::DGModel,
    subsdg::NTuple{NumModels, DGModel};
    numerical_flux_first_order,
    numerical_flux_second_order,
    numerical_flux_gradient,
    state_auxiliary,
    state_gradient_flux,
    states_higher_order,
    diffusion_direction,
    modeldata,
)

Constructs a DGModel from the maindg model and the tuple of subsdg models. The concept of a remainder model is that it computes the contribution of the model after subtracting all of the subcomponents.

By default the numerical fluxes are set to be a tuple of the main models numerical flux and the splitting is done at the PDE level (e.g., the remainder model is calculated prior to discretization). If instead a tuple of numerical fluxes is passed in the main numerical flux is evaluated first and then the subcomponent numerical fluxes are subtracted off. This is discretely different (for the Rusanov / local Lax-Friedrichs flux) than defining a numerical flux for the remainder of the physics model.

The other parameters are set to the value in the maindg component, mainly the data and arrays are aliased to the maindg values.

source
ClimateMachine.DGMethods.auxiliary_field_gradient!Function
auxiliary_field_gradient!(::BalanceLaw, ∇state::MPIStateArray,
                           vars_out, state::MPIStateArray, vars_in, grid;
                           direction = EveryDirection())

Take the gradient of the variables vars_in located in the array state and stores it in the variables vars_out of ∇state. This function computes element wise gradient without accounting for numerical fluxes and hence its primary purpose is to take the gradient of continuous reference fields.

Examples

FT = eltype(state_auxiliary)
grad_Φ = similar(state_auxiliary, vars=@vars(∇Φ::SVector{3, FT}))
auxiliary_field_gradient!(
    model,
    grad_Φ,
    ("∇Φ",),
    state_auxiliary,
    ("orientation.Φ",),
    grid,
)
source
ClimateMachine.DGMethods.courantFunction
courant(local_courant::Function, dg::DGModel, m::BalanceLaw,
        state_prognostic::MPIStateArray, direction=EveryDirection())

Returns the maximum of the evaluation of the function local_courant pointwise throughout the domain. The function local_courant is given an approximation of the local node distance Δx. The direction controls which reference directions are considered when computing the minimum node distance Δx. An example local_courant function is function localcourant(m::AtmosModel, stateprognostic::Vars, state_auxiliary::Vars, diffusive::Vars, Δx) return Δt * cmax / Δx end where Δt is the time step size and cmax is the maximum flow speed in the model.

source
DGMethods.courant(local_cfl, solver_config::SolverConfiguration;
                  Q=solver_config.Q, dt=solver_config.dt)

Returns the maximum of the evaluation of the function local_courant pointwise throughout the domain with the model defined by solver_config. The keyword arguments Q and dt can be used to call the courant method with a different state Q or time step dt than are defined in solver_config.

source

Mathematical Formulation

to be filled

Examples

Attribution

The style of examples we use here is heavily inspired by JuAFEM.jl

to be filled

Custom filter