NAG Library Function Document

nag_1d_spline_interpolant (e01bac)

+− Contents

1 Purpose

2 Specification

3 Description

4 References

5 Arguments

6 Error Indicators and Warnings

7 Accuracy

8 Parallelism and Performance

9 Further Comments

+− 10 Example

10.1 Program Text

10.2 Program Data

10.3 Program Results

1 Purpose

nag_1d_spline_interpolant (e01bac) determines a cubic spline interpolant to a given set of data.

2 Specification

#include <nag.h>

#include <nage01.h>

void	nag_1d_spline_interpolant (Integer m, const double x[], const double y[], Nag_Spline spline, NagError fail)

3 Description

nag_1d_spline_interpolant (e01bac) determines a cubic spline

s (x)

, defined in the range

x_{1} \leq x \leq x_{m}

, which interpolates (passes exactly through) the set of data points

(x_{i}, y_{i})

, for

i = 1, 2, \dots, m

, where

m \geq 4

and

x_{1} < x_{2} < \dots < x_{m}

. Unlike some other spline interpolation algorithms, derivative end conditions are not imposed. The spline interpolant chosen has

m - 4

interior knots

λ_{5}, λ_{6}, \dots, λ_{m}

, which are set to the values of

x_{3}, x_{4}, \dots, x_{m - 2}

respectively. This spline is represented in its B-spline form (see Cox (1975)):

s (x) = \sum_{i = 1}^{m} c_{i} N_{i} (x)

where

N_{i} (x)

denotes the normalized B-spline of degree 3, defined upon the knots

λ_{i}, λ_{i + 1}, \dots, λ_{i + 4}

, and

c_{i}

denotes its coefficient, whose value is to be determined by the function.

The use of B-splines requires eight additional knots

λ_{1}

λ_{2}

λ_{3}

λ_{4}

λ_{m + 1}

λ_{m + 2}

λ_{m + 3}

and

λ_{m + 4}

to be specified; the function sets the first four of these to

x_{1}

and the last four to

x_{m}

The algorithm for determining the coefficients is as described in Cox (1975) except that

Q R

factorization is used instead of

L U

decomposition. The implementation of the algorithm involves setting up appropriate information for the related function nag_1d_spline_fit_knots (e02bac) followed by a call of that function. (For further details of nag_1d_spline_fit_knots (e02bac), see the function document.)

Values of the spline interpolant, or of its derivatives or definite integral, can subsequently be computed as detailed in Section 9.

4 References

Cox M G (1975) An algorithm for spline interpolation J. Inst. Math. Appl. 15 95–108

Cox M G (1977) A survey of numerical methods for data and function approximation The State of the Art in Numerical Analysis (ed D A H Jacobs) 627–668 Academic Press

5 Arguments

1: m – IntegerInput

On entry:

m

, the number of data points.

Constraint:

m \geq 4

2: x[m] – const doubleInput

On entry:

x [i - 1]

must be set to

x_{i}

, the

i

th data value of the independent variable

x

, for

i = 1, 2, \dots, m

Constraint:

x [i] < x [i + 1]

, for

i = 0, 1, \dots, m - 2

3: y[m] – const doubleInput

On entry:

y [i - 1]

must be set to

y_{i}

, the

i

th data value of the dependent variable

y

, for

i = 1, 2, \dots, m

4: spline – Nag_Spline *

Pointer to structure of type Nag_Spline with the following members:

n – IntegerOutput: On exit: the size of the storage internally allocated to $lamda$ . This is set to $m + 4$ .

lamda – double *Output: On exit: the pointer to which storage of size $n$ is internally allocated. $lamda [i - 1]$ contains the $i$ th knot, for $i = 1, 2, \dots, m + 4$ .

c – double *Output: On exit: the pointer to which storage of size $n - 4$ is internally allocated. $c [i - 1]$ contains the coefficient $c_{i}$ of the B-spline $N_{i} (x)$ , for $i = 1, 2, \dots, m$ .

Note that when the information contained in the pointers

lamda

and

c

is no longer of use, or before a new call to nag_1d_spline_interpolant (e01bac) with the same spline, you should free this storage using the NAG macro NAG_FREE. This storage will not have been allocated if this function returns with

fail . code \neq

NE_NOERROR.

5: fail – NagError *Input/Output

The NAG error argument (see Section 3.6 in the Essential Introduction).

6 Error Indicators and Warnings

NE_ALLOC_FAIL: Dynamic memory allocation failed.
NE_INT_ARG_LT: On entry, $m = ⟨value⟩$ .
Constraint: $m \geq 4$ .
NE_NOT_STRICTLY_INCREASING: The sequence x is not strictly increasing: $x [⟨value⟩] = ⟨value⟩$ , $x [⟨value⟩] = ⟨value⟩$ .

7 Accuracy

The rounding errors incurred are such that the computed spline is an exact interpolant for a slightly perturbed set of ordinates

y_{i} + δ y_{i}

. The ratio of the root-mean-square value of the

δ y_{i}

to that of the

y_{i}

is no greater than a small multiple of the relative machine precision.

8 Parallelism and Performance

Not applicable.

9 Further Comments

The time taken by nag_1d_spline_interpolant (e01bac) is approximately proportional to

m

All the

x_{i}

are used as knot positions except

x_{2}

and

x_{m - 1}

. This choice of knots (see Cox (1977)) means that

s (x)

is composed of

m - 3

cubic arcs as follows. If

m = 4

, there is just a single arc space spanning the whole interval

x_{1}

x_{4}

. If

m \geq 5

, the first and last arcs span the intervals

x_{1}

x_{3}

and

x_{m - 2}

x_{m}

respectively. Additionally if

m \geq 6

, the

i

th arc, for

i = 2, 3, \dots, m - 4

, spans the interval

x_{i + 1}

x_{i + 2}

After the call

e01bac(m, x, y, &spline, &fail)

the following operations may be carried out on the interpolant

s (x)

The value of

s (x)

x = xval

can be provided in the variable sval by calling the function

e02bbc(xval, &sval, &spline, &fail)

The values of

s (x)

and its first three derivatives at

x = xval

can be provided in the array sdif of dimension 4, by the call

e02bcc(derivs, xval, sdif, &spline, &fail)

Here derivs must specify whether the left- or right-hand value of the third derivative is required (see nag_1d_spline_deriv (e02bcc) for details). The value of the integral of

s (x)

over the range

x_{1}

x_{m}

can be provided in the variable sint by

e02bdc(&spline, &sint, &fail)

10 Example

The following example program sets up data from 7 values of the exponential function in the interval 0 to 1. nag_1d_spline_interpolant (e01bac) is then called to compute a spline interpolant to these data.

The spline is evaluated by nag_1d_spline_evaluate (e02bbc), at the data points and at points halfway between each adjacent pair of data points, and the spline values and the values of