next up previous index
Next: E211 Cubic Splines Up: CERNLIB Previous: E208 Least Squares

E210 Polynomial Splines / Normalized B-Splines

Routine ID: E210
Author(s): W. Mönch, B. SchorrLibrary: MATHLIB
Submitter: W. MönchSubmitted: 15.03.1993
Language: FortranRevised:

NORBAS ( NORmalized BAsis Splines) is a portable collection of subprograms for various calculations with polynomial splines in one dimension ( 1D) and in two dimensions ( 2D). The polynomial splines are represented as linear combinations of normalized basis splines (B-splines).

On computers other than CDC or Cray, only the double-precision versions DSPKN1, etc. are available. On CDC and Cray computers, only the single-precision versions RSPKN1, etc. are available.

The following outline provides the background material and the notations needed for describing the subprograms and their parameters. For further information about splines and their applications see References, in particular Ref. 7.
Case (1D):

k
Degree (order -1) of the B-spline (0 ≤k ≤25) .
m
Number of spline-knots ( m ≥2  k+2 ) .
i
Index of the B-spline (1 ≤i ≤m-k-1) .
τ

Set of m spline-knots τ= t1,t2, ...,tm , in non-decreasing order, with multiplicity ≤k+1

(i.e. no more than k+1 knots coincide).

[a,b]
Interval, defined by a=tk+1 , b=tm-k .
Bi(x)

Normalized B-spline of degree k over τ with index i. The value of Bi(x) is identically zero outside the interval ti≤x ≤ti+k+1 , and the normalization of Bi(x)

is such that

- &inf;+ &inf;Bi(x) dx = {ti+k+1-tik+1}(i = 1, ..., m-k-1).

s(x)

Polynomial spline at x ∈[a,b] in B-spline representation

y = s(x) = ∑i=1m-k-1 ci Bi(x)  .

Spline interpolation to a data set:

Given a data set xl,yll=1, ...,n ; determine coefficients cii=1, ...,n of a polynomial interpolation spline y = s(x) in B-spline representation with degree k over a set τ of m=n+k+1 knots, such that the following relations (interpolation conditions) hold:

s(xl) = yl(l=1, ...,n).

The existence of a solution of this interpolation problem depends on an appropriate choice of the spline-knots (Ref. 7, Theorem XIII.1 (Schoenberg-Whitney)).

Least squares spline approximation to a data set:

Given a data set xl,yll=1, ...,n ; determine coefficients cii=1, ...,m-k-1 of a polynomial approximation spline y = s(x) in B-spline representation with degree k over a set τ of m ≤n+k+1 knots, such that following least squares problem is solved:

φ(c1, ...,cm-k-1) = ∑l=1n (yl- s(xl))2 = min  !

Variation diminishing spline approximation to a function (Schoenberg):

For a given function y = f(x) on [a,b] this spline approximation is defined by y = s(x) , with
(Ref. 7, p. 158-162)

ci = f(xi); xi = (ti+1+...+ti+k)/k&quad;(i=1,...,m-k-1; k ≥1).

Case (2D):

kx
Degree of one-dimensional B-splines in x-direction (0 ≤kx ≤25) .
ky
Degree of one-dimensional B-splines in y-direction (0 ≤ky ≤25) .
mx
Number of spline-knots in x-direction (mx ≥2  kx+2) .
my
Number of spline-knots in y-direction (my ≥2  ky+2) .
i
Index of B-spline (1 ≤i ≤mx-kx-1) in x-direction.
j
Index of B-spline (1 ≤j ≤my-ky-1) in y-direction.
τx

Set of mx spline-knots τx= tx,1,tx,2, ...,tx,mx in x-direction, in non-decreasing order, with multiplicity ≤kx+1

(i.e. no more than kx+1 knots coincide).

τy

Set of my spline-knots τy= ty,1,ty,2, ...,ty,my in y-direction, in non-decreasing order, with multiplicity ≤ky+1

(i.e. no more than ky+1 knots coincide).

