Exported Functions

RadialBasisFunctions.PHS Method

function PHS(n::T=3; poly_deg::T=2) where {T<:Int}

Convienience contructor for polyharmonic splines.

Arguments

• n: Order of the spline (1, 3, 5, or 7). Higher = smoother.

• poly_deg: Polynomial augmentation degree (default: 2 for quadratic).

See also: IMQ, Gaussian

source

RadialBasisFunctions.classify_stencil Method

classify_stencil(is_boundary, boundary_conditions, eval_idx,
                neighbors, global_to_boundary)

Classify stencil type for dispatch in kernel execution.

source

RadialBasisFunctions.custom Method

custom(data, ℒ; basis=PHS(3; poly_deg=2), eval_points=data, k, adjl, hermite)

Build a RadialBasisOperator with a custom operator function.

Arguments

• data: Vector of data points

• ℒ: Custom function that accepts a basis and returns a callable (x, xᵢ) -> value

Keyword Arguments

• basis: RBF basis (default: PHS(3; poly_deg=2))

• eval_points: Evaluation points (default: data)

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation

Examples

# Custom operator that returns the basis function itself
op = custom(data, basis -> (x, xᵢ) -> basis(x, xᵢ))

# Custom biharmonic operator (∇⁴)
op = custom(data, basis -> ∇²(basis) ∘ ∇²(basis))

source

RadialBasisFunctions.directional Method

directional(data, v; basis=PHS(3; poly_deg=2), eval_points=data, k, adjl, hermite)

Build a RadialBasisOperator for the directional derivative (∇f⋅v).

Arguments

• data: Vector of data points

• v: Direction vector. Can be:

  • A single vector of length Dim (constant direction)

  • A vector of vectors matching length(data) (spatially-varying direction)

Keyword Arguments

• basis: RBF basis (default: PHS(3; poly_deg=2))

• eval_points: Evaluation points (default: data)

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation

Examples

# Constant direction
∂_x = directional(data, [1.0, 0.0])

# Spatially-varying direction (e.g., radial)
normals = [normalize(p) for p in data]
∂_n = directional(data, normals)

See also: gradient, partial, laplacian

source

RadialBasisFunctions.gradient Method

gradient(data, x; basis=PHS(3; poly_deg=2), k, adjl)

One-shot convenience function that creates a gradient operator and applies it to scalar field x.

For repeated evaluations, prefer creating the operator once with gradient(data).

Examples

points = [SVector{2}(rand(2)) for _ in 1:1000]
u = sin.(getindex.(points, 1))
∇u = gradient(points, u)  # One-shot gradient computation

source

RadialBasisFunctions.gradient Method

gradient(data; basis=PHS(3; poly_deg=2), eval_points=data, k, adjl, hermite)

Build a RadialBasisOperator for computing gradients of scalar fields.

This is a convenience alias for jacobian. The gradient of a scalar field is mathematically the Jacobian (a 1×D row vector, returned as a length-D vector per point).

Arguments

• data: Vector of points (e.g., Vector{SVector{2,Float64}})

Keyword Arguments

• basis: RBF basis function (default: PHS(3; poly_deg=2))

• eval_points: Evaluation points (default: data)

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation

Examples

points = [SVector{2}(rand(2)) for _ in 1:1000]
op = gradient(points)

u = sin.(getindex.(points, 1))
∇u = op(u)  # Matrix (1000 × 2)
∂u_∂x = ∇u[:, 1]
∂u_∂y = ∇u[:, 2]

See also: jacobian

source

RadialBasisFunctions.jacobian Method

jacobian(data, x; basis=PHS(3; poly_deg=2), k, adjl)

One-shot convenience function that creates a Jacobian operator and applies it to field x.

For repeated evaluations on the same points, prefer creating the operator once with jacobian(data) and calling it via functor syntax op(x).

Arguments

• data: Vector of points

• x: Field values to differentiate

Keyword Arguments

• basis: RBF basis function (default: PHS(3; poly_deg=2))

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

Examples

points = [SVector{2}(rand(2)) for _ in 1:1000]
u = sin.(getindex.(points, 1))
∇u = jacobian(points, u)  # One-shot gradient computation

source

RadialBasisFunctions.jacobian Method

jacobian(data; basis=PHS(3; poly_deg=2), eval_points=data, k, adjl, hermite)

Build a RadialBasisOperator for computing Jacobians (or gradients for scalar fields).

The Jacobian is the fundamental differential operator. For a scalar field, it computes the gradient. For a vector field, it computes the full Jacobian matrix. The spatial dimension is automatically inferred from the data.

Arguments

• data: Vector of points (e.g., Vector{SVector{2,Float64}})

Keyword Arguments

• basis: RBF basis function (default: PHS(3; poly_deg=2))

• eval_points: Evaluation points (default: data)

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation

Examples

points = [SVector{2}(rand(2)) for _ in 1:1000]
op = jacobian(points)

# Scalar field → gradient
u = sin.(getindex.(points, 1))
∇u = op(u)  # Matrix (1000 × 2)

# Vector field → Jacobian matrix
v = hcat(u, cos.(getindex.(points, 2)))
J = op(v)  # Array (1000 × 2 × 2)

See also: gradient

source

RadialBasisFunctions.laplacian Method

laplacian(data, eval_points, basis, is_boundary, boundary_conditions, normals; k, adjl)

Build a Hermite-compatible RadialBasisOperator for the Laplacian. Maintains backward compatibility with the positional argument API.

source

RadialBasisFunctions.laplacian Method

