KernelInterpolation.jl API

KernelInterpolation

KernelInterpolation.jl is a Julia package that implements methods for multivariate interpolation in arbitrary dimension based on symmetric (conditionally) positive-definite kernels with a focus on radial basis functions. It can be used for classical interpolation of scattered data, as well as for generalized (Hermite-Birkhoff) interpolation by using a meshfree collocation approach. This can be used to solve partial differential equations both stationary ones and time-dependent ones by using some time integration method from OrdinaryDiffEq.jl.

See also: KernelInterpolation.jl

AbstractKernel

An abstract supertype of kernels.

dim(kernel)

Return the dimension of a kernel, i.e. the size of the input vector.

get_name(kernel::AbstractKernel)

Returns the canonical, human-readable name for the given system of equations.

GaussKernel{Dim}(; shape_parameter = 1.0)

Gaussian kernel function with

\[ \phi(r) = \exp(-(\varepsilon r)^2),\]

where $\varepsilon$ is the shape parameter. The Gaussian kernel is always positive definite. See Wendland (2004), p. 74.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539

InverseMultiquadricKernel{Dim}(beta = 0.5; shape_parameter = 1.0)

Inverse multiquadric kernel function with

\[ \phi(r) = (1 + (\varepsilon r)^2)^{-\beta},\]

where $\varepsilon$ is the shape parameter. The inverse multiquadric kernel is always positive definite. See Wendland (2004), p. 76 and p. 95.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539

Matern12Kernel{Dim}(; shape_parameter = 1.0)

Matern kernel with $\nu = 1/2$, i.e.,

\[ \phi(r) = \exp(-\varepsilon r),\]

where $\varepsilon$ is the shape parameter. The Matern kernel is positive definite.

See Wikipedia and Fasshauer (2007), p. 41.

See also MaternKernel, RadialSymmetricKernel.

• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437

Matern32Kernel{Dim}(; shape_parameter = 1.0)

Matern kernel with $\nu = 3/2$, i.e.,

\[ \phi(r) = (1 + \sqrt{3}\varepsilon r)\exp(-\sqrt{3}\varepsilon r),\]

where $\varepsilon$ is the shape parameter. The Matern kernel is positive definite.

See Wikipedia and Fasshauer (2007), p. 41.

See also MaternKernel, RadialSymmetricKernel.

• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437

Matern52Kernel{Dim}(; shape_parameter = 1.0)

Matern kernel with $\nu = 5/2$, i.e.,

\[ \phi(r) = (1 + \sqrt{5}\varepsilon r + 5\cdot(\varepsilon r)^2/3)\exp(-\sqrt{5}\varepsilon r),\]

where $\varepsilon$ is the shape parameter. The Matern kernel is positive definite.

See Wikipedia and Fasshauer (2007), p. 41.

See also MaternKernel, RadialSymmetricKernel.

• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437

Matern72Kernel{Dim}(; shape_parameter = 1.0)

Matern kernel with $\nu = 7/2$, i.e.,

\[ \phi(r) = (1 + \sqrt{7}\varepsilon r + 12\cdot(\varepsilon r)^2/5 + 7\cdot(\varepsilon r)^3/15)\exp(-\sqrt{7}\varepsilon r),\]

where $\varepsilon$ is the shape parameter. The Matern kernel is positive definite.

See Wikipedia and Fasshauer (2007), p. 41.

See also MaternKernel, RadialSymmetricKernel.

• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437

MaternKernel{Dim}(nu = 1.5; shape_parameter = 1.0)

Matern kernel with

\[ \phi_\nu(r) = \frac{2^{1-\nu}}{\Gamma(\nu)}\big(\sqrt{2\nu}\varepsilon r\big)^\nu K_\nu\big(\sqrt{2\nu}\varepsilon r\big),\]

where $\varepsilon$ is the shape parameter. The Matern kernel is positive definite.

See Wikipedia and Fasshauer (2007), p. 41.

See also RadialSymmetricKernel.

• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437

MultiquadricKernel{Dim}(beta = 0.5; shape_parameter = 1.0)

Multiquadric kernel function with

\[ \phi(r) = (1 + (\varepsilon r)^2)^\beta,\]

