Docstrings

HAdaptiveIntegration.Domain.OrthotopeType
struct Orthotope{D,T} <: AbstractDomain{D,T}

Axes-aligned Orthotope in D dimensions, with element type T, given by two points low_corner and high_corner. Note that, we must have low_corner .≤ high_corner.

Fields:

  • corners::SVector{2,SVector{D,T}}: corners[1] is the low corner and corners[2] is the high corner.

Invariants (not check at construction):

  • corners[1] .≤ corners[2]

Constructors:

  • Orthotope(low_corner, high_corner)
  • Orthotope{T}(low_corner, high_corner)
  • Orthotope(corners::SVector{2,SVector{D,T}})
source
HAdaptiveIntegration.Domain.SegmentType
struct Segment{T} <: AbstractDomain{1,T}

Segment in 1 dimension, with element type T, given by xmin and xmax. Note that, we must have xmin ≤ xmax.

Fields:

  • xmin::T: the low corner.
  • xmax::T: the high corner.

Invariants (not check at construction):

  • xmin ≤ xmax

Constructors:

  • Segment(xmin, xmax)
source
HAdaptiveIntegration.Domain.SimplexType
struct Simplex{D,T,N} <: AbstractDomain{D,T}

A simplex in D dimensions with N=D+1 vertices of element type T.

Fields:

  • vertices::SVector{N,SVector{D,T}}: vertices of the simplex.

Invariants (not check at construction):

  • N = D+1

Constructors:

  • Simplex(vertices...)
  • Simplex{T}(vertices...)
  • Simplex(vertices::SVector{N,SVector{D,T}})
source
HAdaptiveIntegration.Domain.map_from_referenceFunction
Φ, μ = map_from_reference(domain::DOM) where {DOM<:AbstractDomain}

Return Φ and μ where Φ is an anonymous function that maps the reference domain to the physical domain domain and μ is the absolute value of the Jacobian's determinant of the map Φ.

source
HAdaptiveIntegration.Domain.map_to_referenceFunction
map_to_reference(domain::DOM) where {DOM<:AbstractDomain}

Return an anonymous function that maps the physical domain to the reference domain.

Constraints:

  • For Segment, must have xmax > xmin.
  • For Simplex{D}, the vertices must form a valid D-dimensional simplex with non-zero volume).
  • For Orthotope, must have high_corner .> low_corner.
source
HAdaptiveIntegration.Domain.reference_simplexMethod
reference_simplex(D::Int, T=float(Int))

Return the reference D-dimensional simplex with element type T, which is the convex hull of the N=D+1 points (0,...,0), (1,0,...,0), (0,1,0,...,0), ..., (0,...,0,1).

source
HAdaptiveIntegration.Domain.subdivide_reference_orthotopeMethod
subdivide_reference_orthotope(::Val{D}, ::Type{T}=float(Int))

Like subdivide_orthotope, but operates on the reference orthotope. Since the output depends only on the dimension D, and the type T used to represent coordinates, this function is generated for each combination of D and T.

source
HAdaptiveIntegration.Domain.subdivide_reference_simplexMethod
subdivide_reference_simplex(::Val{D}, ::Type{T}=float(Int)) where {D,T}

Like subdivide_simplex, but operates on the reference simplex. Since the output depends only on the dimension D, and the type T used to represent coordinates, this function is generated for each combination of D and T.

source
HAdaptiveIntegration.Rule.EmbeddedCubatureType
struct EmbeddedCubature{D,T}

An embedded cubature rule consisting of a high order cubature rule with H nodes and a low order cubature rule with L nodes. Note that the low order cubature uses nodes[1:L] as its nodes. The cubature nodes and weights are assumed to be for the reference domain.

Fields:

  • nodes::Vector{SVector{D,T}}: the cubature nodes.
  • weights_high::Vector{T}: the cubature weights for the high order cubature.
  • weights_low::Vector{T}: the cubature weights for the low order cubature.

Invariants (check at construction):

  • length(nodes) == length(weights_high)
  • length(weights_high) ≥ length(weights_low)
source
HAdaptiveIntegration.Rule.EmbeddedCubatureMethod
(ec::EmbeddedCubature{D,T})(
    fct, domain::AbstractDomain{D,T}, norm=LinearAlgebra.norm
) where {D,T}

Return I_high and norm(I_high - I_low) where I_high and I_low are the result of the high order cubature and the low order cubature on domain. The function fct must take a SVector{D,T} to a return type K, and K must support the multiplication by a scalar and the addition. Note that there is no check, beyond compatibility of dimension and type, that the embedded cubature is for the right domain.

source
HAdaptiveIntegration.Rule.GenzMalikType
struct GenzMalik{D} <: AbstractRule{Orthotope{D}}

Embedded cubature rule for a D-orthotope of high order 7 and low order 5.

Type Parameters:

  • D: The dimension of the orthotope.
source
HAdaptiveIntegration.Rule.GrundmannMoellerType
struct GrundmannMoeller{D} <: AbstractRule{Simplex{D}}

Embedded cubature rule for a D-simplex.

Type Parameters:

  • D: The dimension of the simplex.