[ax,bx]
Interval in x-direction, defined by ax=tx,kx+1 , bx=tx,mx-kx .
[ay,by]
Interval in y-direction, defined by ay=ty,ky+1 , by=ty,my-ky .
Bi(x)

B-spline of degree kx over τx with index i.
Bj(y)

B-spline of degree ky over τy with index j.
Bi,j(x,y)

Product Bi,j(x,y) = Bi(x)  Bj(y) of two one-dimensional B-splines.
s(x,y)

Two-dimensional polynomial spline at (x,y) ∈[ax,bx]x[ay,by]

in B-spline representation

z = s(x,y) = ∑i=1mx-kx-1 ∑j=1my-ky-1 ci,j Bi,j(x,y).

Spline interpolation to a data set:
Given a data set xlx,yly,zlx,lylx=1, ...,nx ;ly=1, ..., ny ; determine coefficients ci,ji=1, ...,nx ; j=1, ...,ny of a two-dimensional polynomial interpolation spline z = s(x,y) in B-spline representation with degrees kx , ky over the sets τx

of mx=nx+kx+1 knots in x-direction and τy of my=ny+ky+1 knots in y-direction, such that following relations (interpolation conditions) hold:

s(xlx,yly) = zlx,ly(lx=1, ...,nx  ;  ly=1, ...,ny) .

The existence of a solution of this interpolation problem depends on an appropriate choice of the spline-knots τx , τy in the two-dimensional interpolation area [ax,bx]x[ay,by] .

Least squares spline approximation to a data set:

Given a data set xlx,yly,zlx,lylx=1, ...,nx ;ly=1, ..., ny ; determine coefficients ci,ji=1, ...,nx ; j=1, ...,ny of a two-dimensional polynomial approximation spline z = s(x,y) in B-spline representation with degrees kx, ky over the sets τx

of mx ≤nx+kx+1 knots in x-direction and τy of my ≤ny+ky+1 knots in y-direction, such that following least squares problem is solved:

φ(c1,1, ...,cmx-kx-1,my-ky-1) =  ∑lx=1nx ∑ly=1ny (zlx,ly-s(xlx,yly))2 = min !

Variation diminishing spline approximation to a function:

For a given function z = f(x,y) on [ax,by]x[ay,by] this two-dimensional spline approximation is defined by z = s(x,y)

on [ax,bx]x[ay,by] , with

ci,j f(xi,yj); xi (tx,i+1+...+tx,i+kx)/kx (i=1,...,mx-kx-1; kx ≥1),

2c yj (ty,j+1+...+ty,j+ky)/ky (j=1,...,my-ky-1; ky ≥1).

The package NORBAS contains FUNCTION and SUBROUTINE subprograms for solving the following problems. To calculate:

(K)
A set τ of m spline-knots in the interval [a,b] for normalized B-splines Bi(x) of degree k, use subprogram tSPKN1 ( 1D). The knots are in non-decreasing order and determined in such a way that the first k+1 knots coincide with a, the last k+1 knots coincide with b, and the remaining (m-2 k-2) knots are equidistant in (a,b) .

Two sets τx , τy of spline-knots in [ax,bx] and [ay,by] for B-splines Bi,j(x,y)

of degrees kx and ky in x- and y-direction, use subprogram tSPKN2 ( 2D). τx and τy , are calculated by the same formulae in x-, and y-direction, as in case ( 1D).

(B)
The function Bi(x) , the n-th derivative {dnBi(x)dxn} , or the integral -&inf;xBi(ξ) dξ

of a normalized B-spline Bi(x) , with fixed degree k and index i over a set τ of spline-knots, use subprogram: tSPNB1 ( 1D).

The function Bi,j(x,y) , the partial derivative {nxny Bi,j(x,y)xnx yny} , or the integral -&inf;x ∫-&inf;y Bi,j(ξ,η) dη dξ

of a two-dimensional B-spline Bi,j(x,y) , with fixed degrees kx, ky and indices i, j over the sets τx , τy of spline-knots, use subprogram tSPNB2 ( 2D).