where $\varepsilon$ is the shape parameter. The multiquadric kernel is conditionally positive definite of order $m = \lceil\beta \rceil$. See Wendland (2004), p. 109.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539

PolyharmonicSplineKernel{Dim}(k)

Polyharmonic spline kernel function with

\[ \phi_k(r) = \begin{cases} r^k, &\text{ if } k \text{ odd}\\ r^k\log(r), &\text{ if } k \text{ even} \end{cases}.\]

The polyharmonic spline is conditionally positive definite of order $m = \lceil k/2\rceil$ for odd k and order $m = k/2 + 1$ for even k. See Wendland (2004), pp. 111–112.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539

RadialCharacteristicKernel{Dim}(beta = 2.0; shape_parameter = 1.0)

Radial characteristic function (or also called truncated power or Askey) kernel function with

\[ \phi(r) = (1 - \varepsilon r)^\beta_+,\]

where $\varepsilon$ is the shape parameter. The radial characteristic function is positive definite if $\beta\ge (d + 1)/2$. It is compactly supported. See Wendland (2004), p. 80, Iske (2018), p. 281.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539
• Armin Iske (2018) Approximation Theory and Algorithms for Data Analysis Texts in Applied Mathematics (Springer) DOI: 10.1007/978-3-030-05228-7

RadialSymmetricKernel

An abstract supertype of radial symmetric kernels. Radial symmetric kernels are generated by an even and continuous function $\Phi: \mathbb{R}^d\to\mathbb{R}$, which is radial-symmetric meaning that there exists a $\phi:[0,\infty]\to\mathbb{R}$ such that

\[ \Phi(x) = \phi(\Vert x\Vert).\]

The kernel is then defined by

\[ K(x, y) = \Phi(x - y).\]

A RadialSymmetricKernel can be evaluated at two points x and y by calling kernel(x, y) or at a single point x by calling kernel(x), which implicitly sets y to zero.

RieszKernel{Dim}(beta; shape_parameter = 1.0)

Riesz kernel with

\[ \phi(r) = -(\varepsilon r)^\beta,\]

where $\varepsilon$ is the shape parameter and $\beta\in (0,2)$. The Riesz kernel is conditionally positive definite of order 1. See Hertrich et al. (2023).

See also RadialSymmetricKernel.

• Johannes Hertrich, Christian Wald, Fabian Altekrüger, Paul Hagemann (2023) Generative Sliced MMD Flows with Riesz Kernels ArXiv: 2305.11463

ThinPlateSplineKernel{Dim}()

Thin plate spline kernel function with

\[ \phi(r) = r^2\log(r),\]

i.e., PolyharmonicSplineKernel with $k = 2$. The thin plate spline is conditionally positive definite of order $m = 2$. See Wendland (2004), p. 112.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539

WendlandKernel{Dim}(k; shape_parameter = 1.0, d = Dim)

Wendland kernel with

\[ \phi_{d,k}(r) = \begin{cases} p_{d,k}(\varepsilon r), \text{ if } 0\le \varepsilon r\le 1\\ 0, \text{ if } \varepsilon r > 1 \end{cases},\]

where $\varepsilon$ is the shape parameter and $p$ is a polynomial with minimal degree. The Wendland kernel is positive definite for d\le Dim and compactly supported. See Wendland (2004), p. 129 or Fasshauer (2007), pp. 87–88.

See also RadialSymmetricKernel.

• Holger Wendland (2004) Scattered Data Approximation Cambridge University Press DOI: 10.1017/CBO9780511617539
• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437

WuKernel{Dim}(l, k; shape_parameter = 1.0)

Wu kernel with

\[ \phi_{l,k}(r) = \begin{cases} p_{l,k}(\varepsilon r), \text{ if } 0\le \varepsilon r\le 1\\ 0, \text{ if } \varepsilon r > 1 \end{cases},\]

where $\varepsilon$ is the shape parameter, $k\le l$, and $p$ is a polynomial of degree $4l - 2k + 1$. The Wu kernel is positive definite for $Dim\le 2k + 1$ and compactly supported. See Fasshauer (2007), pp. 88–90 and Wu (1995).

See also RadialSymmetricKernel.

