/* nag_2d_triang_interp (e01sjc) Example Program. * * Copyright 2005 Numerical Algorithms Group. * * Mark 8, 2004. */ #include #include #include #include int main(void) { /* Scalars */ double xhi, xlo, yhi, ylo; Integer exit_status, i, j, m, nx, ny; /* Arrays */ double *f=0, *grads=0, *pf=0, *px=0, *py=0, *x=0, *y=0; Integer *triang=0; /* Nag Types */ NagError fail; exit_status = 0; INIT_FAIL(fail); Vprintf( "nag_2d_triang_interp (e01sjc) Example Program Results\n"); /* Skip heading in data file */ Vscanf("%*[^\n] "); /* Input the number of nodes. */ Vscanf("%ld%*[^\n] ", &m); if (m >= 3) { /* Allocate memory */ if ( !(f = NAG_ALLOC(m, double)) || !(grads = NAG_ALLOC(2 * m, double)) || !(pf = NAG_ALLOC(25, double)) || !(px = NAG_ALLOC(25, double)) || !(py = NAG_ALLOC(25, double)) || !(x = NAG_ALLOC(m, double)) || !(y = NAG_ALLOC(m, double)) || !(triang = NAG_ALLOC(7*m, Integer)) ) { Vprintf("Allocation failure\n"); exit_status = -1; goto END; } } else { Vprintf("Invalid m.\n"); exit_status = 1; return exit_status; } /* Input the nodes (X,Y) and heights, F. */ for (i = 1; i <= m; ++i) { Vscanf("%lf%lf%lf%*[^\n] ", &x[i - 1], &y[i - 1], &f[i - 1]); } /* Generate the triangulation and gradients. */ /* nag_2d_triang_interp (e01sjc). * A function to generate a two-dimensional surface * interpolating a set of data points, using either the * method of Renka and Cline or the modified Shepard's * method */ nag_2d_triang_interp(m, x, y, f, triang, grads, &fail); if (fail.code != NE_NOERROR) { Vprintf("Error from nag_2d_triang_interp (e01sjc).\n%s\n", fail.message); exit_status = 1; goto END; } /* Evaluate the interpolant on a rectangular grid at NX*NY points */ /* over the domain (XLO to XHI) x (YLO to YHI). */ Vscanf("%ld%lf%lf%*[^\n] ", &nx, &xlo, &xhi); Vscanf("%ld%lf%lf%*[^\n] ", &ny, &ylo, &yhi); if (nx > 0 && ny >0 ) { /* Allocate memory */ if ( !(pf = NAG_ALLOC(nx, double)) || !(px = NAG_ALLOC(nx, double)) || !(py = NAG_ALLOC(ny, double)) ) { Vprintf("Allocation failure\n"); exit_status = -1; goto END; } } else { Vprintf("Invalid nx or ny.\n"); exit_status = 1; return exit_status; } for (i = 1; i <= nx; ++i) { px[i - 1] = (double)(nx - i) / (nx - 1) * xlo + (double)(i - 1) / (nx - 1) * xhi; } for (i = 1; i <= ny; ++i) { py[i - 1] = (double)(ny - i) / (ny - 1) * ylo + (double)(i - 1) / (ny - 1) * yhi; } Vprintf("\n"); Vprintf("%s", " X"); for (i = 1; i <= nx; ++i) { Vprintf("%8.2f", px[i - 1]); } Vprintf("\n"); Vprintf("%s", " Y"); Vprintf("\n"); for (i = ny; i >= 1; --i) { for (j = 1; j <= nx; ++j) { /* nag_2d_triang_eval (e01skc). * A function to evaluate, at a set of points, the * two-dimensional interpolant function generated by * nag_2d_scat_interpolant (e01sac) */ nag_2d_triang_eval(m, x, y, f, triang, grads, px[j - 1], py[i - 1], &pf[j - 1], &fail); if (fail.code != NE_NOERROR) { Vprintf("Error from nag_2d_triang_eval (e01skc).\n%s\n", fail.message); exit_status = 1; goto END; } } Vprintf("%8.2f", py[i - 1]); Vprintf("%3s", ""); for (j = 1; j <= nx; ++j) { Vprintf("%8.2f", pf[j - 1]); } Vprintf("\n"); } END: if (f) NAG_FREE(f); if (grads) NAG_FREE(grads); if (pf) NAG_FREE(pf); if (px) NAG_FREE(px); if (py) NAG_FREE(py); if (x) NAG_FREE(x); if (y) NAG_FREE(y); if (triang) NAG_FREE(triang); return exit_status; }