(P)
The function s(x) , the n-th derivative {dns(x)dxn} , or the integral -&inf;xs(ξ) dξ

of a polynomial spline y = s(x) in B-spline representation with given coefficients ci , use subprogram tSPPS1 ( 1D).

The function s(x,y) , the partial derivative {nxny s(x,y)xnx yny} , or the integral -&inf;x ∫-&inf;y s(ξ,η) dη dξ

of a two-dimensional polynomial spline z = s(x,y)

in B-spline representation with given coefficients ci,j , use subprogram tSPPS2 ( 2D).

(I)
The coefficients ci of a one-dimensional polynomial interpolation spline y = s(x) in B-spline representation to a user-supplied data set xl,yl , use subprogram tSPIN1 ( 1D).

The coefficients ci,j of a two-dimensional polynomial interpolation spline z = s(x,y) in B-spline representation to a user-supplied data set xlx,yly,zlx,ly , use subprogram tSPIN2 ( 2D).

(A)
The coefficients ci of a one-dimensional polynomial least squares approximation spline y=s(x) in B-spline representation to a user-supplied data set xl,yl , use subprogram tSPAP1 ( 1D).

The coefficients ci,j of a two-dimensional polynomial least squares approximation spline z=s(x,y) in B-spline representation to a user-supplied data set xlx,yly,zlx,ly , use subprogram tSPAP2 ( 2D).

(V)
The coefficients ci of a one-dimensional polynomial variation diminishing spline approximation y = s(x)

in B-spline representation to a user-supplied function y = f(x) , use subprogram tSPVD1 ( 1D).

The coefficients ci,j of a two-dimensional polynomial variation diminishing spline approximation z = s(x,y)

in B-spline representation to a user-supplied function z = f(x,y) , use subprogram tSPVD2 ( 2D).

(D)
From given coefficients ci of a one-dimensional polynomial spline y = s(x) in B-spline representation, the corresponding coefficients di of its n-th derivative dns(x)/dxn in B-spline representation, use subprogram tSPCD1 ( 1D).

From given coefficients ci,j of a two-dimensional polynomial spline z = s(x,y) in B-spline representation, the corresponding coefficients di,j of its partial derivative nxny s(x,y)/xnx yny in B-spline representation, use subprogram tSPCD2 ( 2D).

Structure:

SUBROUTINE and FUNCTION subprograms
User Entry Names:
RSPKN1, RSPKN2, RSPNB1, RSPNB2, RSPPS1, RSPPS2, RSPIN1, RSPIN2,
RSPAP1, RSPAP2, RSPVD1, RSPVD2, RSPCD1, RSPCD2,
DSPKN1, DSPKN2, DSPNB1, DSPNB2, DSPPS1, DSPPS2, DSPIN1, DSPIN2,
DSPAP1, DSPAP2, DSPVD1, DSPVD2, DSPCD1, DSPCD2,

Internal Entry Names: RSPAS1, RSPAS2, RSPLKK, DSPAS1 DSPAS2, DSPLKK
Files Referenced: Unit 6
External References:
RGBTRF, RGBTRS, RGESVD, DGBTRF, DGBTRS, DGESVD(F001)
RVSET, RVSUM, RVCPY, RVMPY,
DVSET, DVSUM, DVCPY, DVMPY(F002)
RMCPY, RMMPY, DMCPY, DVMPY(F003)
8lMTLMTR(N002), ABEND(Z035);User supplied FUNCTIONsubprogram

Usage:

For t = R (type REAL), t = D (type DOUBLE PRECISION):

(K) Knots

    CALL tSPKN1(K,M,A,B,T,NERR)
    CALL tSPKN2(KX,KY,MX,MY,AX,BX,AY,BY,TX,TY,NERR)
(B) Normalized B-splines
    tSPNB1(K,M,I,NDER,X,T,NERR)
    tSPNB2(KX,KY,MX,MY,I,J,NDERX,NDERY,X,Y,TX,TY,NERR)