• Gregory Fasshauer (2007) Meshfree Approximation Methods with MATLAB World Scientific DOI: 10.1142/6437
• Zongmin Wu (1995) Compactly supported positive definite radial functions Advances in Computational Mathematics DOI: 10.1007/BF03177517

Phi(kernel, x)

For a RadialSymmetricKernel kernel return value of the multivariate function $\Phi$ defined by $\Phi(x) = \phi(\Vert x\Vert)$.

See also RadialSymmetricKernel, phi.

order(kernel)

Return order of kernel.

phi(kernel, r)

For a RadialSymmetricKernel kernel return value of the univariate function $\phi$ defining the kernel.

See also RadialSymmetricKernel, Phi.

ProductKernel{Dim}(kernels)

Given a vector of kernels, construct a new kernel that multiplies the results of the component kernels, i.e., the new kernel $K$ is given by

\[ K(x, y) = \prod_{i = 1}^n K_i(x, y),\]

where $K_i$ are the component kernels and $n$ the number of kernels. Note that all component kernels need to have the same dim. A ProductKernel can also be constructed using the * operator.

SumKernel{Dim}(kernels)

Given a vector of kernels, construct a new kernel that sums the results of the component kernels, i.e., the new kernel $K$ is given by

\[ K(x, y) = \sum_{i = 1}^n K_i(x, y),\]

where $K_i$ are the component kernels and $n$ the number of kernels. Note that all component kernels need to have the same dim. A SumKernel can also be constructed using the + operator.

TransformationKernel{Dim}(kernel, transformation)

Given a base kernel and a bijective transformation function, construct a new kernel that applies the transformation to both arguments $x$ and $y$, i.e., the new kernel $K_T$ is given by

\[ K_T(x, y) = K(Tx, Ty),\]

where $K$ is the base kernel and $T$ the transformation, i.e. if $K$ is a kernel of dimension $d$, $T$ is a function from dimension Dim to $d$, where Dim is the dimension of the new kernel.

NodeSet(nodes)

Set of interpolation nodes.

distance_matrix(nodeset1::NodeSet, nodeset2::NodeSet)

Compute the distance matrix between two NodeSets, which is a matrix $D$ with $D_{ij} = \|x_i - \xi_j\|$ for all $i$ and $j$, where $x_i$ are the nodes in nodeset1 and $\xi_j$ are the nodes on nodeset2.

empty_nodeset(Dim, RealT = Float64)

Create an empty NodeSet.

homogeneous_hypercube(n, x_min = ntuple(_ -> 0.0, dim), x_max = ntuple(_ -> 1.0, dim); [dim])

If n is integer, create a NodeSet with n homogeneously distributed nodes in every dimension each of dimension dim inside a hypercube defined by the bounds x_min and x_max. If n is a Tuple of length dim, then use as many nodes in each dimension as described by n. The resulting NodeSet will have $n^{\textrm{dim}}$ respectively $\prod_{j = 1}^{\textrm{dim}}n_j$ points. If the bounds are given as single values, they are applied for each dimension. If they are Tuples of size dim, the hypercube has the according bounds. If dim is not given explicitly, it is inferred by the lengths of n, x_min and x_max if possible.

homogeneous_hypercube_boundary(n, x_min = ntuple(_ -> 0.0, dim), x_max = ntuple(_ -> 1.0, dim); [dim])

If n is integer, create a NodeSet with n homogeneously distributed nodes in every dimension each of dimension dim on the boundary of a hypercube defined by the bounds x_min and x_max. If n is a Tuple of length dim, then use as many nodes in each dimension as described by n. If the bounds are given as single values, they are applied for each dimension. If they are Tuples of size dim, the hypercube has the according bounds. If dim is not given explicitly, it is inferred by the lengths of n, x_min and x_max if possible.

random_hypercube([rng], n, x_min = ntuple(_ -> 0.0, dim), x_max = ntuple(_ -> 1.0, dim); [dim])

Create a NodeSet with n random nodes each of dimension dim inside a hypercube defined by the bounds x_min and x_max. If the bounds are given as single values, they are applied for each dimension. If they are Tuples of size dim the hypercube has the according bounds. If dim is not given explicitly, it is inferred by the lengths of x_min and x_max if possible. Optionally, pass a random number generator rng.

