Index

@suitesplines_reexport(pkgs...)

Convenience macro to reexport submodules of SuiteSplines.

• If no arguments are passed, it reexports all submodules in SuiteSplines.SUITESPLINES_PKGS.
• If one or more Symbols are passed (e.g., UnivariateSplines), only those specific submodules are reexported.

This macro expands to a block of @reexport using SuiteSplines.PkgName statements.

Vector of package names included in the SuiteSplines.jl bundle (without .jl).

nullspace(A::AbstractMatrix{T<:Real})

Construct a left basis for the null space of matrix A.

nullspace(a::AbstractVector{T}) where {T<:Real}

Find a sparse basis for the null space of a vector. The entries are computed such that null space operator is maximally sparse and column sum of absolute values equals 1.0 for improved conditioning.

boundary(s; direction=1, component=1)
∂(s; direction=1, component=1)

Return boundary component in direction of the object s.

replace_ref_begin_end!(ex)

Recursively replace occurrences of the symbols :begin and :end in a "ref" expression (i.e. A[...]) ex with the appropriate function calls (firstindex or lastindex). Replacement uses the closest enclosing ref, so A[B[end]] should transform to A[B[lastindex(B)]]

AbstractMapping{Dim,S1,S2}

Abstract type representing a mapping with an Dim dimensional domain, taking values in an S1 x S2 codimensional space.

AbstractSpline{Dim} <: ScalarMapping{Dim}

Abstract subtype of ScalarMapping that represents a AbstractSpline function.

Boundary{S<:AbstractMapping}

Iterator over the boundary components of a mapping of type S.

Degree

Polynomial degree (typeallias of Integer).

Dimension

Dimension (typeallias of Integer).

Regularity

Regularity of the B-spline basis (typeallias of Integer).

ScalarMapping{Dim} <: AbstractMapping{Dim,1,1}

