```/* nag_2d_spline_deriv_rect (e02dhc) Example Program.
*
* NAGPRODCODE Version.
*
* Copyright 2016 Numerical Algorithms Group.
*
* Mark 26, 2016.
*/

#include <math.h>
#include <string.h>
#include <nag.h>
#include <nag_stdlib.h>
#include <nage02.h>
#include <nagx04.h>

#ifdef __cplusplus
extern "C"
{
#endif

static void NAG_CALL print_spline(Integer *ngx, double *gridx, Integer *ngy,
double *gridy, double *z, double *zder,
Integer *exit_status);

#ifdef __cplusplus
}
#endif

#define F(I, J) f[my*(I)+(J)]

int main(void)
{
/* Scalars */
Integer exit_status = 0;
Integer i, j, mx, my, ngx, ngy, nux, nuy, nxest, nyest;
double delta, fp, s, xhi, xlo, yhi, ylo;
/* Arrays */
double *f = 0, *gridx = 0, *gridy = 0, *x = 0, *y = 0, *z = 0, *zder = 0;
/* NAG types */
Nag_2dSpline spline;
Nag_Comm warmstartinf;
Nag_Start startc;
NagError fail;

INIT_FAIL(fail);

printf("nag_2d_spline_deriv_rect (e02dhc) Example Program Results\n");
fflush(stdout);

/* Skip heading in data file */
scanf("%*[^\n] ");

/* Input the number of X, Y co-ordinates MX, MY. */
scanf("%" NAG_IFMT " %" NAG_IFMT "%*[^\n]", &mx, &my);
nxest = mx + 4;
nyest = my + 4;
spline.nx = 4;
spline.ny = 4;

/* Alocations for spline fitting */
if (!(x = NAG_ALLOC((mx), double)) ||
!(y = NAG_ALLOC((my), double)) || !(f = NAG_ALLOC((mx * my), double))
)
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
for (i = 0; i < mx; i++)
scanf("%lf", &x[i]);
scanf("%*[^\n]");
for (i = 0; i < my; i++)
scanf("%lf", &y[i]);
scanf("%*[^\n]");

/* Input the MX*MY function values F at grid points and smoothing factor. */
for (i = 0; i < mx; i++)
for (j = 0; j < my; j++)
scanf("%lf", &F(i, j));
scanf("%*[^\n]");
scanf("%lf%*[^\n]", &s);

/* nag_2d_spline_fit_grid (e02dcc).
* Least squares bicubic spline fit with automatic knot placement,
* two variables (rectangular grid)
*/
startc = Nag_Cold;
nag_2d_spline_fit_grid(startc, mx, x, my, y, f, s, nxest, nyest, &fp,
&warmstartinf, &spline, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_2d_spline_fit_grid (e02dcc)\n%s\n", fail.message);
exit_status = 1;
goto END;
}

printf("\nSpline fit used smoothing factor s = %13.4e.\n", s);
printf("Number of knots in each direction  = %5" NAG_IFMT ",%5" NAG_IFMT
".\n\n", spline.nx, spline.ny);
printf("Sum of squared residuals           = %13.4e.\n", fp);
fflush(stdout);

/* Spline and its derivative to be evaluated on rectangular grid with
* ngx*ngy points on the domain [xlo,xhi] by [ylo,yhi].
*/
scanf("%" NAG_IFMT "%lf%lf%*[^\n]", &ngx, &xlo, &xhi);
scanf("%" NAG_IFMT "%lf%lf%*[^\n]", &ngy, &ylo, &yhi);

/* Allocations for spline evaluation. */
if (!(gridx = NAG_ALLOC((ngx), double)) ||
!(gridy = NAG_ALLOC((ngy), double)) ||
!(z = NAG_ALLOC((ngx * ngy), double)) ||
!(zder = NAG_ALLOC((ngx * ngy), double))
)
{
printf("Allocation failure\n");
exit_status = -2;
goto END;
}

delta = (xhi - xlo) / (double) (ngx - 1);
gridx = xlo;
for (i = 1; i < ngx - 1; i++)
gridx[i] = gridx[i - 1] + delta;
gridx[ngx - 1] = xhi;

delta = (yhi - ylo) / (double) (ngy - 1);
gridy = ylo;
for (i = 1; i < ngy - 1; i++)
gridy[i] = gridy[i - 1] + delta;
gridy[ngy - 1] = yhi;

/* Evaluate spline (nux=nuy=0) using
* nag_2d_spline_deriv_rect (e02dhc).
* Evaluation of spline surface at mesh of points with derivatives
*/
nux = 0;
nuy = 0;
nag_2d_spline_deriv_rect(ngx, ngy, gridx, gridy, nux, nuy, z, &spline,
&fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_2d_spline_deriv_rect (e02dhc))\n%s\n",
fail.message);
exit_status = 2;
goto END;
}

/* Evaluate spline partial derivative of order (nux,nuy) */
scanf("%" NAG_IFMT "%" NAG_IFMT "%*[^\n]", &nux, &nuy);
printf("\nDerivative of spline has order nux, nuy =%5" NAG_IFMT ", %5"
NAG_IFMT ".\n", nux, nuy);
fflush(stdout);
/* nag_2d_spline_deriv_rect (e02dhc).
* Evaluation of spline surface at mesh of points with derivatives
*/
nag_2d_spline_deriv_rect(ngx, ngy, gridx, gridy, nux, nuy, zder, &spline,
&fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_2d_spline_deriv_rect (e02dhc)\n%s\n",
fail.message);
exit_status = 3;
goto END;
}
fflush(stdout);

/* Print tabulated spline and derivative evaluations. */
print_spline(&ngx, gridx, &ngy, gridy, z, zder, &exit_status);

END:
NAG_FREE(f);
NAG_FREE(gridx);
NAG_FREE(gridy);
NAG_FREE(x);
NAG_FREE(y);
NAG_FREE(z);
NAG_FREE(zder);
NAG_FREE(spline.lamda);
NAG_FREE(spline.mu);
NAG_FREE(spline.c);
NAG_FREE(warmstartinf.nag_w);
NAG_FREE(warmstartinf.nag_iw);
return exit_status;
}

static void NAG_CALL print_spline(Integer *ngx, double *gridx, Integer *ngy,
double *gridy, double *z, double *zder,
Integer *exit_status)
{
/* Print spline function and spline derivative evaluation */
Integer indent = 0, ncols = 80;
char formc[] = "%8.3f";
Integer i;
char title;
char *outfile = 0;
char **clabsc = 0, **rlabsc = 0;
Nag_OrderType order;
Nag_MatrixType matrixc;
Nag_DiagType diagc;
Nag_LabelType chlabelc;
NagError fail;

INIT_FAIL(fail);

/* Allocate for row and column label */
if (!(clabsc = NAG_ALLOC(*ngx, char *)) ||
!(rlabsc = NAG_ALLOC(*ngy, char *))
)
{
printf("Allocation failure\n");
*exit_status = -3;
goto END;
}

/* Allocate memory to clabsc and rlabsc elements and generate
* column and row labels to print the results with.
*/
for (i = 0; i < *ngx; i++) {
clabsc[i] = NAG_ALLOC(11, char);
sprintf(clabsc[i], "%5.2f%5s", gridx[i], "");
}
for (i = 0; i < *ngy; i++) {
rlabsc[i] = NAG_ALLOC(11, char);
sprintf(rlabsc[i], "%5.2f%5s", gridy[i], "");
}

order = Nag_ColMajor;
matrixc = Nag_GeneralMatrix;
diagc = Nag_NonUnitDiag;
chlabelc = Nag_CharacterLabels;

/* Print the spline evaluations, z. */
strcpy(title, "Spline evaluated on X-Y grid (X across, Y down):");
printf("\n");
fflush(stdout);

/* nag_gen_real_mat_print_comp (x04cbc).
* Print real general matrix (comprehensive)
*/
nag_gen_real_mat_print_comp(order, matrixc, diagc, *ngy, *ngx, z, *ngy,
formc, title, chlabelc, (const char **) rlabsc,
chlabelc, (const char **) clabsc, ncols,
indent, outfile, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_gen_real_mat_print_comp (x04cbc)\n%s\n",
fail.message);
*exit_status = 4;
goto END;
}

/* Print the spline derivative evaluations, zder. */
strcpy(title, "Spline derivative evaluated on X-Y grid:");
printf("\n");
fflush(stdout);

nag_gen_real_mat_print_comp(order, matrixc, diagc, *ngy, *ngx, zder, *ngy,
formc, title, chlabelc, (const char **) rlabsc,
chlabelc, (const char **) clabsc, ncols, indent,
outfile, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_gen_real_mat_print_comp (x04cbc)\n%s\n",
fail.message);
*exit_status = 5;
goto END;
}

END:
for (i = 0; i < *ngy; i++)
NAG_FREE(rlabsc[i]);
NAG_FREE(rlabsc);
for (i = 0; i < *ngx; i++)
NAG_FREE(clabsc[i]);
NAG_FREE(clabsc);
}
```