hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_interp_1d_spline (e01ba)

Purpose

nag_interp_1d_spline (e01ba) determines a cubic spline interpolant to a given set of data.

Syntax

[lamda, c, ifail] = e01ba(x, y, 'm', m)
[lamda, c, ifail] = nag_interp_1d_spline(x, y, 'm', m)

Description

nag_interp_1d_spline (e01ba) determines a cubic spline s(x)s(x), defined in the range x1xxmx1xxm, which interpolates (passes exactly through) the set of data points (xi,yi)(xi,yi), for i = 1,2,,mi=1,2,,m, where m4m4 and x1 < x2 < < xmx1<x2<<xm. Unlike some other spline interpolation algorithms, derivative end conditions are not imposed. The spline interpolant chosen has m4m-4 interior knots λ5,λ6,,λmλ5,λ6,,λm, which are set to the values of x3,x4,,xm2x3,x4,,xm-2 respectively. This spline is represented in its B-spline form (see Cox (1975)):
m
s(x) = ciNi(x),
i = 1
s(x)=i=1mciNi(x),
where Ni(x)Ni(x) denotes the normalized B-spline of degree 33, defined upon the knots λi,λi + 1,,λi + 4λi,λi+1,,λi+4, and cici denotes its coefficient, whose value is to be determined by the function.
The use of B-splines requires eight additional knots λ1λ1, λ2λ2, λ3λ3, λ4λ4, λm + 1λm+1, λm + 2λm+2, λm + 3λm+3 and λm + 4λm+4 to be specified; nag_interp_1d_spline (e01ba) sets the first four of these to x1x1 and the last four to xmxm.
The algorithm for determining the coefficients is as described in Cox (1975) except that QRQR factorization is used instead of LULU decomposition. The implementation of the algorithm involves setting up appropriate information for the related function nag_fit_1dspline_knots (e02ba) followed by a call of that function. (See nag_fit_1dspline_knots (e02ba) for further details.)
Values of the spline interpolant, or of its derivatives or definite integral, can subsequently be computed as detailed in Section [Further Comments].

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

Parameters

Compulsory Input Parameters

1:     x(m) – double array
m, the dimension of the array, must satisfy the constraint m4m4.
x(i)xi must be set to xixi, the iith data value of the independent variable xx, for i = 1,2,,mi=1,2,,m.
Constraint: x(i) < x(i + 1)xi<xi+1, for i = 1,2,,m1i=1,2,,m-1.
2:     y(m) – double array
m, the dimension of the array, must satisfy the constraint m4m4.
y(i)yi must be set to yiyi, the iith data value of the dependent variable yy, for i = 1,2,,mi=1,2,,m.

Optional Input Parameters

1:     m – int64int32nag_int scalar
Default: The dimension of the arrays x, y. (An error is raised if these dimensions are not equal.)
mm, the number of data points.
Constraint: m4m4.

Input Parameters Omitted from the MATLAB Interface

lck wrk lwrk

Output Parameters

1:     lamda(lck) – double array
lckm + 4lckm+4.
The value of λiλi, the iith knot, for i = 1,2,,m + 4i=1,2,,m+4.
2:     c(lck) – double array
lckm + 4lckm+4.
The coefficient cici of the B-spline Ni(x)Ni(x), for i = 1,2,,mi=1,2,,m. The remaining elements of the array are not used.
3:     ifail – int64int32nag_int scalar
ifail = 0ifail=0 unless the function detects an error (see [Error Indicators and Warnings]).

Error Indicators and Warnings

Errors or warnings detected by the function:
  ifail = 1ifail=1
On entry,m < 4m<4,
orlck < m + 4lck<m+4,
orlwrk < 6 × m + 16lwrk<6×m+16.
  ifail = 2ifail=2
The x-values fail to satisfy the condition
x(1) < x(2) < x(3) < < x(m)x1<x2<x3<<xm.

Accuracy

The rounding errors incurred are such that the computed spline is an exact interpolant for a slightly perturbed set of ordinates yi + δyiyi+δyi. The ratio of the root-mean-square value of the δyiδyi to that of the yiyi is no greater than a small multiple of the relative machine precision.

Further Comments

The time taken by nag_interp_1d_spline (e01ba) is approximately proportional to mm.
All the xixi are used as knot positions except x2x2 and xm1xm-1. This choice of knots (see Cox (1977)) means that s(x)s(x) is composed of m3m-3 cubic arcs as follows. If m = 4m=4, there is just a single arc space spanning the whole interval x1x1 to x4x4. If m5m5, the first and last arcs span the intervals x1x1 to x3x3 and xm2xm-2 to xmxm respectively. Additionally if m6m6, the iith arc, for i = 2,3,,m4i=2,3,,m-4, spans the interval xi + 1xi+1 to xi + 2xi+2.
After the call
[lamda, c, ifail] = e01ba(x, y, lck);
the following operations may be carried out on the interpolant s(x)s(x).
The value of s(x)s(x) at x = xx=x can be provided in the double variable s by the call
[s, ifail] = e02bb(lamda, c, x);
(see nag_fit_1dspline_eval (e02bb)).
The values of s(x)s(x) and its first three derivatives at x = xx=x can be provided in the double array s of dimension 44, by the call
[s, ifail] = e02bc(lamda, c, x, left);
(see nag_fit_1dspline_deriv (e02bc)).
Here left must specify whether the left- or right-hand value of the third derivative is required (see nag_fit_1dspline_deriv (e02bc) for details).
The value of the integral of s(x)s(x) over the range x1x1 to xmxm can be provided in the double variable dint by
[dint, ifail] = e02bd(lamda, c);
(see nag_fit_1dspline_integ (e02bd)).

Example

function nag_interp_1d_spline_example
x = [0;
     0.2;
     0.4;
     0.6;
     0.75;
     0.9;
     1];
y = [1;
     1.22140275816017;
     1.49182469764127;
     1.822118800390509;
     2.117000016612675;
     2.45960311115695;
     2.718281828459045];
[lamda, c, ifail] = nag_interp_1d_spline(x, y)
 

lamda =

         0
         0
         0
         0
    0.4000
    0.6000
    0.7500
    1.0000
    1.0000
    1.0000
    1.0000


c =

    1.0000
    1.1336
    1.3726
    1.7827
    2.1744
    2.4918
    2.7183
         0
         0
         0
         0


ifail =

                    0


function e01ba_example
x = [0;
     0.2;
     0.4;
     0.6;
     0.75;
     0.9;
     1];
y = [1;
     1.22140275816017;
     1.49182469764127;
     1.822118800390509;
     2.117000016612675;
     2.45960311115695;
     2.718281828459045];
[lamda, c, ifail] = e01ba(x, y)
 

lamda =

         0
         0
         0
         0
    0.4000
    0.6000
    0.7500
    1.0000
    1.0000
    1.0000
    1.0000


c =

    1.0000
    1.1336
    1.3726
    1.7827
    2.1744
    2.4918
    2.7183
         0
         0
         0
         0


ifail =

                    0



PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2013