(P) Polynomial spline
    tSPPS1(K,M,NDER,X,T,C,NERR)
    tSPPS2(KX,KY,MX,MY,NDERX,NDERY,X,Y,TX,TY,C,NDIMC,W,NERR)
(I) Spline interpolation
    CALL tSPIN1(K,N,XI,YI,KNOT,T,C,W,IW,NERR)
    CALL tSPIN2(KX,KY,NX,NY,XI,YI,ZI,NDIMZ,KNOT,TX,TY,C,NDIMC,W,IW,NERR)
(A) Least squares spline approximation
    CALL tSPAP1(K,M,N,XI,YI,KNOT,T,C,W,NW,NERR)
    CALL tSPAP2(KX,KY,MX,MY,NX,NY,XI,YI,ZI,NDIMZ,KNOT,TX,TY,C,NDIMC,W,NW,NERR)
(V) Variation diminishing spline approximation
    CALL tSPVD1(F,K,M,T,C,NERR)
    CALL tSPVD2(F,KX,KY,MX,MY,TX,TY,C,NDIMC,NERR)
(D) Coefficients of derivatives
    CALL tSPCD1(K,M,NDER,T,C,D,NERR)
    CALL tSPCD2(KX,KY,MX,MY,NDERX,NDERY,TX,TY,C,NDIMC,D,NERR)

Case (1D):

F
Name of a user-upplied FUNCTION subprogram, declared EXTERNAL in the calling program. This subprogram must provide the value of the function y = f(x)

for variation diminishing spline approximation.

K
( INTEGER) Degree of B-splines (1 ≤K ≤25 for tSPVD1, 0 ≤K ≤25 otherwise).
M
( INTEGER) Number of knots (≥2*K+2 ).
I
( INTEGER) Index of B-splines (1 ≤I ≤M-K-1 ).
N
( INTEGER) Number of data points xl,yl

(≥K+1 for spline interpolation ( tSPIN1); ≥M-K-1 for spline approximation ( tSPAP1)).

NDER
( INTEGER) Selects one out of three cases: (i) integral, (ii) function value, or (iii) derivative.
(-1 ≤NDER ≤K for tSPNB1 and tSPPS1; 1 ≤NDER ≤K for tSPCD1).
=-1: Calculation of the integral of Bi(x)

( tSPNB1), or the integral of s(x) ( tSPPS1).

=0: Calculation of the function value Bi(x)

( tSPNB1), or the function value s(x) ( tSPPS1).

≥1: Calculation of the NDER-th derivative of Bi(x) ( tSPNB1), or the NDER-th derivative of s(x) ( tSPPS1).
X
(Type according to t) Independent variable x of polynomial spline s(x) or B-spline Bi(x) .
XI
(Type according to t) One-dimensional array of length ≥N . On entry, XI(L) must contain the l-th data point xl for spline interpolation ( tSPIN1) or spline approximation ( tSPAP1). The data points must be in ascending order.
YI
(Type according to t) One-dimensional array of length ≥N . On entry, YI(L) must contain the l-th data point yl for spline interpolation ( tSPIN1) or spline approximation ( tSPAP1).
KNOT
( INTEGER) Controls the mode of supplying the knots for spline interpolation or approximation.
=1,2: The knots are computed by the subprograms tSPIN1 and tSPAP1. At the left and right end-point of the interpolation (approximation) interval [x1,xn] arise multiple knots. The remaining knots are either equidistant (KNOT=1 ) or are computed by using the data points xl of interpolation (approximation) (KNOT=2) .
≠1,2: The knots must be supplied by the user.
A,B
(Type according to t) On entry, A and B must contain the left (right) end-point of the interval [a,b] for the calculation of a set of spline knots ( tSPKN1).
T
(Type according to t) One-dimensional array of length ≥M .
For tSPKN1 and for tSPINT1, tSPAP1 if KNOT = 1,2:

On exit, T(J) contains the j-th knot. In the other cases, on entry, T(J) must contain the j-th knot. The knots must be in non-decreasing order with multiplicity ≤K+1 .

