Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

Chapter Contents
Chapter Introduction
NAG Toolbox

# NAG Toolbox: nag_fit_2dspline_derivm (e02dh)

## Purpose

nag_fit_2dspline_derivm (e02dh) computes the partial derivative (of order νx${\nu }_{x}$, νy${\nu }_{y}$), of a bicubic spline approximation to a set of data values, from its B-spline representation, at points on a rectangular grid in the x$x$-y$y$ plane. This function may be used to calculate derivatives of a bicubic spline given in the form produced by nag_interp_2d_spline_grid (e01da), nag_fit_2dspline_panel (e02da), nag_fit_2dspline_grid (e02dc) and nag_fit_2dspline_sctr (e02dd).

## Syntax

[z, ifail] = e02dh(x, y, lamda, mu, c, nux, nuy, 'mx', mx, 'my', my, 'px', px, 'py', py)
[z, ifail] = nag_fit_2dspline_derivm(x, y, lamda, mu, c, nux, nuy, 'mx', mx, 'my', my, 'px', px, 'py', py)

## Description

nag_fit_2dspline_derivm (e02dh) determines the partial derivative ( νx + νy )/( xνx yνy ) $\frac{{\partial }^{{\nu }_{x}+{\nu }_{y}}}{\partial {x}^{{\nu }_{x}}\partial {y}^{{\nu }_{y}}}$ of a smooth bicubic spline approximation s(x,y)$s\left(x,y\right)$ at the set of data points (xq,yr)$\left({x}_{q},{y}_{r}\right)$.
The spline is given in the B-spline representation
 nx − 4 ny − 4 s(x,y) = ∑ ∑ cijMi(x)Nj(y), i = 1 j = 1
$s(x,y) = ∑ i=1 nx-4 ∑ j=1 ny-4 cij Mi(x) Nj(y) ,$
(1)
where Mi(x)${M}_{i}\left(x\right)$ and Nj(y)${N}_{j}\left(y\right)$ denote normalized cubic B-splines, the former defined on the knots λi${\lambda }_{i}$ to λi + 4${\lambda }_{i+4}$ and the latter on the knots μj${\mu }_{j}$ to μj + 4${\mu }_{j+4}$, with nx${n}_{x}$ and ny${n}_{y}$ the total numbers of knots of the computed spline with respect to the x$x$ and y$y$ variables respectively. For further details, see Hayes and Halliday (1974) for bicubic splines and de Boor (1972) for normalized B-splines. This function is suitable for B-spline representations returned by nag_interp_2d_spline_grid (e01da), nag_fit_2dspline_panel (e02da), nag_fit_2dspline_grid (e02dc) and nag_fit_2dspline_sctr (e02dd).
The partial derivatives can be up to order 2$2$ in each direction; thus the highest mixed derivative available is (4)/( x2 y2 ) $\frac{{\partial }^{4}}{\partial {x}^{2}\partial {y}^{2}}$.
The points in the grid are defined by coordinates xq${x}_{\mathit{q}}$, for q = 1,2,,mx$\mathit{q}=1,2,\dots ,{m}_{x}$, along the x$x$ axis, and coordinates yr${y}_{\mathit{r}}$, for r = 1,2,,my$\mathit{r}=1,2,\dots ,{m}_{y}$, along the y$y$ axis.

## References

de Boor C (1972) On calculating with B-splines J. Approx. Theory 6 50–62
Dierckx P (1981) An improved algorithm for curve fitting with spline functions Report TW54 Department of Computer Science, Katholieke Univerciteit Leuven
Dierckx P (1982) A fast algorithm for smoothing data on a rectangular grid while using spline functions SIAM J. Numer. Anal. 19 1286–1304
Hayes J G and Halliday J (1974) The least squares fitting of cubic spline surfaces to general data sets J. Inst. Math. Appl. 14 89–103
Reinsch C H (1967) Smoothing by spline functions Numer. Math. 10 177–183

## Parameters

### Compulsory Input Parameters

