CombinatorialSpaces.jl provides discrete differential operators defined on simplicial sets via the Discrete Exterior Calculus (DEC).
There are two modules for these DEC operators. The first, CombinatorialSpaces.DiscreteExteriorCalculus, serves as our reference implementation. It is suitable for any manifold-like delta dual complex. The second, CombinatorialSpaces.FastDEC, offers more efficient operators, both in construction, and in execution time. The operators offered by both modules agree up to differences introduced by re-ordering floating-point operations.
Certain operators in FastDEC are made more efficient by assuming that the delta dual complex has not been altered in any way after it is created and volumes have been assigned with subdivide_duals!. So, such operators should not be called on such a complex if it has been manually edited.
The discrete exterior calculus (DEC) for simplicial sets.
This module provides the dual complex associated with a delta set (the primal complex), which is a discrete incarnation of Hodge duality, as well as the many operators of the DEC that depend on it, such as the Hodge star, codifferential, wedge product, interior product, and Lie derivative. The main reference for this module is Hirani's 2003 PhD thesis.
Abstract type for dual complex of a 1D delta set.
Abstract type for dual complex of a 2D delta set.
Abstract type for dual complex of a 3D delta set.
Barycenter, aka centroid, of a simplex.
Circumcenter, or center of circumscribed circle, of a simplex.
The circumcenter is calculated by inverting the Cayley-Menger matrix, as explained by Westdendorp. This method of calculation is also used in the package AlphaShapes.jl.
Dual complex of a one-dimensional delta set.
The data structure includes both the primal complex and the dual complex, as well as the mapping between them.
Dual complex of a two-dimensional delta set.
Dual complex of a three-dimensional delta set.
Wrapper for chain of dual cells of dimension n.
In an $N$-dimensional complex, the elementary dual simplices of each $n$-simplex together comprise the dual $(N-n)$-cell of the simplex. Using this correspondence, a basis for primal $n$-chains defines the basis for dual $(N-n)$-chains.
In (Hirani 2003, Definition 3.4.1), the duality operator assigns a certain sign to each elementary dual simplex. For us, all of these signs should be regarded as positive because we have already incorporated them into the orientation of the dual simplices.
Edge in simplicial set: alias for Simplex{1}.
Wrapper for form, aka cochain, on dual cells of dimension n.
Tetrahedron in simplicial set: alias for Simplex{3}.
Triangle in simplicial set: alias for Simplex{2}.
Vertex in simplicial set: alias for Simplex{0}.
Wrapper for vector field on dual vertices.
Embedded dual complex of an embedded 1D delta set.
Although they are redundant information, the lengths of the primal and dual edges are precomputed and stored.
Embedded dual complex of an embedded 2D delta set.
Although they are redundant information, the lengths and areas of the primal/dual edges and triangles are precomputed and stored.
Embedded dual complex of an embedded 3D delta set.
Incenter, or center of inscribed circle, of a simplex.
Oriented dual complex of an oriented 1D delta set.
Oriented dual complex of an oriented 2D delta set.
Oriented dual complex of an oriented 3D delta set.
Wrapper for vector field on primal vertices.
Wedge product of discrete forms.
The wedge product of a $k$-form and an $l$-form is a $(k+l)$-form.
The DEC and related systems have several flavors of wedge product. This one is the discrete primal-primal wedge product introduced in (Hirani, 2003, Chapter 7) and (Desbrun et al 2005, Section 8). It depends on the geometric embedding and requires the dual complex. Note that we diverge from Hirani in that his formulation explicitly divides by (k+1)!. We do not do so in this computation.
Laplace-de Rham operator on discrete forms.
This linear operator on primal $n$-forms is defined by $Δ := δ d + d δ$. Restricted to 0-forms, it reduces to the negative of the Laplace-Beltrami operator ∇²: $Δ f = -∇² f$.
Hodge star operator from primal 1-forms to dual 1-forms.
This specific hodge star implementation is based on the hodge star presented in (Ayoub et al 2020), which generalizes the operator presented in (Hirani 2003). This reproduces the diagonal hodge for a dual mesh generated under circumcentric subdivision and provides off-diagonal correction factors for meshes generated under other subdivision schemes (e.g. barycentric).
Hodge star operator from primal $n$-forms to dual $N-n$-forms.
Some authors, such as (Hirani 2003) and (Desbrun 2005), use the symbol $⋆$ for the duality operator on chains and the symbol $*$ for the Hodge star operator on cochains. We do not explicitly define the duality operator and we use the symbol $⋆$ for the Hodge star.
Alias for the codifferential operator δ.
Discrete exterior derivative of dual form.
Transpose of dual_boundary. For more info, see (Desbrun, Kanso, Tong, 2008: Discrete differential forms for computational modeling, §4.5).
Boundary dual vertices of a dual edge.
This accessor assumes that the simplicial identities for the dual hold.
Point associated with dual vertex of complex.
Boundary dual vertices of a dual tetrahedron.
This accessor assumes that the simplicial identities for the dual hold.
Boundary dual vertices of a dual triangle.
This accessor assumes that the simplicial identities for the dual hold.
Calls the constructor for d's dual type on d, including parameters. Does not call subdivide_duals! on the result. Should work out of the box on new DeltaSet types if (1) their dual type has the same name as their primal type with "Set" substituted by "DualComplex" and (2) their dual type has the same parameter set as their primal type. At the time of writing (PR 117) only "Embedded" types fail criterion (2) and get special treatment.
Examples
s = EmbeddedDeltaSet2D{Bool,SVector{2,Float64}}() dualize(s)::EmbeddedDeltaDualComplex2D{Bool,Float64,SVector{2,Float64}}
Dual vertex corresponding to center of primal edge.
List of elementary dual simplices corresponding to primal simplex.
In general, in an $n$-dimensional complex, the elementary duals of primal $k$-simplices are dual $(n-k)$-simplices. Thus, in 1D dual complexes, the elementary duals of...
- primal vertices are dual edges
- primal edges are (single) dual vertices
In 2D dual complexes, the elementary duals of...
- primal vertices are dual triangles
- primal edges are dual edges
- primal triangles are (single) dual vertices
In 3D dual complexes, the elementary duals of...
- primal vertices are dual tetrahedra
- primal edges are dual triangles
- primal triangles are dual edges
- primal tetrahedra are (single) dual vertices
Alias for the flat operator ♭.
Alias for the flat-sharp dual-to-primal interpolation operator ♭♯.
Alias for the flat-sharp dual-to-primal interpolation matrix ♭♯_mat.
Calculate the center of simplex spanned by given points.
The first argument is a list of points and the second specifies the notion of "center", via an instance of SimplexCenter.
Alias for the Hodge star operator ⋆.
Interior product of a vector field (or 1-form) and a $n$-form.
Specifically, this operation is the primal-dual interior product defined in (Hirani 2003, Section 8.2) and (Desbrun et al 2005, Section 10). Thus it takes a primal vector field (or primal 1-form) and a dual $n$-forms and then returns a dual $(n-1)$-form.
Interior product of a 1-form and a $n$-form, yielding an $(n-1)$-form.
Usually, the interior product is defined for vector fields; this function assumes that the flat operator ♭ (not yet implemented for primal vector fields) has already been applied to yield a 1-form.
Inverse Hodge star operator from dual $N-n$-forms to primal $n$-forms.
Confusingly, this is not the operator inverse of the Hodge star ⋆ because it carries an extra global sign, in analogy to the smooth case (Gillette, 2009, Notes on the DEC, Definition 2.27).
Alias for the Laplace-Beltrami operator ∇².
Alias for the Laplace-de Rham operator Δ.
Alias for Lie derivative operator ℒ.
Lie derivative of $n$-form with respect to a 1-form.
Assumes that the flat operator ♭ has already been applied to the vector field.
p2_d2_interpolation(sd::HasDeltaSet2D)
Generates a sparse matrix that converts data on primal 2-forms into data on dual 2-forms.
Primal vertex associated with a dual simplex.
Alias for the sharp operator ♯.
Compute geometric subdivision for embedded dual complex.
Supports different methods of subdivision through the choice of geometric center, as defined by geometric_center. In particular, barycentric subdivision and circumcentric subdivision are supported.
List of dual simplices comprising the subdivision of a primal simplex.
A primal $n$-simplex is always subdivided into $(n+1)!$ dual $n$-simplices, not be confused with the elementary_duals which have complementary dimension.
The returned list is ordered such that subsimplices with the same primal vertex appear consecutively.
Dual vertex corresponding to center of primal tetrahedron.
Dual vertex corresponding to center of primal triangle.
Dual vertex corresponding to center of primal vertex.
Alias for the wedge product operator ∧.
Codifferential operator from primal $n$ forms to primal $n-1$-forms.
Lie derivative of $n$-form with respect to a vector field (or 1-form).
Specifically, this is the primal-dual Lie derivative defined in (Hirani 2003, Section 8.4) and (Desbrun et al 2005, Section 10).
Laplace-Beltrami operator on discrete forms.
This linear operator on primal $n$-forms defined by $∇² α := -δ d α$, where δ is the codifferential and d is the exterior derivative.
For following texts such as Abraham-Marsden-Ratiu, we take the sign convention that makes the Laplace-Beltrami operator consistent with the Euclidean Laplace operator (the divergence of the gradient). Other authors, such as (Hirani 2003), take the opposite convention, which has the advantage of being consistent with the Laplace-de Rham operator Δ.
Flat operator converting vector fields to 1-forms.
A generic function for discrete flat operators. Currently the DPP-flat from (Hirani 2003, Definition 5.5.2) and (Desbrun et al 2005, Definition 7.3) is implemented, as well as a primal-to-primal flat, which assumes linear-interpolation of the vector field across an edge, determined solely by the values at the endpoints.
See also: the sharp operator ♯.
♭♯(s::HasDeltaSet2D, α::SimplexForm{1})
Make a dual 1-form primal by chaining ♭ᵈᵖ♯ᵈᵈ.
This returns the given dual 1-form as a primal 1-form. See also ♭♯_mat.
♭♯_mat(s::HasDeltaSet2D)
Make a dual 1-form primal by chaining ♭ᵈᵖ♯ᵈᵈ.
This returns a matrix which can be multiplied by a dual 1-form. See also ♭♯.
Sharp operator for converting dual 1-forms to dual vector fields.
This dual-dual sharp uses a method of local linear least squares to provide a tangent vector field.
See also: ♯_mat, which returns a matrix that encodes this operator.
Sharp operator for converting primal 1-forms to primal vector fields.
This the primal-primal sharp from Hirani 2003, Definition 5.8.1 and Remark 2.7.2.
A PP-flat is also defined in (Desbrun et al 2005, Definition 7.4) but differs in two ways: Desbrun et al's notation suggests a unit normal vector, whereas the gradient of Hirani's primal-primal interpolation function is not necessarily a unit vector. More importantly, Hirani's vector is a normal to a different face than Desbrun et al's, with further confusion created by the fact that Hirani's Figure 5.7 agrees with Desbrun et al's description rather than his own. That being said, to the best of our knowledge, our implementation is the correct one and agrees with Hirani's description, if not his figure.
See also: ♭ and ♯_mat, which returns a matrix that encodes this operator.
function ♯_mat(s::AbstractDeltaDualComplex2D, DS::DiscreteSharp)
Sharpen a 1-form into a vector field.
3 primal-primal methods are supported. See ♯_denominator for the distinction between Hirani's method and and an "Alternative" method. Desbrun's definition is selected with DesbrunSharp, and is like Hirani's, save for dividing by the norm twice.
A dual-dual method which uses linear least squares to estimate a vector field is selected with LLSDDSharp.
function ♯_mat(s::AbstractDeltaDualComplex2D, ::LLSDDSharp)
Sharpen a dual 1-form into a DualVectorField, using linear least squares.
Up to floating point error, this method perfectly produces fields which are constant over any triangle in the domain. Assume that the contribution of each half-edge to the value stored on the entire dual edge is proportional to their lengths. Since this least squares method does not perform pre-normalization, the contribution of each half-edge value is proportional to its length on the given triangle. Satisfying the continuous exterior calculus, sharpened vectors are constrained to lie on their triangle (i.e. they are indeed tangent).
It is not known whether this method has been exploited previously in the DEC literature, or defined in code elsewhere.
The discrete exterior calculus (DEC) with high performance in mind.
This module provides similar functionality to the DiscreteExteriorCalculus module but uses assumptions about the ACSet mesh structure to greatly improve performance. Some operators, like the exterior derivative are returned as sparse matrices while others, like the wedge product, are instead returned as functions that will compute the product.
∧(s::HasDeltaSet, α::DualForm{1}, β::SimplexForm{1})
Wedge product of a dual 1-form and a primal 1-form.
Chain the musical isomorphisms to interpolate the dual 1-form to a primal 1-form. Then use the CombinatorialSpaces version of the Hirani primal-primal weddge (without explicitly dividing by 2.)
∧(s::HasDeltaSet, α::SimplexForm{1}, β::DualForm{1})
Wedge product of a primal 1-form and a dual 1-form.
Chain the musical isomorphisms to interpolate the dual 1-form to a primal 1-form, using the linear least squares ♯. Then use the CombinatorialSpaces version of the Hirani primal-primal weddge.
Alias for the averaging operator [`avg₀₁`](@ref).
Alias for the averaging matrix [`avg₀₁_mat`](@ref).
avg₀₁(s::HasDeltaSet, α::SimplexForm{0})
Turn a 0-form into a 1-form by averaging data stored on the face of an edge.
See also avg₀₁_mat.
Averaging matrix from 0-forms to 1-forms.
Given a 0-form, this matrix computes a 1-form by taking the mean of value stored on the faces of each edge.
This matrix can be used to implement a wedge product: (avg₀₁(s)*X) .* Y where X is a 0-form and Y a 1-form, assuming the center of an edge is halfway between its endpoints.
See also avg₀₁.
d0p0interpolation(sd::HasDeltaSet2D; hodge=GeometricHodge())
Generates a sparse matrix that converts data on dual 0-forms into data on primal 0-forms. This uses the p2_d2_interpolation function as an intermediate step.
dec_boundary(n::Int, sd::HasDeltaSet)
Return the boundary operator (as a matrix) for (n+1)-simplices to (n)-simplices
dec_differential(n::Int, sd::HasDeltaSet)
Return the exterior derivative (as a matrix) between n-simplices and (n+1)-simplices
dec_dual_derivative(n::Int, sd::HasDeltaSet)
Return the dual exterior derivative (as a matrix) between dual n-simplices and dual (n+1)-simplices
dec_hodge_star(n::Int, sd::HasDeltaSet; hodge=GeometricHodge())
Return the hodge matrix between n-simplices and dual 'n'-simplices.
dec_inv_hodge_star(n::Int, sd::HasDeltaSet; hodge=GeometricHodge())
Return the inverse hodge matrix between dual n-simplices and 'n'-simplices.
dec_wedge_product(::Type{Tuple{m,n}}, sd::HasDeltaSet, backend=Val{:CPU}, arr_cons=identity, cast_float=nothing) where {m,n}
Return a function that computes the wedge product between a primal m-form and a primal n-form, assuming special properties of the mesh.
It is assumed... ... for the 0-1 wedge product, that the dual vertex on an edge is at the midpoint. ... for the 1-1 wedge product, that the dual mesh simplices are in the default order as returned by the dual complex constructor.
Arguments:
Tuple{m,n}: the degrees of the differential forms. sd: the simplicial complex. backend=Val{:CPU}: a value-type to select special backend logic, if implemented. arr_cons=identity: a constructor of the desired array type on the appropriate backend e.g. MtlArray. cast_float=nothing: a specific Float type to use e.g. Float32. Otherwise, the type of the first differential form will be used.
dec_wedge_product_dd(::Type{Tuple{m,n}}, sd::HasDeltaSet) where {m,n}
Return a function that computes the wedge product between a dual m-form and a dual n-form.
The currently supported dual-dual wedges are 0-1 and 1-0.
dec_wedge_product_dp(::Type{Tuple{m,n}}, sd::HasDeltaSet) where {m,n}
Return a function that computes the wedge product between a dual m-form and a primal n-form.
It is assumed... ... for the 1-0 and 0-1 wedge product, that means are barycentric, and performs bilinear interpolation. It is not known if this definition has appeared in the literature or any code.
The currently supported dual-primal wedges are 0-1, 1-0, and 1-1.
dec_wedge_product_pd(::Type{Tuple{m,n}}, sd::HasDeltaSet) where {m,n}
Return a function that computes the wedge product between a primal m-form and a dual n-form.
See dec_wedge_product_dp for assumptions.
dec_Δ⁻¹(::Type{Val{0}}, s::AbstractGeometricMapSeries; steps = 3, cycles = 5, alg = cg, μ = 2)
Return a function that solves the inverse Laplacian problem.
function interior_product_dd(::Type{Tuple{1,1}}, s::SimplicialSets.HasDeltaSet)
Given a dual 1-form and a dual 1-form, return their interior product as a dual 0-form.
function interior_product_dd(::Type{Tuple{1,1}}, s::SimplicialSets.HasDeltaSet)
Given a dual 1-form and a dual 2-form, return their interior product as a dual 1-form.
function Δᵈ_mat(::Type{Val{0}}, s::SimplicialSets.HasDeltaSet)
Return a function matrix encoding the dual 0-form Laplacian.
function Δᵈ_mat(::Type{Val{2}}, s::SimplicialSets.HasDeltaSet)
Return a function matrix encoding the dual 1-form Laplacian.
function ℒ_dd(::Type{Tuple{1,1}}, s::SimplicialSets.HasDeltaSet)
Given a dual 1-form and a dual 1-form, return their lie derivative as a dual 1-form.
