/* nag_specfun_mathieu_ang_periodic_real (s22cac) Example Program.
 *
 * Copyright 2023 Numerical Algorithms Group.
 *
 * Mark 29.0, 2023.
 */

#include <math.h>
#include <nag.h>
#include <stdio.h>

int main(void) {
  Integer n = 2, mode, ordval, parity, exit_status = 0, n_order;
  double x[2], f[2], f_deriv[2];
  double a, q;
  /* Nag Types */
  NagError fail;

  INIT_FAIL(fail);

  printf("nag_specfun_mathieu_ang_periodic_real (s22cac) ");
  printf("Example Program Results");
  printf("\n      parity      ordval           ");
  printf("a           x   ce_m(x,q)  ce_m'(x,q)\n");

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

  /* Get values */
  scanf("%lf%*[^\n] ", &q);
  scanf("%" NAG_IFMT "%*[^\n] ", &n_order);
  scanf("%lf%lf%*[^\n] ", &x[0], &x[1]);

  mode = 2;

  /* Loop over even orders */
  parity = 0;
  for (ordval = 0; ordval < n_order; ordval++) {
    nag_specfun_mathieu_ang_periodic_real(ordval, q, parity, mode, n, x, f,
                                          f_deriv, &a, &fail);
    if (fail.code != NE_NOERROR) {
      printf("\nError from nag_specfun_mathieu_ang_periodic_real (s22cac).");
      printf("\n%s\n", fail.message);
      exit_status = 1;
    } else
      /* Print values at x=0.0 */
      printf("%12" NAG_IFMT "%12" NAG_IFMT "%12.4f%12.4f%12.4f%12.4f\n", parity,
             ordval, a, x[0], f[0], f_deriv[0]);
  }

  printf("      parity      ordval           ");
  printf("a           x   se_m(x,q)  se_m'(x,q)\n");

  /* Loop over odd orders */
  parity = 1;
  for (ordval = 1; ordval <= n_order; ordval++) {
    nag_specfun_mathieu_ang_periodic_real(ordval, q, parity, mode, n, x, f,
                                          f_deriv, &a, &fail);
    if (fail.code != NE_NOERROR) {
      printf("\nError from nag_specfun_mathieu_ang_periodic_real (s22cac).");
      printf("\n%s\n", fail.message);
      exit_status = 2;
    } else
      /* Print values at x=0.0 */
      printf("%12" NAG_IFMT "%12" NAG_IFMT "%12.4f%12.4f%12.4f%12.4f\n", parity,
             ordval, a, x[0], f[0], f_deriv[0]);
  }

  printf("\n");

  return exit_status;
}