C
(Type according to t) One-dimensional array of length ≥M-K-1 .
For tSPPS1 and tSPCD1: On entry, C(I) must contain the i-th coefficient ci of the polynomial spline s(x) in B-spline representation.
For tSPIN1, tSPAP1 and tSPVD1: On exit, C(I) contains the i-th coefficient ci of the polynomial interpolation or approximation spline s(x) in B-spline representation.
D
(Type according to t) One-dimensional array of length ≥M-K-1 .
On exit, D(I) contains the coefficient di of the NDER-th derivative of a polynomial spline s(x) in B-spline representation.
W
(Type according to t) One-dimensional array of length ≥(3*K+1)*N ( tSPIN1), and of length ≥NW ( tSPAP1); used as working space.
NW
( INTEGER) Length of working array W, (NW ≥N*(n0+5)+n0*(n0+1); n0 = M-K-1 ).
IW
( INTEGER) One-dimensional array of length ≥N , used as working space.
NERR
( INTEGER) Error indicator. On exit:
=0: No error or warning detected.
=1:

At least one of the parameters I, K, M, N, NDER is not in range or A ≤B is not true.

=2: The subprograms tGEQPF, tORMQR, tTRTRS in the Linear Algebra package LAPACK (F001) were unable to solve the linear system of equations for calculating the coefficients of the spline interpolation to a given data set.
Case (2D):
F
Name of a user-upplied FUNCTION subprogram, declared EXTERNAL in the calling program. This subprogram must provide the value of the function z = f(x,y) for variation diminishing spline approximation.
KX,KY
( INTEGER) Degree of one-dimensional B-splines in x-(y-)direction (1 ≤KX ≤25, 1 ≤KY ≤25

for tSPVD2; 0 ≤KX ≤25, 0 ≤KY ≤25 otherwise).

MX,MY
( INTEGER) Number of knots in x-(y-)direction (MX ≥2*KX+2, MY ≥2*KY+2) .
I,J
( INTEGER) Indices of B-splines (1 ≤I ≤MX-KX-1, 1 ≤J ≤MY-KY-1 ).
NX,NY
( INTEGER) Number of data points xlx (yly)

in x-(y-)direction (NX ≥KX+1, NY ≥KY+1

for spline interpolation tSPIN2; NX ≥MX-KX-1, NY ≥MY-KY-1

for spline approximation tSPAP2).

NDERX,
( INTEGER) Selects one out of three cases: (i) integral, (ii) function value, or (iii) derivative.
NDERY
(-1 ≤NDERX ≤KX, -1 ≤NDERY ≤KY

for tSPNB2 and tSPPS2;
1 ≤NDERX ≤KX , 1 ≤NDERY ≤KY for tSPCD2).

=-1:

Calculation of the integral of Bi,j(x,y)

( tSPNB2), or the integral of s(x,y) ( tSPPS2).

=0:

Calculation of the function value Bi.j(x,y)

( tSPNB2), or the function value s(x,y) ( tSPPS2).

≥1: Calculation of the NDERX-th partial derivative with respect to x and the NDERY-th partial derivative with respect to y of Bi,j(x,y) ( tSPNB2), or the calcultion of these derivatives of s(x,y) ( tSPPS2).
Note that in the first two cases NDERX=NDERY=-1, NDERX=NDERY=0 , respectively.
X,Y
(Type according to t) Independent variables x,y of polynomial spline s(x,y) or B-spline Bi,j(x,y) .
XI,YI
(Type according to t) One-dimensional arrays of length ≥NX and ≥NY , respectively. On entry, XI(LX) and YI(LY) must contain the lx-th data point xlx and the ly-th data point yly

for spline interpolation ( tSPIN2) or spline approximation ( tSPAP2). The data points must be in ascending order.