random_hypercube_boundary([rng], n, x_min = ntuple(_ -> 0.0, dim), x_max = ntuple(_ -> 1.0, dim); [dim])

Create a NodeSet with n random nodes each of dimension dim on the boundary of a hypercube defined by the bounds x_min and x_max. If the bounds are given as single values, they are applied for each dimension. If they are Tuples of size dim the hypercube has the according bounds. If dim is not given explicitly, it is inferred by the lengths of x_min and x_max if possible. Optionally, pass a random number generator rng.

random_hypersphere([rng], n, r = 1.0, center = Tuple(zeros(dim)); [dim])

Create a NodeSet with n random nodes each of dimension dim inside a hypersphere with radius r around the center center, which is given as a tuple. If dim is not given explicitly, it is inferred by the length of center if possible. Optionally, pass a random number generator rng.

random_hypersphere_boundary([rng], n, r = 1.0, center = Tuple(zeros(dim)); [dim])

Create a NodeSet with n random nodes each of dimension dim at the boundary of a hypersphere with radius r around the center center, which is given as a tuple. If dim is not given explicitly, it is inferred by the length of center if possible. Optionally, pass a random number generator rng.

separation_distance(nodeset::NodeSet)

Return the separation distance of a NodeSet $X = \{x_1,\ldots, x_n\}$ defined by

\[ q_X = \frac{1}{2}\min_{x_i\neq x_j}\|x_i - x_j\|.\]

AbstractBasis

Abstract type for a basis of a kernel function space. Every basis represents a set of functions, which can be obtained by indexing the basis object. Every basis object holds a kernel function and a NodeSet of centers and potentially more fields depending on the concrete basis type.

LagrangeBasis(centers, kernel, m = order(kernel))

The Lagrange (or cardinal) basis with respect to a kernel and a NodeSet of centers. This basis already includes polynomial augmentation of degree m defaulting to order(kernel). The basis functions are given such that

\[ b_j(x_i) = \delta_{ij},\]

which means that the kernel_matrix of this basis is the identity matrix making it suitable if multiple interpolations with the same centers of the basis and the same kernel, but with different right-hand sides or nodesets are performed. Since the basis already includes polynomials no additional polynomial augmentation is needed for interpolation with this basis.

StandardBasis(centers, kernel)

The standard basis for a function space defined by a kernel and a NodeSet of centers. The basis functions are given by

\[ b_j(x) = K(x, x_j)\]

where K is the kernel and x_j are the nodes in centers.

centers(basis)

Return the centers from a basis object.

interpolation_kernel(basis)

Return the kernel from a basis.

order(basis)

Return the order $m$ of the polynomial, which is needed by this basis for the interpolation, i.e., the polynomial degree plus 1. If $m = 0$, no polynomial is added.

Interpolation

Interpolation object that can be evaluated at a node and represents a kernel interpolation of the form

\[ s(x) = \sum_{j = 1}^N c_jb_j(x) + \sum_{k = 1}^Q d_kp_k(x),\]

where $b_j$ are the basis functions and $p_k$ is a basis of the Q-dimensional space of multivariate polynomials of order order. The additional conditions

\[ \sum_{j = 1}^N c_jp_k(x_j) = 0, \quad k = 1,\ldots, Q\]

are enforced.

See also interpolate.

TemporalInterpolation(ode_sol::ODESolution)

Temporal interpolation of an ODE solution. The result can be evaluated at a time t and a spatial point x. Evaluating the interpolation at a time t returns an Interpolation object that can be evaluated at a spatial point x.

basis(itp)

Return the basis from an interpolation object.

centers(itp::Interpolation)

Return the centers from the basis of an interpolation object.

coefficients(itp::Interpolation)

Obtain all the coefficients of the linear combination for the interpolant, i.e., both the coefficients for the kernel part and for the polynomial part.

See also kernel_coefficients and polynomial_coefficients.

dim(itp::Interpolation)

Return the dimension of the input variables of the interpolation.

interpolate(basis, values, nodeset = centers(basis); m = order(basis),
            regularization = NoRegularization(), factorization_method = nothing)
interpolate(centers, [nodeset,] values, kernel = GaussKernel{dim(nodeset)}();
            m = order(kernel), regularization = NoRegularization(),
            factorization_method = nothing)