laplacian(data; basis=PHS(3; poly_deg=2), eval_points=data, k, adjl, hermite)

Build a RadialBasisOperator for the Laplacian operator (∇²f).

Arguments

• data: Vector of data points

Keyword Arguments

• basis: RBF basis (default: PHS(3; poly_deg=2))

• eval_points: Evaluation points (default: data)

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation

Examples

# Basic usage
op = laplacian(data)

# With custom basis
op = laplacian(data; basis=PHS(5; poly_deg=3))

# With different evaluation points
op = laplacian(data; eval_points=eval_pts)

See also: partial, gradient, directional

source

RadialBasisFunctions.partial Method

partial(data, order, dim; basis=PHS(3; poly_deg=2), eval_points=data, k, adjl, hermite)

Build a RadialBasisOperator for a partial derivative.

Arguments

• data: Vector of data points

• order: Derivative order (1 or 2)

• dim: Dimension index to differentiate

Keyword Arguments

• basis: RBF basis (default: PHS(3; poly_deg=2))

• eval_points: Evaluation points (default: data)

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation

Examples

# First derivative in x-direction
∂x = partial(data, 1, 1)

# Second derivative in y-direction
∂²y = partial(data, 2, 2; basis=PHS(5; poly_deg=4))

See also: laplacian, gradient, directional

source

RadialBasisFunctions.regrid Method

regrid(data, eval_points; basis=PHS(3; poly_deg=2), k, adjl)

Build a RadialBasisOperator for interpolating from data points to eval_points.

Arguments

• data: Source data points

• eval_points: Target evaluation points

Keyword Arguments

• basis: RBF basis (default: PHS(3; poly_deg=2))

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

Examples

# Interpolate from coarse grid to fine grid
coarse = [SVector{2}(rand(2)) for _ in 1:100]
fine = [SVector{2}(rand(2)) for _ in 1:1000]
op = regrid(coarse, fine)

# Apply to field values
u_coarse = sin.(getindex.(coarse, 1))
u_fine = op(u_coarse)

See also: Interpolator

source

RadialBasisFunctions.update_hermite_stencil_data! Method

update_hermite_stencil_data!(hermite_data, global_data, neighbors,
                             is_boundary, boundary_conditions, normals,
                             global_to_boundary)

Populate local Hermite stencil data from global arrays. Used within kernels to extract boundary info for specific neighborhoods.

source

RadialBasisFunctions.∂virtual Method

function ∂virtual(data, eval_points, dim, Δ, basis; k=autoselect_k(data, basis))

Builds a virtual RadialBasisOperator whichi will be evaluated at eval_points where the operator is the partial derivative with respect to dim. Virtual operators interpolate the data to structured points at a distance Δ for which standard finite difference formulas can be applied.

source

RadialBasisFunctions.∂virtual Method

function ∂virtual(data, dim, Δ, basis; k=autoselect_k(data, basis))

Builds a virtual RadialBasisOperator whichi will be evaluated at the input points (data) where the operator is the partial derivative with respect to dim. Virtual operators interpolate the data to structured points at a distance Δ for which standard finite difference formulas can be applied.

source

RadialBasisFunctions.AbstractPHS Type

abstract type AbstractPHS <: AbstractRadialBasis

Supertype of all Polyharmonic Splines.

source

RadialBasisFunctions.AbstractRadialBasis Type

abstract type AbstractRadialBasis <: AbstractBasis end

source

RadialBasisFunctions.BoundaryCondition Type

BoundaryCondition{T}

Unified boundary condition representation: Bu = α_u + β_∂ₙu

Special cases:

• Dirichlet: α=1, β=0

• Neumann: α=0, β=1

• Robin: α≠0, β≠0

• Internal: α=0, β=0 (sentinel for interior points)

source

RadialBasisFunctions.Custom Type

Custom{F<:Function} <: AbstractOperator

Custom operator that applies a user-defined function to basis functions. The function ℒ should accept a basis and return a callable (x, xᵢ) -> value.

source

RadialBasisFunctions.Directional Type

Directional{Dim,T} <: ScalarValuedOperator

Operator for the directional derivative (∇f⋅v), the inner product of the gradient and a direction vector.

source

RadialBasisFunctions.Gaussian Type

Gaussian(ε=1; poly_deg=2)

Gaussian radial basis function: ϕ(r) = exp(-(εr)²)

Arguments

• ε: Shape parameter (must be > 0). Smaller values = wider basis.

• poly_deg: Polynomial augmentation degree (default: 2).

See also: PHS, IMQ

source

RadialBasisFunctions.HermiteStencilData Type

HermiteStencilData{T}

Local stencil data for Hermite interpolation with boundary conditions.

Fields:

• data: Coordinates of k stencil points

• is_boundary: Boolean flags for each point

• boundary_conditions: BC for each point (use Internal() for interior)

• normals: Normal vectors (zero for interior points)

• poly_workspace: Pre-allocated buffer for polynomial operations (avoids allocations in hot path)

Note: For interior points (is_boundary[i] == false), boundary_conditions[i] and normals[i] contain sentinel values and should not be accessed.

source

RadialBasisFunctions.HermiteStencilData Method

Pre-allocation constructor for HermiteStencilData

source

RadialBasisFunctions.IMQ Type

IMQ(ε=1; poly_deg=2)

Inverse Multiquadric radial basis function: ϕ(r) = 1/√((εr)² + 1)

Arguments

• ε: Shape parameter (must be > 0). Smaller values = flatter basis.

• poly_deg: Polynomial augmentation degree (default: 2).

See also: PHS, Gaussian

source