1:     x(mx) – double array
mx, the dimension of the array, must satisfy the constraint mx1${\mathbf{mx}}\ge 1$.
x(q)${\mathbf{x}}\left(q\right)$ must be set to xq${x}_{\mathit{q}}$, the x$x$ coordinate of the q$\mathit{q}$th grid point along the x$x$ axis, for q = 1,2,,mx$\mathit{q}=1,2,\dots ,{m}_{x}$, on which values of the partial derivative are sought.
Constraint: x1 < x2 < < xmx${x}_{1}<{x}_{2}<\cdots <{x}_{{m}_{x}}$.
2:     y(my) – double array
my, the dimension of the array, must satisfy the constraint my1${\mathbf{my}}\ge 1$.
y(r)${\mathbf{y}}\left(\mathit{r}\right)$ must be set to yr${y}_{\mathit{r}}$, the y$y$ coordinate of the r$\mathit{r}$th grid point along the y$y$ axis, for r = 1,2,,my$\mathit{r}=1,2,\dots ,{m}_{y}$ on which values of the partial derivative are sought.
Constraint: y1 < y2 < < ymy${y}_{1}<{y}_{2}<\cdots <{y}_{{m}_{y}}$.
3:     lamda(px) – double array
Contains the position of the knots in the x$x$-direction of the bicubic spline approximation to be differentiated, e.g., lamda as returned by nag_fit_2dspline_grid (e02dc).
4:     mu(py) – double array
Contains the position of the knots in the y$y$-direction of the bicubic spline approximation to be differentiated, e.g., mu as returned by nag_fit_2dspline_grid (e02dc).
5:     c((px4) × (py4)$\left({\mathbf{px}}-4\right)×\left({\mathbf{py}}-4\right)$) – double array
The coefficients of the bicubic spline approximation to be differentiated, e.g., c as returned by nag_fit_2dspline_grid (e02dc).
6:     nux – int64int32nag_int scalar
Specifies the order, νx${\nu }_{x}$ of the partial derivative in the x$x$-direction.
Constraint: 0nux2$0\le {\mathbf{nux}}\le 2$.
7:     nuy – int64int32nag_int scalar
Specifies the order, νy${\nu }_{y}$ of the partial derivative in the y$y$-direction.
Constraint: 0nuy2$0\le {\mathbf{nuy}}\le 2$.

### Optional Input Parameters

1:     mx – int64int32nag_int scalar
Default: The dimension of the array x.
mx${m}_{x}$, the number of grid points along the x$x$ axis.
Constraint: mx1${\mathbf{mx}}\ge 1$.
2:     my – int64int32nag_int scalar
Default: The dimension of the array y.
my${m}_{y}$, the number of grid points along the y$y$ axis.
Constraint: my1${\mathbf{my}}\ge 1$.
3:     px – int64int32nag_int scalar
Default: The dimension of the array lamda.
The total number of knots in the x$x$-direction of the bicubic spline approximation, e.g., the value nx as returned by nag_fit_2dspline_grid (e02dc).
4:     py – int64int32nag_int scalar
Default: The dimension of the array mu.
The total number of knots in the y$y$-direction of the bicubic spline approximation, e.g., the value ny as returned by nag_fit_2dspline_grid (e02dc).

None.

### Output Parameters

1:     z(mx × my${\mathbf{mx}}×{\mathbf{my}}$) – double array
z(my × (q1) + r)${\mathbf{z}}\left({m}_{y}×\left(\mathit{q}-1\right)+\mathit{r}\right)$ contains the derivative (νx + νy)/( xνx yνy ) s (xq,yr) $\frac{{\partial }^{{\nu }_{x}+{\nu }_{y}}}{{\partial x}^{{\nu }_{x}}{\partial y}^{{\nu }_{y}}}s\left({x}_{q},{y}_{r}\right)$, for q = 1,2,,mx$\mathit{q}=1,2,\dots ,{m}_{x}$ and r = 1,2,,my$\mathit{r}=1,2,\dots ,{m}_{y}$.
2:     ifail – int64int32nag_int scalar
${\mathrm{ifail}}={\mathbf{0}}$ unless the function detects an error (see [Error Indicators and Warnings]).

## Error Indicators and Warnings

