Helmholtz decomposition
For vector fields on the grid, we often prefer to work with potentials of that field and compose the vector field from these potentials. This is the so-called Helmholtz decomposition
\[\mathbf{v} = \nabla\phi + \nabla\times\psi\]
where $\psi$ is a vector potential (or streamfunction in two dimensions) and $\phi$ is a scalar potential. Note that, in an unbounded domain, $\phi$ satisfies a Poisson equation:
\[\nabla^2\phi = \nabla\cdot\mathbf{v}\]
and, assuming that $\psi$ has zero 'gauge' (i.e., its divergence is zero), then $\psi$ also satisfies a Poisson equation:
\[\nabla^2\psi = -\nabla\times\mathbf{v}\]
So we can recover the velocity field from the divergence and curl of the vector field (in fluid dynamics, these are the rate of dilatation and the vorticity, respectively) by solving these equations for the potentials and then evaluating their gradient and curl.
With immersed layers, these are extendable to domains with surfaces and their associated jump in the vector field $[\mathbf{v}] = \mathbf{v}^+ - \mathbf{v}^-$. The crucial relationships are
\[\nabla\cdot{\overline{\mathbf{v}}} = \overline{\nabla\cdot\mathbf{v}} + \delta(\chi)\mathbf{n}\cdot[\mathbf{v}]\]
and
\[\nabla\times{\overline{\mathbf{v}}} = \overline{\nabla\times\mathbf{v}} + \delta(\chi)\mathbf{n}\times[\mathbf{v}]\]
These are, respectively, the divergence and curl of the masked vector field, $\overline{\mathbf{v}}$. In ImmersedLayers.jl
, we are always working with the masked form of a field, whether we acknowledge it or not. In contrast, the first terms on the right-hand sides of these are the masked divergence and curl of the vector field in each side of the surface. The difference between these is a type of single-layer surface term–-a 'source sheet' or a 'vortex sheet', respectively–- due to the jump in vector field. If we omit the overbar notation, then the Helmholtz decomposition and the Poisson equations all still hold.
This package has a number of routines that can be used to perform the operations described here. They rely on two caches, ScalarPotentialCache
and VectorPotentialCache
, each of which supports one of the respective potentials. The function vecfield_helmholtz!
carries out all of these operations in order to assemble the vector field $\mathbf{v}$ from the masked divergence and curl fields and the jump in the vector field. (We can also add an additional irrotational, divergence-free vector field during this assembly.
Helmholtz decomposition functions
ImmersedLayers.ScalarPotentialCache
— TypeScalarPotentialCache(base_cache::BasicILMCache)
Create a cache for calculations involving the scalar potential field that induces a vector field on the grid. The base cache must be in support of vector field data.
ImmersedLayers.VectorPotentialCache
— TypeVectorPotentialCache(base_cache::BasicILMCache)
Create a cache for calculations involving the vector potential field (i.e. streamfunction) that induces a vector field on the grid. The base cache must be in support of vector field data.
ImmersedLayers.VectorFieldCache
— TypeVectorFieldCache(base_cache::BasicILMCache)
Create a cache for calculations involving the vector field derived from the vector and scalar potentials. The base cache must be in support of vector field data.
ImmersedLayers.vectorpotential_from_masked_curlv!
— Functionvectorpotential_from_masked_curlv!(ψ::Nodes{Dual},masked_curlv::Nodes{Dual},dv::VectorData,base_cache::BasicILMCache,wcache::VectorPotentialCache)
Return the vector potential field ψ
from the masked curl of the vector field, masked_curlv
($\overline{\nabla\times\mathbf{v}}$), the jump in the vector field across immersed surface dv
. It solves
$L_n\psi = -\overline{\nabla\times\mathbf{v}} - R_n\mathbf{n}\times[\mathbf{v}]$
and returns $\psi$.
ImmersedLayers.scalarpotential_from_masked_divv!
— Functionscalarpotential_from_masked_divv!(ϕ::Nodes{Primal},masked_divv::Nodes{Primal},dv::VectorData,base_cache::BasicILMCache,dcache::ScalarPotentialCache)
Return the scalar potential field ϕ
from the masked divergence of the vector field, masked_divv
($\overline{\nabla\cdot\mathbf{v}}$) the jump in the vector field across immersed surface dv
. It solves
$L_c\phi = \overline{\nabla\cdot\mathbf{v}} + R_c\mathbf{n}\cdot[\mathbf{v}]$
and returns $\phi$.
ImmersedLayers.vectorpotential_from_curlv!
— Functionvectorpotential_from_curlv!(ψ::Nodes{Dual},curlv::Nodes{Dual},base_cache::BasicILMCache)
Return the vector potential field ψ
from the curl of the masked vector field, curlv
($\nabla\times\overline{\mathbf{v}}$), the jump in the vector field across immersed surface dv
. It solves
$L_n\psi = -\nabla\times\overline{\mathbf{v}}$
and returns $\psi$.
ImmersedLayers.vecfield_from_vectorpotential!
— Functionvecfield_from_vectorpotential!(v::Edges{Primal},ψ::Nodes{Dual},base_cache::BasicILMCache)
Return the vector field v
associated with vector potential field ψ
(in 2-d, a scalar streamfunction), via the curl
$\mathbf{v} = \nabla\times\psi \mathbf{e}_z$
ImmersedLayers.masked_curlv_from_curlv_masked!
— Functionmasked_curlv_from_curlv_masked!(masked_curlv::Nodes{Dual},curlv::Nodes{Dual},dv::VectorData,base_cache::BasicILMCache,wcache::VectorPotentialCache)
Return the masked curl of the vector field $\overline{\nabla\times\mathbf{v}}$ (masked_curlv
) from the curl of the masked vector field $\nabla\times\overline{\mathbf{v}}$ (curlv
). It obtains this from
$\overline{\nabla\times\mathbf{v}} = \nabla\times\overline{\mathbf{v}} - R_n\mathbf{n}\times[\mathbf{v}]$
ImmersedLayers.curlv_masked_from_masked_curlv!
— Functioncurlv_masked_from_masked_curlv!(curlv::Nodes{Dual},masked_curlv::Nodes{Dual},dv::VectorData,base_cache::BasicILMCache,wcache::VectorPotentialCache)
Return the curl of the masked vector field $\nabla\times\overline{\mathbf{v}}$ (curlv
) from the masked curl of the vector field $\overline{\nabla\times\mathbf{v}}$ (masked_curlv
). It obtains this from
$\nabla\times\overline{\mathbf{v}} = \overline{\nabla\times\mathbf{v}} + R_n\mathbf{n}\times[\mathbf{v}]$
ImmersedLayers.scalarpotential_from_divv!
— Functionscalarpotential_from_divv!(ϕ::Nodes{Primal},divv::Nodes{Primal},base_cache::BasicILMCache,dcache::ScalarPotentialCache)
Return the scalar potential field ϕ
from the divergence of the masked vector field, divv
($\nabla\cdot\overline{\mathbf{v}}$) the jump in the vector field across immersed surface dv
. It solves
$L_c\phi = \nabla\cdot\overline{\mathbf{v}}$
and returns $\phi$.
ImmersedLayers.masked_divv_from_divv_masked!
— Functionmasked_divv_from_divv_masked!(masked_divv::Nodes{Primal},divv::Nodes{Primal},dv::VectorData,base_cache::BasicILMCache,dcache::ScalarPotentialCache)
Return the masked divergence of the vector field $\overline{\nabla\cdot\mathbf{v}}$ (masked_divv
) from the divergence of the masked vector field $\nabla\cdot\overline{\mathbf{v}}$ (divv
). It obtains this from
$\overline{\nabla\cdot\mathbf{v}} = \nabla\cdot\overline{\mathbf{v}} - R_c\mathbf{n}\cdot[\mathbf{v}]$
ImmersedLayers.divv_masked_from_masked_divv!
— Functiondivv_masked_from_masked_divv!(divv::Nodes{Primal},masked_divv::Nodes{Primal},dv::VectorData,base_cache::BasicILMCache,dcache::ScalarPotentialCache)
Return the divergence of the masked vector field $\nabla\cdot\overline{\mathbf{v}}$ (divv
) from the masked divergence of the vector field $\overline{\nabla\cdot\mathbf{v}}$ (masked_divv
). It obtains this from
$\nabla\cdot\overline{\mathbf{v}} = \overline{\nabla\cdot\mathbf{v}} + R_c\mathbf{n}\cdot[\mathbf{v}]$
ImmersedLayers.vecfield_from_scalarpotential!
— Functionvecfield_from_scalarpotential!(v::Edges{Primal},ϕ::Nodes{Primal},base_cache::BasicILMCache)
Return the vector field v
associated with scalar potential field ϕ
, via the gradient
$\mathbf{v} = \nabla\phi$
ImmersedLayers.vecfield_helmholtz!
— Functionvecfield_helmholtz!(v::Edges{Primal},curlv::Nodes{Dual},divv::Nodes{Primal},dv::VectorData,vp::Union{Edges{Primal},Nothing},base_cache::BasicILMCache,veccache::VectorFieldCache)
Recover the vector field v
from the masked curl field curlv
($\overline{\nabla\times\mathbf{v}}$), divergence field divv
($\overline{\nabla\cdot\mathbf{v}}$), surface vector jump dv
($[\mathbf{v}]$), and additional irrotational, divergence-free vector field vp
. It obtains this from the Helmholtz decomposition,
$\mathbf{v} = \nabla\times\psi + \nabla\phi + \mathbf{v}_p$
in which $\psi$ is the solution of
$L_n\overline{\psi} = -\overline{\nabla\times\mathbf{v}} - R_n\mathbf{n}\times[\mathbf{v}]$
in which $\phi$ is the solution of
$L_c\overline{\phi} = \overline{\nabla\cdot\mathbf{v}} + R_c\mathbf{n}\cdot[\mathbf{v}]$
It is important to stress that curlv
is the masked curl of the vector field, not the curl of the masked vector field. To get this the former from the latter, use masked_curlv_from_curlv_masked!
. Similarly for divv
, which is the masked divergence of the vector field, not the divergence of the masked vector field. One can use masked_divv_from_divv_masked!
for this.
To specify the irrotational, divergence-free vector field vp
, one can also simply provide a tuple, e.g., (1.0,0.0)
, to specify a uniform vector field.
ImmersedLayers.vectorpotential_uniformvecfield!
— Functionvectorpotential_uniformvecfield!(ψ::Nodes{Dual},Vx::Real,Vy::Real,base_cache::BasicILMCache)
Return a vector potential field ψ
associated with a uniform vector field with the components Vx
and Vy
.
ImmersedLayers.scalarpotential_uniformvecfield!
— Functionscalarpotential_uniformvecfield!(ϕ::Nodes{Primal},Vx::Real,Vy::Real,base_cache::BasicILMCache)
Return a vector potential field ϕ
associated with a uniform vector field with the components Vx
and Vy
.
ImmersedLayers.vecfield_uniformvecfield!
— Functionvecfield_uniformvecfield!(v::Edges{Primal},Vx::Real,Vy::Real,base_cache::BasicILMCache)
Return a vector field v
filled with the uniform components Vx
and Vy
, respectively.
This page was generated using Literate.jl.