PP1GD
Evaluates the derivative of a piecewise polynomial on a grid.
Required Arguments
XVEC — Array of length N containing the points at which the piecewise polynomial is to be evaluated. (Input)
The points in XVEC should be strictly increasing.
BREAK — Array of length NINTV + 1 containing the breakpoints for the piecewise polynomial representation. (Input)
BREAK must be strictly increasing.
PPCOEF — Matrix of size KORDER by NINTV containing the local coefficients of the polynomial pieces. (Input)
VALUE — Array of length N containing the values of the IDERIV-th derivative of the piecewise polynomial at the points in XVEC. (Output)
Optional Arguments
IDERIV — Order of the derivative to be evaluated. (Input)
In particular, IDERIV = 0 returns the values of the piecewise polynomial.
Default: IDERIV = 1.
N — Length of vector XVEC. (Input)
Default: N = size (XVEC,1).
KORDER — Order of the polynomial. (Input)
Default: KORDER = size (PPCOEF,1).
NINTV — Number of polynomial pieces. (Input)
Default: NINTV = size (PPCOEF,2).
FORTRAN 90 Interface
Generic: CALL PP1GD (XVEC, BREAK, PPCOEF, VALUE [, …])
Specific: The specific interface names are S_PP1GD and D_PP1GD.
FORTRAN 77 Interface
Single: CALL PP1GD (IDERIV, N, XVEC, KORDER, NINTV, BREAK, PPCOEF, VALUE)
Double: The double precision name is DPP1GD.
Description
The routine PP1GD evaluates a piecewise polynomial function f (or its derivative) at a vector of points. That is, given a vector x of length n satisfying xi < xi + 1 for i = 1, …, n − 1, a derivative value j, and a piecewise polynomial function f that is represented by a breakpoint sequence and coefficient matrix this routine returns the values
f(j) (xi) i = 1,…, n
in the array
VALUE. The functionality of this routine is the same as that of
PPDER called in a loop, however
PP1GD is much more efficient.
Comments
1. Workspace may be explicitly provided, if desired, by use of P21GD/DP21GD. The reference is:
CALL P21GD (IDERIV, N, XVEC, KORDER, NINTV, BREAK, PPCOEF, VALUE, IWK, WORK1, WORK2)
The additional arguments are as follows:
IWK — Array of length N.
WORK1 — Array of length N.
WORK2 — Array of length N.
2. Informational error
Type | Code | Description |
|---|
4 | 4 | The points in XVEC must be strictly increasing. |
Example
To illustrate the use of
PP1GD, we modify the example program for
PPDER. In this example, a piecewise polynomial interpolant to
F is computed. The values of this polynomial are then compared with the exact function values. The routine
PP1GD is based on the routine
PPVALU in de Boor (1978, page 89).
USE IMSL_LIBRARIES
IMPLICIT NONE
INTEGER KORDER, N, NCOEF, NDATA, NKNOT
PARAMETER (KORDER=4, N=20, NCOEF=20, NDATA=20,&
NKNOT=NDATA+KORDER)
!
INTEGER I, NINTV, NOUT, NPPCF
REAL BREAK(NCOEF), BSCOEF(NCOEF), DF, EXP, F,&
FDATA(NDATA), FLOAT, PPCOEF(KORDER,NCOEF), VALUE1(N),&
VALUE2(N), X, XDATA(NDATA), XKNOT(NKNOT), XVEC(N)
INTRINSIC EXP, FLOAT
!
F(X) = X*EXP(X)
DF(X) = (X+1.)*EXP(X)
! Set up interpolation points
DO 10 I=1, NDATA
XDATA(I) = FLOAT(I-1)/FLOAT(NDATA-1)
FDATA(I) = F(XDATA(I))
10 CONTINUE
! Generate knot sequence
CALL BSNAK (NDATA, XDATA, KORDER, XKNOT)
! Compute the B-spline interpolant
CALL BSINT (NCOEF, XDATA, FDATA, KORDER, XKNOT, BSCOEF)
! Convert to piecewise polynomial
CALL BSCPP (KORDER, XKNOT, NCOEF, BSCOEF, NPPCF, BREAK, PPCOEF)
! Compute evaluation points
DO 20 I=1, N
XVEC(I) = FLOAT(I-1)/FLOAT(N-1)
20 CONTINUE
! Compute values of the piecewise
! polynomial
NINTV = NPPCF
CALL PP1GD (XVEC, BREAK, PPCOEF, VALUE1, IDERIV=0, NINTV=NINTV)
! Compute the values of the first
! derivative of the piecewise
! polynomial
CALL PP1GD (XVEC, BREAK, PPCOEF, VALUE2, IDERIV=1, NINTV=NINTV)
! Get output unit number
CALL UMACH (2, NOUT)
! Write heading
WRITE (NOUT,99998)
! Print the results on a uniform
! grid
DO 30 I=1, N
WRITE (NOUT,99999) XVEC(I), VALUE1(I), F(XVEC(I)) - VALUE1(I)&
, VALUE2(I), DF(XVEC(I)) - VALUE2(I)
30 CONTINUE
99998 FORMAT (11X, 'X', 8X, 'S(X)', 7X, 'Error', 7X, 'S''(X)', 7X,&
'Error')
99999 FORMAT (' ', 2F12.3, F12.6, F12.3, F12.6)
END
Output
X S(X) Error S’(X) Error
0.000 0.000 0.000000 1.000 -0.000112
0.053 0.055 0.000000 1.109 0.000030
0.105 0.117 0.000000 1.228 -0.000008
0.158 0.185 0.000000 1.356 0.000002
0.211 0.260 0.000000 1.494 0.000000
0.263 0.342 0.000000 1.643 0.000000
0.316 0.433 0.000000 1.804 -0.000001
0.368 0.533 0.000000 1.978 0.000002
0.421 0.642 0.000000 2.165 0.000001
0.474 0.761 0.000000 2.367 0.000000
0.526 0.891 0.000000 2.584 -0.000001
0.579 1.033 0.000000 2.817 0.000001
0.632 1.188 0.000000 3.068 0.000001
0.684 1.356 0.000000 3.338 0.000001
0.737 1.540 0.000000 3.629 0.000001
0.789 1.739 0.000000 3.941 0.000000
0.842 1.955 0.000000 4.276 -0.000006
0.895 2.189 0.000000 4.636 0.000024
0.947 2.443 0.000000 5.022 -0.000090
1.000 2.718 0.000000 5.436 0.000341