RadialBasisFunctions.Interpolator Type

struct Interpolator

Construct a radial basis interpolation.

source

RadialBasisFunctions.Interpolator Method

function Interpolator(x, y, basis::B=PHS())

Construct a radial basis interpolator.

See also: regrid for local stencil-based interpolation between point sets.

source

RadialBasisFunctions.Jacobian Type

Jacobian{Dim} <: VectorValuedOperator

Operator type for computing Jacobians (and gradients as a special case).

The Jacobian is the fundamental differential operator that computes all partial derivatives. When applied to a scalar field, it produces the gradient. When applied to a vector field, it produces the full Jacobian matrix.

Differentiation increases tensor rank by 1. The output gains a trailing dimension of size D (the spatial dimension).

Input/Output Shapes

• Scalar field Vector{T} (N,) → Gradient Matrix{T} (N_eval × D)

• Vector field Matrix{T} (N × D) → Jacobian Array{T,3} (N_eval × D × D)

• Matrix field Array{T,3} (N × D × D) → 3-tensor Array{T,4} (N_eval × D × D × D)

• General: input shape (N, dims...) → output shape (N_eval, dims..., D)

source

RadialBasisFunctions.Laplacian Type

Laplacian <: ScalarValuedOperator

Operator for the sum of the second derivatives w.r.t. each independent variable (∇²f).

source

RadialBasisFunctions.MonomialBasis Type

struct MonomialBasis{Dim,Deg} <: AbstractBasis

Dim dimensional monomial basis of order Deg.

source

RadialBasisFunctions.PHS1 Type

struct PHS1{T<:Int} <: AbstractPHS

Polyharmonic spline radial basis function: 

source

RadialBasisFunctions.PHS3 Type

struct PHS3{T<:Int} <: AbstractPHS

Polyharmonic spline radial basis function: 

source

RadialBasisFunctions.PHS5 Type

struct PHS5{T<:Int} <: AbstractPHS

Polyharmonic spline radial basis function: 

source

RadialBasisFunctions.PHS7 Type

struct PHS7{T<:Int} <: AbstractPHS

Polyharmonic spline radial basis function: 

source

RadialBasisFunctions.Partial Type

Partial{T<:Int} <: ScalarValuedOperator

Operator for a partial derivative of specified order with respect to a dimension.

source

RadialBasisFunctions.RadialBasisOperator Type

struct RadialBasisOperator

Operator of data using a radial basis with potential monomial augmentation.

source

RadialBasisFunctions.RadialBasisOperator Method

RadialBasisOperator(ℒ, data; eval_points, basis, k, adjl, hermite, device)

Unified constructor with keyword arguments.

Arguments

• ℒ: The operator type (e.g., Laplacian(), Partial(1, 2))

• data: Vector of data points

Keyword Arguments

• eval_points: Evaluation points (default: data)

• basis: RBF basis (default: PHS(3; poly_deg=2))

• k: Stencil size (default: autoselect_k(data, basis))

• adjl: Adjacency list (default: computed via find_neighbors)

• hermite: Optional NamedTuple for Hermite interpolation with fields:

  • is_boundary::Vector{Bool}

  • bc::Vector{<:BoundaryCondition}

  • normals::Vector{<:AbstractVector}

• device: KernelAbstractions backend for weight computation (default: auto-detected from data via get_backend)

Examples

# Basic usage
op = RadialBasisOperator(Laplacian(), data)

# With custom basis and stencil size
op = RadialBasisOperator(Laplacian(), data; basis=PHS(5; poly_deg=3), k=40)

# With different evaluation points
op = RadialBasisOperator(Laplacian(), data; eval_points=eval_pts)

# With Hermite boundary conditions
op = RadialBasisOperator(Laplacian(), data;
    hermite=(is_boundary=is_bound, bc=boundary_conds, normals=normal_vecs))

# With explicit device
using KernelAbstractions
op = RadialBasisOperator(Laplacian(), data; device=CPU())

source

RadialBasisFunctions.RadialBasisOperator Method

function (op::RadialBasisOperator)(y, x)

Evaluate the operator at x in-place and store the result in y.

source

RadialBasisFunctions.RadialBasisOperator Method

function (op::RadialBasisOperator)(x)

Evaluate the operator at x.

source

RadialBasisFunctions.Regrid Type

Regrid

Operator for interpolating from one set of points to another.

source

Private

RadialBasisFunctions._backward_partial_poly_1d! Method

Backward pass for polynomial section in 1D.

source

RadialBasisFunctions._backward_partial_poly_2d! Method

Backward pass for polynomial section in 2D.

source

RadialBasisFunctions._backward_partial_poly_3d! Method

Backward pass for polynomial section in 3D.

source

RadialBasisFunctions._backward_partial_polynomial_section! Method

_backward_partial_polynomial_section!(Δeval_point, Δb, k, nmon, dim, num_ops)

Backward pass through the polynomial section of the RHS for Partial operator.

For monomials in 2D with poly_deg=2 (1, x, y, xy, x², y²): ∂/∂x gives: 0, 1, 0, y, 2x, 0

The gradients of these w.r.t. eval_point are: ∂(0)/∂(x,y) = (0, 0) ∂(1)/∂(x,y) = (0, 0) ∂(0)/∂(x,y) = (0, 0) ∂(y)/∂(x,y) = (0, 1) -> b[4] contributes to Δeval_point[2] ∂(2x)/∂(x,y) = (2, 0) -> b[5] contributes 2 to Δeval_point[1] ∂(0)/∂(x,y) = (0, 0)

