@muladd begin
@inline nvariables(::AbstractEquations{NDIMS, NVARS}) where {NDIMS, NVARS} = NVARS
"""
eachvariable(equations::AbstractEquations)
Return an iterator over the indices that specify the location in relevant data structures
for the variables in `equations`. In particular, not the variables themselves are returned.
"""
@inline eachvariable(equations::AbstractEquations) = Base.OneTo(nvariables(equations))
"""
get_name(equations::AbstractEquations)
Returns the canonical, human-readable name for the given system of equations.
# Examples
```jldoctest
julia> Trixi.get_name(CompressibleEulerEquations1D(1.4))
"CompressibleEulerEquations1D"
```
"""
get_name(equations::AbstractEquations) = equations |> typeof |> nameof |> string
"""
varnames(conversion_function, equations)
Return the list of variable names when applying `conversion_function` to the
conserved variables associated to `equations`. This function is mainly used
internally to determine output to screen and for IO, e.g., for the
[`AnalysisCallback`](@ref) and the [`SaveSolutionCallback`](@ref).
Common choices of the `conversion_function` are [`cons2cons`](@ref) and
[`cons2prim`](@ref).
"""
function varnames end
function get_variable_index(varname, equations;
solution_variables = cons2cons)
index = findfirst(==(varname), varnames(solution_variables, equations))
if isnothing(index)
throw(ArgumentError("$varname is no valid variable."))
end
return index
end
function Base.show(io::IO, equations::AbstractEquations)
@nospecialize equations
print(io, get_name(equations), " with ")
if nvariables(equations) == 1
print(io, "one variable")
else
print(io, nvariables(equations), " variables")
end
end
function Base.show(io::IO, ::MIME"text/plain", equations::AbstractEquations)
@nospecialize equations
if get(io, :compact, false)
show(io, equations)
else
summary_header(io, get_name(equations))
summary_line(io, "#variables", nvariables(equations))
for variable in eachvariable(equations)
summary_line(increment_indent(io),
"variable " * string(variable),
varnames(cons2cons, equations)[variable])
end
summary_footer(io)
end
end
@inline Base.ndims(::AbstractEquations{NDIMS}) where {NDIMS} = NDIMS
Base.broadcastable(equations::AbstractEquations) = (equations,)
"""
flux(u, orientation_or_normal, equations)
Given the conservative variables `u`, calculate the (physical) flux in Cartesian
direction `orientation::Integer` or in arbitrary direction `normal::AbstractVector`
for the corresponding set of governing `equations`.
`orientation` is `1`, `2`, and `3` for the x-, y-, and z-directions, respectively.
"""
function flux end
"""
flux(u, normal_direction::AbstractVector, equations::AbstractEquations{1})
Enables calling `flux` with a non-integer argument `normal_direction` for one-dimensional
equations. Returns the value of `flux(u, 1, equations)` scaled by `normal_direction[1]`.
"""
@inline function flux(u, normal_direction::AbstractVector,
equations::AbstractEquations{1})
return normal_direction[1] * flux(u, 1, equations)
end
"""
rotate_to_x(u, normal, equations)
Apply the rotation that maps `normal` onto the x-axis to the convservative variables `u`.
This is used by [`FluxRotated`](@ref) to calculate the numerical flux of rotationally
invariant equations in arbitrary normal directions.
See also: [`rotate_from_x`](@ref)
"""
function rotate_to_x end
"""
rotate_from_x(u, normal, equations)
Apply the rotation that maps the x-axis onto `normal` to the convservative variables `u`.
This is used by [`FluxRotated`](@ref) to calculate the numerical flux of rotationally
invariant equations in arbitrary normal directions.
See also: [`rotate_to_x`](@ref)
"""
function rotate_from_x end
"""
BoundaryConditionDirichlet(boundary_value_function)
Create a Dirichlet boundary condition that uses the function `boundary_value_function`
to specify the values at the boundary.
This can be used to create a boundary condition that specifies exact boundary values
by passing the exact solution of the equation.
The passed boundary value function will be called with the same arguments as an initial condition function is called, i.e., as
```julia
boundary_value_function(x, t, equations)
```
where `x` specifies the coordinates, `t` is the current time, and `equation` is the corresponding system of equations.
# Examples
```julia
julia> BoundaryConditionDirichlet(initial_condition_convergence_test)
```
"""
struct BoundaryConditionDirichlet{B}
boundary_value_function::B
end
@inline function (boundary_condition::BoundaryConditionDirichlet)(u_inner,
orientation_or_normal,
direction,
x, t,
surface_flux_function,
equations)
u_boundary = boundary_condition.boundary_value_function(x, t, equations)
if iseven(direction)
flux = surface_flux_function(u_inner, u_boundary, orientation_or_normal,
equations)
else
flux = surface_flux_function(u_boundary, u_inner, orientation_or_normal,
equations)
end
return flux
end
@inline function (boundary_condition::BoundaryConditionDirichlet)(u_inner,
orientation_or_normal,
direction,
x, t,
surface_flux_functions::Tuple,
equations)
surface_flux_function, nonconservative_flux_function = surface_flux_functions
u_boundary = boundary_condition.boundary_value_function(x, t, equations)
if iseven(direction)
flux = surface_flux_function(u_inner, u_boundary, orientation_or_normal,
equations)
noncons_flux = nonconservative_flux_function(u_inner, u_boundary,
orientation_or_normal,
equations)
else
flux = surface_flux_function(u_boundary, u_inner, orientation_or_normal,
equations)
noncons_flux = nonconservative_flux_function(u_inner, u_boundary,
orientation_or_normal,
equations)
end
return flux, noncons_flux
end
@inline function (boundary_condition::BoundaryConditionDirichlet)(u_inner,
normal_direction::AbstractVector,
x, t,
surface_flux_function,
equations)
u_boundary = boundary_condition.boundary_value_function(x, t, equations)
flux = surface_flux_function(u_inner, u_boundary, normal_direction, equations)
return flux
end
@inline function (boundary_condition::BoundaryConditionDirichlet)(u_inner,
normal_direction::AbstractVector,
x, t,
surface_flux_functions::Tuple,
equations)
surface_flux_function, nonconservative_flux_function = surface_flux_functions
u_boundary = boundary_condition.boundary_value_function(x, t, equations)
flux = surface_flux_function(u_inner, u_boundary, normal_direction, equations)
noncons_flux = nonconservative_flux_function(u_inner, u_boundary, normal_direction,
equations)
return flux, noncons_flux
end
struct Gradient end
struct Divergence end
"""
BoundaryConditionNeumann(boundary_normal_flux_function)
Similar to `BoundaryConditionDirichlet`, but creates a Neumann boundary condition for parabolic
equations that uses the function `boundary_normal_flux_function` to specify the values of the normal
flux at the boundary.
The passed boundary value function will be called with the same arguments as an initial condition function is called, i.e., as
```julia
boundary_normal_flux_function(x, t, equations)
```
where `x` specifies the coordinates, `t` is the current time, and `equation` is the corresponding system of equations.
"""
struct BoundaryConditionNeumann{B}
boundary_normal_flux_function::B
end
"""
NonConservativeLocal()
Struct used for multiple dispatch on non-conservative flux functions in the format of "local * symmetric".
When the argument `nonconservative_type` is of type `NonConservativeLocal`,
the function returns the local part of the non-conservative term.
"""
struct NonConservativeLocal end
"""
NonConservativeSymmetric()
Struct used for multiple dispatch on non-conservative flux functions in the format of "local * symmetric".
When the argument `nonconservative_type` is of type `NonConservativeSymmetric`,
the function returns the symmetric part of the non-conservative term.
"""
struct NonConservativeSymmetric end
"""
NonConservativeJump()
Struct used for multiple dispatch on non-conservative flux functions in the format of "local * jump".
When the argument `nonconservative_type` is of type `NonConservativeJump`,
the function returns the jump part of the non-conservative term.
"""
struct NonConservativeJump end
"""
FluxNonConservative{STRUCTURE}
Abstract type for non-conservative fluxes that are composed of a local term and a structured two-point
term. The `STRUCTURE` type parameter should be set to [`NonConservativeSymmetric`](@ref) or
[`NonConservativeJump`](@ref), depending on the structure of the non-conservative term.
The abstract type is required for dispatch on the non-conservative type (symmetric / jump)
for the staggered volume flux computation in `calcflux_fhat!`.
"""
abstract type FluxNonConservative{STRUCTURE} end
"""
have_nonconservative_terms(equations)
Trait function determining whether `equations` represent a conservation law
with or without nonconservative terms. Classical conservation laws such as the
[`CompressibleEulerEquations2D`](@ref) do not have nonconservative terms. The
[`IdealGlmMhdEquations2D`](@ref) are an example of equations with nonconservative terms.
The return value will be `True()` or `False()` to allow dispatching on the return type.
"""
have_nonconservative_terms(::AbstractEquations) = False()
"""
n_nonconservative_terms(volume_flux_noncons)
Number of nonconservative terms for a particular nonconservative flux. Even for a specific equation,
this number may vary between different nonconservative fluxes.
This function needs to be specialized only if equations with nonconservative terms are
combined with certain solvers (e.g., subcell limiting).
"""
function n_nonconservative_terms end
"""
Trixi.combine_conservative_and_nonconservative_fluxes(flux, equations)
Trait function indicating whether the given `flux` and `equations` support
fusing the computation of conservative fluxes with nonconservative fluxes.
This is purely a performance optimization for equations with nonconservative
terms (i.e., where `have_nonconservative_terms(equations)` is `Trixi.True()`).
The default value is `Trixi.False()`, i.e., you have to pass a tuple of
numerical fluxes for the conservative and the nonconservative terms, e.g.,
to compute surface terms or the [`VolumeIntegralFluxDifferencing`](@ref).
For some systems and flux implementations, it is cheaper to compute
flux_noncons(u_ll, u_rr, orientation_or_normal_direction, equations)
and
flux_noncons(u_rr, u_ll, orientation_or_normal_direction, equations)
together, or to compute conservative and nonconservative flux contributions in
a single fused kernel. In this case, you should set this trait to be `Trixi.True()`
to take advantage of a more efficient implementation. In this case, you have to
define a single method that computes
flux_cons(u_ll, u_rr, n, equations) + 0.5f0 * flux_noncons(u_ll, u_rr, n, equations)
and
flux_cons(u_ll, u_rr, n, equations) + 0.5f0 * flux_noncons(u_rr, u_ll, n, equations)
together and returns them as a tuple.
See also the test section P4estMesh2D with combine_conservative_and_nonconservative_fluxes in
[Test Performance](https://github.com/trixi-framework/Trixi.jl/blob/main/test/test_performance_specializations_2d.jl).
"""
combine_conservative_and_nonconservative_fluxes(flux, ::AbstractEquations) = False()
"""
have_constant_speed(::AbstractEquations)
Indicates whether the characteristic speeds are constant, i.e., independent of the solution.
Queried in the timestep computation [`StepsizeCallback`](@ref).
This is the default fallback for nonlinear equations.
# Returns
- `False()`
"""
have_constant_speed(::AbstractEquations) = False()
"""
default_analysis_errors(equations)
Default analysis errors (`:l2_error` and `:linf_error`) used by the
[`AnalysisCallback`](@ref).
"""
default_analysis_errors(::AbstractEquations) = (:l2_error, :linf_error)
"""
default_analysis_integrals(equations)
Default analysis integrals used by the [`AnalysisCallback`](@ref).
"""
default_analysis_integrals(::AbstractEquations) = (entropy_timederivative,)
"""
cons2cons(u, equations)
Return the conserved variables `u`. While this function is as trivial as `identity`,
it is also as useful.
"""
@inline cons2cons(u, ::AbstractEquations) = u
@inline Base.first(u, ::AbstractEquations) = first(u)
"""
cons2prim(u, equations)
Convert the conserved variables `u` to the primitive variables for a given set of
`equations`. `u` is a vector type of the correct length `nvariables(equations)`.
Notice the function doesn't include any error checks for the purpose of efficiency,
so please make sure your input is correct.
The inverse conversion is performed by [`prim2cons`](@ref).
"""
function cons2prim end
"""
prim2cons(u, equations)
Convert the primitive variables `u` to the conserved variables for a given set of
`equations`. `u` is a vector type of the correct length `nvariables(equations)`.
Notice the function doesn't include any error checks for the purpose of efficiency,
so please make sure your input is correct.
The inverse conversion is performed by [`cons2prim`](@ref).
"""
function prim2cons end
"""
velocity(u, equations)
Return the velocity vector corresponding to the equations, e.g., fluid velocity for
Euler's equations. The velocity in certain orientation or normal direction (scalar) can be computed
with `velocity(u, orientation, equations)` or `velocity(u, normal_direction, equations)`
respectively. The `velocity(u, normal_direction, equations)` function calls
`velocity(u, equations)` to compute the velocity vector and then normal vector, thus allowing
a general function to be written for the AbstractEquations type. However, the
`velocity(u, orientation, equations)` is written for each equation separately to ensure
only the velocity in the desired direction (orientation) is computed.
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function velocity end
@inline function velocity(u, normal_direction::AbstractVector,
equations::AbstractEquations{2})
vel = velocity(u, equations)
v = vel[1] * normal_direction[1] + vel[2] * normal_direction[2]
return v
end
@inline function velocity(u, normal_direction::AbstractVector,
equations::AbstractEquations{3})
vel = velocity(u, equations)
v = vel[1] * normal_direction[1] + vel[2] * normal_direction[2] +
vel[3] * normal_direction[3]
return v
end
"""
entropy(u, equations)
Return the chosen entropy of the conserved variables `u` for a given set of
`equations`.
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function entropy end
"""
cons2entropy(u, equations)
Convert the conserved variables `u` to the entropy variables for a given set of
`equations` with chosen standard [`entropy`](@ref).
`u` is a vector type of the correct length `nvariables(equations)`.
Notice the function doesn't include any error checks for the purpose of efficiency,
so please make sure your input is correct.
The inverse conversion is performed by [`entropy2cons`](@ref).
"""
function cons2entropy end
"""
entropy2cons(w, equations)
Convert the entropy variables `w` based on a standard [`entropy`](@ref) to the
conserved variables for a given set of `equations`.
`u` is a vector type of the correct length `nvariables(equations)`.
Notice the function doesn't include any error checks for the purpose of efficiency,
so please make sure your input is correct.
The inverse conversion is performed by [`cons2entropy`](@ref).
"""
function entropy2cons end
"""
energy_total(u, equations)
Return the total energy of the conserved variables `u` for a given set of
`equations`, e.g., the [`CompressibleEulerEquations2D`](@ref).
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function energy_total end
"""
energy_kinetic(u, equations)
Return the kinetic energy of the conserved variables `u` for a given set of
`equations`, e.g., the [`CompressibleEulerEquations2D`](@ref).
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function energy_kinetic end
"""
energy_internal(u, equations)
Return the internal energy of the conserved variables `u` for a given set of
`equations`, e.g., the [`CompressibleEulerEquations2D`](@ref).
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function energy_internal end
"""
density(u, equations)
Return the density associated to the conserved variables `u` for a given set of
`equations`, e.g., the [`CompressibleEulerEquations2D`](@ref).
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function density end
"""
pressure(u, equations)
Return the pressure associated to the conserved variables `u` for a given set of
`equations`, e.g., the [`CompressibleEulerEquations2D`](@ref).
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function pressure end
"""
density_pressure(u, equations)
Return the product of the [`density`](@ref) and the [`pressure`](@ref)
associated to the conserved variables `u` for a given set of
`equations`, e.g., the [`CompressibleEulerEquations2D`](@ref).
This can be useful, e.g., as a variable for (shock-cappturing or AMR)
indicators as it combines two variables which must stay positive into one.
Furthermore, this implementation is for media which are described by an
ideal gas law alike equation of state more efficient than
computing [`pressure(u, equations)`](@ref) first and then multiplying with the density.
This is due to the fact that in computation of the pressure,
the kinetic energy needs to be computed, which usually involves
**division** of the squared momenta by the density.
This operation can be avoided!
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
"""
function density_pressure end
"""
waterheight(u, equations)
Return the water height associated to the conserved variables `u` for a given set of
`equations`.
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
!!! note
This function is defined in Trixi.jl to have a common interface for the
methods implemented in the subpackages [TrixiAtmo.jl](https://github.com/trixi-framework/TrixiAtmo.jl)
and [TrixiShallowWater.jl](https://github.com/trixi-framework/TrixiShallowWater.jl).
"""
function waterheight end
"""
waterheight_pressure(u, equations)
Return the product of the [`waterheight`](@ref) and the [`pressure`](@ref)
associated to the conserved variables `u` for a given set of
`equations`.
`u` is a vector of the conserved variables at a single node, i.e., a vector
of the correct length `nvariables(equations)`.
!!! note
This function is defined in Trixi.jl to have a common interface for the
methods implemented in the subpackages [TrixiAtmo.jl](https://github.com/trixi-framework/TrixiAtmo.jl)
and [TrixiShallowWater.jl](https://github.com/trixi-framework/TrixiShallowWater.jl).
"""
function waterheight_pressure end
"""
lake_at_rest_error(u, equations)
Calculate the point-wise error for the lake-at-rest steady state solution.
!!! note
This function is defined in Trixi.jl to have a common interface for the
methods implemented in the subpackages [TrixiAtmo.jl](https://github.com/trixi-framework/TrixiAtmo.jl)
and [TrixiShallowWater.jl](https://github.com/trixi-framework/TrixiShallowWater.jl).
"""
function lake_at_rest_error end
@inline function gradient_conservative(variable, u, equations)
return ForwardDiff.gradient(x -> variable(x, equations), u)
end
include("numerical_fluxes.jl")
abstract type AbstractLinearScalarAdvectionEquation{NDIMS} <:
AbstractEquations{NDIMS, 1} end
include("linear_scalar_advection_1d.jl")
include("linear_scalar_advection_2d.jl")
include("linear_scalar_advection_3d.jl")
abstract type AbstractInviscidBurgersEquation{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("inviscid_burgers_1d.jl")
abstract type AbstractShallowWaterEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
abstract type AbstractCompressibleEulerEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("compressible_euler_1d.jl")
include("compressible_euler_2d.jl")
include("compressible_euler_3d.jl")
include("compressible_euler_quasi_1d.jl")
abstract type AbstractNonIdealCompressibleEulerEquations{NDIMS, NVARS} <:
AbstractCompressibleEulerEquations{NDIMS, NVARS} end
include("equations_of_state.jl")
include("nonideal_compressible_euler.jl")
include("nonideal_compressible_euler_1d.jl")
include("nonideal_compressible_euler_2d.jl")
abstract type AbstractCompressibleEulerMulticomponentEquations{NDIMS, NVARS, NCOMP} <:
AbstractEquations{NDIMS, NVARS} end
include("compressible_euler_multicomponent_1d.jl")
include("compressible_euler_multicomponent_2d.jl")
abstract type AbstractPolytropicEulerEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("polytropic_euler_2d.jl")
include("passive_tracers.jl")
@inline function ncomponents(::AbstractCompressibleEulerMulticomponentEquations{NDIMS,
NVARS,
NCOMP}) where {
NDIMS,
NVARS,
NCOMP
}
return NCOMP
end
"""
eachcomponent(equations::AbstractCompressibleEulerMulticomponentEquations)
Return an iterator over the indices that specify the location in relevant data structures
for the components in `AbstractCompressibleEulerMulticomponentEquations`.
In particular, not the components themselves are returned.
"""
@inline function eachcomponent(equations::AbstractCompressibleEulerMulticomponentEquations)
return Base.OneTo(ncomponents(equations))
end
abstract type AbstractIdealGlmMhdEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("ideal_glm_mhd_1d.jl")
include("ideal_glm_mhd_2d.jl")
include("ideal_glm_mhd_3d.jl")
abstract type AbstractIdealGlmMhdMulticomponentEquations{NDIMS, NVARS, NCOMP} <:
AbstractEquations{NDIMS, NVARS} end
include("ideal_glm_mhd_multicomponent_1d.jl")
include("ideal_glm_mhd_multicomponent_2d.jl")
abstract type AbstractIdealGlmMhdMultiIonEquations{NDIMS, NVARS, NCOMP} <:
AbstractEquations{NDIMS, NVARS} end
include("ideal_glm_mhd_multiion.jl")
include("ideal_glm_mhd_multiion_2d.jl")
include("ideal_glm_mhd_multiion_3d.jl")
@inline function ncomponents(::AbstractIdealGlmMhdMulticomponentEquations{NDIMS, NVARS,
NCOMP}) where {
NDIMS,
NVARS,
NCOMP
}
return NCOMP
end
"""
eachcomponent(equations::AbstractIdealGlmMhdMulticomponentEquations)
Return an iterator over the indices that specify the location in relevant data structures
for the components in `AbstractIdealGlmMhdMulticomponentEquations`.
In particular, not the components themselves are returned.
"""
@inline function eachcomponent(equations::AbstractIdealGlmMhdMulticomponentEquations)
return Base.OneTo(ncomponents(equations))
end
@inline function ncomponents(::AbstractIdealGlmMhdMultiIonEquations{NDIMS, NVARS,
NCOMP}) where {
NDIMS,
NVARS,
NCOMP
}
return NCOMP
end
"""
eachcomponent(equations::AbstractIdealGlmMhdMultiIonEquations)
Return an iterator over the indices that specify the location in relevant data structures
for the components in `AbstractIdealGlmMhdMultiIonEquations`.
In particular, not the components themselves are returned.
"""
@inline function eachcomponent(equations::AbstractIdealGlmMhdMultiIonEquations)
return Base.OneTo(ncomponents(equations))
end
abstract type AbstractHyperbolicDiffusionEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("hyperbolic_diffusion_1d.jl")
include("hyperbolic_diffusion_2d.jl")
include("hyperbolic_diffusion_3d.jl")
abstract type AbstractLatticeBoltzmannEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("lattice_boltzmann_2d.jl")
include("lattice_boltzmann_3d.jl")
abstract type AbstractAcousticPerturbationEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("acoustic_perturbation_2d.jl")
abstract type AbstractLinearizedEulerEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("linearized_euler_1d.jl")
include("linearized_euler_2d.jl")
include("linearized_euler_3d.jl")
abstract type AbstractEquationsParabolic{NDIMS, NVARS, GradientVariables} <:
AbstractEquations{NDIMS, NVARS} end
abstract type AbstractTrafficFlowLWREquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("traffic_flow_lwr_1d.jl")
abstract type AbstractMaxwellEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("maxwell_1d.jl")
abstract type AbstractLinearElasticityEquations{NDIMS, NVARS} <:
AbstractEquations{NDIMS, NVARS} end
include("linear_elasticity_1d.jl")
end