Interpolate the values evaluated at the nodes in the nodeset to a function using the kernel kernel and polynomials up to a order m (i.e. degree - 1), i.e., determine the coefficients $c_j$ and $d_k$ in the expansion

\[ s(x) = \sum_{j = 1}^N c_jb_j(x) + \sum_{k = 1}^Q d_kp_k(x),\]

where $b_j$ are the basis functions in the basis and $s(x)$ the interpolant $s(x_j) = f(x_j)$, where $f(x_j)$ are given by values, $x_j$ are the nodes in the nodeset, and $p_k$ is a basis of the $Q$-dimensional space of multivariate polynomials with maximum degree of m - 1. If m = 0, no polynomial is added. The additional conditions

\[ \sum_{j = 1}^N c_jp_k(x_j) = 0, \quad k = 1,\ldots, Q = \begin{pmatrix}m - 1 + d\\d\end{pmatrix}\]

are enforced. Returns an Interpolation object.

If nodeset is provided, the interpolant is a least squares approximation with a different set of nodes as the centers used for the basis. Otherwise, nodeset is set to centers(basis) or centers.

A regularization can be applied to the kernel matrix using the regularization argument, cf. regularize!. In addition, the factorization_method can be specified to determine how the system matrix is factorized. By default, the system matrix is just wrapped as a Symmetric matrix for interpolation and no factorization is applied for a least squares solution, but you can, e.g., also explicitly use cholesky, lu, or qr factorization.

interpolation_kernel(itp)

Return the kernel from an interpolation object.

kernel_coefficients(itp::Interpolation)

Obtain the coefficients of the kernel part of the linear combination for the interpolant.

See also coefficients and polynomial_coefficients.

kernel_inner_product(itp1, itp2)

Inner product of the native space for two interpolants itp1 and itp2 with the same kernel. The inner product is defined as

\[ \langle f, g\rangle_K = \sum_{i = 1}^N\sum_{j = 1}^Mc_i^fc_j^gK(x_i, \xi_j)\]

for the interpolants $f(x) = \sum_{i = 1}^Nc_i^fK(x, x_i)$ and $g(x) = \sum_{j = 1}^Mc_j^gK(x, \xi_j)$.

See also kernel_norm.

kernel_norm(itp)

Norm of the native space defined by the kernel of the interpolant itp. The norm is defined as

\[ \|f\|_K^2 = \sum_{i,j=1}^Nc_ic_jK(x_i, x_j)\]

for the interpolant $f(x) = \sum_{j = 1}^nc_jK(x, x_j)$.

See also kernel_inner_product.

nodeset(itp)

Return the node set from an interpolation object.

order(itp)

Return the order $m$ of the polynomial used for the interpolation, i.e., the polynomial degree plus 1. If $m = 0$, no polynomial is added.

polynomial_basis(itp::Interpolation)

Return a vector of the polynomial basis functions used for the interpolation.

See also polyvars.

polynomial_coefficients(itp::Interpolation)

Obtain the coefficients of the polynomial part of the linear combination for the interpolant.

See also coefficients and kernel_coefficients.

polyvars(itp::Interpolation)

Return a vector of the polynomial variables.

See also polynomial_basis.

system_matrix(itp::Interpolation)

Return the system matrix, i.e., the matrix $A$ in the linear system

\[ Ac = f,\]

where $c$ are the coefficients of the kernel interpolant and $f$ the vector of known values. The exact form of $A$ differs depending on which method is used.

AbstractRegularization

An abstract supertype of regularizations. A regularization implements a function regularize! that takes a matrix and returns a regularized version of it.

L2Regularization(regularization_parameter::Real)

A regularization that adds a multiple of the identity matrix to the input matrix.

NoRegularization()

A regularization that does nothing.

regularize!(A, reg::AbstractRegularization)

Apply the regularization reg to the matrix A in place.

EllipticOperator(A, b, c)

Linear second-order elliptic operator with matrix $A(x)\in\mathbb{R}^{d\times d}$, vector $b(x)\in\mathbb{R}^d$, and scalar $c(x)$. The operator is defined as