This is equivalent to computing the mixed second derivatives ∂²pⱼ/∂x[dim]∂x[d].

source

RadialBasisFunctions._build_collocation_matrix! Method

_build_collocation_matrix!(A, data, basis, mon, k)

Build RBF collocation matrix. Works for both interior stencils (AbstractVector data) and Hermite stencils (HermiteStencilData) via dispatch helpers.

Matrix structure:

┌─────────────────┬─────────┐
│  Φ(xᵢ, xⱼ)      │ P(xᵢ)   │  k×k RBF + k×nmon polynomial
├─────────────────┼─────────┤
│  P(xⱼ)ᵀ         │   0     │  nmon×k poly + nmon×nmon zero
└─────────────────┴─────────┘

For Hermite stencils with Neumann/Robin conditions, basis functions are modified via _rbf_entry and _poly_entry! dispatch to maintain matrix symmetry.

source

RadialBasisFunctions._build_rhs! Method

_build_rhs!(b, ℒrbf::Tuple, ℒmon::Tuple, data, eval_point, basis, mon, k)

Build RHS vector for multiple operators. Works for both interior stencils (AbstractVector) and Hermite stencils (HermiteStencilData) via dispatch helpers.

source

RadialBasisFunctions._build_rhs! Method

_build_rhs!(b, ℒrbf, ℒmon, data, eval_point, basis, mon, k)

Build RHS vector for single operator. Works for both interior stencils (AbstractVector) and Hermite stencils (HermiteStencilData) via dispatch helpers.

source

RadialBasisFunctions._build_stencil! Method

_build_stencil!(A, b, ℒrbf, ℒmon, data, eval_point, basis, mon, k)

Assemble complete stencil: build collocation matrix, build RHS, solve for weights. Works for both interior stencils (AbstractVector) and Hermite stencils (HermiteStencilData) via dispatch helpers in _build_collocation_matrix! and _build_rhs!.

Returns: weights (first k rows of solution, size k×num_ops)

source

RadialBasisFunctions._build_stencil! Method

_build_stencil!(λ, A, b, ℒrbf, ℒmon, data, eval_point, basis, mon, k)

In-place variant: writes solution into pre-allocated λ buffer, returns view of first k rows. Avoids allocating the solution vector and the slice on every call.

source

RadialBasisFunctions._build_weights Method

_build_weights(ℒ, data, eval_points, adjl, basis)

Apply operator to basis functions and route to weight computation.

source

RadialBasisFunctions._build_weights Method

_build_weights(data, eval_points, adjl, basis, ℒrbf, ℒmon, mon;
              batch_size=10, device=CPU())

Build weights for interior-only problems (no boundary conditions). Creates empty BoundaryData to indicate all interior points.

source

RadialBasisFunctions._build_weights Method

_build_weights(ℒ, data, eval_points, adjl, basis,
              is_boundary, boundary_conditions, normals)

Generic Hermite dispatcher for operators. Applies operator to basis and routes to Hermite weight computation.

This eliminates repetitive _build_weights methods across operator files. Note: Type constraint removed to avoid circular dependency with operators.jl

source

RadialBasisFunctions._build_weights Method

_build_weights(data, eval_points, adjl, basis, ℒrbf, ℒmon, mon,
              is_boundary, boundary_conditions, normals;
              batch_size=10, device=CPU())

Build weights for problems with boundary conditions using Hermite interpolation. Exact allocation: Dirichlet points get single entry, others get full stencil.

source

RadialBasisFunctions._build_weights Method

_build_weights(ℒ, op)

Entry point from operator construction. Extracts configuration from operator and routes to appropriate implementation.

source

RadialBasisFunctions._construct_sparse Method

_construct_sparse(I, J, V, N_eval, N_data, num_ops)

Construct sparse matrix/vector from COO arrays.