ZI
(Type according to t) Two-dimensional array of dimension (NDIMZ, ≥NY) . On entry, ZI(LX,LY) must contain the (lx,ly) -th data point zlx,ly for spline interpolation ( tSPIN2) or spline approximation ( tSPAP2).
NDIMZ
( INTEGER) Declared first dimension of a two-dimensional array ZI in the calling program (≥NX ).
KNOT
( INTEGER) Controls the mode of supplying the knots for spline interpolation or approximation.
=1,2: The set of knots are computed by subprograms tSPIN2 and tSPAP2. At the left and right end-points of the interpolation (approximation) intervals [x1,xnx], [y1,yny] arise multiple knots. The remaining knots are either equidistant (KNOT=1 ) or are computed by using the data points xlx,yly of interpolation (approximation) (KNOT=2 ).
≠1,2:

The knots must be supplied by the user.

AX,BX;
(Type according to t) On entry, AX, BX; AY, BY must contain the left (right) end-points of the
AY,BY
intervals [ax,bx]; [ay,by] for the calculation of a set of spline knots in x-(y-)direction, respectively, by tSPKN2.
TX,TY
(Type according to t) One-dimensional arrays of length ≥MX and ≥MY , repectively.
For tSPKN2 and for tSPIN2, tSPAP2 if KNOT = 1,2:

On exit, TX(J) and TY(J) contain the j-th knot in x-(y-)direction. In the other cases, on entry, TX(J) and TY(J) must contain the j-th knot in x-(y-)direction. The knots must be in non-decreasing order with multiplicity ≤KX+1 and ≤KY+1 , respectively.

C
(Type according to t) Two-dimensional array of dimension (NDIMC , ≥MY-KY-1) .
For tSPPS2, tSPCD2: On entry, C(I,J) must contain the (i,j) -th coefficient ci,j of the polynomial spline s(x,y)

in B-spline representation.
For tSPIN2, tSPAP2, tSPVD2: On exit, C(I,J) contains the (i,j) -th coefficient ci,j of the polynomial interpolation or approximation spline s(x,y) in B-spline representation.

NDIMC
( INTEGER) Declared first dimension of a two-dimensional array C in the calling program
(≥MX-KX-1 ).
D
(Type according to t) Two-dimensional array of dimension (NDIMC , ≥MY-KY-1) .
On exit, D(I,J) contains the coefficient di,j of the partial derivative of order nx, ny with respect to x and y of a polynomial spline s(x,y) in B-spline representation.
W
(Type according to t) One-dimensional array of length ≥MY-KY-1 ( tSPPS2),
≥(3*KX*NY+2)*NX*NY ( tSPIN2), and of length ≥NW ( tSPAP2), used as working space.
NW
( INTEGER) Length of a one-dimensional array W, used as working space
(≥NX*NY*(n0+6)+n0*(n0+1); n0 = (MX-KX-1)*(MY-KY-1) ).
IW
( INTEGER) One-dimensional array of length ≥NX*NY , used as working space.
NERR
( INTEGER) Error indicator. On exit:
= 0: No error or warning detected.
= 1: At least one of the parameters I, J, KX, KY, MX, MY, NX, NY, NDERX, NDERY is not in range or at least one of the relations AX ≤BX , AY ≤BY is not true.
= 2: The routines tGEQPF, tORMQR, tTRTRS in the Linear Algebra package LAPACK (F001) were unable to solve the linear system of equations for calculation coefficients of spline interpolation to a given data set.

Examples:

