## NormalHermiteSplines.jl

Multivariate Normal Hermite-Birkhoff Interpolating Splines in Julia
Author IgorKohan
Popularity
11 Stars
Updated Last
2 Years Ago
Started In
July 2020

# Multivariate Normal Hermite-Birkhoff Interpolating Splines in Julia

## Instalation

Start Julia and run the following commands:

``````julia> using Pkg
``````

## Example usage

To use this package, begin your code with

``````using NormalHermiteSplines
``````

The first example is the function `φ(x,y)=sin(4.0*sqrt(x^2+y^2))` defined in the region `Ω=[-1,1]x[-1,1]`.

We'll construct an interpolating normal spline using this function and its gradient values sampled on set of 1000 Halton nodes ([1]) distributed in the [-1,1]x[-1,1] square.

``````    using NormalHermiteSplines

nodes = get_2D_halton_nodes(1000)         # generates Halton data set in [0, 1] x [0, 1]
n_1 = size(nodes, 2)
u = Vector{Float64}(undef, n_1)           # function values
d_nodes = Matrix{Float64}(undef, 2, n_1)  # directional derivative nodes
es = Matrix{Float64}(undef, 2, n_1)       # derivative directions
du = Vector{Float64}(undef, n_1)          # directional derivative values

for i = 1:n_1
nodes[1, i] = nodes[1, i] * 2.0 - 1.0   # transforming Halton nodes to [-1, 1] x [-1, 1]
nodes[2, i] = nodes[2, i] * 2.0 - 1.0
d_nodes[1,i] = nodes[1,i]
d_nodes[2,i] = nodes[2,i]
d = sqrt(nodes[1, i]^2 + nodes[2, i]^2) # here d > 0
u[i] = sin(4.0 * d)
grad[1] = 4.0 * nodes[1, i] * cos(4.0 * d) / d
grad[2] = 4.0 * nodes[2, i] * cos(4.0 * d) / d
es[1,i] = grad[1]                       # no need to normalize 'es' vectors
end

# Hermite spline must be constructed with RK_H1 or RK_H2 kernel.
# Here value of the 'scaling parameter' ε is estimated in the interpolate procedure.
rk = RK_H1()
#
spline = interpolate(nodes, u, d_nodes, es, du, rk)

grid = get_2D_grid2(100)     # regular grid covering [-1, 1] x [-1, 1]
σ = evaluate(spline, grid)   # spline values on the grid

σ1 = evaluate_one(spline, [0.5; 0.5])
#  ≈ 0.308071
#
#  exact function value at [0.5; 0.5]:
#  ≈ 0.308072

#  ≈ -2.69080
#  ≈ -2.69090
#
# exact function gradient at [0.5; 0.5]:
#  ≈ -2.69086
#  ≈ -2.69086
``````

The spline surface and filled 2-D contour plots:

Approximation error plots:

Spline was evaluated on a uniform Cartesian grid of size 101x101. Accuracy of the interpolation was measured by calculating the Root Mean Square Error (RMSE) and the Maximum Absolute Error (MAE). For this case `RMSE`: 1.6E-03, `MAE`: 1.1E-01, estimated value of the scaling parameter `ε` is 8.8, estimation of the Gram matrix condition number is 1.0E+11.

The second example is the function `Ψ(x,y,z)=cos(π*x)*cos(y-0.5)*sin(π*(z-0.5))` defined in the region `Ω=[0,1]x[0,1]x[0,1]`.

We'll construct an interpolating normal spline using function `Ψ` values sampled on set of 1000 non-uniform random nodes distributed in the unit cube [0,1]x[0,1]x[0,1].

``````    using NormalHermiteSplines

nodes = get_3D_random_grid(9)       # generates 1000 non-uniform random nodes
n_1 = size(nodes, 2)
u = Vector{Float64}(undef, n_1)     # function values
for i = 1:n_1
x = nodes[1,i]
y = nodes[2,i]
z = nodes[3,i]
u[i] = cos(π*x)*cos(y-0.5)*sin(π*(z-0.5))
end

# Here spline is being constructed with RK_H2 kernel,
# the 'scaling parameter' ε is defined explicitly.
rk = RK_H2(5.0)
#
spline = interpolate(nodes, u, rk)
grid = get_3D_grid(50) # creates the uniform Cartesian grid of size 51x51x51 in [0, 1]x[0, 1]x[0, 1]
σ = evaluate(spline, grid)

σ1 = evaluate_one(spline, [0.8; 0.6; 0.8])
#  ≈ -0.65122
#
#  exact function value at [0.8; 0.6; 0.8]:
#  ≈ -0.65124

g1 = evaluate_gradient(spline, [0.8; 0.6; 0.8])
#  ≈ -1.4862
#  ≈  0.0653
#  ≈ -1.4863
#
# exact function gradient at [0.8; 0.6; 0.8]:
#  ≈ -1.4865
#  ≈  0.0653
#  ≈ -1.4865
``````

The spline plots:

Approximation error plots:

Spline was evaluated on a uniform Cartesian grid of size 51x51x51. For this case `RMSE`: 2.9E-03, `MAE`: 5.1E-02, value of the scaling parameter `ε` is 5.0, estimation of the Gram matrix condition number is 1.0E+14.

Further examples are given in documentation.

The normal splines method for one-dimensional function interpolation and linear ordinary differential and integral equations was proposed in [2]. An idea of the multivariate splines in Sobolev space was initially formulated in [8], however it was not well-suited to solving real-world problems. Using that idea the multivariate generalization of the normal splines method was developed for two-dimensional problem of low-range computerized tomography in [3] and applied for solving a mathematical economics problem in [4]. At the same time an interpolation scheme with Matérn kernels was developed in [9], this scheme coincides with interpolating normal splines method. Further results related to applications of the normal splines method were reported at the seminars and conferences [5,6,7].

## Documentation

References

[2] V. Gorbunov, The method of normal spline collocation. USSR Computational Mathematics and Mathematical Physics, Vol. 29, No. 1, 1989

[3] I. Kohanovsky, Normal Splines in Computing Tomography (Нормальные сплайны в вычислительной томографии). Avtometriya, No.2, 1995

[4] V. Gorbunov, I. Kohanovsky, K. Makedonsky, Normal splines in reconstruction of multi-dimensional dependencies. Papers of WSEAS International Conference on Applied Mathematics, Numerical Analysis Symposium, Corfu, 2004

[5] I. Kohanovsky, Multidimensional Normal Splines and Problem of Physical Field Approximation, International Conference on Fourier Analysis and its Applications, Kuwait, 1998.

[6] I. Kohanovsky, Inequality-Constrained Multivariate Normal Splines with Some Applications in Finance. 27th GAMM-Seminar on Approximation of Multiparametric functions, Max-Planck-Institute for Mathematics in the Sciences, Leipzig, Germany, 2011.

[7] V. Gorbunov, I. Kohanovsky, Heterogeneous Parallel Method for the Construction of Multi-dimensional Smoothing Splines. ESCO 2014 4th European Seminar on Computing, 2014

[8] A. Imamov, M. Dzhurabaev, Splines in S.L. Sobolev spaces (Сплайны в пространствах С.Л.Соболева). Deposited manuscript. Dep. UzNIINTI, No 880, 1989.

[9] J. Dix, R. Ogden, An Interpolation Scheme with Radial Basis in Sobolev Spaces H^s(R^n). Rocky Mountain J. Math. Vol. 24, No.4, 1994.