\[ \mathcal{L}u = -\sum_{i,j = 1}^d a_{ij}(x)\partial_{x_i,x_j}^2u + \sum_{i = 1}^db_i(x)\partial_{x_i}u + c(x)u.\]

A, b and c are space-dependent functions returning a matrix, a vector, and a scalar, respectively. The matrix A should be symmetric and positive definite for any input x. The operator can be called with a RadialSymmetricKernel and points x and y to evaluate the operator of the kernel at x - y. It can also be called with an Interpolation object and a point x to evaluate the elliptic operator of the interpolation at x. Note that this is only supported for the kernel part of the interpolation, i.e. the polynomial part, if existent, is ignored.

Gradient()

The gradient operator. It can be called with a RadialSymmetricKernel and points x and y to evaluate the gradient of the kernel at x - y. It can also be called with an Interpolation object and a point x to evaluate the gradient of the interpolation at x. Note that this is only supported for the kernel part of the interpolation, i.e. the polynomial part, if existent, is ignored.

Laplacian()

The Laplacian operator. It can be called with a RadialSymmetricKernel and points x and y to evaluate the Laplacian of the kernel at x - y. It can also be called with an Interpolation object and a point x to evaluate the Laplacian of the interpolation at x. Note that this is only supported for the kernel part of the interpolation, i.e. the polynomial part, if existent, is ignored.

PartialDerivative(i)

Partial derivative operator with respect to the i-th component. The operator can be called with a RadialSymmetricKernel and points x and y to evaluate the derivative of the kernel at x - y. It can also be called with an Interpolation object and a point x to evaluate the first partial derivative of the interpolation at x in the i-th direction. Note that this is only supported for the kernel part of the interpolation, i.e. the polynomial part, if existent, is ignored.

AdvectionDiffusionEquation(diffusivity, advection_velocity, f)

Advection-diffusion equation with diffusivity diffusivity and advection velocity advection_velocity. The advection-diffusion equation is defined as

\[ \partial_t u + \mathbf{a}\cdot\nabla u = \kappa\Delta u + f,\]

where $\mathbf{a}$ is the advection velocity, $\kappa$ is the diffusivity, and $f$ is the right-hand side, which can be a time- and space-dependent function or a vector.

AdvectionEquation(advection_velocity)

Advection equation with advection velocity advection_velocity. The advection equation is defined as

\[ \partial_t u + \mathbf{a}\cdot\nabla u = f,\]

where $\mathbf{a}$ is the advection velocity and $f$ a source term.

EllipticEquation(A, b, c, f)

Libear second-order elliptic equation with matrix A, vector b, and scalar c and right-hand side f. The elliptic equation is defined as

\[ \mathcal{L}u = -\sum_{i,j = 1}^d a_{ij}(x)\partial_{x_i,x_j}^2u + \sum_{i = 1}^db_i(x)\partial_{x_i}u + c(x)u = f,\]

where A, b and c are space-dependent functions returning a matrix, a vector, and a scalar, respectively.

See also EllipticOperator.

HeatEquation(diffusivity, f)

Heat equation with thermal diffusivity diffusivity. The heat equation is defined as

\[ \partial_t u = \kappa\Delta u + f,\]

where $\kappa$ is the thermal diffusivity and $f$ is the right-hand side, which can be a time- and space-dependent function or a vector.

PoissonEquation(f)

Poisson equation with right-hand side f, which can be a space-dependent function or a vector. The Poisson equation is defined as

\[ -\Delta u = f\]

See also Laplacian.

Semidiscretization(spatial_discretization, initial_condition)
Semidiscretization(equations, nodeset_inner, boundary_condition, nodeset_boundary, [centers,]
                   initial_condition, kernel = GaussKernel{dim(nodeset_inner)}())

Semidiscretization of a partial differential equation with Dirichlet boundary conditions and initial condition initial_condition. The boundary_condition function can be time- and space-dependent. The initial_condition function is time- and space-dependent to be able to reuse it as analytical solution if available. If no analytical solution is available, the time variable can be ignored in the initial_condition function. The centers are the centers of the kernel functions. By default, centers is set to merge(nodeset_inner, nodeset_boundary). Note that centers needs to have the center number of nodes as the number of nodes in the domain and on the boundary because OrdinaryDiffEq.jl does not support DAEs with rectangular mass matrices.