Abstract subtype of AbstractMapping` that represents scalar-valued mappings.

insert!(v, index, item)

The Base.Insert! function is extended to operate on NonDecreasingVector and IncreasingVector. An argument-error will be thrown if the item leads to an unsorted sequence.

Examples:

julia> v = NonDecreasingVector([1,2,3,4]);

julia> insert!(v, 3, 2)
5-element Vector{Int64}:
 1
 2
 2
 3
 4

julia> v = IncreasingVector([1,2,4,5]);

julia> insert!(v, 3, 3)
5-element Vector{Int64}:
 1
 2
 3
 4
 5

construct_vector(u, m)

Construct a vector by repeating the elements in u m times.

Examples:

julia> u = [0.0,1.0,2.5,3.0]
4-element Vector{Float64}:
 0.0
 1.0
 2.5
 3.0

julia> m = [2,1,2,3]
4-element Vector{Int64}:
 2
 1
 2
 3

julia> construct_vector(u, m)
8-element Vector{Float64}:
 0.0
 0.0
 1.0
 2.5
 2.5
 3.0
 3.0
 3.0

deconstruct_vector(v)

Decompose a sequence into a new sequence, a multiplicity vector, and vectors specifying indexing into the sequence.

Examples:

julia> v = [0.0,0.0,1.0,2.0,2.0,3.0,3.0,3.0];

julia> u, m, ia, ic = deconstruct_vector(v)
([0.0, 1.0, 2.0, 3.0], [2, 1, 2, 3], [1, 1, 2, 3, 3, 4, 4, 4], [2, 3, 5, 8])

julia> construct_vector(u,m)==v
true

julia> u[ia]==v
true

julia> v[ic]==u
true

global_insert(v, k)

Uniformly subdivide new values into an IncreasingVector or NonDecreasingVector.

Examples:

julia> v = IncreasingVector([0.0,1.0,2.0])
3-element IncreasingVector{Float64}:
 0.0
 1.0
 2.0

julia> global_insert(v, 3)
9-element IncreasingVector{Float64}:
 0.0
 0.25
 0.5
 0.75
 1.0
 1.25
 1.5
 1.75
 2.0

When applied to a NonDecreasingVector only the non-zero length intervals are subdivided

julia> v = NonDecreasingVector([0.0,0.0,1.0,1.0,2.0])
5-element NonDecreasingVector{Float64}:
 0.0
 0.0
 1.0
 1.0
 2.0

julia> global_insert(v, 3)
11-element NonDecreasingVector{Float64}:
 0.0
 0.0
 0.25
 0.5
 0.75
 1.0
 1.0
 1.25
 1.5
 1.75
 2.0

IncreasingRange{T<:Real}

Construct an increasing range of real numbers.

IncreasingSequence{T}

Abstract type for increasing sequences of real numbers.

IncreasingVector{T<:Real}

Construct a vector with an increasing set of real numbers.

Examples:

The element type T is a subtype of Real. Hence it is possible to make IncreasingVector{Int64} as well as IncreasingVector{Float64}

julia> IncreasingVector([0,1,2,3])
4-element IncreasingVector{Int64}:
 0
 1
 2
 3

julia> IncreasingVector([0.0,1.0,2.0,3.0])
4-element IncreasingVector{Float64}:
 0.0
 1.0
 2.0
 3.0

It is also possible to extract the unique values of a NonDecreasingVector into an IncreasingVector.

julia> v = NonDecreasingVector([0.0,0.0,1.0,2.0,2.0,3.0,3.0])
7-element NonDecreasingVector{Float64}:
 0.0
 0.0
 1.0
 2.0
 2.0
 3.0
 3.0

julia> IncreasingVector(v)
4-element IncreasingVector{Float64}:
 0.0
 1.0
 2.0
 3.0

Interval(a::T, b::T)

Create an interval from a point a to a point b.

NonDecreasingSequence{T}

Abstract type for non-decreasing sequences of real numbers.

NonDecreasingVector{T}

Construct a vector with a non-decreasing set of real numbers.

Examples:

julia> v = NonDecreasingVector([0.0,0.0,1.0,2.0,2.0])
5-element NonDecreasingVector{Float64}:
 0.0
 0.0
 1.0
 2.0
 2.0

julia> u, m = deconstruct_vector(v)
([0.0, 1.0, 2.0], [2, 1, 2], [1, 1, 2, 3, 3], [2, 3, 5])

julia> u
3-element IncreasingVector{Float64}:
 0.0
 1.0
 2.0

julia> NonDecreasingVector(u, m) == v
true

It is also possible to extract the unique values of a NonDecreasingVector into an IncreasingVector.

julia> v = NonDecreasingVector([0.0,0.0,1.0,2.0,2.0,3.0,3.0])
7-element NonDecreasingVector{Float64}:
 0.0
 0.0
 1.0
 2.0
 2.0
 3.0
 3.0

julia> IncreasingVector(v)
4-element IncreasingVector{Float64}:
 0.0
 1.0
 2.0
 3.0

SortedSequence{T}

Abstract type for sorted sequences of real numbers.

Unique(v)

Convenience iterator that lazily returns a tuple with the unique consecutive values and their multiplicities.

Examples:

julia> v = [1,1,3,4,5,5,5,6,6,7];

julia> for item in Unique(v)
           @show item
       end
item = (1, 2)
item = (3, 1)
item = (4, 1)
item = (5, 3)
item = (6, 2)
item = (7, 1)

collect(X)

Collect values of CartesianProduct{T,Dim} in Dim-dimensional arrays. Similar to ndgrid in Matlab.

Examples:

julia> x, y = collect([1,2] ⨱ [2,3,4]);

julia> x
2×3 Matrix{Int64}:
 1  1  1
 2  2  2

julia> y
2×3 Matrix{Int64}:
 2  3  4
 2  3  4

adjoint(K::BoxProduct)

Compute the adjoint of a Box product.

adjoint(K::KroneckerProduct)

Compute the adjoint of a Kronecker product.

inv(K::BoxProduct)

Compute the inverse of a Box product.

inv(K::KroneckerProduct)

Compute the inverse of a Kronecker product.

collect(K::KroneckerProduct)

Collects a lazy instance of the KroneckerProduct type into a full, native matrix. Equivalent with Matrix(K::KroneckerProduct).

transpose(K::BoxProduct)

Compute the transpose of a Box product.

transpose(K::KroneckerProduct)

Compute the transpose of a Kronecker product.

det(K::BoxProduct)

Compute the determinant of a Box product.

det(K::KroneckerProduct)

Compute the determinant of a Kronecker product.

eigen(K::KroneckerProduct)

If the matrices of a KroneckerProduct instance are square, performs the eigenvalue decompositon on them and return an Eigen type.

isposdef(K::KroneckerProduct)

Test whether a Kronecker product is positive definite (and Hermitian) by trying to perform a Cholesky factorization of K.

issymmetric(K::KroneckerProduct)

Checks if a Kronecker product is symmetric.

lmul!(a::Number, K::KroneckerProduct)

Scale a KroneckerProduct K inplace by a factor a by rescaling the left matrix.

logdet(K::KroneckerProduct)

Compute the logarithm of the determinant of a Kronecker product.

rmul!(K::KroneckerProduct, a::Number)

Scale a KroneckerProduct K inplace by a factor a by rescaling the right matrix.

tr(K::KroneckerProduct)

Compute the trace of a Kronecker product.

order(M::AbstractMatrix)

Returns the order of a matrix, i.e. how many matrices are involved in the Kronecker product (default to 1 for general matrices).

box(A, B)
box(A, B...)

Compute the box product, which is defined as a perfect shuffle of a kronecker product.

issquare(A::AbstractMatrix)

Checks if an array is a square matrix.

Matrix(K::KroneckerProduct)

Converts a GeneralizedKroneckerProduct instance to a Matrix type.

BoxProduct{M,T,N} <: AbstractKron{M,T,N}

The box product is defined in the following paper.

Olsen, Peder A., Steven J. Rennie, and Vaibhava Goel. "Efficient automatic differentiation of matrix functions." In Recent Advances in Algorithmic Differentiation, pp. 71-81. Springer, Berlin, Heidelberg, 2012.

KroneckerProduct{T,N,S}  <: AbstractVecOrMat{T}

Matrix type that stores the matrices in a Kronecker product and lazily evaluates products, sums, etc.

l2_error(f; to, relative=false, quadrule=standard_quadrature_rule(f,g))

Compute the (relative) L^2 error of f with respect to g. Any f and g will do as long as @cartesian and standard_quadrature_rule(f,g) are implemented.

EvaluationSet{S1,S2, S}

Set of evaluation points that are used as input and output of functions, mappings and fields.

ResizableArray{T,Dim} <: AbstractArray{T, Dim}

Array-type that may be resized. New resources will be allocated if the new array-length exceeds maxlen.

@evaluate! Y = f(x)
@evaluate! Y += f(x)
@evaluate! Y -= f(x)
@evaluate! Y *= f(x)
@evaluate! Y /= f(x)

Fast update routine for evaluating a mapping on a grid $x ...$. array Y needs to be preallocated.

@evaluate Y  = f(x)

Fast evaluation routine for evaluating a mapping on a grid $x ...$. The result is stored in array $Y$. Allocation is done automatically and EvaluationSet's of the correct size and type are cached for each mapping using an LRUCache.

Binomial(Val(n), Val(k)) = binomial(n,k)

@generated function for optimized computation of binomials.

b2c(simplex)

Return the mapping (matrix) that maps barycentric coordinates to Cartesian coordinates.

barycentriccoords(S::AbstractSimplex, x)

Compute the barycentric coordinates of a point, or matrix of points x.

bernstein(p::Degree, control::AbstractVector{T}, λ::BPoint{D,T})

Compute Bezier polynomial at barycentric point λ with the degrees of freedom in control.

bernstein(X::MultiIndex{D}, λ::BPoint{D,T})

Compute single Bernstein basis function corresponding to multi-index X at barycentric point λ.

cartesiancoords(S::AbstractSimplex, λ)

Compute the Cartesian coordinates of a barycentric point or vector of barycentric points λ.

decasteljau!(control::Vector{T}, λ::BPoint{D,T}, MultiIndices)

Classic DeCasteljau algorithm (can compute derivatives too!) on a (D-1)-dimensional simplex.

degree_elevation_operator(p::Degree, d::Dimension)

Return a matrix that maps degrees of freedom of a Bezier polynomial of degree p to the degrees of freedom corresponding to its order elevated Bezier polynomial.

grevillepoint(X::MultiIndex)

Compute Greville point associated with multi-index label X.

jacobian(simplex)

Return the Jacobian mapping of a simplex.

kformbasisconv(n::Dimension, k::Form, r::Degree)

Compute the linear mapping from Bernstein polynomials to the space of k-forms Pᵣ⁻Λᵏ(T). This linear mapping is encoded in Pattern{T}. We refer to the paper. For more information, check out the periodic table of finite elements

linindex(i, j, k...)

Return the linear index of a multi-index label.

mass_matrix(S::AbstractSimplex{D}, p::Degree, q::Degree, T::Type{<:Real} = Type{Float64})

Compute Bernstein mass matrix for degree (p,q) and dimension d.

MultiIndex(::LinearIndex{P<:Degree, D<:Dimension, I<:Integer})

Return the multi-index associated with linear index I and simplicial domain of degree P and dimension D.

step_decasteljau!(control::Vector{T}, λ::BPoint{D,T}, MultiIndices)

Computes a single step of the recursive DeCasteljau algorithm on a (D-1)-dimensional simplex.

subsets(v::NTuple{n,Int}, k::Int)

Generate the binomial(n,k) unique subsets of length k from set v and its complement.

Examples:

julia> α, β = subsets((1,2,3), 1);

julia> α
3-element Array{Tuple{Int64},1}:
 (1,)
 (2,)
 (3,)

julia> β
3-element Array{Tuple{Int64,Int64},1}:
 (2, 3)
 (1, 3)
 (1, 2)

vol(simplex)

Compute the volume of an affine simplex.

dimension(simplex)

Return the dimension of a simplex.

dimsplinespace(p::Degree, d::Dimension)

Compute the dimension of the Bezier spline space on a simplex of degree p and numer of vertices d.

 BezierSimplex{D,T<:Real} <: AbstractSimplex{D}

Simplex of dimension D with Vertices in real number type T.

LinearIndex(P::Degree, D::Dimension, I::Integer)

A linear index I of a simplicial discretization with Degree P and dimension D.

LinearIndex(::MultiIndex{X})

Compute the LinearIndex associated with a MultiIndex.

MultiIndices(P::Degree, D::Dimension)
MultiIndices(D::Dimension)
MultiIndices(s::Simplex)
MultiIndices(s::BezierSimplex)

Returns an iterator corresponding to a Bezier simplex of polynomial degree P and dimension D.

 Simplex{D,T<:Real} <: AbstractSimplex{D}

Simplex of dimension D with Vertices in real number type T.

dimsplinespace(p, U)

Compute the dimension of the spline space defined by degree $p$ and knotsvector $U$.

refine(p::Degree, R::AbstractRefinement)

Change in polynomial degree based on the refinement type.

refine(p::KnotVector, R::AbstractRefinement)

Change in knot vector based on the refinement type.

approximate_collocation_inverse(p, U)
approximate_collocation_inverse(p, U, k1)

Compute an approximate inverse $\mathsf{A} ≈ \mathsf{B}^{-1}$, where $\mathsf{B}$ is the consistent B-spline interpolation matrix. The approximation is of order $k_1 \leq p+1$. The quasi-interpolant is designed to reproduce polynomials

\[ \sum_{i=1}^n \mu_i(x^k) B_{i,p}(x) = x^k  \quad \text{for } k=0,1,...,k_1-1\]

and is based on the following paper

T. Lyche and L. L. Schumaker. Local spline approximation methods. Journal of Approximation Theory, 15(4):294-325, Dec. 1975.

Examples:

julia> p  = Degree(3);

julia> U  = KnotVector([0.0,1.0,1.5,2.5,3.0], [p+1,1,1,2,p+1]);

julia> y  = grevillepoints(p, U);

The approximate collocation inverse, using full approximation order $k_1=p+1$ and the consistent collocation matrix are computed as follows

julia> A = approximate_collocation_inverse(p, U, p+1);

julia> B = bspline_interpolation_matrix(p, U, y, 1)[1];

It can be verified that the quasi-interpolant reproduces constants, linears, quadratic and cubic polynomials

julia> B * A * (y.^0) ≈ y.^0
true

julia> B * A * (y.^1) ≈ y.^1
true

julia> B * A * (y.^2) ≈ y.^2
true

julia> B * A * (y.^3) ≈ y.^3
true

approximate_l2_inverse(p, U, k1)

Compute an approximate inverse $\mathsf{S} ≈ \mathsf{M}^{-1}$, where $\mathsf{M}$ is the consisten B-spline mass matrix. The approximation is of order $k_1 \leq p+1$. The $L^2$ quasi-interpolant is designed to reproduce polynomials

\[ \sum_{i=1}^n \mu_i(x^k) B_{i,p}(x) = x^k  \quad \text{for } k=0,1,...,k_1-1\]

and is based on the following paper

Chui, Charles K., Wenjie He, and Joachim Stöckler. "Nonstationary tight wavelet frames, I: Bounded intervals." Applied and Computational Harmonic Analysis 17.2 (2004): 141-197.

bezier_extraction_operator(p, U)

Compute the Bézier extraction operators corresponding to a B-spline basis of polynomial degree $p$ and knot vector $U$. The output is a 3-dimensional array $\mathsf{C} \in \mathbb{R}^{(p+1) \times (p+1) \times n_{el}}$.

The implementation is based on the following paper

Borden, Michael J., et al. "Isogeometric finite element data structures based on Bézier extraction of NURBS." International Journal for Numerical Methods in Engineering 87.1‐5 (2011): 15-47.

Examples:

julia> p = Degree(2); U = KnotVector([0.0,0.0,0.0,1.0,3.0,3.0,4.0,4.0,4.0]);

julia> bezier_extraction_operator(p, U)
3×3×3 Array{Float64, 3}:
[:, :, 1] =
 1.0  0.0       0.0
 0.0  1.0       0.0
 0.0  0.666667  0.333333

[:, :, 2] =
 0.666667  0.333333  0.0
 0.0       1.0       0.0
 0.0       0.0       1.0

[:, :, 3] =
 1.0  0.0  0.0
 0.0  1.0  0.0
 0.0  0.0  1.0

bspline_integral_value(p, U, i::Int)
bspline_integral_value(p, U)

Compute the integral of a b-spline function defined by degree p and knotvector U.

bspline_interpolation_matrix(p, U, i, x, nout=1)

Given the polynomial degree $p$ and knot-vector $U$ compute the collocation matrices of the B-spline basis functions sampled at a set of points $x$. The output is a set of sparse matrices whose rows correspond to the evaluation points and whose columns are associated with the B-spline degrees of freedom.

Examples:

julia> p = Degree(2);

julia> U = KnotVector([0.0,1.0,1.5,2.5,3.0], [3,1,1,2,3]);

julia> x = IncreasingVector([0.5,1.5,2.75]);

julia> B = bspline_interpolation_matrix(p, U, x, p+1);

The interpolation matrix for quadratic B-splines is

julia> B[1]
3×7 SparseArrays.SparseMatrixCSC{Float64, Int64} with 9 stored entries:
 0.25  0.583333  0.166667   ⋅         ⋅     ⋅    ⋅ 
  ⋅     ⋅        0.666667  0.333333  0.0    ⋅    ⋅ 
  ⋅     ⋅         ⋅         ⋅        0.25  0.5  0.25

The interpolation matrix for linear and constant B-splines is given by

julia> B[2]
3×6 SparseArrays.SparseMatrixCSC{Float64, Int64} with 6 stored entries:
 0.5  0.5   ⋅    ⋅    ⋅    ⋅ 
  ⋅    ⋅   1.0  0.0   ⋅    ⋅ 
  ⋅    ⋅    ⋅    ⋅   0.5  0.5

julia> B[3]
3×5 SparseArrays.SparseMatrixCSC{Float64, Int64} with 3 stored entries:
 1.0   ⋅    ⋅    ⋅    ⋅ 
  ⋅    ⋅   1.0   ⋅    ⋅ 
  ⋅    ⋅    ⋅    ⋅   1.0

It can be easily verified that the B-splines satisfy a partition of unity

julia> sum(B[1],dims=2)
3×1 Matrix{Float64}:
 1.0
 1.0
 1.0

julia> sum(B[2],dims=2)
3×1 Matrix{Float64}:
 1.0
 1.0
 1.0

julia> sum(B[3],dims=2)
3×1 Matrix{Float64}:
 1.0
 1.0
 1.0

The collocation matrices can be used to solve an interpolation problem or a PDE using the collocation method. For example, the interpolation matrix evaluated at the Greville points can be computed as follows

julia> x = grevillepoints(p, U);

julia> B = bspline_interpolation_matrix(p, U, x, p+1)[1]
7×7 SparseArrays.SparseMatrixCSC{Float64, Int64} with 21 stored entries:
 1.0   0.0        0.0        ⋅          ⋅     ⋅    ⋅ 
 0.25  0.583333   0.166667   ⋅          ⋅     ⋅    ⋅ 
  ⋅    0.0833333  0.833333  0.0833333   ⋅     ⋅    ⋅ 
  ⋅     ⋅         0.166667  0.583333   0.25   ⋅    ⋅ 
  ⋅     ⋅          ⋅         ⋅         1.0   0.0  0.0
  ⋅     ⋅          ⋅         ⋅         0.25  0.5  0.25
  ⋅     ⋅          ⋅         ⋅         0.0   0.0  1.0

The matrix $\mathsf{B}$ can then be used to interpolate a univariate function. For example, up to quadratic polynomials are exactly reproduced

# sample grid to test polynomial reproduction
julia> y = global_insert(x, 4);

julia> C = bspline_interpolation_matrix(p, U, y, 1)[1];

julia> C * (B \ x.^0) ≈ y.^0
true

julia> C * (B \ x.^1) ≈ y.^1
true

julia> C * (B \ x.^2) ≈ y.^2
true

Here matrix $\mathsf{C}$ evaluates the B-spline at a refined partition $y$.

bsplinebasisfuns(p, U, i::Integer, x::T, nout=1)
bsplinebasisfuns(p, U, x::T, nout=1)

Compute the $p+1$ non-zero B-spline basis-functions at site $x \in [U_{i},U_{i+1})$. The output is a matrix where column $j$ (rows $1:p+2-j$), for $j=1,...,n_{out}$ correspond to the $p+2-j$ non-zero B-spline basis functions of degree $p+1-j$. The remaining entries correspond to supports. If the active knot-span is not provided then it will be computed.

Examples:

julia> p = Degree(2);

julia> U = KnotVector([0.0,0.0,0.0,1.0,2.0,2.0,2.0]);

julia> x = 1.5
1.5

julia> span = findspan(p, U, x)
4

julia> B = bsplinebasisfuns(p, U, span, x, p+1)
3×3 Matrix{Float64}:
 0.125  0.5  1.0
 0.625  0.5  1.0
 0.25   1.0  2.0

The knot span index can be omited, in which case it is computed on the fly.

julia> bsplinebasisfuns(p, U, x, p+1)
3×3 Matrix{Float64}:
 0.125  0.5  1.0
 0.625  0.5  1.0
 0.25   1.0  2.0

By prescribing a vector of $m$ sites $x \in \mathbb{R}^m$ we can compute the B-splines at more than one point. The function can be called with or without prescribing the knot-span index

bsplinebasisfuns(p, U, span::Vector{Integer}, u::Vector{T}, nout=1)
bsplinebasisfuns(p, U, u::SortedSequence{T}, nout=1)

The output is a 3-dimensional array $B \in \mathbb{R}^{(p+1) \times m \times n_{out}}$.

Examples:

julia> x = IncreasingVector([0.5, 1.5]);

julia> span = findspan(p, U, x)
2-element NonDecreasingVector{Int64}:
 3
 4

julia> bsplinebasisfuns(p, U, span, x, p+1)
3×2×3 Array{Float64, 3}:
[:, :, 1] =
 0.25   0.125
 0.625  0.625
 0.125  0.25

[:, :, 2] =
 0.5  0.5
 0.5  0.5
 2.0  1.0

[:, :, 3] =
 1.0  1.0
 1.0  1.0
 1.0  2.0

Alternatively, the non-vanishing B-splines up to and including degree $p$ can be computed at $n_{el}$ non-zero knot-spans (elements) with $m$ sites each, that is, $x \in \mathbb{R}^{m \times n_{el}}$. Again, the function can be called with or without prescribing the knot-span index

bsplinebasisfuns(p, U, span::Vector{Int64}, u::Matrix{T}, nout=1)
bsplinebasisfuns(p, U, u::Matrix{T}, nout=1)

The output is a 4-dimensional array $B \in \mathbb{R}^{(p+1) \times m \times n_{el} \times n_{out}}$.

Examples:

julia> x = [0.25 1.25;
            0.75 1.75];

julia> bsplinebasisfuns(p, U, x, p+1)
3×2×2×3 Array{Float64, 4}:
[:, :, 1, 1] =
 0.5625   0.0625
 0.40625  0.65625
 0.03125  0.28125

[:, :, 2, 1] =
 0.28125  0.03125
 0.65625  0.40625
 0.0625   0.5625

[:, :, 1, 2] =
 0.75  0.25
 0.25  0.75
 2.0   2.0

[:, :, 2, 2] =
 0.75  0.25
 0.25  0.75
 1.0   1.0

[:, :, 1, 3] =
 1.0  1.0
 1.0  1.0
 1.0  1.0

[:, :, 2, 3] =
 1.0  1.0
 1.0  1.0
 2.0  2.0

ders_bspline_interpolation_matrix(p, U, i, x, nout=1)

Given the polynomial degree $p$ and knot-vector $U$ compute the collocation matrices of the first $n_{out}$ derivatives of the B-spline basis functions sampled at a set of points $x$. The output is a set of sparse matrices whose rows correspond to the evaluation points and whose columns are associated with the B-spline degrees of freedom.

Examples:

julia> p = Degree(2);

julia> U = KnotVector([0.0,1.0,1.5,2.5,3.0], [3,1,1,2,3]);

julia> x = IncreasingVector([0.5,1.5,2.75]);

julia> B = ders_bspline_interpolation_matrix(p, U, x, p+1);

The output is a set of SparseMatrixCSC matrices. For example the B-splines evaluate at points $x$ are

julia> B[1]
3×7 SparseArrays.SparseMatrixCSC{Float64, Int64} with 9 stored entries:                       
 0.25  0.583333  0.166667   ⋅         ⋅     ⋅    ⋅                                            
  ⋅     ⋅        0.666667  0.333333  0.0    ⋅    ⋅                                            
  ⋅     ⋅         ⋅         ⋅        0.25  0.5  0.25                                          

Derivatives of quadratic B-splines evaluated at $x$

julia> B[2]
3×7 SparseArrays.SparseMatrixCSC{Float64, Int64} with 9 stored entries:
 -1.0  0.333333   0.666667   ⋅         ⋅    ⋅    ⋅ 
   ⋅    ⋅        -1.33333   1.33333   0.0   ⋅    ⋅ 
   ⋅    ⋅          ⋅         ⋅       -2.0  0.0  2.0

Second derivatives of quadratic B-splines evaluated at $x$

julia> B[3]
3×7 SparseArrays.SparseMatrixCSC{Float64, Int64} with 9 stored entries:
 2.0  -3.33333  1.33333    ⋅        ⋅      ⋅    ⋅ 
  ⋅     ⋅       1.33333  -3.33333  2.0     ⋅    ⋅ 
  ⋅     ⋅        ⋅         ⋅       8.0  -16.0  8.0

dersbsplinebasisfuns(p, U, span::Int, u::T, nout=1)
dersbsplinebasisfuns(p, U, u::T, nout=1)

Compute the non-vanishing B-spline basis-functions of degree $p$ and their $n_{out}-1$'th order derivatives at the site $x$. The output is a matrix $B \in \mathbb{R}^{(p+1) \times n_{out}}$ where the $i$'th row corresponds to the $i$'th non-vanishing B-spline at $x$ and the $j$'th column corresponds to the $j-1$'th derivative.

The implementation of this function is based on

Piegl, Les, and Wayne Tiller. The NURBS book. Springer Science & Business Media, 2012.

Examples:

julia> p = Degree(2);

julia> U = KnotVector([0.0,0.0,0.0,1.0,2.0,2.0,2.0]);

julia> x = 1.5;

julia> span = findspan(p, U, x)
4

julia> dersbsplinebasisfuns(p, U, span, x, p+1)
3×3 Matrix{Float64}:
 0.125  -0.5   1.0
 0.625  -0.5  -3.0
 0.25    1.0   2.0

dersbsplinebasisfuns can also be called without providing the span-indices, in which case they will be calculated on the fly.

julia> dersbsplinebasisfuns(p, U, x, p+1)
3×3 Matrix{Float64}:
 0.125  -0.5   1.0
 0.625  -0.5  -3.0
 0.25    1.0   2.0

By prescribing a vector of $m$ sites $x \in \mathbb{R}^m$ we can compute the B-splines and their derivatives at more than one point. The function can be called with or without prescribing the knot-span index

dersbsplinebasisfuns(p, U, span::Vector{Integer}, u::Vector{T}, nout=1)
dersbsplinebasisfuns(p, U, u::SortedSequence{T}, nout=1)

The output is a 3-dimensional array $B \in \mathbb{R}^{(p+1) \times m \times n_{out}}$.

Examples:

julia> x = IncreasingVector([0.5, 1.5]);

julia> span = findspan(p, U, x)
2-element NonDecreasingVector{Int64}:
 3
 4

julia> dersbsplinebasisfuns(p, U, span, x, p+1)
3×2×3 Array{Float64, 3}:
[:, :, 1] =
 0.25   0.125
 0.625  0.625
 0.125  0.25

[:, :, 2] =
 -1.0  -0.5
  0.5  -0.5
  0.5   1.0

[:, :, 3] =
  2.0   1.0
 -3.0  -3.0
  1.0   2.0

Alternatively, the non-vanishing derivatives of the B-splines can be computed at $n_{el}$ non-zero knot-spans (elements) with $m$ sites each, that is, $x \in \mathbb{R}^{m \times n_{el}}$. Again, the function can be called with or without prescribing the knot-span index

dersbsplinebasisfuns(p, U, span::Vector{Int64}, u::Matrix{T}, nout=1)
dersbsplinebasisfuns(p, U, u::Matrix{T}, nout=1)

The output is a 4-dimensional array $B \in \mathbb{R}^{(p+1) \times m \times n_{el} \times n_{out}}$.

Examples:

julia> x = [0.25 1.25;
            0.75 1.75];

julia> dersbsplinebasisfuns(p, U, x, p+1)
3×2×2×3 Array{Float64, 4}:
[:, :, 1, 1] =
 0.5625   0.0625
 0.40625  0.65625
 0.03125  0.28125

[:, :, 2, 1] =
 0.28125  0.03125
 0.65625  0.40625
 0.0625   0.5625

[:, :, 1, 2] =
 -1.5   -0.5
  1.25  -0.25
  0.25   0.75

[:, :, 2, 2] =
 -0.75  -0.25
  0.25  -1.25
  0.5    1.5

[:, :, 1, 3] =
  2.0   2.0
 -3.0  -3.0
  1.0   1.0

[:, :, 2, 3] =
  1.0   1.0
 -3.0  -3.0
  2.0   2.0

distribute_points(S::SplineSpace, V::SplineSpace; add_boundary_points::Bool=false)

Compute distribution of quadrature points in weigthed quadrature with a test-space $V$ and target space for quadrature $S$. The nodes are distribute in such a way that all conditions for exact quadrature are satisfied with a minimum number of points.

The distribution of quadrature points in weighted quadrature is based on the following paper.

Hiemstra, René R., et al. "Fast formation and assembly of finite element matrices with application to isogeometric linear elasticity." Computer Methods in Applied Mechanics and Engineering 355 (2019): 234-260.

Examples:

julia> S = SplineSpace(2, IncreasingVector([0.0,1.0,2.0,3.0]), [3,2,2,3]);

julia> V = SplineSpace(2, IncreasingVector([0.0,1.0,2.0,3.0]), [3,1,1,3]);

julia> UnivariateSplines.distribute_points(S, V)
8-element IncreasingVector{Float64}:
 0.16666666666666666
 0.5
 0.8333333333333334
 1.25
 1.75
 2.1666666666666665
 2.5
 2.8333333333333335

extraction_operator(p::Degree, U::KnotVector; [cperiodic, cleft, cright])

Compute an extraction operator that creates a subspace of the SplineSpace with for example peroidicity built-in.

findspan(p, U, u)

Given polynomial degree $p$ and knotvector $U$ determine the knot span index $i$ of a point $x$ such that $u \in [U_i, U_{i+1})$.

Examples:

julia> p = Degree(2)
2

julia> U = KnotVector([0.0,1.0,1.5,2.5,3.0],[3,1,1,2,3])
10-element NonDecreasingVector{Float64}:
 0.0
 0.0
 0.0
 1.0
 1.5
 2.5
 2.5
 3.0
 3.0
 3.0

julia> x = 2.0
2.0

julia> span = findspan(p, U, x)
5

julia> U[span] ≤ x < U[span+1]
true

It is also possible to compute the knot-spans of all values in an 'IncreasingVector{T}'

julia> x = IncreasingVector([0.25, 0.75, 1.0, 1.25, 2.5, 2.75]);

julia> span = findspan(p, U, x)
6-element NonDecreasingVector{Int64}:
 3
 3
 4
 4
 7
 7

grevillepoints(p, U, i)
grevillepoints(p, U)

Compute the $i$'th Greville Absissa corresponding to the $i$'th B-spline defined by degree 'p' and knotvector 'U' of spline space.

Examples:

julia> p = Degree(2);

julia> U = KnotVector([0.0,1.0,2.5,3.0], [3,1,2,1])
7-element NonDecreasingVector{Float64}:
 0.0
 0.0
 0.0
 1.0
 2.5
 2.5
 3.0

julia> grevillepoints(p, U, 2)
0.5

julia> grevillepoints(p, U)
4-element IncreasingVector{Float64}:
 0.0
 0.5
 1.75
 2.5

h_refinement_operator!(p::Degree, U::KnotVector{T}, u::AbstractVector{T}) where {T<:Real}
h_refinement_operator!(p::Degree, U::KnotVector{T}, u::T) where {T<:Real}
h_refinement_operator!(p::Degree, U::KnotVector{T}, span::Integer, u::T) where {T<:Real}

Insert one or more knots $u$ into knotvector $U$ and output the transformation operator from the coarse to the refined space. The knotvector is updated in-place. The output is represented as a SparseMatrixCSC matrix.

Examples:

julia> p = Degree(2)
2

julia> U = KnotVector([0.0,1.0,2.0,3.0,4], [3,1,1,2,3])
10-element NonDecreasingVector{Float64}:
 0.0
 0.0
 0.0
 1.0
 2.0
 3.0
 3.0
 4.0
 4.0
 4.0

julia> C = h_refinement_operator!(p, U, [0.5, 1.5, 2.5, 3.5]); # sparse matrix

julia> Matrix(C)
11×7 Matrix{Float64}:
 1.0  0.0   0.0   0.0   0.0  0.0  0.0
 0.5  0.5   0.0   0.0   0.0  0.0  0.0
 0.0  0.75  0.25  0.0   0.0  0.0  0.0
 0.0  0.25  0.75  0.0   0.0  0.0  0.0
 0.0  0.0   0.75  0.25  0.0  0.0  0.0
 0.0  0.0   0.25  0.75  0.0  0.0  0.0
 0.0  0.0   0.0   0.5   0.5  0.0  0.0
 0.0  0.0   0.0   0.0   1.0  0.0  0.0
 0.0  0.0   0.0   0.0   0.5  0.5  0.0
 0.0  0.0   0.0   0.0   0.0  0.5  0.5
 0.0  0.0   0.0   0.0   0.0  0.0  1.0

julia> @show U;
U = [0.0, 0.0, 0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.0, 3.5, 4.0, 4.0, 4.0]

normalization_weights(X)

Determine a set of normalization weights that result in a more uniform set oftype( weighted quadrature weights when solving for the least norm solution. )

nquadpoints(S::SplineSpace, V::SplineSpace; add_boundary_points=false, add_additional_points=0, min_points=1)

Determine the minimum number of quadrature points corresponding to a weighted quadrature rule with test-space $V$ and target-space $S$.

Examples:

julia> S = SplineSpace(2, IncreasingVector([0.0,1.0,2.0,3.0]), [3,2,2,3]);

julia> V = SplineSpace(2, IncreasingVector([0.0,1.0,2.0,3.0]), [3,1,1,3]);

julia> A = UnivariateSplines.nquadpoints(S, V; add_boundary_points=false)
3-element Vector{Int64}:
 3
 2
 3

num_elements(U::KnotVector)

Count the number of non-zero knot spans or elements in the knotvector.

onebasisfuneval(U, u)

Compute value of a single B-spline basis-function defined by local knotvector 'U' evaluated at point 'x'.

The implementation of this function is taken from

Piegl, Les, and Wayne Tiller. The NURBS book. Springer Science & Business Media, 2012.

refinement_operator(p::Degree, U::KnotVector, method::AbstractRefinement)

Compute the two_scale_operator given p and U and the applicable refinement strategy: IgaBase.hRefinement(), IgaBase.pRefinement(), IgaBase.kRefinement(). See IgaBase package for definition of refinement methods.

solve_leastnorm_without_constraints(A, b, S=BigFloat)

Solve underdetermined system of equation $A x = b$ using QR factorization with arbitrary precision.

table_required_points(S, V)

Determine the number of exactness conditions in subsets of elements. The output is a matrix $\mathsf{A}$ where $\mathsf{A}_{ij}$ denotes the minimum number of points required in interval $i:j$ to perform exact weighted quadrature.

Examples:

julia> S = SplineSpace(2, IncreasingVector([0.0,1.0,2.0,3.0]), [3,2,2,3]);

julia> V = SplineSpace(2, IncreasingVector([0.0,1.0,2.0,3.0]), [3,1,1,3]);

julia> A = UnivariateSplines.table_required_points(S, V)
3×3×1 Array{Int64, 3}:
[:, :, 1] =
 3  5  7
 2  5  0
 3  0  0

two_scale_operator(p::Degree, U::KnotVector{T}, q::Degree, V::KnotVector{T}) where {T<:Real}

The two_scale_operator encodes the main logic for all spline refinemement operations.

The arguments p and U correspond to the current spline space and q and V are the "target" degree knot vector.

unit_integral_rescaling(p, U, i::Int)
unit_integral_rescaling(p, U)

Compute the scaling factors that normalize the B-spline basis functions to have unit integral.

wquadweights!(M, B, W, C)

Compute the quadrature weights in weighted quadrature such that M ≈ W' * B.

Bspline{T} <: ScalarMapping{1}

A univariate B-spline function of field T.

BsplineBasis{T<:Real} <: Basis{T,3}

Special array type for storing evaluating B-spline and their derivatives.

BsplineEvaluationCache <: EvaluationCache{1}

Cache that caches basis function evaluation at a set of evaluation points.

BsplineSupport(S::SplineSpace)

Iterator that outputs the supporting elements of each B-spline basis-function as a UnitRange{Int64}.

Example:

julia> S = SplineSpace(2, [0.0,2.0,3.0], [3,1,3]);

julia> for α in BsplineSupport(S)
           @show α
       end
α = 1:1
α = 1:2
α = 1:2
α = 2:2

GeneralizedGaussrule(p, r, nel, a, b)

Compute a generalized Gaussian quadrature rule for a target spline-space $\mathbb{S}^p_{r}(a,b)$ defined in a uniform partition with $n_{el}$ elements.

These rules are made available in the following paper

Hiemstra, René R., et al. "Optimal and reduced quadrature rules for tensor product and hierarchically refined splines in isogeometric analysis." Computer Methods in Applied Mechanics and Engineering 316 (2017): 966-1004.

KnotSpanIndices

A sorted sequence of integers that refer to the corresponding non-zero knot span of a KnotVector.

KnotVector{T<:Real}

A knot-vector is a non-decreasing sequence of real numbers. (typeallias of NonDecreasingVector{T})

SpanIndex(S::SplineSpace)

Iterator that outputs the span index corresponding to each element in the partition

Example:

julia> S = SplineSpace(2, [0.0,1.0,2.0,3.0,4.0], [3,1,2,1,3]);

julia> for s in SpanIndex(S)
           @show s
       end
s = 3
s = 4
s = 6
s = 7

SplineSpace(p, U)
SplineSpace(p, x, m)
SplineSpace(p, Interval(a, b), num_elements)

Definition of a spline space by means of the polynomial degree $p$, a sequence of break-points $x$ and the knot-multiplicity $m$ or by prescribing the knot-vector $U$.

Examples:

julia> S = SplineSpace(2, [0.0,2.0,3.0], [3,1,3])
SplineSpace(degree = 2, interval = [0.0, 3.0], dimension = 4)

julia> dimsplinespace(S)
4

julia> Degree(S)
2

julia> KnotVector(S)
7-element NonDecreasingVector{Float64}:
 0.0
 0.0
 0.0
 2.0
 3.0
 3.0
 3.0

BsplineSupport(S::SplineSpace)

Iterator that outputs the supporting elements of each B-spline basis-function as a UnitRange{Int64}.

WeightedQuadrule{T<:Real} <: AbstractQuadrule{1}

Weighted quadrature rules are test function specific quadrature rules that exactly integrate all functions in a targetspace. The rules are computed in high precision and are accurate up to 16 digits.

QuadratureRule(Q::TensorProduct{Dim,<:AbstractQuadrule{1}})

Generate a tensor product quadrature rule where the quadrature points are stored as a CartesianProduct and the weights as a KroneckerProduct vector allowing for easy evaluation and numerical integration based on Kronecker products.

TPEvaluationCache{Dim}

Type that caches spline basis functions as a KroneckerProduct{Dim} and a grid of points as a CartesianProduct{Dim}. The type is mutable such that the cached variables can be updated at runtime.

TensorProductBspline{Dim,Ders,T}(space, coeffs)
TensorProductBspline{Dim,Ders,T}(space; codim=1)

Construct a Dim dimensional TensorProductBspline patch that is mapped to Codim dimensional Euclidean space.

QuadraturePoints(acc::ElementAccessor{Dim}, element::Element{Dim}) where Dim

Get the quadrature points on this element.

QuadraturePoints(acc::PatchAccessor{Dim}) where Dim

Get the quadrature points on this element.

QuadratureWeights(acc::ElementAccessor{Dim}, element::Element{Dim}) where Dim

Get the quadrature weights on this element.

QuadratureWeights(acc::PatchAccessor{Dim}) where Dim

Get the quadrature weights on this element.

TestFunctions(acc::ElementAccessor{Dim}, element::Element{Dim}) where Dim

Get the test functions on an element as a KroneckerProduct matrix.

TestFunctions(acc::PatchAccessor{Dim}) where Dim

Get the test functions on the patch as a KroneckerProduct matrix.

TestIndices(acc::PatchtAccessor)

Get global indices of the test functions.

TestIndices(acc::ElementAccessor{Dim}, element::Element{Dim})

Get global indices of the test functions on an element.

TrialFunctions(acc::ElementAccessor{Dim}, element::Element{Dim})

Get the trial functions on an element as a KroneckerProduct matrix.

TrialFunctions(acc::PatchAccessor)

Get the trial functions on the patch as a KroneckerProduct matrix.

TrialIndices(acc::PatchAccessor)

Get global indices of the trial functions.

TrialIndices(acc::ElementAccessor{Dim}, element::Element{Dim})

Get global indices of the trial functions on an element.

coface(element::Element)

If element is on the boundary then element.parent is a view and this function returns its single co-face.

pre_allocate_arrays(T, nquad::NTuple{3,Int}, m::NTuple{3,Int}, n::NTuple{3,Int})

Preallocate arrays of element type==T for use in sum factorization of matrices.

pre_allocate_arrays(T, nquad::NTuple{3,Int}, m::NTuple{3,Int})

Preallocate arrays of element type==T for use in sum factorization of vectors.

pre_allocate_arrays(T, nquad::NTuple{2,Int}, m::NTuple{2,Int}, n::NTuple{2,Int})

Preallocate arrays of element type==T for use in sum factorization of matrices.

pre_allocate_arrays(T, nquad::NTuple{2,Int}, m::NTuple{2,Int})

Preallocate arrays of element type==T for use in sum factorization of vectors.

sumfact!(A::Tuple{<:Array,<:Array,<:Array}, C, testfuns)

Apply sum factorization. A and B are allocated using pre_allocate_arrays(...). testfuns and trialfuns are tuples of arrays that represent the test and trial functions evaluated at the quadrature nodes. It is assumed that the quadrature weights have been incorporated into the test functions.

sumfact!(A::Tuple{<:Array,<:Array,<:Array}, B, C, testfuns, trialfuns)

Apply sum factorization. A and B are allocated using pre_allocate_arrays(...). testfuns and trialfuns are tuples of arrays that represent the test and trial functions evaluated at the quadrature nodes. It is assumed that the quadrature weights have been incorporated into the test functions

sumfact!(A::Tuple{<:Array,<:Array}, C, testfuns)

Apply sum factorization. A is pre-allocated using pre_allocate_rhs_arrays(...). testfuns is a tuple of matrices that represent the test functions evaluated at the quadrature nodes. It is assumed that the quadrature weights have been incorporated into the test functions.

sumfact!(A::Tuple{<:Array,<:Array}, B, C, testfuns, trialfuns)

Apply sum factorization. A and B are allocated using pre_allocate_arrays(...). testfuns and trialfuns are tuples of arays that represent the test and trial functions evaluated at the quadrature nodes. It is assumed that the quadrature weights have been incorporated into the test functions.

QuadratureRule(acc::ElementAccessor{Dim}, element::Element{Dim}) where Dim

Get the quadrature rule on this element.

QuadratureRule(acc::PatchAccessor{Dim}) where Dim

Get the quadrature rule on the patch.

Element{Dim, T}

Dim-dimensional element storing its domain of definition as a Cartesian{Dim, Interval{T}} and its element number as a CartesianIndex{Dim}

ElementAccessor{Dim, X<:IncreasingSequence, Data}

Datastructure that provides access to element data such as the element domain, element trialfunctions, element testfunctions and element quadrature rule.

Elements(partition::CartesianProduct{Dim, <:IncreasingSequence{<:Real}})

Return an iterator over all elements in a partition.

PatchAccessor{Dim, X, ..., Data}

Datastructure that provides access to element data such as the element domain, element trialfunctions, element testfunctions and element quadrature rule.

UnivariatePatchData{T<:Real, Q<:AbstractQuadrule{1}}

Datastructure that stores all the univariate patch data. This data is conveniently accessed by an element accessor.

active_splines(U::TensorProduct{2, <:SplineSpace}, E::Array)

Determine the active B-splines and boundary B-splines given a distance function ϕ.

A[k] == 0 => no physical elements in support of the function A[k] == 1 => boundary function with a physical cut element in its support A[k] == 2 => at least one physical element is fully inside support A[k] == 3 => all physical elements are fully inside support

active_splines(U::TensorProduct{Dim, <:SplineSpace}, E::Array)

Determine the active B-splines and boundary B-splines given a distance function ϕ.

A[k] == 0 => no physical elements in support of the function A[k] == 1 => boundary function with a physical cut element in its support A[k] == 2 => at least one physical element is fully inside support A[k] == 3 => all physical elements are fully inside support

active_splines(U::TensorProduct{2, <:SplineSpace}, E::Array)

Determine the active B-splines and boundary B-splines given a distance function ϕ.

A[k] == 0 => no physical elements in support of the function A[k] == 1 => boundary function with a physical cut element in its support A[k] == 2 => at least one physical element is fully inside support A[k] == 3 => all physical elements are fully inside support

compute_extension_coefficients(p, kts, span, k)

Compute the extension coefficients using dual functionals. This works for general non-uniform knot vectors and is based on the paper [Höllig, Klaus, and Ulrich Reif. "Nonuniform web-splines." Computer Aided Geometric Design 20, no. 5 (2003): 277-294.]

elementisinside(Y::Matrix{Bool})

Return a boolean array that specifies whether an element is inside or outside, respectively

find_extension_indices(active_functions, gperm, mi)

Find the extension indices that are closest to the boundary function with CartesianIndex mi.

find_extension_indices(active_functions, gperm, mi)

Find the extension indices that are closest to the boundary function with CartesianIndex mi.

pointisinside(phi, partition::CartesianProduct)

Return a boolean array that specifies whether a point is inside or outside, respectively

spline_extension_operator(U::TensorProduct{2, <:SplineSpace}, F)

Compute an extension operator that stabilizes the splinespace according to the definition of [Höllig, Klaus, Ulrich Reif, and Joachim Wipper. "Weighted extended B-spline approximation of Dirichlet problems." SIAM Journal on Numerical Analysis 39, no. 2 (2001): 442-462.]

ClosestExtensionArray(F, U::TensorProduct{Dim,<:SplineSpace})

Defines a closest index array, according to the definition of [Höllig, Klaus, Ulrich Reif, and Joachim Wipper. "Weighted extended B-spline approximation of Dirichlet problems." SIAM Journal on Numerical Analysis 39, no. 2 (2001): 442-462.]

ClosestExtensionArray(F, U::TensorProduct{Dim,<:SplineSpace})

Computes a closest index array, according to the definition of [Höllig, Klaus, Ulrich Reif, and Joachim Wipper. "Weighted extended B-spline approximation of Dirichlet problems." SIAM Journal on Numerical Analysis 39, no. 2 (2001): 442-462.]

GeneratePerimeter{Dim,Indices}

Datastructure that allows iteration over active neighborhoods of boundary functions

IPerm{Dim,Indices}

Iterator that returns active neighborhoods around a B-spline.

ImmersedQuadRule(map::AlgoimMapping, xa::Real, ya::Real, xb::Real, yb::Real, qo::Int64)

Compute a algoim quadrature rule in bounding box [xa, ya] × [xb, yb] based on a Gauss-Legendre rule of qo points.