Future GPU support: convert to device-sparse format here (see #88)

source

RadialBasisFunctions._forward_with_cache Method

_forward_with_cache(data, eval_points, adjl, basis, ℒrbf, ℒmon, mon, ℒType)

Forward pass that builds weights while caching intermediate results for backward pass.

Returns: (W, cache) where W is the sparse weight matrix and cache contains per-stencil factorizations and solutions needed for the pullback.

source

RadialBasisFunctions._get_grad_funcs Method

_get_grad_funcs(OpType, basis, ℒ)

Get gradient functions for the given operator type and basis. Returns (grad_Lφ_x, grad_Lφ_xi) tuple.

source

RadialBasisFunctions._get_rhs_closures Method

_get_rhs_closures(OpType, ℒ, basis)

Get operator-specific closures for backward_stencil_with_ε!. Returns (poly_backward!, ∂Lφ_∂ε_fn) keyword arguments.

• Partial: polynomial section backward + partial ε derivative

• Laplacian: no polynomial backward + laplacian ε derivative

source

RadialBasisFunctions._interpolator_constructor_backward Method

_interpolator_constructor_backward(Δrbf_weights, Δmon_weights, A, k)

Backward pass for the Interpolator constructor w.r.t. y (the data values).

Given cotangents of rbf_weights and monomial_weights, computes the cotangent of y using the implicit function theorem. Since w = A⁻¹ [y; 0] and A is constant w.r.t. y:

Δy = (A⁻¹ [Δrbf_weights; Δmon_weights])[1:k]

Used by both Mooncake and potentially Enzyme extensions.

source

RadialBasisFunctions._interpolator_point_gradient! Method

_interpolator_point_gradient!(Δx, interp::Interpolator, x, Δy)

Accumulate the gradient of interp(x) * Δy into Δx.

RBF contribution: Σᵢ wᵢ * Δy * ∇φ(x, xᵢ) Polynomial contribution: Σⱼ wⱼ * Δy * ∇pⱼ(x)

source

RadialBasisFunctions._num_ops Method

Get number of operators (type-stable)

source

RadialBasisFunctions._optype Method

_optype(ℒ)

Map operator instance to its abstract type for dispatch in AD rules.

source

RadialBasisFunctions._prepare_buffer Method

Prepare RHS buffer with correct type (type-stable)

source

RadialBasisFunctions.allocate_sparse_arrays Method

allocate_sparse_arrays(TD, k, N_eval, num_ops, adjl, boundary_data)

Allocate sparse matrix arrays for COO format sparse matrix construction. Exactly counts non-zeros: interior points get k entries, Dirichlet points get 1 entry.

source

RadialBasisFunctions.autoselect_k Method

autoselect_k(data::Vector, basis<:AbstractRadialBasis)

See Bayona, 2017 - https://doi.org/10.1016/j.jcp.2016.12.008

source

RadialBasisFunctions.backward_collocation! Method

backward_collocation!(Δdata, ΔA, neighbors, data, basis, mon, k)

Chain rule through collocation matrix construction.

The collocation matrix has structure: A[i,j] = φ(xi, xj) for i,j ≤ k (RBF block) A[i,k+j] = pⱼ(xi) for i ≤ k (polynomial block)

For RBF block (using ∇φ from existing basis_rules): Δxi += ΔA[i,j] * ∇φ(xi, xj) Δxj -= ΔA[i,j] * ∇φ(xi, xj) (by symmetry of φ(x-y))

For polynomial block: Δxi += ΔA[i,k+j] * ∇pⱼ(xi)

Note: A is symmetric, so we need to handle both triangles.

source

RadialBasisFunctions.backward_collocation_ε! Method

backward_collocation_ε!(Δε_acc, ΔA, neighbors, data, basis, k)

Compute gradient contribution to shape parameter ε from collocation matrix.

Uses implicit differentiation: Δε += Σᵢⱼ ΔA[i,j] * ∂A[i,j]/∂ε where A[i,j] = φ(xi, xj) for the RBF block.

source

RadialBasisFunctions.backward_linear_solve! Method

backward_linear_solve!(ΔA, Δb, Δw, cache)

Compute cotangents of collocation matrix A and RHS vector b from cotangent of weights Δw.

Given: Aλ = b, w = λ[1:k] We have: Δλ = [Δw; 0] (padded with zeros for monomial part)

Using implicit function theorem: η = A⁻ᵀ Δλ ΔA = -η λᵀ Δb = η

source

RadialBasisFunctions.backward_rhs! Method

backward_rhs!(Δdata, Δeval_point, Δb, neighbors, eval_point, data, basis, k, grad_Lφ_x, grad_Lφ_xi; poly_backward!=nothing)

Chain rule through RHS vector construction for any operator.

RBF section (shared by all operators): b[i] = ℒφ(eval_point, xi) → accumulate ∂/∂eval_point and ∂/∂xi

Polynomial section (Partial only — Laplacian gives constants, no gradient): b[k+j] = ℒpⱼ(eval_point) → passed as optional poly_backward! closure

source

RadialBasisFunctions.backward_rhs_ε! Method

backward_rhs_ε!(Δε_acc, Δb, neighbors, eval_point, data, basis, k, ∂Lφ_∂ε_fn)

Compute gradient contribution to shape parameter ε from RHS.

∂Lφ_∂ε_fn(x, xi) returns ∂(ℒφ)/∂ε for the specific operator:

• Laplacian: (x, xi) -> ∂Laplacian_φ_∂ε(basis, x, xi)

• Partial: (x, xi) -> ∂Partial_φ_∂ε(basis, dim, x, xi)

source

RadialBasisFunctions.backward_stencil_with_ε! Method

backward_stencil_with_ε!(Δdata, Δeval_point, Δε_acc, Δw, cache, neighbors, eval_point, data, basis, mon, k, grad_Lφ_x, grad_Lφ_xi; poly_backward!=nothing, ∂Lφ_∂ε_fn=nothing)

Complete backward pass for a single stencil including shape parameter gradient.

Combines:

1. backward_linear_solve! → compute ΔA, Δb from Δw

2. backward_collocation! → chain ΔA to Δdata

3. backward_collocation_ε! → chain ΔA to Δε

4. backward_rhs! → chain Δb to Δdata and Δeval_point

5. backward_rhs_ε! → chain Δb to Δε

Optional closures:

• poly_backward!: polynomial section gradient (Partial only)

• ∂Lφ_∂ε_fn: shape parameter derivative function (IMQ/Gaussian only)

source

RadialBasisFunctions.build_weights_kernel Method

build_weights_kernel(data, eval_points, adjl, basis, ℒrbf, ℒmon, mon,
                    boundary_data; batch_size, device)

Main orchestrator for weight computation. Currently CPU-only. GPU stencil solve is not yet supported — see GitHub issue #88.

source

RadialBasisFunctions.build_weights_pullback_loop! Method

build_weights_pullback_loop!(Δdata, Δeval, Δε_acc, ΔW_extractor, cache, adjl,
    eval_points, data, basis, mon, ℒ, OpType, grad_Lφ_x, grad_Lφ_xi)

Shared stencil iteration loop for _build_weights pullback across all AD backends.

ΔW_extractor(eval_idx, neighbors, k) is a callable that returns the stencil cotangent matrix Δw given the eval index, neighbor list, and stencil size. This abstracts over the different ways each AD framework stores cotangents (dense matrix, nzval vector, etc.).

source

RadialBasisFunctions.calculate_batch_range Method

Calculate batch index range for kernel execution

source

RadialBasisFunctions.compute_hermite_poly_entry! Method

compute_hermite_poly_entry!(a, i, data, mon)

Compute polynomial entries for Hermite interpolation. Dispatches based on boundary point type for type stability.

source

RadialBasisFunctions.compute_hermite_rbf_entry Method

compute_hermite_rbf_entry(i, j, data, ops)

Compute single RBF matrix entry with Hermite boundary modifications. Dispatches based on point types (Interior/Dirichlet/NeumannRobin).

Uses BasisOperators for efficient evaluation (avoids functor construction in hot loop).

source

RadialBasisFunctions.construct_global_to_boundary Method

construct_global_to_boundary(is_boundary)

Construct mapping from global indices to boundary-only indices. For boundary points: global_to_boundary[i] = boundary array index For interior points: global_to_boundary[i] = 0 (sentinel)

source

RadialBasisFunctions.count_nonzeros Method

count_nonzeros(adjl, is_boundary, boundary_conditions)

Count exact number of non-zero entries for optimized allocation. Returns (total_nnz, row_offsets) where row_offsets[i] is the starting position for row i.

source

RadialBasisFunctions.extract_stencil_cotangent Method

extract_stencil_cotangent(ΔW, eval_idx, neighbors, k, num_ops)

Extract cotangent values for a single stencil from a dense/sparse matrix cotangent. Used by the Enzyme extension.

source

RadialBasisFunctions.extract_stencil_cotangent_from_nzval Method

extract_stencil_cotangent_from_nzval(ΔW_nzval, W, eval_idx, neighbors, k)

Extract cotangent values for a single stencil from sparse matrix nzval gradient. Used by Mooncake extension where gradients are stored in fdata.nzval.

source

RadialBasisFunctions.fill_dirichlet_entry! Method

Fill Dirichlet identity row for optimized allocation

source

RadialBasisFunctions.fill_entries! Method

Fill sparse matrix entries using indexed storage (row_offsets)

source

RadialBasisFunctions.grad_applied_laplacian_wrt_x Method

grad_applied_laplacian_wrt_x(basis)

Get gradient of applied Laplacian operator w.r.t. evaluation point.

source

RadialBasisFunctions.grad_applied_laplacian_wrt_xi Method

grad_applied_laplacian_wrt_xi(basis)

Get gradient of applied Laplacian operator w.r.t. data point. By symmetry, always the negation of the _wrt_x version.

source

RadialBasisFunctions.grad_applied_partial_wrt_x Method

grad_applied_partial_wrt_x(basis, dim)

Get gradient of applied partial derivative operator w.r.t. evaluation point.

source

RadialBasisFunctions.grad_applied_partial_wrt_xi Method

grad_applied_partial_wrt_xi(basis, dim)

Get gradient of applied partial derivative operator w.r.t. data point. By symmetry, always the negation of the _wrt_x version.

source

RadialBasisFunctions.grad_laplacian_gaussian_wrt_x Method

grad_laplacian_gaussian_wrt_x(ε)

Returns a function computing ∂/∂x[j] of [∇²φ] for Gaussian.

Mathematical derivation: ∇²φ = (4ε⁴r² - 2ε²D) * φ

∂(∇²φ)/∂x[j] = φ * δ_j * 4ε⁴ * [2 + D - 2ε²r²]

source

RadialBasisFunctions.grad_laplacian_imq_wrt_x Method

grad_laplacian_imq_wrt_x(ε)

Returns a function computing ∂/∂x[j] of [∇²φ] for IMQ.

Mathematical derivation: ∇²φ = sum_i [∂²φ/∂x[i]²] = 3ε⁴r²/s^(5/2) - D*ε²/s^(3/2)

∂(∇²φ)/∂x[j] = δ_j * [3(D+2)ε⁴/s^(5/2) - 15ε⁶r²/s^(7/2)]

source

RadialBasisFunctions.grad_laplacian_phs1_wrt_x Method

grad_laplacian_phs1_wrt_x()

Returns a function computing ∂/∂x[j] of [∇²φ] for PHS1.

Mathematical derivation: ∇²φ = 2/r ∂/∂x[j] [2/r] = -2 * δ_j / r³

Note: At r=0, we return 0 to avoid singularity.

source

RadialBasisFunctions.grad_laplacian_phs3_wrt_x Method

grad_laplacian_phs3_wrt_x()

Returns a function computing ∂/∂x[j] of [∇²φ] for PHS3.

Mathematical derivation: ∇²φ = 12r ∂/∂x[j] [12r] = 12 * δ_j / r

source

RadialBasisFunctions.grad_laplacian_phs5_wrt_x Method

grad_laplacian_phs5_wrt_x()

Returns a function computing ∂/∂x[j] of [∇²φ] for PHS5.

Mathematical derivation: ∇²φ = 30r³ ∂/∂x[j] [30r³] = 30 * 3r² * δ_j / r = 90 * r * δ_j

source

RadialBasisFunctions.grad_laplacian_phs7_wrt_x Method

grad_laplacian_phs7_wrt_x()

Returns a function computing ∂/∂x[j] of [∇²φ] for PHS7.

Mathematical derivation: ∇²φ = 56r⁵ ∂/∂x[j] [56r⁵] = 56 * 5r⁴ * δ_j / r = 280 * r³ * δ_j

source

RadialBasisFunctions.grad_partial_gaussian_wrt_x Method

grad_partial_gaussian_wrt_x(ε, dim)

Returns a function computing ∂/∂x[j] of [∂φ/∂x[dim]] for Gaussian.

Mathematical derivation: φ = exp(-ε²r²) ∂φ/∂x[dim] = -2ε² * δ_d * φ

∂²φ/∂x[j]∂x[dim] = φ * [-2ε² * δ_{j,dim} + 4ε⁴ * δ_d * δ_j]

For j == dim: φ * (4ε⁴ * δ_d² - 2ε²) For j != dim: φ * 4ε⁴ * δ_d * δ_j

source

RadialBasisFunctions.grad_partial_imq_wrt_x Method

grad_partial_imq_wrt_x(ε, dim)

Returns a function computing ∂/∂x[j] of [∂φ/∂x[dim]] for IMQ.

Mathematical derivation: Let s = ε²r² + 1, δ_d = x[dim] - xi[dim] ∂φ/∂x[dim] = -ε² * δ_d / s^(3/2)

∂²φ/∂x[j]∂x[dim] = -ε² * [δ_{j,dim} / s^(3/2) - δ_d * (3/2) * s^(-5/2) * 2ε² * δ_j] = -ε² * δ_{j,dim} / s^(3/2) + 3ε⁴ * δ_d * δ_j / s^(5/2)

For j == dim: -ε² / s^(3/2) + 3ε⁴ * δ_d² / s^(5/2) For j != dim: 3ε⁴ * δ_d * δ_j / s^(5/2)

source

RadialBasisFunctions.grad_partial_phs1_wrt_x Method

grad_partial_phs1_wrt_x(dim)

Returns a function computing ∂/∂x[j] of [∂φ/∂x[dim]] for PHS1.

Mathematical derivation: ∂φ/∂x[dim] = δ_d / r

∂²φ/∂x[j]∂x[dim] = (δ_{j,dim} * r - δ_d * δ_j / r) / r² = δ_{j,dim} / r - δ_d * δ_j / r³

Note: At r=0, the derivative is singular but we return 0 since the RBF value itself is 0 at r=0, so this term doesn't contribute to the gradient.

source

RadialBasisFunctions.grad_partial_phs3_wrt_x Method

grad_partial_phs3_wrt_x(dim)

Returns a function computing ∂/∂x[j] of [∂φ/∂x[dim]] for PHS3.

Mathematical derivation: ∂φ/∂x[dim] = 3 * δ_d * r where δ_d = x[dim] - xi[dim], r = ||x - xi||

∂²φ/∂x[j]∂x[dim] = 3 * (δ_{j,dim} * r + δ_d * δ_j / r)

For j == dim: 3 * (r + δ_d² / r) For j != dim: 3 * δ_d * δ_j / r

source

RadialBasisFunctions.grad_partial_phs5_wrt_x Method

grad_partial_phs5_wrt_x(dim)

Returns a function computing ∂/∂x[j] of [∂φ/∂x[dim]] for PHS5.

Mathematical derivation: ∂φ/∂x[dim] = 5 * δ_d * r³

∂²φ/∂x[j]∂x[dim] = 5 * (δ_{j,dim} * r³ + δ_d * 3r * δ_j) = 5 * (δ_{j,dim} * r³ + 3 * δ_d * δ_j * r)

source

RadialBasisFunctions.grad_partial_phs7_wrt_x Method

grad_partial_phs7_wrt_x(dim)

Returns a function computing ∂/∂x[j] of [∂φ/∂x[dim]] for PHS7.

Mathematical derivation: ∂φ/∂x[dim] = 7 * δ_d * r⁵

∂²φ/∂x[j]∂x[dim] = 7 * (δ_{j,dim} * r⁵ + δ_d * 5r³ * δ_j)

source

RadialBasisFunctions.hermite_mono_rhs! Method

hermite_mono_rhs!(bmono, ℒmon, mon, eval_point, is_bound, bc, normal, workspace)

Apply boundary conditions to monomial operator evaluation. Dispatches based on boundary point type for type stability.

source

RadialBasisFunctions.hermite_rbf_rhs Method

hermite_rbf_rhs(ℒrbf, eval_point, data_point, is_bound, bc, normal)

Apply boundary conditions to RBF operator evaluation.

• Interior/Dirichlet: standard evaluation ℒΦ(x_eval, x_data)

• Neumann/Robin: α_ℒΦ + β_ℒ(∂ₙΦ)

source

RadialBasisFunctions.launch_kernel! Method

launch_kernel!(...)

Launch parallel CPU kernel for weight computation. Handles Dirichlet/Interior/Hermite stencil classification via dispatch.

source

RadialBasisFunctions.negate_grad Method

negate_grad(grad_fn)

Given a gradient function grad_fn(x, xi), returns (x, xi) -> -grad_fn(x, xi). All _wrt_xi functions are the negation of their _wrt_x counterparts by symmetry.

source

RadialBasisFunctions.operator_arity Method

Extract operator arity at compile time

source

RadialBasisFunctions.point_type Method

Determine boundary type of a single point

source

RadialBasisFunctions.∂Laplacian_φ_∂ε Method

∂Laplacian_φ_∂ε(basis::Gaussian, x, xi)

Derivative of Laplacian of Gaussian basis w.r.t. shape parameter ε.

∇²φ = (-2ε²D + 4ε⁴r²) exp(-ε²r²), where D = dimension ∂(∇²φ)/∂ε = exp(-ε²r²) [-4εD + 16ε³r² + 4ε³r²D - 8ε⁵r⁴]

source

RadialBasisFunctions.∂Laplacian_φ_∂ε Method

∂Laplacian_φ_∂ε(basis::IMQ, x, xi)

Derivative of Laplacian of IMQ basis w.r.t. shape parameter ε.

Let s = ε²r² + 1, D = dimension ∇²φ = -ε²D/s^(3/2) + 3ε⁴r²/s^(5/2) ∂(∇²φ)/∂ε = ∂/∂ε[-ε²D s^(-3/2) + 3ε⁴r² s^(-5/2)]

source

RadialBasisFunctions.∂Partial_φ_∂ε Method

∂Partial_φ_∂ε(basis::Gaussian, dim::Int, x, xi)

Derivative of first partial derivative of Gaussian basis w.r.t. shape parameter ε.

∂φ/∂x_dim = -2ε²(x_dim - xi_dim) exp(-ε²r²) ∂/∂ε[∂φ/∂x_dim] = 4ε(x_dim - xi_dim)(ε²r² - 1) exp(-ε²r²)

source

RadialBasisFunctions.∂Partial_φ_∂ε Method

∂Partial_φ_∂ε(basis::IMQ, dim::Int, x, xi)

Derivative of first partial derivative of IMQ basis w.r.t. shape parameter ε.

∂φ/∂x_dim = ε²(xi_dim - x_dim) s^(-3/2) ∂/∂ε[∂φ/∂x_dim] = 2ε(xi_dim - x_dim) s^(-3/2) + ε²(xi_dim - x_dim)(-3/2)s^(-5/2) · 2εr² = (xi_dim - x_dim)[2ε s^(-3/2) - 3ε³r² s^(-5/2)]

source

RadialBasisFunctions.∂φ_∂ε Method

∂φ_∂ε(basis::Gaussian, x, xi)

Derivative of Gaussian basis function w.r.t. shape parameter ε.

φ(r) = exp(-ε²r²) ∂φ/∂ε = -2εr² exp(-ε²r²)

source

RadialBasisFunctions.∂φ_∂ε Method

∂φ_∂ε(basis::IMQ, x, xi)

Derivative of IMQ basis function w.r.t. shape parameter ε.

φ(r) = (ε²r² + 1)^(-1/2) ∂φ/∂ε = -εr² (ε²r² + 1)^(-3/2)

source

RadialBasisFunctions.AbstractBasis Type

abstract type AbstractBasis end

source

RadialBasisFunctions.BasisOperators Type

BasisOperators{B,G,Hess}

Bundle of pre-constructed basis operators for efficient evaluation in hot loops. Avoids repeated functor construction inside hermite_rbf_dispatch.

Fields:

• φ: The basis function itself

• ∇φ: Gradient operator (pre-constructed ∇(basis))

• Hφ: Hessian operator (pre-constructed H(basis))

Usage:

ops = BasisOperators(basis)
# In hot loop:
φ_val = ops.φ(x, xᵢ)
grad = ops.∇φ(x, xᵢ)      # Returns vector
hess = ops.Hφ(x, xᵢ)      # Returns matrix
Dφ = dot(n, grad)         # Directional derivative
D²φ = dot(ni, hess * nj)  # Second directional derivative

source

RadialBasisFunctions.BasisOperators Method

Construct BasisOperators from a basis function.

source

RadialBasisFunctions.BoundaryData Type

BoundaryData{T}

Wrapper for global boundary information (replaces fragile tuples).

source

RadialBasisFunctions.BoundaryPointType Type

Trait types for individual point boundary classification

source

RadialBasisFunctions.D Type

D{B<:AbstractRadialBasis,V}

Directional derivative operator functor. Construct with D(basis, v). Computes the derivative of the basis function in direction v.

source

RadialBasisFunctions.D² Type

D²{B<:AbstractRadialBasis,V1,V2}

Directional second derivative operator functor. Construct with D²(basis, v1, v2).

source

RadialBasisFunctions.H Type

H{B<:AbstractRadialBasis}

Hessian operator functor. Construct with H(basis). Returns the Hessian matrix of the basis function.

source

RadialBasisFunctions.OperatorArity Type

Operator arity traits for compile-time dispatch

source

RadialBasisFunctions.StencilForwardCache Type

StencilForwardCache{T}

Per-stencil storage from forward pass needed for backward pass.

• lambda: Full solution vector (k+nmon) × num_ops from solving Aλ = b

• A_mat: The symmetric collocation matrix (stored for backprop)

• k: Number of RBF neighbors in stencil

• nmon: Number of monomial basis functions

source

RadialBasisFunctions.StencilType Type

Trait types for stencil classification

source

RadialBasisFunctions.WeightsBuildForwardCache Type

WeightsBuildForwardCache{T}

Global cache storing all stencil results and references to inputs.

• stencil_caches: Vector of StencilForwardCache, one per evaluation point

• k: Stencil size (number of neighbors)

• nmon: Number of monomial basis functions

• num_ops: Number of operators (1 for scalar, D for gradient)

source

RadialBasisFunctions.∂ Type

∂{B<:AbstractRadialBasis}

Partial derivative operator functor. Construct with ∂(basis, dim).

source

RadialBasisFunctions.∂² Type

∂²{B<:AbstractRadialBasis}

Second partial derivative operator functor. Construct with ∂²(basis, dim).

source

RadialBasisFunctions.∇ Type

∇{B<:AbstractRadialBasis}

Gradient operator functor. Construct with ∇(basis).

source

RadialBasisFunctions.∇² Type

∇²{B<:AbstractRadialBasis}

Laplacian operator functor. Construct with ∇²(basis).

source