Docstrings
Inti.Inti
— Modulemodule Inti
Library for solving integral equations using Nyström methods.
Inti.COMPRESSION_METHODS
— Constantconst COMPRESSION_METHODS = [:none, :hmatrix, :fmm]
Available compression methods for the dense linear operators in Inti
.
Inti.CORRECTION_METHODS
— Constantconst CORRECTION_METHODS = [:none, :dim, :hcubature]
Available correction methods for the singular and nearly-singular integrals in Inti
.
Inti.ENTITIES
— Constantconst ENTITIES
Dictionary mapping EntityKey
to GeometricEntity
. Contains all entities created in a given session.
Inti.SAME_POINT_TOLERANCE
— ConstantSAME_POINTS_TOLERANCE
Two points x
and y
are considerd the same if norm(x-y) ≤ SAME_POINT_TOLERANCE
.
Inti.AbstractDifferentialOperator
— Typeabstract type AbstractDifferentialOperator{N}
A partial differential operator in dimension N
.
AbstractDifferentialOperator
types are used to define AbstractKernel
s related to fundamental solutions of differential operators.
Inti.AbstractKernel
— Typeabstract type AbstractKernel{T}
A kernel functions K
with the signature K(target,source)::T
.
See also: SingleLayerKernel
, DoubleLayerKernel
, AdjointDoubleLayerKernel
, HyperSingularKernel
Inti.AbstractMesh
— Typeabstract type AbstractMesh{N,T}
An abstract mesh structure in dimension N
with primite data of type T
(e.g. Float64
for double precision representation).
Concrete subtypes of AbstractMesh
should implement ElementIterator
for accessing the mesh elements.
See also: Mesh
Inti.AdjointDoubleLayerKernel
— Typestruct AdjointDoubleLayerKernel{T,Op} <: AbstractKernel{T}
Given an operator Op
, construct its free-space adjoint double-layer kernel. This corresponds to the transpose(γ₁,ₓ[G])
, where G
is the SingleLayerKernel
. For operators such as Laplace
or Helmholtz
, this is simply the normal derivative of the fundamental solution respect to the target variable.
Inti.BlockArray
— Typestruct BlockArray{T<:StaticArray,N,S} <: AbstractMatrix{T,N}
A struct which behaves like an Array{T,N}
, but with the underlying data
stored as a Matrix{S}
, where S::Number = eltype(T)
is the scalar type associated with T
. This allows for the use of blas
routines under-the-hood, while providing a convenient interface for handling matrices over StaticArray
s.
using StaticArrays
+Docstrings · Inti.jl Docstrings
Inti.Inti
— Modulemodule Inti
Library for solving integral equations using Nyström methods.
sourceInti.COMPRESSION_METHODS
— Constantconst COMPRESSION_METHODS = [:none, :hmatrix, :fmm]
Available compression methods for the dense linear operators in Inti
.
sourceInti.CORRECTION_METHODS
— Constantconst CORRECTION_METHODS = [:none, :dim, :hcubature]
Available correction methods for the singular and nearly-singular integrals in Inti
.
sourceInti.ENTITIES
— Constantconst ENTITIES
Dictionary mapping EntityKey
to GeometricEntity
. Contains all entities created in a given session.
sourceInti.SAME_POINT_TOLERANCE
— ConstantSAME_POINTS_TOLERANCE
Two points x
and y
are considerd the same if norm(x-y) ≤ SAME_POINT_TOLERANCE
.
sourceInti.AbstractDifferentialOperator
— Typeabstract type AbstractDifferentialOperator{N}
A partial differential operator in dimension N
.
AbstractDifferentialOperator
types are used to define AbstractKernel
s related to fundamental solutions of differential operators.
sourceInti.AbstractKernel
— Typeabstract type AbstractKernel{T}
A kernel functions K
with the signature K(target,source)::T
.
See also: SingleLayerKernel
, DoubleLayerKernel
, AdjointDoubleLayerKernel
, HyperSingularKernel
sourceInti.AbstractMesh
— Typeabstract type AbstractMesh{N,T}
An abstract mesh structure in dimension N
with primite data of type T
(e.g. Float64
for double precision representation).
Concrete subtypes of AbstractMesh
should implement ElementIterator
for accessing the mesh elements.
See also: Mesh
sourceInti.AdjointDoubleLayerKernel
— Typestruct AdjointDoubleLayerKernel{T,Op} <: AbstractKernel{T}
Given an operator Op
, construct its free-space adjoint double-layer kernel. This corresponds to the transpose(γ₁,ₓ[G])
, where G
is the SingleLayerKernel
. For operators such as Laplace
or Helmholtz
, this is simply the normal derivative of the fundamental solution respect to the target variable.
sourceInti.BlockArray
— Typestruct BlockArray{T<:StaticArray,N,S} <: AbstractMatrix{T,N}
A struct which behaves like an Array{T,N}
, but with the underlying data
stored as a Matrix{S}
, where S::Number = eltype(T)
is the scalar type associated with T
. This allows for the use of blas
routines under-the-hood, while providing a convenient interface for handling matrices over StaticArray
s.
using StaticArrays
T = SMatrix{2,2,Int,4}
B = Inti.BlockArray{T}([i*j for i in 1:4, j in 1:4])
@@ -8,7 +8,7 @@
2×2 Inti.BlockArray{SMatrix{2, 2, Int64, 4}, 2, Int64}:
[1 2; 2 4] [3 4; 6 8]
[3 6; 4 8] [9 12; 12 16]
-
sourceInti.DimParameters
— Typestruct DimParameters
Parameters associated with the density interpolation method used in bdim_correction
.
sourceInti.Domain
— Typestruct Domain
Representation of a geometrical domain formed by a set of entities with the same geometric dimension. For basic set operations on domains are supported (union, intersection, difference, etc), and they all return a new Domain
object.
Calling keys(Ω)
returns the set of EntityKey
s that make up the domain; given a key, the underlying entities can be accessed with global_get_entity(key)
.
sourceInti.Domain
— MethodDomain([f::Function,] keys)
Create a domain from a set of EntityKey
s. Optionally, a filter function f
can be passed to filter the entities.
Note that all entities in a domain must have the same geometric dimension.
sourceInti.Domain
— MethodDomain(f::Function, msh::AbstractMesh)
Call Domain(f, ents)
on ents = entities(msh).
sourceInti.DoubleLayerKernel
— Typestruct DoubleLayerKernel{T,Op} <: AbstractKernel{T}
Given an operator Op
, construct its free-space double-layer kernel. This corresponds to the γ₁
trace of the SingleLayerKernel
. For operators such as Laplace
or Helmholtz
, this is simply the normal derivative of the fundamental solution respect to the source variable.
sourceInti.Elastostatic
— Typestruct Elastostatic{N,T} <: AbstractDifferentialOperator{N}
Elastostatic operator in N
dimensions: -μΔu - (μ+λ)∇(∇⋅u)
Note that the displacement $u$ is a vector of length N
since this is a vectorial problem.
sourceInti.ElementIterator
— Typestruct ElementIterator{E,M} <: AbstractVector{E}
Structure to lazily access elements of type E
in a mesh of type M
. This is particularly useful for LagrangeElement
s, where the information to reconstruct the element is stored in the mesh connectivity matrix.
sourceInti.EmbeddedQuadrature
— Typestruct EmbeddedQuadrature{L,H,D} <: ReferenceQuadrature{D}
A quadrature rule for the reference shape D
based on a high-order quadrature of type H
and a low-order quadrature of type L
. The low-order quadrature rule is embedded in the sense that its n
nodes are exactly the first n
nodes of the high-order quadrature rule.
sourceInti.EntityKey
— TypeEntityKey
Used to represent the key of a GeometricEntity
, comprised of a dim
and a tag
field, where dim
is the geometrical dimension of the entity, and tag
is a unique integer identifying the entity.
The sign of the tag
field is used to distinguish the orientation of the entity, and is ignored when comparing two EntityKey
s for equality.
sourceInti.Fejer
— Typestruct Fejer{N}
N
-point Fejer's first quadrature rule for integrating a function over [0,1]
. Exactly integrates all polynomials of degree ≤ N-1
.
using Inti
+
sourceInti.DimParameters
— Typestruct DimParameters
Parameters associated with the density interpolation method used in bdim_correction
.
sourceInti.Domain
— Typestruct Domain
Representation of a geometrical domain formed by a set of entities with the same geometric dimension. For basic set operations on domains are supported (union, intersection, difference, etc), and they all return a new Domain
object.
Calling keys(Ω)
returns the set of EntityKey
s that make up the domain; given a key, the underlying entities can be accessed with global_get_entity(key)
.
sourceInti.Domain
— MethodDomain([f::Function,] keys)
Create a domain from a set of EntityKey
s. Optionally, a filter function f
can be passed to filter the entities.
Note that all entities in a domain must have the same geometric dimension.
sourceInti.Domain
— MethodDomain(f::Function, msh::AbstractMesh)
Call Domain(f, ents)
on ents = entities(msh).
sourceInti.DoubleLayerKernel
— Typestruct DoubleLayerKernel{T,Op} <: AbstractKernel{T}
Given an operator Op
, construct its free-space double-layer kernel. This corresponds to the γ₁
trace of the SingleLayerKernel
. For operators such as Laplace
or Helmholtz
, this is simply the normal derivative of the fundamental solution respect to the source variable.
sourceInti.Elastostatic
— Typestruct Elastostatic{N,T} <: AbstractDifferentialOperator{N}
Elastostatic operator in N
dimensions: -μΔu - (μ+λ)∇(∇⋅u)
Note that the displacement $u$ is a vector of length N
since this is a vectorial problem.
sourceInti.ElementIterator
— Typestruct ElementIterator{E,M} <: AbstractVector{E}
Structure to lazily access elements of type E
in a mesh of type M
. This is particularly useful for LagrangeElement
s, where the information to reconstruct the element is stored in the mesh connectivity matrix.
sourceInti.EmbeddedQuadrature
— Typestruct EmbeddedQuadrature{L,H,D} <: ReferenceQuadrature{D}
A quadrature rule for the reference shape D
based on a high-order quadrature of type H
and a low-order quadrature of type L
. The low-order quadrature rule is embedded in the sense that its n
nodes are exactly the first n
nodes of the high-order quadrature rule.
sourceInti.EntityKey
— TypeEntityKey
Used to represent the key of a GeometricEntity
, comprised of a dim
and a tag
field, where dim
is the geometrical dimension of the entity, and tag
is a unique integer identifying the entity.
The sign of the tag
field is used to distinguish the orientation of the entity, and is ignored when comparing two EntityKey
s for equality.
sourceInti.Fejer
— Typestruct Fejer{N}
N
-point Fejer's first quadrature rule for integrating a function over [0,1]
. Exactly integrates all polynomials of degree ≤ N-1
.
using Inti
q = Inti.Fejer(;order=10)
@@ -16,7 +16,7 @@
# output
-true
sourceInti.Gauss
— Typestruct Gauss{D,N} <: ReferenceQuadrature{D}
Tabulated N
-point symmetric Gauss quadrature rule for integration over D
.
sourceInti.GaussLegendre
— Typestruct GaussLegendre{N,T}
N
-point Gauss-Legendre quadrature rule for integrating a function over [0,1]
. Exactly integrates all polynomials of degree ≤ 2N-1
.
using Inti
+true
sourceInti.Gauss
— Typestruct Gauss{D,N} <: ReferenceQuadrature{D}
Tabulated N
-point symmetric Gauss quadrature rule for integration over D
.
sourceInti.GaussLegendre
— Typestruct GaussLegendre{N,T}
N
-point Gauss-Legendre quadrature rule for integrating a function over [0,1]
. Exactly integrates all polynomials of degree ≤ 2N-1
.
using Inti
q = Inti.GaussLegendre(;order=10)
@@ -24,14 +24,14 @@
# output
-true
sourceInti.GeometricEntity
— Typestruct GeometricEntity
Geometrical objects such as lines, surfaces, and volumes.
Geometrical entities are stored in a global ENTITIES
dictionary mapping EntityKey
to the corresponding GeometricEntity
, and usually entities are manipulated through their keys.
A GeometricEntity
can also contain a pushforward
field used to parametrically represent the entry as the image of a reference domain (pushforward.domain
) under some function (pushforward.parametrization
).
Note that entities are manipulated through their keys, and the GeometricEntity
constructor returns the key of the created entity; to retrieve the entity, use the global_get_entity
function.
sourceInti.GeometricEntity
— MethodGeometricEntity(shape::String [; translation, rotation, scaling, kwargs...])
Constructs a geometric entity with the specified shape and optional parameters, and returns its key
.
Arguments
shape::String
: The shape of the geometric entity.translation
: The translation vector of the geometric entity. Default is SVector(0, 0, 0)
.rotation
: The rotation vector of the geometric entity. Default is SVector(0, 0, 0)
.scaling
: The scaling vector of the geometric entity. Default is SVector(1, 1, 1)
.kwargs...
: Additional keyword arguments to be passed to the shape constructor.
Supported shapes
sourceInti.Helmholtz
— MethodHelmholtz(; k, dim)
Helmholtz operator in dim
dimensions: -Δu - k²u
.
The parameter k
can be a real or complex number. For purely imaginary wavenumbers, consider using the Yukawa
kernel.
sourceInti.HyperRectangle
— Typestruct HyperRectangle{N,T} <: ReferenceInterpolant{ReferenceHyperCube{N},T}
Axis-aligned hyperrectangle in N
dimensions given by low_corner::SVector{N,T}
and high_corner::SVector{N,T}
.
sourceInti.HyperSingularKernel
— Typestruct HyperSingularKernel{T,Op} <: AbstractKernel{T}
Given an operator Op
, construct its free-space hypersingular kernel. This corresponds to the transpose(γ₁,ₓγ₁[G])
, where G
is the SingleLayerKernel
. For operators such as Laplace
or Helmholtz
, this is simply the normal derivative respect to the target variable of the DoubleLayerKernel
.
sourceInti.IntegralOperator
— Typestruct IntegralOperator{T} <: AbstractMatrix{T}
A discrete linear integral operator given by
\[I[u](x) = \int_{\Gamma\_s} K(x,y)u(y) ds_y, x \in \Gamma_{t}\]
where $\Gamma_s$ and $\Gamma_t$ are the source and target domains, respectively.
sourceInti.IntegralPotential
— Typestruct IntegralPotential
Represent a potential given by a kernel
and a quadrature
over which integration is performed.
IntegralPotential
s are created using IntegralPotential(kernel, quadrature)
.
Evaluating an integral potential requires a density σ
(defined over the quadrature nodes of the source mesh) and a point x
at which to evaluate the integral
\[\int_{\Gamma} K(oldsymbol{x},oldsymbol{y})\sigma(y) ds_y, x \not \in \Gamma\]
Assuming 𝒮
is an integral potential and σ
is a vector of values defined on quadrature
, calling 𝒮[σ]
creates an anonymous function that can be evaluated at any point x
.
sourceInti.Kronrod
— Typestruct Kronrod{D,N} <: ReferenceQuadrature{D}
N
-point Kronrod rule obtained by adding n+1
points to a Gauss quadrature containing n
points. The order is either 3n + 1
for n
even or 3n + 2
for n
odd.
sourceInti.LagrangeCube
— Typeconst LagrangeSquare = LagrangeElement{ReferenceSquare}
sourceInti.LagrangeElement
— Typestruct LagrangeElement{D,Np,T} <: ReferenceInterpolant{D,T}
A polynomial p : D → T
uniquely defined by its Np
values on the Np
reference nodes of D
.
The return type T
should be a vector space (i.e. support addition and multiplication by scalars). For istance, T
could be a number or a vector, but not a Tuple
.
sourceInti.LagrangeLine
— Typeconst LagrangeLine = LagrangeElement{ReferenceLine}
sourceInti.LagrangeSquare
— Typeconst LagrangeSquare = LagrangeElement{ReferenceSquare}
sourceInti.LagrangeTetrahedron
— Typeconst LagrangeTetrahedron = LagrangeElement{ReferenceTetrahedron}
sourceInti.LagrangeTriangle
— Typeconst LagrangeTriangle = LagrangeElement{ReferenceTriangle}
sourceInti.Laplace
— MethodLaplace(; dim)
Laplace's differential operator in dim
dimension: $-Δu$. ```
Note the negative sign in the definition.
sourceInti.Mesh
— Typestruct Mesh{N,T} <: AbstractMesh{N,T}
Unstructured mesh defined by a set of nodes
(of type
SVector{N,T}`), and a dictionary mapping element types to connectivity matrices. Each columns of a given connectivity matrix stores the integer tags of the nodes in the mesh comprising the element.
Additionally, the mesh contains a mapping from EntityKey
s to the tags of the elements composing the entity. This can be used to extract submeshes from a given mesh using e.g. view(msh,Γ)
or msh[Γ]
, where Γ
is a Domain
.
See elements
for a way to iterate over the elements of a mesh.
sourceInti.ModifiedHelmholtz
— Typeconst ModifiedHelmholtz
Type alias for the Yukawa
operator.
sourceInti.MultiIndex
— TypeMultiIndex{N}
Wrapper around NTuple{N,Int}
mimicking a multi-index in ℤ₀ᴺ
.
sourceInti.ParametricElement
— TypeParametricElement{D,T,F} <: ReferenceInterpolant{D,T}
An element represented through a explicit function f
mapping D
into the element. For performance reasons, f
should take as input a StaticVector
and return a StaticVector
or StaticArray
.
See also: ReferenceInterpolant
, LagrangeElement
sourceInti.ParametricElement
— MethodParametricElement(f, d::HyperRectangle)
Construct the element defined as the image of f
over d
.
sourceInti.PolynomialSpace
— Typestruct PolynomialSpace{D,K}
The space of all polynomials of degree ≤K
, commonly referred to as ℙₖ
.
The type parameter D
, of singleton type, is used to determine the reference domain of the polynomial basis. In particular, when D
is a hypercube in d
dimensions, the precise definition is ℙₖ = span{𝐱ᶿ : 0≤max(θ)≤ K}
; when D
is a d
-dimensional simplex, the space is ℙₖ = span{𝐱ᶿ : 0≤sum(θ)≤ K}
, where θ ∈ 𝐍ᵈ
is a multi-index.
See also: monomial_basis
, lagrange_basis
sourceInti.Quadrature
— Typestruct Quadrature{N,T} <: AbstractVector{QuadratureNode{N,T}}
A collection of QuadratureNode
s used to integrate over an AbstractMesh
.
sourceInti.Quadrature
— MethodQuadrature(Ω::Domain; meshsize, qorder)
Construct a Quadrature
over the domain Ω
with a mesh of size meshsize
and quadrature order qorder
.
sourceInti.Quadrature
— MethodQuadrature(msh::AbstractMesh, etype2qrule::Dict)
+true
sourceInti.GeometricEntity
— Typestruct GeometricEntity
Geometrical objects such as lines, surfaces, and volumes.
Geometrical entities are stored in a global ENTITIES
dictionary mapping EntityKey
to the corresponding GeometricEntity
, and usually entities are manipulated through their keys.
A GeometricEntity
can also contain a pushforward
field used to parametrically represent the entry as the image of a reference domain (pushforward.domain
) under some function (pushforward.parametrization
).
Note that entities are manipulated through their keys, and the GeometricEntity
constructor returns the key of the created entity; to retrieve the entity, use the global_get_entity
function.
sourceInti.GeometricEntity
— MethodGeometricEntity(shape::String [; translation, rotation, scaling, kwargs...])
Constructs a geometric entity with the specified shape and optional parameters, and returns its key
.
Arguments
shape::String
: The shape of the geometric entity.translation
: The translation vector of the geometric entity. Default is SVector(0, 0, 0)
.rotation
: The rotation vector of the geometric entity. Default is SVector(0, 0, 0)
.scaling
: The scaling vector of the geometric entity. Default is SVector(1, 1, 1)
.kwargs...
: Additional keyword arguments to be passed to the shape constructor.
Supported shapes
sourceInti.Helmholtz
— MethodHelmholtz(; k, dim)
Helmholtz operator in dim
dimensions: -Δu - k²u
.
The parameter k
can be a real or complex number. For purely imaginary wavenumbers, consider using the Yukawa
kernel.
sourceInti.HyperRectangle
— Typestruct HyperRectangle{N,T} <: ReferenceInterpolant{ReferenceHyperCube{N},T}
Axis-aligned hyperrectangle in N
dimensions given by low_corner::SVector{N,T}
and high_corner::SVector{N,T}
.
sourceInti.HyperSingularKernel
— Typestruct HyperSingularKernel{T,Op} <: AbstractKernel{T}
Given an operator Op
, construct its free-space hypersingular kernel. This corresponds to the transpose(γ₁,ₓγ₁[G])
, where G
is the SingleLayerKernel
. For operators such as Laplace
or Helmholtz
, this is simply the normal derivative respect to the target variable of the DoubleLayerKernel
.
sourceInti.IntegralOperator
— Typestruct IntegralOperator{T} <: AbstractMatrix{T}
A discrete linear integral operator given by
\[I[u](x) = \int_{\Gamma\_s} K(x,y)u(y) ds_y, x \in \Gamma_{t}\]
where $\Gamma_s$ and $\Gamma_t$ are the source and target domains, respectively.
sourceInti.IntegralPotential
— Typestruct IntegralPotential
Represent a potential given by a kernel
and a quadrature
over which integration is performed.
IntegralPotential
s are created using IntegralPotential(kernel, quadrature)
.
Evaluating an integral potential requires a density σ
(defined over the quadrature nodes of the source mesh) and a point x
at which to evaluate the integral
\[\int_{\Gamma} K(oldsymbol{x},oldsymbol{y})\sigma(y) ds_y, x \not \in \Gamma\]
Assuming 𝒮
is an integral potential and σ
is a vector of values defined on quadrature
, calling 𝒮[σ]
creates an anonymous function that can be evaluated at any point x
.
sourceInti.Kronrod
— Typestruct Kronrod{D,N} <: ReferenceQuadrature{D}
N
-point Kronrod rule obtained by adding n+1
points to a Gauss quadrature containing n
points. The order is either 3n + 1
for n
even or 3n + 2
for n
odd.
sourceInti.LagrangeCube
— Typeconst LagrangeSquare = LagrangeElement{ReferenceSquare}
sourceInti.LagrangeElement
— Typestruct LagrangeElement{D,Np,T} <: ReferenceInterpolant{D,T}
A polynomial p : D → T
uniquely defined by its Np
values on the Np
reference nodes of D
.
The return type T
should be a vector space (i.e. support addition and multiplication by scalars). For istance, T
could be a number or a vector, but not a Tuple
.
sourceInti.LagrangeLine
— Typeconst LagrangeLine = LagrangeElement{ReferenceLine}
sourceInti.LagrangeSquare
— Typeconst LagrangeSquare = LagrangeElement{ReferenceSquare}
sourceInti.LagrangeTetrahedron
— Typeconst LagrangeTetrahedron = LagrangeElement{ReferenceTetrahedron}
sourceInti.LagrangeTriangle
— Typeconst LagrangeTriangle = LagrangeElement{ReferenceTriangle}
sourceInti.Laplace
— MethodLaplace(; dim)
Laplace's differential operator in dim
dimension: $-Δu$. ```
Note the negative sign in the definition.
sourceInti.Mesh
— Typestruct Mesh{N,T} <: AbstractMesh{N,T}
Unstructured mesh defined by a set of nodes
(of type
SVector{N,T}`), and a dictionary mapping element types to connectivity matrices. Each columns of a given connectivity matrix stores the integer tags of the nodes in the mesh comprising the element.
Additionally, the mesh contains a mapping from EntityKey
s to the tags of the elements composing the entity. This can be used to extract submeshes from a given mesh using e.g. view(msh,Γ)
or msh[Γ]
, where Γ
is a Domain
.
See elements
for a way to iterate over the elements of a mesh.
sourceInti.ModifiedHelmholtz
— Typeconst ModifiedHelmholtz
Type alias for the Yukawa
operator.
sourceInti.MultiIndex
— TypeMultiIndex{N}
Wrapper around NTuple{N,Int}
mimicking a multi-index in ℤ₀ᴺ
.
sourceInti.ParametricElement
— TypeParametricElement{D,T,F} <: ReferenceInterpolant{D,T}
An element represented through a explicit function f
mapping D
into the element. For performance reasons, f
should take as input a StaticVector
and return a StaticVector
or StaticArray
.
See also: ReferenceInterpolant
, LagrangeElement
sourceInti.ParametricElement
— MethodParametricElement(f, d::HyperRectangle)
Construct the element defined as the image of f
over d
.
sourceInti.PolynomialSpace
— Typestruct PolynomialSpace{D,K}
The space of all polynomials of degree ≤K
, commonly referred to as ℙₖ
.
The type parameter D
, of singleton type, is used to determine the reference domain of the polynomial basis. In particular, when D
is a hypercube in d
dimensions, the precise definition is ℙₖ = span{𝐱ᶿ : 0≤max(θ)≤ K}
; when D
is a d
-dimensional simplex, the space is ℙₖ = span{𝐱ᶿ : 0≤sum(θ)≤ K}
, where θ ∈ 𝐍ᵈ
is a multi-index.
See also: monomial_basis
, lagrange_basis
sourceInti.Quadrature
— Typestruct Quadrature{N,T} <: AbstractVector{QuadratureNode{N,T}}
A collection of QuadratureNode
s used to integrate over an AbstractMesh
.
sourceInti.Quadrature
— MethodQuadrature(Ω::Domain; meshsize, qorder)
Construct a Quadrature
over the domain Ω
with a mesh of size meshsize
and quadrature order qorder
.
sourceInti.Quadrature
— MethodQuadrature(msh::AbstractMesh, etype2qrule::Dict)
Quadrature(msh::AbstractMesh, qrule::ReferenceQuadrature)
-Quadrature(msh::AbstractMesh; qorder)
Construct a Quadrature
for msh
, where for each element type E
in msh
the reference quadrature q = etype2qrule[E]
is used. When a single qrule
is passed, it is used for all element types in msh
.
If an order
keyword is passed, a default quadrature of the desired order is used for each element type usig _qrule_for_reference_shape
.
For co-dimension one elements, the normal vector is also computed and stored in the QuadratureNode
s.
sourceInti.QuadratureNode
— TypeQuadratureNode{N,T<:Real}
A point in ℝᴺ
with a weight
for performing numerical integration. A QuadratureNode
can optionally store a normal
vector.
sourceInti.ReferenceCube
— Typeconst ReferenceCube = ReferenceHyperCube{3}
Singleton type representing the unit cube [0,1]³
.
sourceInti.ReferenceHyperCube
— Typestruct ReferenceHyperCube{N} <: ReferenceShape{N}
Singleton type representing the axis-aligned hypercube in N
dimensions with the lower corner at the origin and the upper corner at (1,1,…,1)
.
sourceInti.ReferenceInterpolant
— Typeabstract type ReferenceInterpolant{D,T}
Interpolanting function mapping points on the domain D<:ReferenceShape
(of singleton type) to a value of type T
.
Instances el
of ReferenceInterpolant
are expected to implement:
el(x̂)
: evaluate the interpolation scheme at the (reference) coordinate x̂ ∈ D
.jacobian(el,x̂)
: evaluate the jacobian matrix of the interpolation at the (reference) coordinate x ∈ D
.
Note For performance reasons, both el(x̂)
and jacobian(el,x̂)
should take as input a StaticVector
and output a static vector or static array.
sourceInti.ReferenceLine
— Typeconst ReferenceLine = ReferenceHyperCube{1}
Singleton type representing the [0,1]
segment.
sourceInti.ReferenceQuadrature
— Typeabstract type ReferenceQuadrature{D}
A quadrature rule for integrating a function over the domain D <: ReferenceShape
.
Calling x,w = q()
returns the nodes x
, given as SVector
s, and weights w
, for performing integration over domain(q)
.
sourceInti.ReferenceShape
— Typeabstract type ReferenceShape
A fixed reference domain/shape. Used mostly for defining more complex shapes as transformations mapping an ReferenceShape
to some region of ℜᴹ
.
See e.g. ReferenceLine
or ReferenceTriangle
for some examples of concrete subtypes.
sourceInti.ReferenceSimplex
— Typestruct ReferenceSimplex{N}
Singleton type representing the N-simplex with N+1 vertices (0,...,0),(0,...,0,1),(0,...,0,1,0),(1,0,...,0)
sourceInti.ReferenceSquare
— Typeconst ReferenceSquare = ReferenceHyperCube{2}
Singleton type representing the unit square [0,1]²
.
sourceInti.ReferenceTetrahedron
— Typestruct ReferenceTetrahedron
Singleton type representing the tetrahedron with vertices (0,0,0),(0,0,1),(0,1,0),(1,0,0)
sourceInti.ReferenceTriangle
— Typestruct ReferenceTriangle
Singleton type representing the triangle with vertices (0,0),(1,0),(0,1)
sourceInti.SingleLayerKernel
— Typestruct SingleLayerKernel{T,Op} <: AbstractKernel{T}
The free-space single-layer kernel (i.e. the fundamental solution) of an Op <: AbstractDifferentialOperator
.
sourceInti.Stokes
— MethodStokes(; μ, dim)
Stokes operator in dim
dimensions: $[-μΔu + ∇p, ∇⋅u]$.
sourceInti.SubMesh
— Typestruct SubMesh{N,T} <: AbstractMesh{N,T}
View into a parent
mesh over a given domain
.
A submesh implements the interface for AbstractMesh
; therefore you can iterate over elements of the submesh just like you would with a mesh.
Construct SubMesh
s using view(parent,Ω::Domain)
.
sourceInti.TensorProductQuadrature
— TypeTensorProductQuadrature{N,Q}
A tensor-product of one-dimension quadrature rules. Integrates over [0,1]^N
.
Examples
qx = Inti.Fejer(10)
+Quadrature(msh::AbstractMesh; qorder)
Construct a Quadrature
for msh
, where for each element type E
in msh
the reference quadrature q = etype2qrule[E]
is used. When a single qrule
is passed, it is used for all element types in msh
.
If an order
keyword is passed, a default quadrature of the desired order is used for each element type usig _qrule_for_reference_shape
.
For co-dimension one elements, the normal vector is also computed and stored in the QuadratureNode
s.
sourceInti.QuadratureNode
— TypeQuadratureNode{N,T<:Real}
A point in ℝᴺ
with a weight
for performing numerical integration. A QuadratureNode
can optionally store a normal
vector.
sourceInti.ReferenceCube
— Typeconst ReferenceCube = ReferenceHyperCube{3}
Singleton type representing the unit cube [0,1]³
.
sourceInti.ReferenceHyperCube
— Typestruct ReferenceHyperCube{N} <: ReferenceShape{N}
Singleton type representing the axis-aligned hypercube in N
dimensions with the lower corner at the origin and the upper corner at (1,1,…,1)
.
sourceInti.ReferenceInterpolant
— Typeabstract type ReferenceInterpolant{D,T}
Interpolanting function mapping points on the domain D<:ReferenceShape
(of singleton type) to a value of type T
.
Instances el
of ReferenceInterpolant
are expected to implement:
el(x̂)
: evaluate the interpolation scheme at the (reference) coordinate x̂ ∈ D
.jacobian(el,x̂)
: evaluate the jacobian matrix of the interpolation at the (reference) coordinate x ∈ D
.
Note For performance reasons, both el(x̂)
and jacobian(el,x̂)
should take as input a StaticVector
and output a static vector or static array.
sourceInti.ReferenceLine
— Typeconst ReferenceLine = ReferenceHyperCube{1}
Singleton type representing the [0,1]
segment.
sourceInti.ReferenceQuadrature
— Typeabstract type ReferenceQuadrature{D}
A quadrature rule for integrating a function over the domain D <: ReferenceShape
.
Calling x,w = q()
returns the nodes x
, given as SVector
s, and weights w
, for performing integration over domain(q)
.
sourceInti.ReferenceShape
— Typeabstract type ReferenceShape
A fixed reference domain/shape. Used mostly for defining more complex shapes as transformations mapping an ReferenceShape
to some region of ℜᴹ
.
See e.g. ReferenceLine
or ReferenceTriangle
for some examples of concrete subtypes.
sourceInti.ReferenceSimplex
— Typestruct ReferenceSimplex{N}
Singleton type representing the N-simplex with N+1 vertices (0,...,0),(0,...,0,1),(0,...,0,1,0),(1,0,...,0)
sourceInti.ReferenceSquare
— Typeconst ReferenceSquare = ReferenceHyperCube{2}
Singleton type representing the unit square [0,1]²
.
sourceInti.ReferenceTetrahedron
— Typestruct ReferenceTetrahedron
Singleton type representing the tetrahedron with vertices (0,0,0),(0,0,1),(0,1,0),(1,0,0)
sourceInti.ReferenceTriangle
— Typestruct ReferenceTriangle
Singleton type representing the triangle with vertices (0,0),(1,0),(0,1)
sourceInti.SingleLayerKernel
— Typestruct SingleLayerKernel{T,Op} <: AbstractKernel{T}
The free-space single-layer kernel (i.e. the fundamental solution) of an Op <: AbstractDifferentialOperator
.
sourceInti.Stokes
— MethodStokes(; μ, dim)
Stokes operator in dim
dimensions: $[-μΔu + ∇p, ∇⋅u]$.
sourceInti.SubMesh
— Typestruct SubMesh{N,T} <: AbstractMesh{N,T}
View into a parent
mesh over a given domain
.
A submesh implements the interface for AbstractMesh
; therefore you can iterate over elements of the submesh just like you would with a mesh.
Construct SubMesh
s using view(parent,Ω::Domain)
.
sourceInti.TensorProductQuadrature
— TypeTensorProductQuadrature{N,Q}
A tensor-product of one-dimension quadrature rules. Integrates over [0,1]^N
.
Examples
qx = Inti.Fejer(10)
qy = Inti.Fejer(15)
-q = Inti.TensorProductQuadrature(qx,qy)
sourceInti.VioreanuRokhlin
— Typestruct VioreanuRokhlin{D,N} <: ReferenceQuadrature{D}
Tabulated N
-point Vioreanu-Rokhlin quadrature rule for integration over D
.
sourceInti.Yukawa
— MethodYukawa(; λ, dim)
Yukawa operator, also known as modified Helmholtz, in dim
dimensions: $-Δu + λ²u$.
The parameter λ
is a positive number. Note the negative sign in front of the Laplacian.
sourceBase.iterate
— Functioniterate(Ω::Domain)
Iterating over a domain means iterating over its entities.
sourceInti._copyto!
— Method_copyto!(target,source)
Defaults to Base.copyto!
, but includes some specialized methods to copy from a Matrix
of SMatrix
to a Matrix
of Number
s and viceversa.
sourceInti._get_gauss_qcoords_and_qweights
— Method_get_gauss_and_qweights(R::Type{<:ReferenceShape{D}}, N) where D
Returns the N
-point symmetric gaussian qnodes and qweights (x, w)
for integration over R
.
sourceInti._get_vioreanurokhlin_qcoords_and_qweights
— Method_get_vioreanurokhlin_qcoords_and_qweights(R::Type{<:ReferenceShape{D}}, N) where D
Returns the N
-point Vioreanu-Rokhlin qnodes and qweights (x, w)
for integration over R
.
sourceInti._green_multiplier
— Method_green_multiplier(s::Symbol)
Return -1.0
if s == :inside
, 0.0
if s == :outside
, and -0.5
if s == :on
; otherwise, throw an error. The orientation is relative to the normal of the bounding curve/surface.
sourceInti._green_multiplier
— Method_green_multiplier(x, quad)
Helper function to help determine the constant σ in the Green identity S[γ₁u](x)
- D[γ₀u](x) + σ*u(x) = 0. This can be used as a predicate to determine whether a
point is inside a domain or not.
sourceInti._meshgen
— Method_meshgen(f,d::HyperRectangle,sz)
Create prod(sz)
elements of ParametricElement
type representing the push forward of f
on each of the subdomains defined by a uniform cartesian mesh of d
of size sz
.
sourceInti._normal
— Method_normal(jac::SMatrix{M,N})
Given a an M
by N
matrix representing the jacobian of a codimension one object, compute the normal vector.
sourceInti._qrule_for_reference_shape
— Method_qrule_for_reference_shape(ref,order)
Given a ref
erence shape and a desired quadrature order
, return an appropiate quadrature rule.
sourceInti.acorn
— Methodacorn(; translation, rotation, scaling, labels)
Create an acorn entity in 3D, and apply optional transformations. Returns the key.
sourceInti.adaptive_correction
— Methodadaptive_correction(iop::IntegralOperator; tol, maxdist = farfield_distance(iop; atol), maxsplit = 1000])
Given an integral operator iop
, this function provides a sparse correction to iop
for the entries i,j
such that the distance between the i
-th target and the j
-th source is less than maxdist
.
Choosing maxdist
is a trade-off between accuracy and efficiency. The smaller the value, the fewer corrections are needed, but this may compromise the accuracy. For a fixed quadrature, the size of maxdist
has to grow as the tolerance tol
decreases. The default [farfield_distance(iop; tol)](@ref) provides a heuristic to determine a suitable
maxdist`.
The correction is computed by using the adaptive_integration
routine, with a tolerance atol
and a maximum number of subdivisions maxsplit
; see adaptive_integration
for more details.
sourceInti.adaptive_integration
— Methodadaptive_integration(f, τ̂::RefernceShape; kwargs...)
-adaptive_integration(f, qrule::EmbeddedQuadrature; kwargs...)
Use an adaptive procedure to estimate the integral of f
over τ̂ = domain(qrule)
. The following optional keyword arguments are available:
atol::Real=0.0
: absolute tolerance for the integral estimatertol::Real=0.0
: relative tolerance for the integral estimatemaxsplit::Int=1000
: maximum number of times to split the domainnorm::Function=LinearAlgebra.norm
: norm to use for error estimatesbuffer::BinaryHeap
: a pre-allocated buffer to use for the adaptive procedure (see allocate_buffer)
sourceInti.adaptive_integration_singular
— Methodadaptive_integration_singular(f, τ̂, x̂ₛ; kwargs...)
Similar to adaptive_integration
, but indicates that f
has an isolated (integrable) singularity at x̂ₛ ∈ x̂ₛ
.
The integration is performed by splitting τ̂
so that x̂ₛ
is a fixed vertex, guaranteeing that f
is never evaluated at x̂ₛ
. Aditionally, a suitable change of variables may be applied to alleviate the singularity and improve the rate of convergence.
sourceInti.adj_double_layer_hypersingular
— Methodadj_double_layer_hypersingular(; op, target, source, compression,
-correction)
Similar to single_double_layer
, but for the adjoint double-layer and hypersingular operators. See the documentation of [single_double_layer
] for a description of the arguments.
sourceInti.allocate_buffer
— Methodallocate_buffer(f, quad::EmbeddedQuadrature)
Create the buffer
needed for the call adaptive_integration(f, τ̂; buffer, ...)
.
sourceInti.ambient_dimension
— Functionambient_dimension(x)
Dimension of the ambient space where x
lives. For geometrical objects this can differ from its geometric_dimension
; for example a triangle in ℝ³
has ambient dimension 3
but geometric dimension 2
, while a curve in ℝ³
has ambient dimension 3 but geometric dimension 1.
sourceInti.assemble_fmm
— Methodassemble_fmm(iop; atol)
Set up a 2D or 3D FMM for evaluating the discretized integral operator iop
associated with the op
. In 2D the FMM2D
or FMMLIB2D
library is used (whichever was most recently loaded) while in 3D FMM3D
is used.
FMMLIB2D FMMLIB2D does no checking for if the targets and sources coincide, and will return Inf
values if iop.target !== iop.source
, but there is a point x ∈ iop.target
such that x ∈ iop.source
.
sourceInti.assemble_hmatrix
— Methodassemble_hmatrix(iop[; atol, rank, rtol, eta])
Assemble an H-matrix representation of the discretized integral operator iop
using the HMatrices.jl
library.
See the documentation of HMatrices
for more details on usage and other keyword arguments.
sourceInti.assemble_matrix
— Methodassemble_matrix(iop::IntegralOperator; threads = true)
Assemble a dense matrix representation of an IntegralOperator
.
sourceInti.bdim_correction
— Methodbdim_correction(op,X,Y,S,D; green_multiplier, kwargs...)
Given a op
and a (possibly innacurate) discretizations of its single and double-layer operators S
and D
(taking a vector of values on Y
and returning a vector on of values on X
), compute corrections δS
and δD
such that S + δS
and D + δD
are more accurate approximations of the underlying single- and double-layer integral operators.
See [7] for more details on the method.
Arguments
Required:
op
must be an AbstractDifferentialOperator
Y
must be a Quadrature
object of a closed surfaceX
is either inside, outside, or on Y
S
and D
are approximations to the single- and double-layer operators for op
taking densities in Y
and returning densities in X
.green_multiplier
(keyword argument) is a vector with the same length as X
storing the value of μ(x)
for x ∈ X
in the Green identity S\[γ₁u\](x) - D\[γ₀u\](x) + μ*u(x) = 0
. See _green_multiplier
.
Optional kwargs
:
parameters::DimParameters
: parameters associated with the density interpolation methodderivative
: if true, compute the correction to the adjoint double-layer and hypersingular operators instead. In this case, S
and D
should be replaced by a (possibly innacurate) discretization of adjoint double-layer and hypersingular operators, respectively.maxdist
: distance beyond which interactions are considered sufficiently far so that no correction is needed. This is used to determine a threshold for nearly-singular corrections when X
and Y
are different surfaces. When X === Y
, this is not needed.
sourceInti.bean
— Methodbean(; translation, rotation, scaling, labels)
Create a bean entity in 3D, and apply optional transformations. Returns the key.
sourceInti.boundary
— Methodboundary(Ω::Domain)
Return the external boundaries of a domain.
See also: external_boundary
, internal_boundary
, skeleton
.
sourceInti.boundary_idxs
— Methodboundary_idxs(el::LagrangeElement)
The indices of the nodes in el
that define the boundary of the element.
sourceInti.cart2sph
— Methodcart2sph(x,y,z)
Map cartesian coordinates x,y,z
to spherical ones r, θ, φ
representing the radius, elevation, and azimuthal angle respectively. The convention followed is that 0 ≤ θ ≤ π
and -π < φ ≤ π
. Same as the cart2sph
function in MATLAB.
sourceInti.connectivity
— Methodconnectivity(msh::AbstractMesh,E::DataType)
Return the connectivity matrix for elements of type E
in msh
. The integer tags in the matrix refer to the points in nodes(msh)
sourceInti.coords
— Methodcoords(q)
Return the spatial coordinates of q
.
sourceInti.cushion
— Methodcushion(; translation, rotation, scaling, labels)
Create a cushion entity in 3D, and apply optional transformations. Returns the key.
sourceInti.decompose
— Functiondecompose(s::ReferenceShape,x)
Decompose an ReferenceShape
into LagrangeElement
s so that x
is a fixed vertex of the children elements.
The decomposed elements may be oriented differently than the parent, and thus care has to be taken regarding e.g. normal vectors.
sourceInti.degree
— Methoddegree(el::LagrangeElement)
-degree(el::Type{<:LagrangeElement})
The polynomial degree el
.
sourceInti.dimension
— Methoddimension(space)
The length of a basis for space
; i.e. the number of linearly independent elements required to span space
.
sourceInti.dom2elt
— Methoddom2elt(m::Mesh,Ω,E)::Vector{Int}
Compute the element indices idxs
of the elements of type E
composing Ω
.
sourceInti.dom2qtags
— Methoddom2qtags(Q::Quadrature, dom::Domain)
Given a domain, return the indices of the quadratures nodes in Q
associated to its quadrature.
sourceInti.domain
— Functiondomain(f)
Given a function-like object f: Ω → R
, return Ω
.
sourceInti.domain
— Methoddomain(msh::AbstractMesh)
Return a [Domain
] containing of all entities covered by the mesh.
sourceInti.domain
— Methoddomain(Q::Quadrature)
The Domain
over which Q
performs integration.
sourceInti.domain
— Methoddomain(q::ReferenceQuadrature)
The domain of integratino for quadrature rule q
.
sourceInti.element_types
— Functionelement_types(msh::AbstractMesh)
Return the element types present in the msh
.
sourceInti.elements
— Functionelements(msh::AbstractMesh [, E::DataType])
Return the elements of a msh
. Passing and element type E
will restricts to elements of that type.
A common pattern to avoid type-instabilies in performance critical parts of the code is to use a function barrier, as illustrated below:
for E in element_types(msh)
+q = Inti.TensorProductQuadrature(qx,qy)
sourceInti.VioreanuRokhlin
— Typestruct VioreanuRokhlin{D,N} <: ReferenceQuadrature{D}
Tabulated N
-point Vioreanu-Rokhlin quadrature rule for integration over D
.
sourceInti.Yukawa
— MethodYukawa(; λ, dim)
Yukawa operator, also known as modified Helmholtz, in dim
dimensions: $-Δu + λ²u$.
The parameter λ
is a positive number. Note the negative sign in front of the Laplacian.
sourceBase.iterate
— Functioniterate(Ω::Domain)
Iterating over a domain means iterating over its entities.
sourceInti._copyto!
— Method_copyto!(target,source)
Defaults to Base.copyto!
, but includes some specialized methods to copy from a Matrix
of SMatrix
to a Matrix
of Number
s and viceversa.
sourceInti._get_gauss_qcoords_and_qweights
— Method_get_gauss_and_qweights(R::Type{<:ReferenceShape{D}}, N) where D
Returns the N
-point symmetric gaussian qnodes and qweights (x, w)
for integration over R
.
sourceInti._get_vioreanurokhlin_qcoords_and_qweights
— Method_get_vioreanurokhlin_qcoords_and_qweights(R::Type{<:ReferenceShape{D}}, N) where D
Returns the N
-point Vioreanu-Rokhlin qnodes and qweights (x, w)
for integration over R
.
sourceInti._green_multiplier
— Method_green_multiplier(s::Symbol)
Return -1.0
if s == :inside
, 0.0
if s == :outside
, and -0.5
if s == :on
; otherwise, throw an error. The orientation is relative to the normal of the bounding curve/surface.
sourceInti._green_multiplier
— Method_green_multiplier(x, quad)
Helper function to help determine the constant σ in the Green identity S[γ₁u](x)
- D[γ₀u](x) + σ*u(x) = 0. This can be used as a predicate to determine whether a
point is inside a domain or not.
sourceInti._meshgen
— Method_meshgen(f,d::HyperRectangle,sz)
Create prod(sz)
elements of ParametricElement
type representing the push forward of f
on each of the subdomains defined by a uniform cartesian mesh of d
of size sz
.
sourceInti._normal
— Method_normal(jac::SMatrix{M,N})
Given a an M
by N
matrix representing the jacobian of a codimension one object, compute the normal vector.
sourceInti._qrule_for_reference_shape
— Method_qrule_for_reference_shape(ref,order)
Given a ref
erence shape and a desired quadrature order
, return an appropiate quadrature rule.
sourceInti.acorn
— Methodacorn(; translation, rotation, scaling, labels)
Create an acorn entity in 3D, and apply optional transformations. Returns the key.
sourceInti.adaptive_correction
— Methodadaptive_correction(iop::IntegralOperator; tol, maxdist = farfield_distance(iop; atol), maxsplit = 1000])
Given an integral operator iop
, this function provides a sparse correction to iop
for the entries i,j
such that the distance between the i
-th target and the j
-th source is less than maxdist
.
Choosing maxdist
is a trade-off between accuracy and efficiency. The smaller the value, the fewer corrections are needed, but this may compromise the accuracy. For a fixed quadrature, the size of maxdist
has to grow as the tolerance tol
decreases. The default [farfield_distance(iop; tol)](@ref) provides a heuristic to determine a suitable
maxdist`.
The correction is computed by using the adaptive_integration
routine, with a tolerance atol
and a maximum number of subdivisions maxsplit
; see adaptive_integration
for more details.
sourceInti.adaptive_integration
— Methodadaptive_integration(f, τ̂::RefernceShape; kwargs...)
+adaptive_integration(f, qrule::EmbeddedQuadrature; kwargs...)
Use an adaptive procedure to estimate the integral of f
over τ̂ = domain(qrule)
. The following optional keyword arguments are available:
atol::Real=0.0
: absolute tolerance for the integral estimatertol::Real=0.0
: relative tolerance for the integral estimatemaxsplit::Int=1000
: maximum number of times to split the domainnorm::Function=LinearAlgebra.norm
: norm to use for error estimatesbuffer::BinaryHeap
: a pre-allocated buffer to use for the adaptive procedure (see allocate_buffer)
sourceInti.adaptive_integration_singular
— Methodadaptive_integration_singular(f, τ̂, x̂ₛ; kwargs...)
Similar to adaptive_integration
, but indicates that f
has an isolated (integrable) singularity at x̂ₛ ∈ x̂ₛ
.
The integration is performed by splitting τ̂
so that x̂ₛ
is a fixed vertex, guaranteeing that f
is never evaluated at x̂ₛ
. Aditionally, a suitable change of variables may be applied to alleviate the singularity and improve the rate of convergence.
sourceInti.adj_double_layer_hypersingular
— Methodadj_double_layer_hypersingular(; op, target, source, compression,
+correction)
Similar to single_double_layer
, but for the adjoint double-layer and hypersingular operators. See the documentation of [single_double_layer
] for a description of the arguments.
sourceInti.allocate_buffer
— Methodallocate_buffer(f, quad::EmbeddedQuadrature)
Create the buffer
needed for the call adaptive_integration(f, τ̂; buffer, ...)
.
sourceInti.ambient_dimension
— Functionambient_dimension(x)
Dimension of the ambient space where x
lives. For geometrical objects this can differ from its geometric_dimension
; for example a triangle in ℝ³
has ambient dimension 3
but geometric dimension 2
, while a curve in ℝ³
has ambient dimension 3 but geometric dimension 1.
sourceInti.assemble_fmm
— Methodassemble_fmm(iop; atol)
Set up a 2D or 3D FMM for evaluating the discretized integral operator iop
associated with the op
. In 2D the FMM2D
or FMMLIB2D
library is used (whichever was most recently loaded) while in 3D FMM3D
is used.
FMMLIB2D FMMLIB2D does no checking for if the targets and sources coincide, and will return Inf
values if iop.target !== iop.source
, but there is a point x ∈ iop.target
such that x ∈ iop.source
.
sourceInti.assemble_hmatrix
— Methodassemble_hmatrix(iop[; atol, rank, rtol, eta])
Assemble an H-matrix representation of the discretized integral operator iop
using the HMatrices.jl
library.
See the documentation of HMatrices
for more details on usage and other keyword arguments.
sourceInti.assemble_matrix
— Methodassemble_matrix(iop::IntegralOperator; threads = true)
Assemble a dense matrix representation of an IntegralOperator
.
sourceInti.bdim_correction
— Methodbdim_correction(op,X,Y,S,D; green_multiplier, kwargs...)
Given a op
and a (possibly innacurate) discretizations of its single and double-layer operators S
and D
(taking a vector of values on Y
and returning a vector on of values on X
), compute corrections δS
and δD
such that S + δS
and D + δD
are more accurate approximations of the underlying single- and double-layer integral operators.
See [7] for more details on the method.
Arguments
Required:
op
must be an AbstractDifferentialOperator
Y
must be a Quadrature
object of a closed surfaceX
is either inside, outside, or on Y
S
and D
are approximations to the single- and double-layer operators for op
taking densities in Y
and returning densities in X
.green_multiplier
(keyword argument) is a vector with the same length as X
storing the value of μ(x)
for x ∈ X
in the Green identity S\[γ₁u\](x) - D\[γ₀u\](x) + μ*u(x) = 0
. See _green_multiplier
.
Optional kwargs
:
parameters::DimParameters
: parameters associated with the density interpolation methodderivative
: if true, compute the correction to the adjoint double-layer and hypersingular operators instead. In this case, S
and D
should be replaced by a (possibly innacurate) discretization of adjoint double-layer and hypersingular operators, respectively.maxdist
: distance beyond which interactions are considered sufficiently far so that no correction is needed. This is used to determine a threshold for nearly-singular corrections when X
and Y
are different surfaces. When X === Y
, this is not needed.
sourceInti.bean
— Methodbean(; translation, rotation, scaling, labels)
Create a bean entity in 3D, and apply optional transformations. Returns the key.
sourceInti.boundary
— Methodboundary(Ω::Domain)
Return the external boundaries of a domain.
See also: external_boundary
, internal_boundary
, skeleton
.
sourceInti.boundary_idxs
— Methodboundary_idxs(el::LagrangeElement)
The indices of the nodes in el
that define the boundary of the element.
sourceInti.cart2sph
— Methodcart2sph(x,y,z)
Map cartesian coordinates x,y,z
to spherical ones r, θ, φ
representing the radius, elevation, and azimuthal angle respectively. The convention followed is that 0 ≤ θ ≤ π
and -π < φ ≤ π
. Same as the cart2sph
function in MATLAB.
sourceInti.connectivity
— Methodconnectivity(msh::AbstractMesh,E::DataType)
Return the connectivity matrix for elements of type E
in msh
. The integer tags in the matrix refer to the points in nodes(msh)
sourceInti.coords
— Methodcoords(q)
Return the spatial coordinates of q
.
sourceInti.cushion
— Methodcushion(; translation, rotation, scaling, labels)
Create a cushion entity in 3D, and apply optional transformations. Returns the key.
sourceInti.decompose
— Functiondecompose(s::ReferenceShape,x)
Decompose an ReferenceShape
into LagrangeElement
s so that x
is a fixed vertex of the children elements.
The decomposed elements may be oriented differently than the parent, and thus care has to be taken regarding e.g. normal vectors.
sourceInti.degree
— Methoddegree(el::LagrangeElement)
+degree(el::Type{<:LagrangeElement})
The polynomial degree el
.
sourceInti.dimension
— Methoddimension(space)
The length of a basis for space
; i.e. the number of linearly independent elements required to span space
.
sourceInti.dom2elt
— Methoddom2elt(m::Mesh,Ω,E)::Vector{Int}
Compute the element indices idxs
of the elements of type E
composing Ω
.
sourceInti.dom2qtags
— Methoddom2qtags(Q::Quadrature, dom::Domain)
Given a domain, return the indices of the quadratures nodes in Q
associated to its quadrature.
sourceInti.domain
— Functiondomain(f)
Given a function-like object f: Ω → R
, return Ω
.
sourceInti.domain
— Methoddomain(msh::AbstractMesh)
Return a [Domain
] containing of all entities covered by the mesh.
sourceInti.domain
— Methoddomain(Q::Quadrature)
The Domain
over which Q
performs integration.
sourceInti.domain
— Methoddomain(q::ReferenceQuadrature)
The domain of integratino for quadrature rule q
.
sourceInti.element_types
— Functionelement_types(msh::AbstractMesh)
Return the element types present in the msh
.
sourceInti.elements
— Functionelements(msh::AbstractMesh [, E::DataType])
Return the elements of a msh
. Passing and element type E
will restricts to elements of that type.
A common pattern to avoid type-instabilies in performance critical parts of the code is to use a function barrier, as illustrated below:
for E in element_types(msh)
_long_computation(elements(msh, E), args...)
end
@@ -39,11 +39,11 @@
for el in iter # the type of el is known at compile time
# do something with el
end
-end
where a dynamic dispatch is performed only on the element types (typically small for a given mesh).
sourceInti.ellipsoid
— Methodellipsoid(; translation, rotation, scaling, labels)
Create an ellipsoid entity in 3D, and apply optional transformations. Returns the key of the created entity.
sourceInti.ent2etags
— Methodent2etags(msh::AbstractMesh)
Return a dictionary mapping entities to a dictionary of element types to element tags.
sourceInti.entities
— Methodentities(Ω::Domain)
Return all entities making up a domain (as a set of EntityKey
s).
sourceInti.etype_to_nearest_points
— Methodetype_to_nearest_points(X,Y::Quadrature; maxdist)
For each element el
in Y.mesh
, return a list with the indices of all points in X
for which el
is the nearest element. Ignore indices for which the distance exceeds maxdist
.
sourceInti.external_boundary
— Methodexternal_boundary(Ω::Domain)
Return the external boundaries inside a domain. These are entities in the skeleton of Ω which are not in the internal boundaries of Ω.
See also: internal_boundary
, skeleton
.
sourceInti.farfield_distance
— Methodfarfield_distance(iop::IntegralOperator; tol, maxiter = 10)
-farfield_distance(K, Q::Quadrature; tol, maxiter = 10)
Return an estimate of the distance d
such that the (absolute) quadrature error of the integrand y -> K(x,y)
is below tol
for x
at a distance d
from the center of the largest element in Q
; when an integral operator is passed, we have Q::Quadrature = source(iop)
and K = kernel(iop)
.
The estimate is computed by finding the first integer n
such that the quadrature error on the largest element τ
lies below tol
for points x
satisfying dist(x,center(τ)) = n*radius(τ)
.
Note that the desired tolerance may not be achievable if the quadrature rule is not accurate enough, or if τ
is not sufficiently small, and therefore a maximum number of iterations maxiter
is provided to avoid an infinite loops. In such cases, it is recommended that you either increase the quadrature order, or decrease the mesh size.
Note: this is obviously a heuristic, and may not be accurate in all cases.
sourceInti.fibonnaci_points_sphere
— Methodfibonnaci_points_sphere(N,r,c)
Return N
points distributed (roughly) in a uniform manner on the sphere of radius r
centered at c
.
sourceInti.flip_normal
— Methodflip_normal(q::QuadratureNode)
Return a new QuadratureNode
with the normal vector flipped.
sourceInti.gauss_curvature
— Methodgauss_curvature(Q::Quadrature)
Compute the gauss_curvature
at each quadrature node in Q
.
sourceInti.gauss_curvature
— Methodgauss_curvature(τ, x̂)
Calculate the Gaussian curvature of the element τ
at the parametric coordinate x̂
.
sourceInti.geometric_dimension
— Functiongeometric_dimension(x)
NNumber of degrees of freedom necessary to locally represent the geometrical object. For example, lines have geometric dimension of 1 (whether in ℝ²
or in ℝ³
), while surfaces have geometric dimension of 2.
sourceInti.global_get_entity
— Methodglobal_get_entity(k::EntityKey)
Retrieve the GeometricEntity
corresponding to the EntityKey
k
from the global ENTITIES
dictionary.
sourceInti.gmsh_curve
— Methodgmsh_curve(f::Function, a, b; npts=100, meshsize = 0, tag=-1)
Create a curve in the current gmsh
model given by {f(t) : t ∈ (a,b) }
where f
is a function from ℝ
to ℝ^3
. The curve is approximated by C² b-splines passing through npts
equispaced in parameter space. If a meshsize
is given, gmsh
will use it when meshing the curve.
sourceInti.hessian
— Methodhesssian(el,x)
Given a (possibly vector-valued) functor f : 𝐑ᵐ → 𝐅ⁿ
, return the n × m × m
matrix Aᵢⱼⱼ = ∂²fᵢ/∂xⱼ∂xⱼ
. By default ForwardDiff
is used to compute the hessian, but you should overload this method for specific f
if better performance and/or precision is required.
Note: both x
and f(x)
are expected to be of SVector
type.
sourceInti.image
— Functionimage(f)
Given a function-like object f: Ω → R
, return f(Ω)
.
sourceInti.import_mesh
— Methodimport_mesh(filename = nothing; dim=3)
Open filename
and create a Mesh
from the gmsh
model in it.
If filename
is nothing
, the current gmsh
model is used. Note that this assumes that the Gmsh API has been initialized through gmsh.initialize
.
Passing dim=2
will create a two-dimensional mesh by projecting the original mesh onto the x,y
plane.
sourceInti.integrate
— Methodintegrate(f,quad::Quadrature)
Compute ∑ᵢ f(qᵢ)wᵢ
, where the qᵢ
are the quadrature nodes of quad
, and wᵢ
are the quadrature weights.
Note that you must define f(::QuadratureNode)
: use q.coords
and q.normal
if you need to access the coordinate or normal vector at que quadrature node.
sourceInti.integrate
— Methodintegrate(f,q::ReferenceQuadrature)
-integrate(f,x,w)
Integrate the function f
using the quadrature rule q
. This is simply sum(f.(x) .* w)
, where x
and w
are the quadrature nodes and weights, respectively.
The function f
should take an SVector
as input.
sourceInti.integrate_with_error_estimate
— Functionintegrate_with_error_estimate(f, quad::EmbeddedQuadrature, norm = LinearAlgebra.norm)
Return I, E
where I
is the estimated integral of f
over domain(quad)
using the high-order quadrature and E
is the error estimate obtained by taking the norm
of the difference between the high and low-order quadratures in quad
.
sourceInti.integration_measure
— Methodintegration_measure(f, x̂)
Given the Jacobian matrix J
of a transformation f : ℝᴹ → ℝᴺ
compute the integration measure √det(JᵀJ)
at the parametric coordinate x̂
sourceInti.interface_method
— Methodinterface_method(x)
A method of an abstract type
for which concrete subtypes are expected to provide an implementation.
sourceInti.internal_boundary
— Methodinternal_boundary(Ω::Domain)
Return the internal boundaries of a Domain
. These are entities in skeleton(Ω)
which appear at least twice as a boundary of entities in Ω
.
sourceInti.interpolation_order
— Methodinterpolation_order(qrule::ReferenceQuadrature)
The interpolation order of a quadrature rule is defined as the the smallest k
such that there exists a unique polynomial in PolynomialSpace{D,k}
that minimizes the error in approximating the function f
at the quadrature nodes.
For an N
-point Gauss quadrature rule on the segment, the interpolation order is N-1
since N
points uniquely determine a polynomial of degree N-1
.
For a triangular reference domain, the interpolation order is more difficult to define. An unisolvent three-node quadrature on the triangular, for example, has an interpolation order k=1
since the three nodes uniquely determine a linear polynomial, but a four-node quadrature may also have an interpolation order k=1
since for k=2
there are multiple polynomials that pass through the four nodes.
sourceInti.jacobian
— Methodjacobian(f,x)
Given a (possibly vector-valued) functor f : 𝐑ᵐ → 𝐅ⁿ
, return the n × m
matrix Aᵢⱼ = ∂fᵢ/∂xⱼ
. By default ForwardDiff
is used to compute the jacobian, but you should overload this method for specific f
if better performance and/or precision is required.
Note: both x
and f(x)
are expected to be of SVector
type.
sourceInti.kress_change_of_variables
— Methodkress_change_of_variables(P)
Return a change of variables mapping [0,1]
to [0,1]
with the property that the first P-1
derivatives of the transformation vanish at x=0
.
sourceInti.kress_change_of_variables_periodic
— Methodkress_change_of_variables_periodic(P)
Like kress_change_of_variables
, this change of variables maps the interval [0,1]
onto itself, but the first P
derivatives of the transformation vanish at both endpoints (thus making it a periodic function).
This change of variables can be used to periodize integrals over the interval [0,1]
by mapping the integrand into a new integrand that vanishes (to order P) at both endpoints.
sourceInti.lagrange_basis
— Methodlagrange_basis(nodes,[sp::AbstractPolynomialSpace])
Return the set of n
polynomials in sp
taking the value of 1
on node i
and 0
on nodes j ≂̸ i
for 1 ≤ i ≤ n
.
sourceInti.lagrange_basis
— Methodlagrange_basis(qrule::ReferenceQuadrature)
Return a function L : ℝᴺ → ℝᵖ
where N
is the dimension of the domain of qrule
, and p
is the number of nodes in qrule
. The function L
is a polynomial in polynomial_space(qrule)
, and L(xⱼ)[i] = δᵢⱼ
(i.e. the i
th component of L
is the i
th Lagrange basis).
sourceInti.lagrange_basis
— Methodlagrange_basis(E::Type{<:LagrangeElement})
Return the Lagrange basis B
for the element E
. Evaluating B(x)
yields the value of each basis function at x
.
sourceInti.line
— Methodline(a,b)
Create a [GeometricEntity
] representing a straight line connecting points a
and b
. The points a
and b
can be either SVector
s or a Tuple
.
The parametrization of the line is given by f(u) = a + u(b - a)
, where 0 ≤ u ≤ 1
.
sourceInti.mean_curvature
— Methodmean_curvature(Q::Quadrature)
Compute the mean_curvature
at each quadrature node in Q
.
sourceInti.mean_curvature
— Methodmean_curvature(τ, x̂)
Calculate the mean curvature of the element τ
at the parametric coordinate x̂
.
sourceInti.measure
— Functionmeasure(k::EntityKey, rtol)
Compute the length/area/volume of the entity k
using an adaptive quadrature with a relative tolerance rtol
. Assumes that the entity has an explicit parametrization.
sourceInti.meshgen!
— Methodmeshgen!(mesh,Ω,sz)
Similar to meshgen
, but append entries to mesh
.
sourceInti.meshgen
— Methodmeshgen(Ω, n)
+end
where a dynamic dispatch is performed only on the element types (typically small for a given mesh).
sourceInti.ellipsoid
— Methodellipsoid(; translation, rotation, scaling, labels)
Create an ellipsoid entity in 3D, and apply optional transformations. Returns the key of the created entity.
sourceInti.ent2etags
— Methodent2etags(msh::AbstractMesh)
Return a dictionary mapping entities to a dictionary of element types to element tags.
sourceInti.entities
— Methodentities(Ω::Domain)
Return all entities making up a domain (as a set of EntityKey
s).
sourceInti.etype_to_nearest_points
— Methodetype_to_nearest_points(X,Y::Quadrature; maxdist)
For each element el
in Y.mesh
, return a list with the indices of all points in X
for which el
is the nearest element. Ignore indices for which the distance exceeds maxdist
.
sourceInti.external_boundary
— Methodexternal_boundary(Ω::Domain)
Return the external boundaries inside a domain. These are entities in the skeleton of Ω which are not in the internal boundaries of Ω.
See also: internal_boundary
, skeleton
.
sourceInti.farfield_distance
— Methodfarfield_distance(iop::IntegralOperator; tol, maxiter = 10)
+farfield_distance(K, Q::Quadrature; tol, maxiter = 10)
Return an estimate of the distance d
such that the (absolute) quadrature error of the integrand y -> K(x,y)
is below tol
for x
at a distance d
from the center of the largest element in Q
; when an integral operator is passed, we have Q::Quadrature = source(iop)
and K = kernel(iop)
.
The estimate is computed by finding the first integer n
such that the quadrature error on the largest element τ
lies below tol
for points x
satisfying dist(x,center(τ)) = n*radius(τ)
.
Note that the desired tolerance may not be achievable if the quadrature rule is not accurate enough, or if τ
is not sufficiently small, and therefore a maximum number of iterations maxiter
is provided to avoid an infinite loops. In such cases, it is recommended that you either increase the quadrature order, or decrease the mesh size.
Note: this is obviously a heuristic, and may not be accurate in all cases.
sourceInti.fibonnaci_points_sphere
— Methodfibonnaci_points_sphere(N,r,c)
Return N
points distributed (roughly) in a uniform manner on the sphere of radius r
centered at c
.
sourceInti.flip_normal
— Methodflip_normal(q::QuadratureNode)
Return a new QuadratureNode
with the normal vector flipped.
sourceInti.gauss_curvature
— Methodgauss_curvature(Q::Quadrature)
Compute the gauss_curvature
at each quadrature node in Q
.
sourceInti.gauss_curvature
— Methodgauss_curvature(τ, x̂)
Calculate the Gaussian curvature of the element τ
at the parametric coordinate x̂
.
sourceInti.geometric_dimension
— Functiongeometric_dimension(x)
NNumber of degrees of freedom necessary to locally represent the geometrical object. For example, lines have geometric dimension of 1 (whether in ℝ²
or in ℝ³
), while surfaces have geometric dimension of 2.
sourceInti.global_get_entity
— Methodglobal_get_entity(k::EntityKey)
Retrieve the GeometricEntity
corresponding to the EntityKey
k
from the global ENTITIES
dictionary.
sourceInti.gmsh_curve
— Methodgmsh_curve(f::Function, a, b; npts=100, meshsize = 0, tag=-1)
Create a curve in the current gmsh
model given by {f(t) : t ∈ (a,b) }
where f
is a function from ℝ
to ℝ^3
. The curve is approximated by C² b-splines passing through npts
equispaced in parameter space. If a meshsize
is given, gmsh
will use it when meshing the curve.
sourceInti.hessian
— Methodhesssian(el,x)
Given a (possibly vector-valued) functor f : 𝐑ᵐ → 𝐅ⁿ
, return the n × m × m
matrix Aᵢⱼⱼ = ∂²fᵢ/∂xⱼ∂xⱼ
. By default ForwardDiff
is used to compute the hessian, but you should overload this method for specific f
if better performance and/or precision is required.
Note: both x
and f(x)
are expected to be of SVector
type.
sourceInti.image
— Functionimage(f)
Given a function-like object f: Ω → R
, return f(Ω)
.
sourceInti.import_mesh
— Methodimport_mesh(filename = nothing; dim=3)
Open filename
and create a Mesh
from the gmsh
model in it.
If filename
is nothing
, the current gmsh
model is used. Note that this assumes that the Gmsh API has been initialized through gmsh.initialize
.
Passing dim=2
will create a two-dimensional mesh by projecting the original mesh onto the x,y
plane.
sourceInti.integrate
— Methodintegrate(f,quad::Quadrature)
Compute ∑ᵢ f(qᵢ)wᵢ
, where the qᵢ
are the quadrature nodes of quad
, and wᵢ
are the quadrature weights.
Note that you must define f(::QuadratureNode)
: use q.coords
and q.normal
if you need to access the coordinate or normal vector at que quadrature node.
sourceInti.integrate
— Methodintegrate(f,q::ReferenceQuadrature)
+integrate(f,x,w)
Integrate the function f
using the quadrature rule q
. This is simply sum(f.(x) .* w)
, where x
and w
are the quadrature nodes and weights, respectively.
The function f
should take an SVector
as input.
sourceInti.integrate_with_error_estimate
— Functionintegrate_with_error_estimate(f, quad::EmbeddedQuadrature, norm = LinearAlgebra.norm)
Return I, E
where I
is the estimated integral of f
over domain(quad)
using the high-order quadrature and E
is the error estimate obtained by taking the norm
of the difference between the high and low-order quadratures in quad
.
sourceInti.integration_measure
— Methodintegration_measure(f, x̂)
Given the Jacobian matrix J
of a transformation f : ℝᴹ → ℝᴺ
compute the integration measure √det(JᵀJ)
at the parametric coordinate x̂
sourceInti.interface_method
— Methodinterface_method(x)
A method of an abstract type
for which concrete subtypes are expected to provide an implementation.
sourceInti.internal_boundary
— Methodinternal_boundary(Ω::Domain)
Return the internal boundaries of a Domain
. These are entities in skeleton(Ω)
which appear at least twice as a boundary of entities in Ω
.
sourceInti.interpolation_order
— Methodinterpolation_order(qrule::ReferenceQuadrature)
The interpolation order of a quadrature rule is defined as the the smallest k
such that there exists a unique polynomial in PolynomialSpace{D,k}
that minimizes the error in approximating the function f
at the quadrature nodes.
For an N
-point Gauss quadrature rule on the segment, the interpolation order is N-1
since N
points uniquely determine a polynomial of degree N-1
.
For a triangular reference domain, the interpolation order is more difficult to define. An unisolvent three-node quadrature on the triangular, for example, has an interpolation order k=1
since the three nodes uniquely determine a linear polynomial, but a four-node quadrature may also have an interpolation order k=1
since for k=2
there are multiple polynomials that pass through the four nodes.
sourceInti.jacobian
— Methodjacobian(f,x)
Given a (possibly vector-valued) functor f : 𝐑ᵐ → 𝐅ⁿ
, return the n × m
matrix Aᵢⱼ = ∂fᵢ/∂xⱼ
. By default ForwardDiff
is used to compute the jacobian, but you should overload this method for specific f
if better performance and/or precision is required.
Note: both x
and f(x)
are expected to be of SVector
type.
sourceInti.kress_change_of_variables
— Methodkress_change_of_variables(P)
Return a change of variables mapping [0,1]
to [0,1]
with the property that the first P-1
derivatives of the transformation vanish at x=0
.
sourceInti.kress_change_of_variables_periodic
— Methodkress_change_of_variables_periodic(P)
Like kress_change_of_variables
, this change of variables maps the interval [0,1]
onto itself, but the first P
derivatives of the transformation vanish at both endpoints (thus making it a periodic function).
This change of variables can be used to periodize integrals over the interval [0,1]
by mapping the integrand into a new integrand that vanishes (to order P) at both endpoints.
sourceInti.lagrange_basis
— Methodlagrange_basis(nodes,[sp::AbstractPolynomialSpace])
Return the set of n
polynomials in sp
taking the value of 1
on node i
and 0
on nodes j ≂̸ i
for 1 ≤ i ≤ n
.
sourceInti.lagrange_basis
— Methodlagrange_basis(qrule::ReferenceQuadrature)
Return a function L : ℝᴺ → ℝᵖ
where N
is the dimension of the domain of qrule
, and p
is the number of nodes in qrule
. The function L
is a polynomial in polynomial_space(qrule)
, and L(xⱼ)[i] = δᵢⱼ
(i.e. the i
th component of L
is the i
th Lagrange basis).
sourceInti.lagrange_basis
— Methodlagrange_basis(E::Type{<:LagrangeElement})
Return the Lagrange basis B
for the element E
. Evaluating B(x)
yields the value of each basis function at x
.
sourceInti.line
— Methodline(a,b)
Create a [GeometricEntity
] representing a straight line connecting points a
and b
. The points a
and b
can be either SVector
s or a Tuple
.
The parametrization of the line is given by f(u) = a + u(b - a)
, where 0 ≤ u ≤ 1
.
sourceInti.mean_curvature
— Methodmean_curvature(Q::Quadrature)
Compute the mean_curvature
at each quadrature node in Q
.
sourceInti.mean_curvature
— Methodmean_curvature(τ, x̂)
Calculate the mean curvature of the element τ
at the parametric coordinate x̂
.
sourceInti.measure
— Functionmeasure(k::EntityKey, rtol)
Compute the length/area/volume of the entity k
using an adaptive quadrature with a relative tolerance rtol
. Assumes that the entity has an explicit parametrization.
sourceInti.meshgen
— Functionmeshgen(Ω, n)
meshgen(Ω, n_dict)
-meshgen(Ω; meshsize)
Generate a Mesh
for the domain Ω
where each curve is meshed using n
elements. Passing a dictionary allows for a finer control; in such cases, n_dict[ent]
should return an integer for each entity ent
in Ω
of geometric_dimension
one.
Alternatively, a meshsize
can be passed, in which case, the number of elements is computed as so as to obtain an average mesh size of meshsize
. Note that the actual mesh size may vary significantly for each element if the parametrization is far from uniform.
This function requires the entities forming Ω
to have an explicit parametrization.
Mesh quality The quality of the generated mesh created using meshgen
depends on the quality of the underlying parametrization. For complex surfaces, you are better off using a proper mesher such as gmsh
.
sourceInti.monomial_basis
— Functionmonomial_basis(sp::PolynomialSpace)
Return a function f : ℝᴺ → ℝᵈ
, where N
is the dimension of the domain of sp
containing a basis of monomials 𝐱ᶿ
spanning the polynomial space PolynomialSpace
.
sourceInti.near_interaction_list
— Methodnear_interaction_list(X,Y::AbstractMesh; tol)
For each element el
of type E
in Y
, return the indices of the points in X
which are closer than tol
to the center
of el
.
This function returns a dictionary where e.g. dict[E][5] --> Vector{Int}
gives the indices of points in X
which are closer than tol
to the center of the fifth element of type E
.
If tol
is a Dict
, then tol[E]
is the tolerance for elements of type E
.
sourceInti.new_tag
— Methodnew_tag(dim)
Return a new tag for an entity of dimension dim
so that EntityKey(dim, tag)
is not already in ENTITIES
.
sourceInti.nodes
— Methodnodes(msh::SubMesh)
A view of the nodes of the parent mesh belonging to the submesh. The ordering is given by the nodetags
function.
sourceInti.nodetags
— Methodnodetags(msh::SubMesh)
Return the tags of the nodes in the parent mesh belonging to the submesh.
sourceInti.normal
— Methodnormal(el, x̂)
Return the normal vector of el
at the parametric coordinate x̂
.
sourceInti.normal
— Methodnormal(q)
Return the normal vector of q
, if it exists.
sourceInti.notimplemented
— Methodnotimplemented()
Things which should probably be implemented at some point.
sourceInti.order
— Methodorder(q::ReferenceQuadrature)
A quadrature of order p
(sometimes called degree of precision) integrates all polynomials of degree ≤ p
but not ≤ p + 1
.
sourceInti.order
— Methodorder(el::LagrangeElement)
The order of the element's interpolating polynomial (e.g. a LagrangeLine
with 2
nodes defines a linear polynomial, and thus has order 1
).
sourceInti.parametric_curve
— Methodparametric_curve(f, a::Real, b::Real)
Create a [GeometricEntity
] representing a parametric curve defined by the {f(t) | a ≤ t ≤ b}
. The function f
should map a scalar to an SVector
.
Flipping the orientation is supported by passing a > b
.
sourceInti.parametric_surface
— Function parametric_surface(f, lc, hc, boundary = nothing; kwargs...)
Create a parametric surface defined by the function f
over the rectangular domain defined by the lower corner lc
and the upper corner hc
. The optional boundary
argument can be used to specify the boundary curves of the surface.
Arguments
f
: A function that takes two arguments x
and y
and returns a tuple (u, v)
representing the parametric coordinates of the surface at (x, y)
.lc
: A 2-element array representing the lower corner of the rectangular domain.hc
: A 2-element array representing the upper corner of the rectangular domain.boundary
: An optional array of boundary curves that define the surface.
Keyword Arguments
kwargs
: Additional keyword arguments that can be passed to the GeometricEntity
constructor.
Returns
- The key of the created
GeometricEntity
.
sourceInti.polynomial_solutions_vdim
— Functionpolynomial_solutions_vdim(op, order[, center])
For every monomial term pₙ
of degree order
, compute a polynomial Pₙ
such that ℒ[Pₙ] = pₙ
, where ℒ
is the differential operator associated with op
. This function returns {pₙ,Pₙ,γ₁Pₙ}
, where γ₁Pₙ
is the generalized Neumann trace of Pₙ
.
Passing a point center
will shift the monomials and solutions accordingly.
sourceInti.polynomial_space
— Methodpolynomial_space(qrule::ReferenceQuadrature)
Return a PolynomialSpace
associated with the interpolation_order
of the quadrature nodes of qrule
.
sourceInti.qcoords
— Methodqcoords(q)
Return the coordinate of the quadrature nodes associated with q
.
sourceInti.quadrature_to_node_vals
— Methodquadrature_to_node_vals(Q::Quadrature, qvals::AbstractVector)
Given a vector qvals
of scalar values at the quadrature nodes of Q
, return a vector ivals
of scalar values at the interpolation nodes of Q.mesh
.
sourceInti.qweights
— Methodqweights(q)
Return the quadrature weights associated with q
.
sourceInti.reference_nodes
— Methodreference_nodes(el::LagrangeElement)
-reference_nodes(::Type{<:LagrangeElement})
Return the reference nodes on domain(el)
used for the polynomial interpolation. The function values on these nodes completely determines the interpolating polynomial.
We use the same convention as gmsh
for defining the reference nodes and their order (see node ordering on gmsh
documentation).
sourceInti.return_type
— Methodreturn_type(f[,args...])
The type returned by f(args...)
, where args
is a tuple of types. Falls back to Base.promote_op
by default.
A functors of type T
with a knonw return type should extend return_type(::T,args...)
to avoid relying on promote_op
.
sourceInti.rotation_matrix
— Methodrotation_matrix(rot)
Constructs a rotation matrix given the rotation angles around the x, y, and z axes.
Arguments
rot
: A tuple or vector containing the rotation angles in radians for each axis.
Returns
R::SMatrix
: The resulting rotation matrix.
sourceInti.single_double_layer
— Methodsingle_double_layer(; op, target, source::Quadrature, compression,
-correction, derivative = false)
Construct a discrete approximation to the single- and double-layer integral operators for op
, mapping values defined on the quadrature nodes of source
to values defined on the nodes of target
. If derivative = true
, return instead the adjoint double-layer and hypersingular operators (which are the derivative of the single- and double-layer, respectively).
You must choose a compression
method and a correction
method, as described below.
Compression
The compression
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the dense linear operators should be compressed. The available options are:
(method = :none, )
: no compression is performed, the resulting matrices are dense.(method =:hmatrix, tol)
: the resulting operators are compressed using hierarchical matrices with an absolute tolerance tol
(defaults to 1e-8
).(method = :fmm, tol)
: the resulting operators are compressed using the fast multipole method with an absolute tolerance tol
(defaults to 1e-8
).
Correction
The correction
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the singular and nearly-singular integrals should be computed. The available options are:
(method = :none, )
: no correction is performed. This is not recommented, as the resulting approximation will be inaccurate if the source and target are not sufficiently far apart.(method = :dim, maxdist, target_location)
: use the density interpolation method to compute the correction. maxdist
specifies the distance between source and target points above which no correction is performed (defaults to Inf
). target_location
should be either :inside
, :outside
, or :on
, and specifies where the target
points lie relative to the to the
sourcecurve/surface (which is assumed to be closed). When
target === source,
target_location` is not needed.
sourceInti.single_double_layer_potential
— Methodsingle_double_layer_potential(; op, source)
Return the single- and double-layer potentials for op
as IntegralPotential
s.
sourceInti.skeleton
— Methodskeleton(Ω::Domain)
Return all the boundaries of the domain, i.e. the domain's skeleton.
sourceInti.stack_weakdeps_env!
— Methodstack_weakdeps_env!(; verbose = false, update = false)
Push to the load stack an environment providing the weak dependencies of Inti.jl. This allows benefiting from additional functionalities of Inti.jl which are powered by weak dependencies without having to manually install them in your environment.
Set update=true
if you want to update the weakdeps
environment.
Warning Calling this function can take quite some time, especially the first time around, if packages have to be installed or precompiled. Run in verbose
mode to see what is happening.
Examples:
Inti.stack_weakdeps_env!()
-using HMatrices
sourceInti.standard_basis_vector
— Methodstandard_basis_vector(k, ::Val{N})
Create an SVector
of length N with a 1 in the kth position and zeros elsewhere.
sourceInti.svector
— Methodsvector(f,n)
Create an SVector
of length n, computing each element as f(i), where i is the index of the element.
sourceInti.torus
— Methodtorus(; r, R, translation, rotation, scaling, labels)
Create a torus entity in 3D, and apply optional transformations. Returns the key. The parameters r
and R
are the inner and outer radii of the torus.
sourceInti.uniform_points_circle
— Methoduniform_points_circle(N,r,c)
Return N
points uniformly distributed on a circle of radius r
centered at c
.
sourceInti.vdim_correction
— Methodvdim_correction(op,X,Y,Y_boundary,S,D,V; green_multiplier, kwargs...)
Compute a correction to the volume potential V : Y → X
such that V + δV
is a more accurate approximation of the underlying volume potential operator. The correction is computed using the (volume) density interpolation method.
This function requires a op::AbstractDifferentialOperator
, a target set X
, a source quadrature Y
, a boundary quadrature Y_boundary
, approximations S : Y_boundary -> X
and D : Y_boundary -> X
to the single- and double-layer potentials (correctly handling nearly-singular integrals), and a naive approximation of the volume potential V
. The green_multiplier
is a vector of the same length as X
storing the value of μ(x)
for x ∈ X
in the Green identity (see _green_multiplier
).
See [8] for more details on the method.
Optional kwargs
:
interpolation_order
: the order of the polynomial interpolation. By default, the maximum order of the quadrature rules is used.maxdist
: distance beyond which interactions are considered sufficiently far so that no correction is needed. This is used to determine a threshold for nearly-singular corrections.center
: the center of the basis functions. By default, the basis functions are centered at the origin.shift
: a boolean indicating whether the basis functions should be shifted and rescaled to each element.
sourceInti.vdim_mesh_center
— Methodvdim_mesh_center(msh)
Point x
which minimizes ∑ (x-xⱼ)²/r²ⱼ, where xⱼ and rⱼ are the circumcenter and circumradius of the elements of msh
, respectively.
sourceInti.vertices
— Methodvertices(el::LagrangeElement)
Coordinates of the vertices of el
.
sourceInti.vertices_idxs
— Methodvertices_idxs(el::LagrangeElement)
The indices of the nodes in el
that define the vertices of the element.
sourceInti.volume_potential
— Methodvolume_potential(; op, target, source::Quadrature, compression, correction)
Compute the volume potential operator for a given PDE.
Arguments
op
: The PDE (Partial Differential Equation) to solve.target
: The target domain where the potential is computed.source
: The source domain where the potential is generated.compression
: The compression method to use for the potential operator.correction
: The correction method to use for the potential operator.
Returns
The volume potential operator V
that represents the interaction between the target and source domains.
Compression
The compression
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the dense linear operators should be compressed. The available options are:
(method = :none, )
: no compression is performed, the resulting matrices are dense.(method =:hmatrix, tol)
: the resulting operators are compressed using hierarchical matrices with an absolute tolerance tol
(defaults to 1e-8
).(method = :fmm, tol)
: the resulting operators are compressed using the fast multipole method with an absolute tolerance tol
(defaults to 1e-8
).
Correction
The correction
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the singular and nearly-singular integrals should be computed. The available options are:
(method = :none, )
: no correction is performed. This is not recommented, as the resulting approximation will be inaccurate if the source and target are not sufficiently far apart.(method = :dim, maxdist, target_location)
: use the density interpolation method to compute the correction. maxdist
specifies the distance between source and target points above which no correction is performed (defaults to Inf
). target_location
should be either :inside
, :outside
, or :on
, and specifies where the target
points lie relative to the to the
source's boundary. When
target === source,
target_location` is not needed.
Details
The volume potential operator is computed by assembling the integral operator V
using the single-layer kernel G
. The operator V
is then compressed using the specified compression method. If no compression is specified, the operator is returned as is. If a correction method is specified, the correction is computed and added to the compressed operator.
sourceSettings
This document was generated with Documenter.jl version 1.7.0 on Saturday 19 October 2024. Using Julia version 1.11.1.
+meshgen(Ω; meshsize)
Generate a Mesh
for the domain Ω
where each curve is meshed using n
elements. Passing a dictionary allows for a finer control; in such cases, n_dict[ent]
should return an integer for each entity ent
in Ω
of geometric_dimension
one.
Alternatively, a meshsize
can be passed, in which case, the number of elements is computed as so as to obtain an average mesh size of meshsize
. Note that the actual mesh size may vary significantly for each element if the parametrization is far from uniform.
This function requires the entities forming Ω
to have an explicit parametrization.
The quality of the generated mesh created using meshgen
depends on the quality of the underlying parametrization. For complex surfaces, you are better off using a proper mesher such as gmsh
.
Inti.meshgen!
— Methodmeshgen!(mesh,Ω,sz)
Similar to meshgen
, but append entries to mesh
.
Inti.monomial_basis
— Functionmonomial_basis(sp::PolynomialSpace)
Return a function f : ℝᴺ → ℝᵈ
, where N
is the dimension of the domain of sp
containing a basis of monomials 𝐱ᶿ
spanning the polynomial space PolynomialSpace
.
Inti.near_interaction_list
— Methodnear_interaction_list(X,Y::AbstractMesh; tol)
For each element el
of type E
in Y
, return the indices of the points in X
which are closer than tol
to the center
of el
.
This function returns a dictionary where e.g. dict[E][5] --> Vector{Int}
gives the indices of points in X
which are closer than tol
to the center of the fifth element of type E
.
If tol
is a Dict
, then tol[E]
is the tolerance for elements of type E
.
Inti.new_tag
— Methodnew_tag(dim)
Return a new tag for an entity of dimension dim
so that EntityKey(dim, tag)
is not already in ENTITIES
.
Inti.nodes
— Methodnodes(msh::SubMesh)
A view of the nodes of the parent mesh belonging to the submesh. The ordering is given by the nodetags
function.
Inti.nodetags
— Methodnodetags(msh::SubMesh)
Return the tags of the nodes in the parent mesh belonging to the submesh.
Inti.normal
— Methodnormal(el, x̂)
Return the normal vector of el
at the parametric coordinate x̂
.
Inti.normal
— Methodnormal(q)
Return the normal vector of q
, if it exists.
Inti.notimplemented
— Methodnotimplemented()
Things which should probably be implemented at some point.
Inti.order
— Methodorder(q::ReferenceQuadrature)
A quadrature of order p
(sometimes called degree of precision) integrates all polynomials of degree ≤ p
but not ≤ p + 1
.
Inti.order
— Methodorder(el::LagrangeElement)
The order of the element's interpolating polynomial (e.g. a LagrangeLine
with 2
nodes defines a linear polynomial, and thus has order 1
).
Inti.parametric_curve
— Methodparametric_curve(f, a::Real, b::Real)
Create a [GeometricEntity
] representing a parametric curve defined by the {f(t) | a ≤ t ≤ b}
. The function f
should map a scalar to an SVector
.
Flipping the orientation is supported by passing a > b
.
Inti.parametric_surface
— Function parametric_surface(f, lc, hc, boundary = nothing; kwargs...)
Create a parametric surface defined by the function f
over the rectangular domain defined by the lower corner lc
and the upper corner hc
. The optional boundary
argument can be used to specify the boundary curves of the surface.
Arguments
f
: A function that takes two argumentsx
andy
and returns a tuple(u, v)
representing the parametric coordinates of the surface at(x, y)
.lc
: A 2-element array representing the lower corner of the rectangular domain.hc
: A 2-element array representing the upper corner of the rectangular domain.boundary
: An optional array of boundary curves that define the surface.
Keyword Arguments
kwargs
: Additional keyword arguments that can be passed to theGeometricEntity
constructor.
Returns
- The key of the created
GeometricEntity
.
Inti.polynomial_solutions_vdim
— Functionpolynomial_solutions_vdim(op, order[, center])
For every monomial term pₙ
of degree order
, compute a polynomial Pₙ
such that ℒ[Pₙ] = pₙ
, where ℒ
is the differential operator associated with op
. This function returns {pₙ,Pₙ,γ₁Pₙ}
, where γ₁Pₙ
is the generalized Neumann trace of Pₙ
.
Passing a point center
will shift the monomials and solutions accordingly.
Inti.polynomial_space
— Methodpolynomial_space(qrule::ReferenceQuadrature)
Return a PolynomialSpace
associated with the interpolation_order
of the quadrature nodes of qrule
.
Inti.qcoords
— Methodqcoords(q)
Return the coordinate of the quadrature nodes associated with q
.
Inti.quadrature_to_node_vals
— Methodquadrature_to_node_vals(Q::Quadrature, qvals::AbstractVector)
Given a vector qvals
of scalar values at the quadrature nodes of Q
, return a vector ivals
of scalar values at the interpolation nodes of Q.mesh
.
Inti.qweights
— Methodqweights(q)
Return the quadrature weights associated with q
.
Inti.reference_nodes
— Methodreference_nodes(el::LagrangeElement)
+reference_nodes(::Type{<:LagrangeElement})
Return the reference nodes on domain(el)
used for the polynomial interpolation. The function values on these nodes completely determines the interpolating polynomial.
We use the same convention as gmsh
for defining the reference nodes and their order (see node ordering on gmsh
documentation).
Inti.return_type
— Methodreturn_type(f[,args...])
The type returned by f(args...)
, where args
is a tuple of types. Falls back to Base.promote_op
by default.
A functors of type T
with a knonw return type should extend return_type(::T,args...)
to avoid relying on promote_op
.
Inti.rotation_matrix
— Methodrotation_matrix(rot)
Constructs a rotation matrix given the rotation angles around the x, y, and z axes.
Arguments
rot
: A tuple or vector containing the rotation angles in radians for each axis.
Returns
R::SMatrix
: The resulting rotation matrix.
Inti.single_double_layer
— Methodsingle_double_layer(; op, target, source::Quadrature, compression,
+correction, derivative = false)
Construct a discrete approximation to the single- and double-layer integral operators for op
, mapping values defined on the quadrature nodes of source
to values defined on the nodes of target
. If derivative = true
, return instead the adjoint double-layer and hypersingular operators (which are the derivative of the single- and double-layer, respectively).
You must choose a compression
method and a correction
method, as described below.
Compression
The compression
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the dense linear operators should be compressed. The available options are:
(method = :none, )
: no compression is performed, the resulting matrices are dense.(method =:hmatrix, tol)
: the resulting operators are compressed using hierarchical matrices with an absolute tolerancetol
(defaults to1e-8
).(method = :fmm, tol)
: the resulting operators are compressed using the fast multipole method with an absolute tolerancetol
(defaults to1e-8
).
Correction
The correction
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the singular and nearly-singular integrals should be computed. The available options are:
(method = :none, )
: no correction is performed. This is not recommented, as the resulting approximation will be inaccurate if the source and target are not sufficiently far apart.(method = :dim, maxdist, target_location)
: use the density interpolation method to compute the correction.maxdist
specifies the distance between source and target points above which no correction is performed (defaults toInf
).target_location
should be either:inside
,:outside
, or:on
, and specifies where thetarget
points lie relative to the to the
sourcecurve/surface (which is assumed to be closed). When
target === source,
target_location` is not needed.
Inti.single_double_layer_potential
— Methodsingle_double_layer_potential(; op, source)
Return the single- and double-layer potentials for op
as IntegralPotential
s.
Inti.skeleton
— Methodskeleton(Ω::Domain)
Return all the boundaries of the domain, i.e. the domain's skeleton.
Inti.stack_weakdeps_env!
— Methodstack_weakdeps_env!(; verbose = false, update = false)
Push to the load stack an environment providing the weak dependencies of Inti.jl. This allows benefiting from additional functionalities of Inti.jl which are powered by weak dependencies without having to manually install them in your environment.
Set update=true
if you want to update the weakdeps
environment.
Calling this function can take quite some time, especially the first time around, if packages have to be installed or precompiled. Run in verbose
mode to see what is happening.
Examples:
Inti.stack_weakdeps_env!()
+using HMatrices
Inti.standard_basis_vector
— Methodstandard_basis_vector(k, ::Val{N})
Create an SVector
of length N with a 1 in the kth position and zeros elsewhere.
Inti.svector
— Methodsvector(f,n)
Create an SVector
of length n, computing each element as f(i), where i is the index of the element.
Inti.torus
— Methodtorus(; r, R, translation, rotation, scaling, labels)
Create a torus entity in 3D, and apply optional transformations. Returns the key. The parameters r
and R
are the inner and outer radii of the torus.
Inti.uniform_points_circle
— Methoduniform_points_circle(N,r,c)
Return N
points uniformly distributed on a circle of radius r
centered at c
.
Inti.vdim_correction
— Methodvdim_correction(op,X,Y,Y_boundary,S,D,V; green_multiplier, kwargs...)
Compute a correction to the volume potential V : Y → X
such that V + δV
is a more accurate approximation of the underlying volume potential operator. The correction is computed using the (volume) density interpolation method.
This function requires a op::AbstractDifferentialOperator
, a target set X
, a source quadrature Y
, a boundary quadrature Y_boundary
, approximations S : Y_boundary -> X
and D : Y_boundary -> X
to the single- and double-layer potentials (correctly handling nearly-singular integrals), and a naive approximation of the volume potential V
. The green_multiplier
is a vector of the same length as X
storing the value of μ(x)
for x ∈ X
in the Green identity (see _green_multiplier
).
See [8] for more details on the method.
Optional kwargs
:
interpolation_order
: the order of the polynomial interpolation. By default, the maximum order of the quadrature rules is used.maxdist
: distance beyond which interactions are considered sufficiently far so that no correction is needed. This is used to determine a threshold for nearly-singular corrections.center
: the center of the basis functions. By default, the basis functions are centered at the origin.shift
: a boolean indicating whether the basis functions should be shifted and rescaled to each element.
Inti.vdim_mesh_center
— Methodvdim_mesh_center(msh)
Point x
which minimizes ∑ (x-xⱼ)²/r²ⱼ, where xⱼ and rⱼ are the circumcenter and circumradius of the elements of msh
, respectively.
Inti.vertices
— Methodvertices(el::LagrangeElement)
Coordinates of the vertices of el
.
Inti.vertices_idxs
— Methodvertices_idxs(el::LagrangeElement)
The indices of the nodes in el
that define the vertices of the element.
Inti.volume_potential
— Methodvolume_potential(; op, target, source::Quadrature, compression, correction)
Compute the volume potential operator for a given PDE.
Arguments
op
: The PDE (Partial Differential Equation) to solve.target
: The target domain where the potential is computed.source
: The source domain where the potential is generated.compression
: The compression method to use for the potential operator.correction
: The correction method to use for the potential operator.
Returns
The volume potential operator V
that represents the interaction between the target and source domains.
Compression
The compression
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the dense linear operators should be compressed. The available options are:
(method = :none, )
: no compression is performed, the resulting matrices are dense.(method =:hmatrix, tol)
: the resulting operators are compressed using hierarchical matrices with an absolute tolerancetol
(defaults to1e-8
).(method = :fmm, tol)
: the resulting operators are compressed using the fast multipole method with an absolute tolerancetol
(defaults to1e-8
).
Correction
The correction
argument is a named tuple with a method
field followed by method-specific fields. It specifies how the singular and nearly-singular integrals should be computed. The available options are:
(method = :none, )
: no correction is performed. This is not recommented, as the resulting approximation will be inaccurate if the source and target are not sufficiently far apart.(method = :dim, maxdist, target_location)
: use the density interpolation method to compute the correction.maxdist
specifies the distance between source and target points above which no correction is performed (defaults toInf
).target_location
should be either:inside
,:outside
, or:on
, and specifies where thetarget
points lie relative to the to the
source's boundary. When
target === source,
target_location` is not needed.
Details
The volume potential operator is computed by assembling the integral operator V
using the single-layer kernel G
. The operator V
is then compressed using the specified compression method. If no compression is specified, the operator is returned as is. If a correction method is specified, the correction is computed and added to the compressed operator.