| |
- Here we have a code for computing cubic interpolating splines.
These may serve for approximating a given function f
by a two times differentiable function, composed from polynomials of degree not larger than three.
While interpolation by polynomials gives useful results on small intervals only, this type
of approximation is not sensitive against the total length of the interval but only to the distance of the
abscissae.
- The function obtained, we name it S here, is piecewise defined, and
each such piece lives between two successive abscissae of the grid
a = x1 < x2 < ... < xn = b
- These pieces, named Si, i=1,...,n-1 here, are determined from the following conditions:
- Si(x) is a polynomial of degree not greater 3.
- S(x) = Si(x) for xi <= x <= xi+1
- Si(xi) = f(xi) (or yi).
- Si(xi+1) = f(xi+1).
- S'i(xi+1) =
S'i+1(xi+1), i=1,..,n-2 (continuity of S').
- S''i(xi+1) =
S''i+1(xi+1), i=1,...,n-2 (continuity of S'').
Observe that the continuity of S itself is obtained by the interpolation condition.
- These conditions lead to linear system of 4(n-2)+2 equations for the
4(n-1) coefficients
which is broken down here in a special way into a tridiagonal linear system and n-1
uncoupled systems of dimension 3 with a triangular structure.
- This construction is always possible, but not unique.
There remain two degrees of freedom, which allow to fix properties of S
at the boundaries of its interval of definition.
Several choices are possible and in use, for example continuity of S''' at
x2 and xn-1 (giving the ordinary cubic interpolation in case n=4)
, periodicity of S' and S'',
which of course is useful only if also f is periodic (i.e. y1=yn)
and the two used here: the so called ''natural'' spline with S''(a) = S''(b) = 0 and
the hermitian spline with S'(a) = f'(a) and S'(b) = f'(b).
The natural spline got its name from the spline of the shipbuilder, a thin elastic ruler,
and the hermitian from Hermites interpolation scheme, where values of f
and f' are interpolated simultaneously.
- The spline has the remarkable property of minimizing the total curvature integral
∫ab (g''(x))2dx
over the set of all two times continuously differentiable functions g which satisfy the interpolation
and the boundary conditions specified.
- The approximation error of those two constructions is of the order
h4, with h = max(xi+1-xi),
provided that f is four times continuously differentiable on [a,b].
This holds true for the hermitian spline on the whole interval and for the natural one
in [x1+ε,xn-ε] with the error constant growing to infinity
for ε->0.
If in case of the hermitian spline one replaces the true values f'(a), f'(b) by the derivative of the cubic interpolation
polynomial for the first resp. last four points of the grid at a resp. b, then
the error constant raises a bit, but the overall quality is unchanged. Hence one can use this
if no derivative is available explicitly.
|