Calculate

  1. The coefficients ci of a polynomial interpolation spline y=s(x) of degree k=2 in B-spline representation to a given data set (xl,yl)l=1,...,6 ;
  2. The corresponding coefficients di of the first derivative y'={ds(x)dx} ;
  3. The values of s(x), &quad;{ds(x)dx}

    and 0xs(ξ)dξ for x=0(0.1)1 .

      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
      DIMENSION XI(6),YI(6),T(9),C(6),D(6),W(42),IW(6)
      DATA K,N,NDER,KNOT / 2,6,1,1 /
      DATA XI / 0D0,0.20D0,0.40D0,0.60D0,0.80D0,1.00D0 /
      DATA YI / 1D0,0.66D0,0.47D0,0.38D0,0.35D0,0.34D0 /
      M=N+K+1
 
      CALL DSPIN1(K,N,XI,YI,KNOT,T,C,W,IW,NERR)
      CALL DSPCD1(K,M,NDER,T,C,D,NERR)
 
      WRITE(6,1000) K,N,(T(I),I=1,M)
      DO 20 J=0,10
      X=J/1D1
 
      SPL0=DSPPS1(K,M, 0,X,T,C,NERR)
      SPL1=DSPPS1(K,M, 1,X,T,D,NERR)
      SPLI=DSPPS1(K,M,-1,X,T,C,NERR)
 
   20 WRITE(6,1010) J,X,SPL0,SPL1,SPLI
 1000 FORMAT(...)
 1010 FORMAT(...)
      END
 DEGREE OF POLYNOMIAL SPLINE:  2          NUMBER OF DATA POINTS:  6
 KNOTS T :        0.00  0.00  0.00  0.25  0.50  0.75  1.00  1.00  1.00
 
     J        X         S(X)          DS(X)         IN(X)
     0      0.00      1.000000     -2.119921      0.000000
     1      0.10      0.809004     -1.700000      0.090100
     2      0.20      0.660000     -1.280079      0.163201
     3      0.30      0.550992     -0.940017      0.223467
     4      0.40      0.470000     -0.679816      0.274299
     5      0.50      0.415028     -0.419615      0.318334
     6      0.60      0.380000     -0.280953      0.357970
     7      0.70      0.358838     -0.142290      0.394796
     8      0.80      0.350000     -0.065306      0.430174
     9      0.90      0.344235     -0.050000      0.464873
    10      1.00      0.340000     -0.034694      0.499072

Error handling:

Error E210.1: K|KX,KYnot in range, Error E210.5: NDER|NDERX,NDERYnot in range,
Error E210.2: M|MX,MYnot in range, Error E210.6: A,B|AX,BX;AY,BYinconsistent,
Error E210.3: I|I,Jnot in range, Error E210.7: NDERX|NDERYinconsistent.
Error E210.4: N|NX,NYnot in range,

In all cases, NERR is set =1 (see above). A message is written on Unit 6, unless subroutine MTLSET (N002) has been called.

References:

  1. J.H. Ahlberg, E.N. Nilson, J.L. Walsh, The Theory of Splines and their Applications, Academic Press, New York, 1967.
  2. M.J. Marsden, An identity for spline functions with applications to variation diminishing spline approximation, J. Appr. Theory 3 (1970), 7-49.
  3. C. de Boor, On calculating with B-splines, J. Appr. Theory 6 (1972), 50-62.
  4. M.G. Cox, The numerical evaluation of B-splines, J. Inst. Maths Applics 10 (1972), 134-149.
  5. J.G. Hayes, J. Halliday, The least-squares fitting of cubic spline surfaces to general data sets, J. Inst. Maths Applics 14 (1974), 89-103.
  6. M.G. Cox, An algorithm for spline interpolation, J. Inst. Maths Applics 15 (1975), 95-108.
  7. C. de Boor, A Practical Guide to Splines, Springer-Verlag, Berlin (1978).
  8. P. Lancester, K. Salkauskas, Curve and Surface Fitting - An Introduction, Academic Press, New York, 1986.
  9. J.C. Mason, M.S. Cox (Eds.), Algorithms for Approximation, Clarendon Press, Oxford, 1987.
  10. J.W. Schmidt, H. Späth (Eds.), Splines in Numerical Analysis, Akademie-Verlag, Berlin, 1989.
  11. H. Späth, Eindimensionale Spline-Interpolations-Algorithmen; Zweidimensionale Spline-Interpolations-Algorithmen, (2 Bd.) R. Oldenbourg, München 1990/1991.

E211



next up previous index
Next: E211 Cubic Splines Up: CERNLIB Previous: E208 Least Squares


Janne Saarela
Mon Apr 3 15:06:23 METDST 1995