Decapodes.Canon.Environment.glenConstant

Glens Law

Source

Nye, J. F. (1957). The Distribution of Stress and Velocity in Glaciers and Ice-Sheets. Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 239(1216), 113–133. http://www.jstor.org/stable/100184

Model

Γ::Form1
              
(A, ρ, g, n)::Constant
              
Γ == (2 / (n + 2)) * A * (ρ * g) ^ n
Decapodes.Canon.Environment.halfar_eq2Constant

Halfar (Eq. 2)

Source

Halfar, P. (1981), On the dynamics of the ice sheets, J. Geophys. Res., 86(C11), 11065–11072, doi:10.1029/JC086iC11p11065

Model

h::Form0
              
Γ::Form1
              
n::Constant
              
∂ₜ(h) == (∘(⋆, d, ⋆))(((Γ * d(h)) ∧ mag(♯(d(h))) ^ (n - 1)) ∧ h ^ (n + 2))
Decapodes.Canon.Environment.tracerConstant

Tracer

Source

Model

(c, C, F, c_up)::Form0
              
(v, V, q)::Form1
              
c_up == (((-1 * (⋆)(L(v, (⋆)(c))) - (⋆)(L(V, (⋆)(c)))) - (⋆)(L(v, (⋆)(C)))) - (∘(⋆, d, ⋆))(q)) + F
Decapodes.Canon.Biology.klausmeier_2aConstant

Klausmeier (Eq. 2a)

Source

Klausmeier, CA. “Regular and irregular patterns in semiarid vegetation.” Science (New York, N.Y.) vol. 284,5421 (1999): 1826-8. doi:10.1126/science.284.5421.1826

Model

(n, w)::DualForm0
              
dX::Form1
              
(a, ν)::Constant
              
∂ₜ(w) == ((a - w) - w * n ^ 2) + ν * ℒ(dX, w)
Decapodes.Canon.Biology.lejeuneConstant

Lejeune

Source

Lejeune, O., & Tlidi, M. (1999). A Model for the Explanation of Vegetation Stripes (Tiger Bush). Journal of Vegetation Science, 10(2), 201–208. https://doi.org/10.2307/3237141

Model

ρ::Form0
              
(μ, Λ, L)::Constant
              
∂ₜ(ρ) == (ρ * (((1 - μ) + (Λ - 1) * ρ) - ρ * ρ) + 0.5 * (L * L - ρ) * Δ(ρ)) - 0.125 * ρ * Δ(ρ) * Δ(ρ)
Decapodes.Canon.Physics.ficks_lawConstant

Ficks Law

Source

Equation for diffusion first stated by Adolf Fick. The diffusion flux is proportional to the concentration gradient.

Model

C::Form0
              
ϕ::Form1
              
ϕ == k(d₀(C))
Decapodes.Canon.Physics.jko_schemeConstant

Jordan-Kinderlehrer-Otto

Source

Jordan, R., Kinderlehrer, D., & Otto, F. (1998). The Variational Formulation of the Fokker–Planck Equation. In SIAM Journal on Mathematical Analysis (Vol. 29, Issue 1, pp. 1–17). Society for Industrial & Applied Mathematics (SIAM). https://doi.org/10.1137/s0036141096303359

Model

(ρ, Ψ)::Form0
              
β⁻¹::Constant
              
∂ₜ(ρ) == (∘(⋆, d, ⋆))(d(Ψ) ∧ ρ) + β⁻¹ * Δ(ρ)
Decapodes.Canon.Physics.mohamed_flowConstant

Mohamed Eq. 10, N2

Source

Model

𝐮::Form1
              
(P, 𝑝ᵈ)::Form0
              
(negone, half, μ)::Constant
              
∂ₜ(𝐮) == 𝐮̇
              
𝑝ᵈ == P + half * i(𝐮, 𝐮)
              
𝐮̇ == μ * (∘(d, ⋆, d, ⋆))(𝐮) + negone * (⋆₁⁻¹)(𝐮 ∧₁₀ₚᵈ (⋆)(d(𝐮))) + d(𝑝ᵈ)
Decapodes.Canon.Physics.momentumConstant

Momentum

Source

Model

(f, b)::Form0
              
(v, V, g, Fᵥ, uˢ, v_up)::Form1
              
τ::Form2
              
U::Parameter
              
uˢ̇ == ∂ₜ(uˢ)
              
v_up == (((((((-1 * L(v, v) - L(V, v)) - L(v, V)) - f ∧ v) - (∘(⋆, d, ⋆))(uˢ) ∧ v) - d(p)) + b ∧ g) - (∘(⋆, d, ⋆))(τ)) + uˢ̇ + Fᵥ
              
uˢ̇ == force(U)
Decapodes.Canon.Physics.navier_stokesConstant

