/* 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); printf("nag_2d_triang_interp (e01sjc) Example Program Results\n"); /* Skip heading in data file */ scanf("%*[^\n] "); /* Input the number of nodes. */ scanf("%ld%*[^\n] ", &m); if (m >= 3) { /* Allocate memory */ if (!(f = NAG_ALLOC(m, double)) || !(grads = NAG_ALLOC(2 * m, double)) || !(x = NAG_ALLOC(m, double)) || !(y = NAG_ALLOC(m, double)) || !(triang = NAG_ALLOC(7*m, Integer))) { printf("Allocation failure\n"); exit_status = -1; goto END; } } else { printf("Invalid m.\n"); exit_status = 1; goto END; } /* Input the nodes (X,Y) and heights, F. */ for (i = 1; i <= m; ++i) { scanf("%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) { printf("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). */ scanf("%ld%lf%lf%*[^\n] ", &nx, &xlo, &xhi); scanf("%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))) { printf("Allocation failure\n"); exit_status = -1; goto END; } } else { printf("Invalid nx or ny.\n"); exit_status = 1; goto END; } 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; } printf("\n"); printf("%s", " X"); for (i = 1; i <= nx; ++i) { printf("%8.2f", px[i - 1]); } printf("\n"); printf("%s", " Y"); printf("\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_triang_interp (e01sjc). */ 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) { printf("Error from nag_2d_triang_eval (e01skc).\n%s\n", fail.message); exit_status = 1; goto END; } } printf("%8.2f", py[i - 1]); printf("%3s", ""); for (j = 1; j <= nx; ++j) { printf("%8.2f", pf[j - 1]); } printf("\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; }