This is a tutorial on how to use the Julia package CubicHermiteSpline.jl, which performs a cubic Hermite spline interpolation on an array of data points, $(x_i, y_i)$, given that their associated gradients at each point, $k_i$, are known in advance.

## Getting started

First, the package can be easily installed from the Julia REPL:

1
2
julia> using Pkg


Or in the package mode (typing ] in REPL):

1


Then you can load the package:

1
2
3
4
using CubicHermiteSpline
import Plots
using Plots: plot, plot!
using LaTeXStrings


Plotting related packages are also loaded. Let’s tweak the appearance of our plots:

1
2
3
4
5
6
7
8
9
10
11
12
Plots.default(size=(600, 370))
fntf = :Helvetica
titlefont = Plots.font(fntf, pointsize=12)
guidefont = Plots.font(fntf, pointsize=12)
tickfont = Plots.font(fntf, pointsize=9)
legendfont = Plots.font(fntf, pointsize=8)
Plots.default(fontfamily=fntf)
Plots.default(titlefont=titlefont, guidefont=guidefont, tickfont=tickfont, legendfont=legendfont)
Plots.default(minorticks=true)
Plots.default(linewidth=1.2)
Plots.default(foreground_color_legend=nothing)
Plots.default(legend=false)


Now we are ready to go.

## Example 1: interpolating cubic polynomial

The cubic polynomial of the form

should be exactly interpolated by the cubic Hermite spline interpolation.

Below we use CubicHermiteSpline.jl to demonstrate this fact. First we define a typical cubic polynomial:

1
f(x) = x^3 - 3x^2 + 2x - 5;


Its gradient are available in an analytical form as

1
df(x) = 3x^2 - 6x + 2;


The cubic polynomial is plotted.

1
2
3
4
5
6
a = 0.0
b = 2.5
x_exact = range(a, b, step=0.01);
y_exact = f.(x_exact);

plot(x_exact, y_exact)


Generate a set of data points to be interpolated in the cubic polynomial curve on a evenly spaced grid:

1
2
3
x_input = range(a, b, step=0.5);
y_input = f.(x_input);
k_input = df.(x_input);


Create an interpolation instance:

1
spl = CubicHermiteSplineInterpolation(x_input, y_input, k_input);


Perform the actual cubic Hermite spline interpolation:

1
2
x_interp = range(a, b, step=0.01);
y_interp = spl(x_interp);


Plotting input data points, exact solution, and the interpolated result in a single plot shows that the interpolation is exact for cubic polynomials.

1
2
3
Plots.scatter(x_input, y_input, label="input data", legend=:topleft)
plot!(x_exact, y_exact, label="exact")
plot!(x_interp, y_interp, label="interpolated")


Do the interpolation again for data points on an irregular grid:

1
2
3
4
5
6
7
8
9
10
11
12
x_input = [a, 0.1, 0.6, 1.75, b];
y_input = f.(x_input);
k_input = df.(x_input);

spl = CubicHermiteSplineInterpolation(x_input, y_input, k_input);

x_interp = range(a, b, step=0.01);
y_interp = spl(x_interp);

Plots.scatter(x_input, y_input, label="input data", legend=:topleft)
plot!(x_exact, y_exact, label="exact")
plot!(x_interp, y_interp, label="interpolated")


## Example 2: interpolating smooth functions

Here we use the spherical bessel function of the first kind:

Below we show that the cubic Hermite spline performs impressively well for interpolating such smooth functions.

1
2
3
4
5
6
7
8
9
10
11
g(x) = (sin(x) - x*cos(x)) / x^2;

# dg(x) = ((x^2 - 2) * sin(x) + 2x*cos(x)) / x^3;
dg(x) = (sin(x) - 2*g(x)) / x;

a = 0.01
b = 8.01
x_exact = range(a, b, step=0.01);
y_exact = g.(x_exact);

plot(x_exact, y_exact)


Perform cubic Hermite interpolation on an even grid:

1
2
3
4
5
6
7
8
9
10
11
12
x_input = range(a, b, step=0.5);
y_input = g.(x_input);
k_input = dg.(x_input);

spl = CubicHermiteSplineInterpolation(x_input, y_input, k_input);

x_interp = range(a, b, step=0.01);
y_interp = spl(x_interp);

Plots.scatter(x_input, y_input, label="input data", legend=:topright)
plot!(x_exact, y_exact, label="exact")
plot!(x_interp, y_interp, label="interpolated")


It can be seen from above plot that the interpolation works perfectly for known data points on a dense even grid.

Interpolate the function on an unevenly spaced grid:

1
2
3
4
5
6
7
8
9
10
11
12
x_input = [a, 0.6, 1.5, 3.4, 4.5, 5.0, 6.2, 7.5, b];
y_input = g.(x_input);
k_input = dg.(x_input);

spl = CubicHermiteSplineInterpolation(x_input, y_input, k_input);

x_interp = range(a, b, step=0.01);
y_interp = spl(x_interp);

Plots.scatter(x_input, y_input, label="input data", legend=:topright)
plot!(x_exact, y_exact, label="exact")
plot!(x_interp, y_interp, label="interpolated")


Note that the accuracy of the interpolation degrades near the extrema of the function if there is no input data point near that region.

## Acknowledgements

• This work is partially supported by the General Program of the National Natural Science Foundation of China (No. 21873021) and Shanghai Pujiang Program (No. 18PJ1401200).
• The captured output text blocks in this post is automatically generated by using Literate.jl.
Updated on