Fields:

  • order_high::Int: the high order of the rule.
  • order_low::Int: the low order of the rule.

Invariants (check at construction):

  • order_high and order_low must be odd.
  • must have order_high > order_low ≥ 1.
source
HAdaptiveIntegration.Rule.TabulatedEmbeddedCubatureType
@kwdef struct TabulatedEmbeddedCubature{DOM<:AbstractDomain} <: AbstractRule{DOM}

An embedded cubature rule consisting of a high order cubature rule and a low order cubature rule. Note that the low order cubature uses nodes[1:L] as its nodes where L is the length of the weights_low. The cubature nodes and weights are assumed to be for the reference domain (use the reference_domain function to get the reference domain).

Fields:

  • description::String: description of the embedded cubature.
  • reference::String: where the values are found.
  • order_high::Int: order of the high order cubature.
  • order_low::Int: order of the low order cubature.
  • precision::Int: number of significant digits on the node and weight values, 10^-precision give the relative precision of the values.
  • nodes::Vector{Vector{String}}: the cubature nodes.
  • weights_high::Vector{String}: the cubature weights for the high order cubature.
  • weights_low::Vector{String}: the cubature weights for the low order cubature.

Invariants (check at construction):

  • length(nodes) == length(weights_high)
  • length(weights_high) ≥ length(weights_low)
  • order_high ≥ order_low ≥ 0
  • precision ≥ 0
source
HAdaptiveIntegration.allocate_bufferMethod
allocate_buffer(
    fct,
    domain::AbstractDomain{D,T};
    embedded_cubature::EmbeddedCubature{D,T}=default_embedded_cubature(domain),
    norm=LinearAlgebra.norm
) where {D,T}

Allocate and return a buffer that can be passed to the integrate function to improve performance by reducing memory allocations when integrate is called multiple times.

source
HAdaptiveIntegration.increase_precisionFunction
increase_precision(
    tec::TabulatedEmbeddedCubature,
    ::Type{T};
    x_atol=10 * eps(T),
    f_atol=10 * eps(T),
    maxiter::Int=16,
)

Increase the precision of a TabulatedEmbeddedCubature to match the precision of the floating point type T. This is an extension, you need using ForwardDiff for using it.

Arguments

  • tec::TabulatedEmbeddedCubature: the tabulated embedded cubature to increase the precision of.
  • ::Type{T}: the floating point type to increase the precision to.

Optional arguments

  • x_atol=10 * eps(T): the absolute tolerance for the change in the variables (nodes and weights) between Newton iterations.
  • f_atol=10 * eps(T): the absolute tolerance for the change in the function values (integral constraints) between Newton iterations.
  • maxiter::Int=16: the maximum number of Newton iterations to perform.
source
HAdaptiveIntegration.integrateMethod
integrate(
    fct,
    domain::AbstractDomain{D,T};
    embedded_cubature::EmbeddedCubature{D,T}=default_embedded_cubature(domain),
    subdiv_algo=default_subdivision(domain),
    buffer=nothing,
    norm=LinearAlgebra.norm,
    atol=zero(T),
    rtol=(atol > zero(T)) ? zero(T) : sqrt(eps(T)),
    maxsubdiv=2^(13 + D),
    callback=(I, E, nb_subdiv, buffer) -> nothing,
) where {D,T}

Return I and E where I is the integral of the function fct over domain and E is an error estimate.

Arguments

  • fct: a function that must take a SVector{D,T} to a return type K, with K must support the multiplication by a scalar of type T and the addition.
  • domain::AbstractDomain{D,T}: the integration domain. Currently, we support Segment, Triangle, Rectangle, Tetrahedron, Cuboid, d-dimensional Simplex, and d-dimensional Orthotope.

Optional arguments

  • embedded_cubature::EmbeddedCubature{D,T}=default_embedded_cubature(domain): the embedded cubature, each supported domain has a default_embedded_cubature.
  • subdiv_algo=default_subdivision(domain): the subdivision algorithm, each domain has a default_subdivision.
  • buffer=nothing: heap use to do the adaptive algorithm, can be allocated using allocate_buffer, which might result in performance gain if multiple call to integrate is perform.
  • norm=LinearAlgebra.norm: norm used to estimate the error.
  • atol=zero(T): absolute tolerance.
  • rtol=(atol > zero(T)) ? zero(T) : sqrt(eps(T)): relative tolerance.
  • maxsubdiv=2^(13 + D): maximum number of subdivision.
  • callback=(I, E, nb_subdiv, buffer) -> nothing: a callback function called for each estimated value of I and E, including the initial estimate (nb_subdiv=0) and after each subdivision. The callback receives the current integral I, error estimate E, number of subdivisions nb_subdiv, and buffer.
source
HAdaptiveIntegration.resumMethod
resum(buffer; norm=LinearAlgebra.norm)

Re-sum the integral and error estimate from a provided buffer. This function is more expensive than a sum because it uses the Kahan-Babuška-Neumaier [1] summation algorithm to reduce numerical error due to floating-point numbers.

[1] Klein, A. A Generalized Kahan-Babuška-Summation-Algorithm. Computing 76, 279-293 (2006). https://doi.org/10.1007/s00607-005-0139-x

source