/* nag_ode_bvp_ps_lin_cgl_deriv (d02udc) Example Program.
*
* Copyright 2022 Numerical Algorithms Group.
*
* Mark 28.7, 2022.
*/
#include <math.h>
#include <nag.h>
#ifdef __cplusplus
extern "C" {
#endif
static double NAG_CALL fcn(double x);
static double NAG_CALL deriv(double x);
#ifdef __cplusplus
}
#endif
int main(void) {
/* Scalars */
Integer exit_status = 0;
Integer i, n;
double a = 0.0, b = 1.5, scale;
double teneps = 100.0 * nag_machine_precision;
double uxerr = 0.0;
/* Arrays */
double *f = 0, *fd = 0, *x = 0;
/* NAG types */
Nag_Boolean reqerr = Nag_FALSE;
NagError fail;
INIT_FAIL(fail);
printf("nag_ode_bvp_ps_lin_cgl_deriv (d02udc) Example Program Results\n\n");
/* Skip heading in data file */
scanf("%*[^\n] ");
scanf("%" NAG_IFMT "%*[^\n] ", &n);
if (!(f = NAG_ALLOC((n + 1), double)) || !(fd = NAG_ALLOC((n + 1), double)) ||
!(x = NAG_ALLOC((n + 1), double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
/* nag_ode_bvp_ps_lin_cgl_grid (d02ucc).
* Generate Chebyshev Gauss-Lobatto solution grid.
*/
nag_ode_bvp_ps_lin_cgl_grid(n, a, b, x, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_ode_bvp_ps_lin_cgl_grid (d02ucc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
/* Evaluate the function on Chebyshev grid. */
for (i = 0; i < n + 1; i++)
f[i] = fcn(x[i]);
/* nag_ode_bvp_ps_lin_cgl_deriv (d02udc).
* Differentiate a function using function values on Chebyshev grid.
*/
nag_ode_bvp_ps_lin_cgl_deriv(n, f, fd, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_ode_bvp_ps_lin_cgl_deriv (d02udc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
scale = 2.0 / (b - a);
for (i = 0; i < n + 1; i++)
fd[i] = scale * fd[i];
/* Print function and its derivative. */
printf("Original function f and numerical derivative fx\n\n");
printf("%8s%11s%11s\n", "x", "f", "fx");
for (i = 0; i < n + 1; i++)
printf("%10.4f %10.4f %10.4f\n", x[i], f[i], fd[i]);
if (reqerr) {
for (i = 0; i < n + 1; i++)
uxerr = MAX(uxerr, fabs(fd[i] - deriv(x[i])));
printf("fx is within a multiple %8" NAG_IFMT " of machine precision.\n",
100 * ((Integer)(uxerr / teneps) + 1));
}
END:
NAG_FREE(f);
NAG_FREE(fd);
NAG_FREE(x);
return exit_status;
}
static double NAG_CALL fcn(double x) { return 2.0 * x + exp(-x); }
static double NAG_CALL deriv(double x) { return 2.0 - exp(-x); }