Navier-Stokes

Source

Partial differential equations which describe the motion of viscous fluid surfaces.

Model

(V, V̇, G)::Form1{X}
              
(ρ, ṗ, p)::Form0{X}
              
V̇ == neg₁(L₁′(V, V)) + div₁(kᵥ(Δ₁(V) + third(d₀(δ₁(V)))), avg₀₁(ρ)) + d₀(half(i₁′(V, V))) + neg₁(div₁(d₀(p), avg₀₁(ρ))) + G
              
∂ₜ(V) == V̇
              
ṗ == neg₀((⋆₀⁻¹)(L₀(V, (⋆₀)(p))))
              
∂ₜ(p) == ṗ
Decapodes.Canon.Physics.oscillatorConstant

Oscillator

Source

Equation governing the motion of an object whose acceleration is negatively-proportional to its position.

Model

X::Form0
              
V::Form0
              
k::Constant
              
∂ₜ(X) == V
              
∂ₜ(V) == -k * X
Decapodes.Canon.Physics.poiseuilleConstant

Poiseuille

Source

A relation between the pressure drop in an incompressible and Newtownian fluid in laminar flow flowing through a long cylindrical pipe.

Model

P::Form0
              
q::Form1
              
(R, μ̃)::Constant
              
Δq == Δ(q)
              
∇P == d(P)
              
∂ₜ(q) == q̇
              
q̇ == μ̃ * ∂q(Δq) + ∇P + R * q
Decapodes.Canon.Physics.poiseuille_densityConstant

Poiseuille Density

Source

Model

q::Form1
              
(P, ρ)::Form0
              
(k, R, μ̃)::Constant
              
∂ₜ(q) == q̇
              
∇P == d(P)
              
q̇ == (μ̃ * ∂q(Δ(q)) - ∇P) + R * q
              
P == k * ρ
              
∂ₜ(ρ) == ρ̇
              
ρ_up == (∘(⋆, d, ⋆))(-1 * (ρ ∧₀₁ q))
              
ρ̇ == ∂ρ(ρ_up)
Decapodes.Canon.Physics.schroedingerConstant

Schoedinger

Source

The evolution of the wave function over time.

Model

(i, h, m)::Constant
              
V::Parameter
              
Ψ::Form0
              
∂ₜ(Ψ) == (((-1 * h ^ 2) / (2m)) * Δ(Ψ) + V * Ψ) / (i * h)
Decapodes.CartesianPointType
CartesianPoint{T}(p)

a point in cartesian coordinates, intended as a wrapper around Point3 from GeometryBasics.

Decapodes.SpherePointType
SpherePoint{T}(p)

a point in spherical coordinates, intended as a wrapper around Point3 from GeometryBasics.

Decapodes.TangentBasisMethod
tb(w)

Take a linear combinations of the tangent vectors at the base point. Use this to get a vector tangent to the sphere in the coordinate system of the base point. If the base point is in spherical coordinates, this is the identity, if the base point is in cartesian coordinates, it returns the tangent vector in cartesian coordinates.

Decapodes.find_unreachable_tvarsMethod
function find_unreachable_tvars(d)

Determine if the given Decapode can be compiled to an explicit time-stepping simulation. Use the simple check that one can traverse the Decapode starting from the state variables and reach all TVars (ignoring ∂ₜ edges).

Decapodes.gensimMethod
function gensim(d::AbstractNamedDecapode; dimension::Int=2)

Generate a simulation function from the given Decapode. The returned function can then be combined with a mesh and a function describing function mappings to return a simulator to be passed to solve.

Decapodes.Canon.Chemistry.GrayScottConstant

Gray-Scott

Source

A model of reaction-diffusion

Model

(U, V)::Form0
              
UV2::Form0
              
(U̇, V̇)::Form0
              
(f, k, rᵤ, rᵥ)::Constant
              
UV2 == U .* (V .* V)
              
U̇ == (rᵤ * Δ(U) - UV2) + f * (1 .- U)
              
V̇ == (rᵥ * Δ(V) + UV2) - (f + k) .* V
              
∂ₜ(U) == U̇
              
∂ₜ(V) == V̇
Decapodes.Canon.Chemistry.brusselatorConstant

Brusselator

Source

A model of reaction-diffusion for an oscillatory chemical system.

Model

(U, V)::Form0{X}
              
(U2V, One)::Form0{X}
              
(U̇, V̇)::Form0{X}
              
α::Constant{X}
              
F::Parameter{X}
              
U2V == (U .* U) .* V
              
U̇ == ((1 + U2V) - 4.4U) + α * Δ(U) + F
              
V̇ == (3.4U - U2V) + α * Δ(U)
              
∂ₜ(U) == U̇
              
∂ₜ(V) == V̇