Errors or warnings detected by the function:
ifail = 1${\mathbf{ifail}}=1$
Constraint: 0nux2$0\le {\mathbf{nux}}\le 2$.
ifail = 2${\mathbf{ifail}}=2$
Constraint: 0nuy2$0\le {\mathbf{nuy}}\le 2$.
ifail = 3${\mathbf{ifail}}=3$
Constraint: mx1${\mathbf{mx}}\ge 1$.
ifail = 4${\mathbf{ifail}}=4$
Constraint: my1${\mathbf{my}}\ge 1$.
ifail = 5${\mathbf{ifail}}=5$
Constraint: x(i1)x(i)${\mathbf{x}}\left(\mathit{i}-1\right)\le {\mathbf{x}}\left(\mathit{i}\right)$, for i = 2,3,,mx$\mathit{i}=2,3,\dots ,{\mathbf{mx}}$.
ifail = 6${\mathbf{ifail}}=6$
Constraint: y(i1)y(i)${\mathbf{y}}\left(\mathit{i}-1\right)\le {\mathbf{y}}\left(\mathit{i}\right)$, for i = 2,3,,my$\mathit{i}=2,3,\dots ,{\mathbf{my}}$.
ifail = 999${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.

## Accuracy

On successful exit, the partial derivatives on the given mesh are accurate to machine precision with respect to the supplied bicubic spline. Please refer to Section [Accuracy] in (e01da), (e02da), (e02dc) and (e02dd) of the function document for the respective function which calculated the spline approximant for details on the accuracy of that approximation.

None.

## Example

```function nag_fit_2dspline_derivm_example
mx = 11;
my = 9;
start = 'C';
x = [0:0.5:5];
y = [0:0.5:4];
f = [ 1.0000e+00, 8.8758e-01, 5.4030e-01, 7.0737e-02,-4.1515e-01, ...
-8.0114e-01,-9.7999e-01,-9.3446e-01,-6.5664e-01, 1.5000e+00, ...
1.3564e+00, 8.2045e-01, 1.0611e-01,-6.2422e-01,-1.2317e+00, ...
-1.4850e+00,-1.3047e+00,-9.8547e-01, 2.0600e+00, 1.7552e+00, ...
1.0806e+00, 1.5147e-01,-8.3229e-01,-1.6023e+00,-1.9700e+00, ...
-1.8729e+00,-1.4073e+00, 2.5700e+00, 2.1240e+00, 1.3508e+00, ...
1.7684e-01,-1.0404e+00,-2.0029e+00,-2.4750e+00,-2.3511e+00, ...
-1.6741e+00, 3.0000e+00, 2.6427e+00, 1.6309e+00, 2.1221e-01, ...
-1.2484e+00,-2.2034e+00,-2.9700e+00,-2.8094e+00,-1.9809e+00, ...
3.5000e+00, 3.1715e+00, 1.8611e+00, 2.4458e-01,-1.4565e+00, ...
-2.8640e+00,-3.2650e+00,-3.2776e+00,-2.2878e+00, 4.0400e+00, ...
3.5103e+00, 2.0612e+00, 2.8595e-01,-1.6946e+00,-3.2046e+00, ...
-3.9600e+00,-3.7958e+00,-2.6146e+00, 4.5000e+00, 3.9391e+00, ...
2.4314e+00, 3.1632e-01,-1.8627e+00,-3.6351e+00,-4.4550e+00, ...
-4.2141e+00,-2.9314e+00, 5.0400e+00, 4.3879e+00, 2.7515e+00, ...
3.5369e-01,-2.0707e+00,-4.0057e+00,-4.9700e+00,-4.6823e+00, ...
-3.2382e+00, 5.5050e+00, 4.8367e+00, 2.9717e+00, 3.8505e-01, ...
-2.2888e+00,-4.4033e+00,-5.4450e+00,-5.1405e+00,-3.5950e+00, ...
6.0000e+00, 5.2755e+00, 3.2418e+00, 4.2442e-01,-2.4769e+00, ...
-4.8169e+00,-5.9300e+00,-5.6387e+00,-3.9319e+00];
s = 0.1;
ngx = 6;
xlo = 0;
xhi = 5;
ngy = 5;
ylo = 0;
yhi = 4;
nux = int64(1);
nuy = int64(0);

nxest = mx + 4;
nyest = my + 4;
lamda = zeros(nxest, 1);
mu    = zeros(nyest, 1);
wrk   = zeros(4*(mx+my) + 11*(nxest+nyest) + nxest*my + max(my,nxest) + 54, 1);
iwrk  = zeros(3 + mx + my + nxest + nyest, 1, 'int64');

% Determine the spline approximation
[nx, lamda, ny, mu, c, fp, wrk, iwrk, ifail] = ...
nag_fit_2dspline_grid(start, x, y, f, s, int64(0), lamda, int64(0), mu, wrk, iwrk);

fprintf('\nSpline fit used smoothing factor S = %13.4e.\n', s);
fprintf('Number of knots in each direction  = %d, %d\n', nx, ny);
fprintf('\nSum of squared residuals           = %13.4e.\n', fp);

% Evaluate the spline derivative on a different rectangular grid with
% ngx*ngy points over the domain (xlo to xhi) x (ylo to yhi).
gridx = [xlo:(xhi-xlo)/(ngx-1):xhi];
gridy = [ylo:(yhi-ylo)/(ngy-1):yhi];

% Evaluate spline (nux=nuy=0)
[z, ifail] = nag_fit_2dspline_derivm(gridx, gridy, lamda(1:double(nx)), mu(1:double(ny)), ...
c, int64(0), int64(0));

% Evaluate spline partial derivative of order (nux, nuy)
[zder, ifail] = nag_fit_2dspline_derivm(gridx, gridy, lamda(1:double(nx)), mu(1:double(ny)), ...
c, nux, nuy);

% Spline partial derivatives to evaluate have order nux,nuy.
fprintf('\nDerivative of spline has order nux, nuy = %d, %d.\n\n', nux, nuy);

% Print the spline evaluations
rlabs = cell(size(gridy));
for i = 1:ngy
rlabs{i} = num2str(gridy(i));
end
clabs = cell(size(gridx));
for i = 1:ngx
clabs{i} = num2str(gridx(i));
end
title = 'Spline evaluated on X-Y grid (X across, Y down):';
[ifail] = nag_file_print_matrix_real_gen_comp('G', 'N', reshape(z, ngy, ngx), 'F8.3', title, ...
'C', rlabs, 'C', clabs, int64(80), int64(0));
fprintf('\n');
title = 'Spline derivative evaluated on X-Y grid:';
[ifail] = nag_file_print_matrix_real_gen_comp('G', 'N', reshape(zder, ngy, ngx), 'F8.3', title, ...
'C', rlabs, 'C', clabs, int64(80), int64(0))
```
```

Spline fit used smoothing factor S =    1.0000e-01.
Number of knots in each direction  = 10, 13

Sum of squared residuals           =    1.0004e-01.

Derivative of spline has order nux, nuy = 1, 0.

Spline evaluated on X-Y grid (X across, Y down):
0       1       2       3       4       5
0    0.992   2.043   3.029   4.014   5.021   5.997
1    0.541   1.088   1.607   2.142   2.705   3.239
2   -0.417  -0.829  -1.241  -1.665  -2.083  -2.485
3   -0.978  -1.975  -2.914  -3.913  -4.965  -5.924
4   -0.648  -1.363  -1.991  -2.606  -3.251  -3.933

Spline derivative evaluated on X-Y grid:
0       1       2       3       4       5
0    1.093   1.013   0.970   1.004   1.001   0.939
1    0.565   0.531   0.515   0.558   0.559   0.499
2   -0.429  -0.404  -0.421  -0.423  -0.412  -0.389
3   -1.060  -0.951  -0.949  -1.048  -1.031  -0.861
4   -0.779  -0.661  -0.608  -0.628  -0.663  -0.701

ifail =

0

```
```function e02dh_example
mx = 11;
my = 9;
start = 'C';
x = [0:0.5:5];
y = [0:0.5:4];
f = [ 1.0000e+00, 8.8758e-01, 5.4030e-01, 7.0737e-02,-4.1515e-01, ...
-8.0114e-01,-9.7999e-01,-9.3446e-01,-6.5664e-01, 1.5000e+00, ...
1.3564e+00, 8.2045e-01, 1.0611e-01,-6.2422e-01,-1.2317e+00, ...
-1.4850e+00,-1.3047e+00,-9.8547e-01, 2.0600e+00, 1.7552e+00, ...
1.0806e+00, 1.5147e-01,-8.3229e-01,-1.6023e+00,-1.9700e+00, ...
-1.8729e+00,-1.4073e+00, 2.5700e+00, 2.1240e+00, 1.3508e+00, ...
1.7684e-01,-1.0404e+00,-2.0029e+00,-2.4750e+00,-2.3511e+00, ...
-1.6741e+00, 3.0000e+00, 2.6427e+00, 1.6309e+00, 2.1221e-01, ...
-1.2484e+00,-2.2034e+00,-2.9700e+00,-2.8094e+00,-1.9809e+00, ...
3.5000e+00, 3.1715e+00, 1.8611e+00, 2.4458e-01,-1.4565e+00, ...
-2.8640e+00,-3.2650e+00,-3.2776e+00,-2.2878e+00, 4.0400e+00, ...
3.5103e+00, 2.0612e+00, 2.8595e-01,-1.6946e+00,-3.2046e+00, ...
-3.9600e+00,-3.7958e+00,-2.6146e+00, 4.5000e+00, 3.9391e+00, ...
2.4314e+00, 3.1632e-01,-1.8627e+00,-3.6351e+00,-4.4550e+00, ...
-4.2141e+00,-2.9314e+00, 5.0400e+00, 4.3879e+00, 2.7515e+00, ...
3.5369e-01,-2.0707e+00,-4.0057e+00,-4.9700e+00,-4.6823e+00, ...
-3.2382e+00, 5.5050e+00, 4.8367e+00, 2.9717e+00, 3.8505e-01, ...
-2.2888e+00,-4.4033e+00,-5.4450e+00,-5.1405e+00,-3.5950e+00, ...
6.0000e+00, 5.2755e+00, 3.2418e+00, 4.2442e-01,-2.4769e+00, ...
-4.8169e+00,-5.9300e+00,-5.6387e+00,-3.9319e+00];
s = 0.1;
ngx = 6;
xlo = 0;
xhi = 5;
ngy = 5;
ylo = 0;
yhi = 4;
nux = int64(1);
nuy = int64(0);

nxest = mx + 4;
nyest = my + 4;
lamda = zeros(nxest, 1);
mu    = zeros(nyest, 1);
wrk   = zeros(4*(mx+my) + 11*(nxest+nyest) + nxest*my + max(my,nxest) + 54, 1);
iwrk  = zeros(3 + mx + my + nxest + nyest, 1, 'int64');

% Determine the spline approximation
[nx, lamda, ny, mu, c, fp, wrk, iwrk, ifail] = ...
e02dc(start, x, y, f, s, int64(0), lamda, int64(0), mu, wrk, iwrk);

fprintf('\nSpline fit used smoothing factor S = %13.4e.\n', s);
fprintf('Number of knots in each direction  = %d, %d\n', nx, ny);
fprintf('\nSum of squared residuals           = %13.4e.\n', fp);

% Evaluate the spline derivative on a different rectangular grid with
% ngx*ngy points over the domain (xlo to xhi) x (ylo to yhi).
gridx = [xlo:(xhi-xlo)/(ngx-1):xhi];
gridy = [ylo:(yhi-ylo)/(ngy-1):yhi];

% Evaluate spline (nux=nuy=0)
[z, ifail] = e02dh(gridx, gridy, lamda(1:double(nx)), mu(1:double(ny)), ...
c, int64(0), int64(0));

% Evaluate spline partial derivative of order (nux, nuy)
[zder, ifail] = e02dh(gridx, gridy, lamda(1:double(nx)), mu(1:double(ny)), ...
c, nux, nuy);

% Spline partial derivatives to evaluate have order nux,nuy.
fprintf('\nDerivative of spline has order nux, nuy = %d, %d.\n\n', nux, nuy);

% Print the spline evaluations
rlabs = cell(size(gridy));
for i = 1:ngy
rlabs{i} = num2str(gridy(i));
end
clabs = cell(size(gridx));
for i = 1:ngx
clabs{i} = num2str(gridx(i));
end
title = 'Spline evaluated on X-Y grid (X across, Y down):';
[ifail] = x04cb('G', 'N', reshape(z, ngy, ngx), 'F8.3', title, ...
'C', rlabs, 'C', clabs, int64(80), int64(0));
fprintf('\n');
title = 'Spline derivative evaluated on X-Y grid:';
[ifail] = x04cb('G', 'N', reshape(zder, ngy, ngx), 'F8.3', title, ...
'C', rlabs, 'C', clabs, int64(80), int64(0))
```
```

Spline fit used smoothing factor S =    1.0000e-01.
Number of knots in each direction  = 10, 13

Sum of squared residuals           =    1.0004e-01.

Derivative of spline has order nux, nuy = 1, 0.

Spline evaluated on X-Y grid (X across, Y down):
0       1       2       3       4       5
0    0.992   2.043   3.029   4.014   5.021   5.997
1    0.541   1.088   1.607   2.142   2.705   3.239
2   -0.417  -0.829  -1.241  -1.665  -2.083  -2.485
3   -0.978  -1.975  -2.914  -3.913  -4.965  -5.924
4   -0.648  -1.363  -1.991  -2.606  -3.251  -3.933

Spline derivative evaluated on X-Y grid:
0       1       2       3       4       5
0    1.093   1.013   0.970   1.004   1.001   0.939
1    0.565   0.531   0.515   0.558   0.559   0.499
2   -0.429  -0.404  -0.421  -0.423  -0.412  -0.389
3   -1.060  -0.951  -0.949  -1.048  -1.031  -0.861
4   -0.779  -0.661  -0.608  -0.628  -0.663  -0.701

ifail =

0

```