See also SpatialDiscretization, semidiscretize.

SpatialDiscretization(equations, nodeset_inner, boundary_condition, nodeset_boundary, basis)
SpatialDiscretization(equations, nodeset_inner, boundary_condition, nodeset_boundary,
                      [centers,] kernel = GaussKernel{dim(nodeset_inner)}())

Spatial discretization of a partial differential equation with Dirichlet boundary conditions. The nodeset_inner are the nodes in the domain and nodeset_boundary are the nodes on the boundary. The boundary_condition is a function describing the Dirichlet boundary conditions. The centers are the centers of the kernel functions. By default, centers is set to merge(nodeset_inner, nodeset_boundary). Otherwise, a least squares problem is solved.

See also Semidiscretization, solve_stationary.

semidiscetize(semi::Semidiscretization, tspan)

Wrap a Semidiscretization object into an ODEProblem object with time span tspan.

solve_stationary(spatial_discretization)

Solve a stationary partial differential equation discretized as spatial_discretization with Dirichlet boundary conditions by non-symmetric collocation (Kansa method). Returns an Interpolation object.

interpolation_matrix(centers, kernel, ps, regularization = NoRegularization(); factorization_method = Symmetric)
interpolation_matrix(basis, ps, regularization; factorization_method = Symmetric)

Return the interpolation matrix for the basis, polynomials ps, and regularization. For the StandardBasis, the interpolation matrix is defined as

\[ A = \begin{pmatrix}K & P\\P^T & 0\end{pmatrix},\]

where $K$ is the regularize!d kernel_matrix and $P$ the polynomial_matrix. If a node set of centers and a kernel are given, the interpolation matrix is computed for the StandardBasis. Additionally, you can specify a factorization_method to use for the system matrix. By default, the system matrix is just wrapped as a Symmetric matrix, but you can, e.g., also explicitly use cholesky, lu, or qr factorization.

kernel_matrix(basis, nodeset = centers(basis))
kernel_matrix(nodeset1[, nodeset2], kernel)

Return the kernel matrix for the nodes and kernel. The kernel matrix is defined as

\[ A_{ij} = b_j(x_i),\]

where $b_i$ are the basis function in the basis and x_i are the nodes in the nodeset. If two nodesets and a kernel are given, the kernel matrix is computed for the StandardBasis meaning

\[ A_{ij} = K(\xi_j, x_i),\]

where $\xi_j$ are the nodes/centers in nodeset1, $x_i$ are the nodes in nodeset2, and K is the kernel. If nodeset2 is not given, it defaults to nodeset1.

least_squares_matrix(basis, nodeset, ps, regularization = NoRegularization();
                     factorization_method = Matrix)
least_squares_matrix(centers, nodeset, kernel, ps, regularization = NoRegularization();
                     factorization_method = Matrix)

Return the least squares matrix for the basis, nodeset, polynomials ps, and regularization. For the StandardBasis, the least squares matrix is defined as

\[ A = \begin{pmatrix}K & P_1\\P_2^T & 0\end{pmatrix},\]

