Evaluates the derivative of a three-dimensional
tensor-product spline, given its tensor-product
B-spline representation on a
grid.
IXDER — Order of the X-derivative. (Input)
IYDER — Order of the Y-derivative. (Input)
IZDER — Order of the Z-derivative. (Input)
XVEC — Array of
length NX
containing the x-coordinates at which the spline is to be
evaluated. (Input)
The points in XVEC should be
strictly increasing.
YVEC — Array of
length NY
containing the y-coordinates at which the spline is to be
evaluated. (Input)
The points in YVEC should be
strictly increasing.
ZVEC — Array of
length NY
containing the y-coordinates at which the spline is to be
evaluated. (Input)
The points in YVEC should be
strictly increasing.
KXORD — Order of the spline in the x-direction. (Input)
KYORD — Order of the spline in the y-direction. (Input)
KZORD — Order of the spline in the z-direction. (Input)
XKNOT — Array of
length NXCOEF +
KXORD containing
the knot sequence in the x-direction. (Input)
XKNOT must be
nondecreasing.
YKNOT — Array of
length NYCOEF +
KYORD containing
the knot sequence in the y-direction. (Input)
YKNOT must be
nondecreasing.
ZKNOT — Array of
length NZCOEF +
KZORD containing
the knot sequence in the z-direction. (Input)
ZKNOT must be
nondecreasing.
BSCOEF — Array of
length NXCOEF
* NYCOEF * NZCOEF containing the
tensor-product
B-spline coefficients. (Input)
BSCOEF is treated
internally as a matrix of size NXCOEF by NYCOEF by NZCOEF.
VALUE — Array of
size NX by NY by NZ containing the
values of the (IXDER, IYDER, IZDER) derivative of
the spline on the NX by NY by NZ grid.
(Output)
VALUE(I, J, K) contains the
derivative of the spline at the point (XVEC(I), YVEC(J), ZVEC(K)).
NX — Number of
grid points in the x-direction. (Input)
Default: NX = size (XVEC,1).
NY — Number of
grid points in the y-direction. (Input)
Default: NY = size (YVEC,1).
NZ — Number of
grid points in the z-direction. (Input)
Default: NZ = size (ZVEC,1).
NXCOEF — Number
of B-spline coefficients in the x-direction.
(Input)
Default: NXCOEF = size (XKNOT,1) – KXORD.
NYCOEF — Number
of B-spline coefficients in the y-direction.
(Input)
Default: NYCOEF = size (YKNOT,1) – KYORD.
NZCOEF — Number
of B-spline coefficients in the z-direction.
(Input)
Default: NZCOEF = size (ZKNOT,1) – KZORD.
LDVALU — Leading
dimension of VALUE exactly as
specified in the dimension statement of the calling program.
(Input)
Default: LDVALU = size
(VALUE,1).
MDVALU — Middle
dimension of VALUE exactly as
specified in the dimension statement of the calling program.
(Input)
Default: MDVALU = size
(VALUE,2).
Generic: CALL BS3GD (IXDER, IYDER, IZDER, XVEC, YVEC, ZVEC, KXORD, KYORD, KZORD, XKNOT, YKNOT, ZKNOT, BSCOEF, VALUE [,…])
Specific: The specific interface names are S_BS3GD and D_BS3GD.
Single: CALL BS3GD (IXDER, IYDER, IZDER, NX, XVEC, NY, YVEC, NZ, ZVEC, KXORD, KYORD, KZORD, XKNOT, YKNOT, ZKNOT, NXCOEF, NYCOEF, NZCOEF, BSCOEF, VALUE, LDVALU, MDVALU)
Double: The double precision name is DBS3GD.
The routine BS3GD evaluates a partial derivative of a trivariate tensor-product spline (represented as a linear combination of tensor-product B-splines) on a grid. For more information, see de Boor (1978, pages 351− 353).
This routine returns the value of the function
s(p,q,r) on the grid
(xi, yj, zk) for
i = 1, …, nx,
j =
1, …, ny, and
k = 1, …,
nz given the coefficients c by computing (for all (x,
y, z) on the grid)
where kx, ky, and kz are the orders of the
splines. (These numbers are passed to the subroutine in KXORD,
KYORD,
and KZORD,
respectively.) Likewise, tx, ty, and tz are the corresponding
knot sequences (XKNOT,
YKNOT,
and ZKNOT).
The grid must be ordered in the sense that
xi < xi +
1, yj < yj +
1, and zk < zk +
1.
1. Workspace may be explicitly provided, if desired, by use of B23GD/DB23GD. The reference is:
CALL B23GD ((IXDER, IYDER, IZDER, NX, XVEC, NY, YVEC, NZ, ZVEC, KXORD, KYORD, KZORD, XKNOT, YKNOT, ZKNOT, NXCOEF, NYCOEF, NZCOEF, BSCOEF, VALUE, LDVALU, MDVALU, LEFTX, LEFTY, LEFTZ, A, B, C, DBIATX, DBIATY, DBIATZ, BX, BY, BZ)
The additional arguments are as follows:
LEFTX — Work array of length NX.
LEFTY — Work array of length NY.
LEFTZ — Work array of length NZ.
A — Work array of length KXORD * KXORD.
B — Work array of length KYORD * KYORD.
C — Work array of length KZORD * KZORD.
DBIATX — Work array of length KXORD * (IXDER + 1).
DBIATY — Work array of length KYORD * (IYDER + 1).
DBIATZ — Work array of length KZORD * (IZDER + 1).
BX — Work array of length KXORD * NX.
BY — Work array of length KYORD * NY.
BZ — Work array of length KZORD * NZ.
2. Informational errors
Type Code
3
1
XVEC(I) does not satisfy
XKNOT(KXORD) ≤ XVEC(I) ≤ XKNOT(NXCOEF + 1).
3
2
YVEC(I) does not satisfy
YKNOT(KYORD) ≤ YVEC(I) ≤ YKNOT(NYCOEF + 1).
3
3
ZVEC(I) does not satisfy
ZKNOT(KZORD) ≤ ZVEC(I) ≤ ZKNOT(NZCOEF + 1).
4 4 XVEC is not strictly increasing.
4 5 YVEC is not strictly increasing.
4 6 ZVEC is not strictly increasing.
In this example, a spline interpolant s to a function f(x, y, z) = x4 + y(xz)3 is constructed using BS3IN. Next, BS3GD is used to compute s(2,0,1)(x, y, z) on the grid. The values of this partial derivative and the error are computed on a 4 × 4 × 2 grid and then displayed.
USE BS3GD_INT
USE BS3IN_INT
USE BSNAK_INT
USE UMACH_INT
IMPLICIT NONE
INTEGER KXORD, KYORD, KZORD, LDF, LDVAL, MDF, MDVAL, NXDATA,&
NXKNOT, NYDATA, NYKNOT, NZ, NZDATA, NZKNOT
PARAMETER (KXORD=5, KYORD=2, KZORD=3, LDVAL=4, MDVAL=4,&
NXDATA=21, NYDATA=6, NZ=2, NZDATA=8, LDF=NXDATA,&
MDF=NYDATA, NXKNOT=NXDATA+KXORD, NYKNOT=NYDATA+KYORD,&
NZKNOT=NZDATA+KZORD)
!
INTEGER I, J, K, L, NOUT, NXCOEF, NYCOEF, NZCOEF
REAL BSCOEF(NXDATA,NYDATA,NZDATA), F, F201,&
FDATA(LDF,MDF,NZDATA), FLOAT, VALUE(LDVAL,MDVAL,NZ),&
X, XDATA(NXDATA), XKNOT(NXKNOT), XVEC(LDVAL), Y,&
YDATA(NYDATA), YKNOT(NYKNOT), YVEC(MDVAL), Z,&
ZDATA(NZDATA), ZKNOT(NZKNOT), ZVEC(NZ)
INTRINSIC FLOAT
!
!
!
F(X,Y,Z) = X*X*X*X + X*X*X*Y*Z*Z*Z
F201(X,Y,Z) = 18.0*X*Y*Z
!
CALL UMACH (2, NOUT)
! Set up X interpolation points
DO 10 I=1, NXDATA
XDATA(I) = 2.0*(FLOAT(I-1)/FLOAT(NXDATA-1)) - 1.0
10 CONTINUE
! Set up Y interpolation points
DO 20 I=1, NYDATA
YDATA(I) = FLOAT(I-1)/FLOAT(NYDATA-1)
20 CONTINUE
! Set up Z interpolation points
DO 30 I=1, NZDATA
ZDATA(I) = FLOAT(I-1)/FLOAT(NZDATA-1)
30 CONTINUE
! Generate knots
CALL BSNAK (NXDATA, XDATA, KXORD, XKNOT)
CALL BSNAK (NYDATA, YDATA, KYORD, YKNOT)
CALL BSNAK (NZDATA, ZDATA, KZORD, ZKNOT)
! Generate FDATA
DO 50 K=1, NZDATA
DO 40 I=1, NYDATA
DO 40 J=1, NXDATA
FDATA(J,I,K) = F(XDATA(J),YDATA(I),ZDATA(K))
40 CONTINUE
50 CONTINUE
! Interpolate
CALL BS3IN (XDATA, YDATA, ZDATA, FDATA, KXORD, KYORD,&
KZORD, XKNOT, YKNOT, ZKNOT, BSCOEF)
!
NXCOEF = NXDATA
NYCOEF = NYDATA
NZCOEF = NZDATA
! Print over a grid of
! [-1.0,1.0] x [0.0,1.0] x [0.0,1.0]
! at 32 points.
DO 60 I=1, 4
XVEC(I) = 2.0*(FLOAT(I-1)/3.0) - 1.0
60 CONTINUE
DO 70 J=1, 4
YVEC(J) = FLOAT(J-1)/3.0
70 CONTINUE
DO 80 L=1, 2
ZVEC(L) = FLOAT(L-1)
80 CONTINUE
CALL BS3GD (2, 0, 1, XVEC, YVEC, ZVEC, KXORD, KYORD,&
KZORD, XKNOT, YKNOT, ZKNOT, BSCOEF, VALUE)
!
!
WRITE (NOUT,99999)
DO 110 I=1, 4
DO 100 J=1, 4
DO 90 L=1, 2
WRITE (NOUT,'(5F13.4)') XVEC(I), YVEC(J), ZVEC(L),&
VALUE(I,J,L),&
F201(XVEC(I),YVEC(J),ZVEC(L)) -&
VALUE(I,J,L)
90 CONTINUE
100 CONTINUE
110 CONTINUE
99999 FORMAT (44X, '(2,0,1)', /, 10X, 'X', 11X, 'Y', 10X, 'Z', 10X,&
'S (X,Y,Z) Error')
STOP
END
(2,0,1)
X
Y
Z
S (X,Y,Z) Error
-1.0000
0.0000 0.0000
-0.0005
0.0005
-1.0000
0.0000
1.0000 0.0002
-0.0002
-1.0000
0.3333
0.0000 0.0641
-0.0641
-1.0000
0.3333 1.0000
-5.9360 -0.0640
-1.0000
0.6667
0.0000 0.1274
-0.1274
-1.0000
0.6667 1.0000
-11.8730 -0.1270
-1.0000
1.0000
0.0000 0.1911
-0.1911
-1.0000
1.0000 1.0000
-17.8086 -0.1914
-0.3333
0.0000
0.0000
0.0000
0.0000
-0.3333
0.0000
1.0000
0.0000
0.0000
-0.3333
0.3333
0.0000 0.0212
-0.0212
-0.3333
0.3333 1.0000
-1.9788 -0.0212
-0.3333
0.6667
0.0000 0.0425
-0.0425
-0.3333
0.6667 1.0000
-3.9575 -0.0425
-0.3333
1.0000
0.0000 0.0637
-0.0637
-0.3333
1.0000 1.0000
-5.9363
-0.0637
0.3333
0.0000
0.0000
0.0000
0.0000
0.3333
0.0000
1.0000
0.0000
0.0000
0.3333
0.3333 0.0000
-0.0212
0.0212
0.3333
0.3333
1.0000
1.9788
0.0212
0.3333
0.6667 0.0000
-0.0425
0.0425
0.3333
0.6667
1.0000
3.9575
0.0425
0.3333
1.0000 0.0000
-0.0637
0.0637
0.3333
1.0000
1.0000
5.9363
0.0637
1.0000
0.0000 0.0000
-0.0005
0.0005
1.0000
0.0000
1.0000
0.0000
0.0000
1.0000
0.3333 0.0000
-0.0637
0.0637
1.0000
0.3333
1.0000
5.9359
0.0641
1.0000
0.6667 0.0000
-0.1273
0.1273
1.0000
0.6667 1.0000
11.8733
0.1267
1.0000
1.0000 0.0000
-0.1912
0.1912
1.0000
1.0000 1.0000
17.8096 0.1904
Visual Numerics, Inc. PHONE: 713.784.3131 FAX:713.781.9260 |