```/* nag_2d_spline_fit_scat (e02ddc) Example Program.
*
* Copyright 2017 Numerical Algorithms Group.
*
* Mark 26.1, 2017.
*
*/

#include <nag.h>
#include <stdio.h>
#include <nag_stdlib.h>
#include <nage02.h>

int main(void)
{

Integer exit_status = 0, i, j, m, npx, npy, nx, nxest, ny, nyest, rank;
NagError fail;
Nag_2dSpline spline;
Nag_Start start;
double delta, *f = 0, *fg = 0, fp, *px = 0, *py = 0, s, warmstartinf;
double *weights = 0, *x = 0, xhi, xlo, *y = 0, yhi, ylo;

INIT_FAIL(fail);

/* Initialize spline */
spline.lamda = 0;
spline.mu = 0;
spline.c = 0;

nxest = 14;
nyest = 14;
printf("nag_2d_spline_fit_scat (e02ddc) Example Program Results\n");
scanf("%*[^\n]"); /* Skip heading in data file */
/* Input the number of data-points m. */
scanf("%" NAG_IFMT "", &m);
if (m >= 16) {
if (!(f = NAG_ALLOC(m, double)) ||
!(weights = NAG_ALLOC(m, double)) ||
!(x = NAG_ALLOC(m, double)) || !(y = NAG_ALLOC(m, double)))
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
}
else {
printf("Invalid m.\n");
exit_status = 1;
return exit_status;
}
/* Input the data-points and the weights. */
for (i = 0; i < m; i++)
scanf("%lf%lf%lf%lf", &x[i], &y[i], &f[i], &weights[i]);
start = Nag_Cold;
if (scanf("%lf", &s) != EOF) {
/* Determine the spline approximation. */

/* nag_2d_spline_fit_scat (e02ddc).
* Least squares bicubic spline fit with automatic knot
* placement, two variables (scattered data)
*/
nag_2d_spline_fit_scat(start, m, x, y, f, weights, s, nxest, nyest, &fp,
&rank, &warmstartinf, &spline, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_2d_spline_fit_scat (e02ddc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
nx = spline.nx;
ny = spline.ny;
printf("\nCalling with smoothing factor s = %13.4e, nx = %1" NAG_IFMT ","
" ny = %1" NAG_IFMT "\n", s, nx, ny);

printf("rank deficiency = %1" NAG_IFMT "\n\n",
(nx - 4) * (ny - 4) - rank);

/* Print the knot sets, lamda and mu. */
printf("Distinct knots in x direction located  at\n");
for (j = 3; j < spline.nx - 3; j++)
printf("%12.4f%s", spline.lamda[j],
((j - 3) % 5 == 4 || j == spline.nx - 4) ? "\n" : " ");
printf("\nDistinct knots in y direction located  at\n");
for (j = 3; j < spline.ny - 3; j++)
printf("%12.4f%s", spline.mu[j],
((j - 3) % 5 == 4 || j == spline.ny - 4) ? "\n" : " ");
printf("\nThe B-spline coefficients:\n\n");
for (i = 0; i < ny - 4; i++) {
for (j = 0; j < nx - 4; j++)
printf("%9.2f", spline.c[i + j * (ny - 4)]);
printf("\n");
}

printf("\n Sum of squared residuals fp = %13.4e\n", fp);
if (nx == 8 && ny == 8)
printf("The spline is the least squares bi-cubic polynomial\n");

/* Evaluate the spline on a rectangular grid at npx*npy points
* over the domain (xlo to xhi) x (ylo to yhi).
*/
scanf("%" NAG_IFMT "%lf%lf", &npx, &xlo, &xhi);
scanf("%" NAG_IFMT "%lf%lf", &npy, &ylo, &yhi);

if (npx >= 1 && npy >= 1) {
if (!(fg = NAG_ALLOC(npx * npy, double)) ||
!(px = NAG_ALLOC(npx, double)) || !(py = NAG_ALLOC(npy, double)))
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
}
else {
printf("Invalid npx or npy.\n");
exit_status = 1;
return exit_status;
}
delta = (xhi - xlo) / (npx - 1);
for (i = 0; i < npx; i++)
px[i] = MIN(xhi, xlo + i * delta);
for (i = 0; i < npy; i++)
py[i] = MIN(yhi, ylo + i * delta);

/* nag_2d_spline_eval_rect (e02dfc).
* Evaluation of bicubic spline, at a mesh of points
*/
nag_2d_spline_eval_rect(npx, npy, px, py, fg, &spline, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_2d_spline_eval_rect (e02dfc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
printf("\nValues of computed spline:\n\n");
printf("          x");
for (i = 0; i < npx; i++)
printf("%8.2f ", px[i]);
printf("\n     y\n");
for (i = npy - 1; i >= 0; i--) {
printf("%8.2f   ", py[i]);
for (j = 0; j < npx; j++)
printf("%8.2f ", fg[npy * j + i]);
printf("\n");
}
/* Free memory used by spline */
NAG_FREE(spline.lamda);
NAG_FREE(spline.mu);
NAG_FREE(spline.c);
NAG_FREE(fg);
NAG_FREE(px);
NAG_FREE(py);
}
END:
NAG_FREE(f);
NAG_FREE(weights);
NAG_FREE(x);
NAG_FREE(y);
return exit_status;
}
```