where $K$ is the regularize!d kernel_matrix, $P_1$ the polynomial_matrix for the nodeset and $P_2$ the polynomial_matrixfor thecenters. If anodesetandkernelare given, the least squares matrix is computed for the [StandardBasis](@ref). Additionally, you can specify afactorization_methodto use for the system matrix. By default, the system matrix is not factorized, but you can, e.g., also explicitly use theqr` factorization.

operator_matrix(diff_op_or_pde, nodeset_inner, nodeset_boundary, kernel)

Compute the operator matrix $L$ discretizing $\mathcal{L}$ for a given kernel. The operator matrix is defined as

\[ L = A_\mathcal{L} A^{-1},\]

where $A_\mathcal{L}$ is the matrix of the differential operator (defined by the equations), and $A$ the kernel matrix.

See also pde_boundary_matrix and kernel_matrix.

pde_boundary_matrix(diff_op_or_pde, nodeset_inner, nodeset_boundary, [centers,] kernel)

Compute the matrix of a partial differential equation (or differential operator) with a given kernel. The matrix is defined as

\[ A_\mathcal{L} = \begin{pmatrix}\tilde A_\mathcal{L}\\\tilde A\end{pmatrix},\]

where $\tilde A_\mathcal{L}$ is the matrix of the differential operator (defined by the equations) for the inner nodes $x_i$:

\[ (\tilde A_\mathcal{L})_{ij} = \mathcal{L}K(x_i, \xi_j),\]

and $\tilde A$ is the kernel matrix for the boundary nodes:

\[ \tilde A_{ij} = K(x_i, \xi_j),\]

where $\mathcal{L}$ is the differential operator (defined by the equations), $K$ the kernel, $x_i$ are the nodes in nodeset_boundary and $\xi_j$ are the centers. By default, centers is set to merge(nodeset_inner, nodeset_boundary).

See also pde_matrix and kernel_matrix.

pde_matrix(diff_op_or_pde, nodeset1, nodeset2, kernel)

Compute the matrix of a partial differential equation (or differential operator) with a given kernel. The matrix is defined as

\[ (\tilde A_\mathcal{L})_{ij} = \mathcal{L}K(x_i, \xi_j),\]

where $\mathcal{L}$ is the differential operator (defined by the equations), $K$ the kernel, $x_i$ are the nodes in nodeset1 and $\xi_j$ are the nodes in nodeset2.

polynomial_matrix(nodeset, ps)

Return the polynomial matrix for the nodeset and polynomials. The polynomial matrix is defined as

\[ A_{ij} = p_j(x_i),\]

where $x_i$ are the nodes in the nodeset and $p_j$ the polynomials.

AliveCallback(io::IO = stdout; interval::Integer=0, dt=nothing)

Inexpensive callback showing that a simulation is still running by printing some information such as the current time to the screen every interval time steps or after a time of dt in terms of integration time by adding additional (shortened) time steps where necessary (note that this may change the solution).

SummaryCallback(io::IO = stdout)

Create and return a callback that resets the timer at the beginning of a simulation and prints the timer values at the end of the simulation.

SaveSolutionCallback(; interval::Integer=0,
                     dt=nothing,
                     save_initial_solution=true,
                     save_final_solution=true,
                     output_directory="out",
                     extra_functions=(),
                     keys=append!(["itp"], "value_" .* string.(eachindex(extra_functions))))

Save the current numerical solution in regular intervals in VTK format as a Paraview Collection (.pvd). Either pass interval to save every interval time steps or pass dt to save in intervals of dt in terms of integration time by adding additional (shortened) time steps where necessary (note that this may change the solution). The interpolation object will always be saved at the inner and boundary nodes of the corresponding Semidiscretization. You can pass extra functions (time- and space-dependent) or vectors to save at these nodes via extra_functions. The corresponding keys in the .vtu files can be specified by keys.

See also add_to_pvd, vtk_save.

add_to_pvd(filename, pvd, time, nodeset::NodeSet, functions_or_vectors...;
           keys = "value_" .* string.(eachindex(functions_or_vectors)))

Same as vtk_save, but appends the data to a Paraview collection file pvd at time time. In contrast to vtk_save, the functions are time- and space-dependent.

vtk_read(filename)

Read a set of nodes from a VTK file and return them as a NodeSet. Note that the data will always be returned as a 3D NodeSet, even if the data is 1D or 2D. The point data will be returned as a dictionary with the keys being the names of the data arrays in the VTK file.

See also vtk_save.

vtk_save(filename, nodeset::NodeSet, functions_or_vectors...;
         keys = "value_" .* string.(eachindex(functions_or_vectors)))

Save a NodeSet to a VTK file. You can optionally pass a list of space-dependent functions or vectors to save the values of the functions at the nodes. The functions can also be passed as Interpolation or directly as vectors. The optional keyword argument keys is used to specify the names of the data arrays in the VTK file.

See also add_to_pvd, vtk_read.

default_example()

Return the path to an example that can be used to quickly see KernelInterpolation.jl in action. See also examples_dir and get_examples.

examples_dir()

Return the directory where the example files provided with KernelInterpolation.jl are located.

Examples

readdir(examples_dir())

get_examples()

Return a list of all examples that are provided by KernelInterpolation.jl. See also examples_dir and default_example.