/* test_aquad3.c adaptive quadrature test */
#include <stdio.h>
#include <math.h>
#include "aquad3.h"

double f(double x)
{
  return 1.0/x;
}

double f1(double x)
{
  return 1.0/((x-0.3)*(x-0.3)+0.01) + 1.0/((x-0.9)*(x-0.9)+0.04) -6;
}

int main(int argc, char * argv[])
{
  double xmin = 0.1;
  double xmax = 2.0;
  double eps  = 0.001;
  double area, err, exact;

  printf("test_aquad3.c testing aquad3.c 1/x  eps=%g \n", eps);
  area = aquad3(f, xmin, xmax, eps);
  exact = log(xmax)-log(xmin);
  err = area-exact;
  printf("xmin=%g, xmax=%g, area=%g, exact=%g, err=%g \n\n",
         xmin, xmax, area, exact, err);

  xmin = 0.01;
  area = aquad3(f, xmin, xmax, eps);
  exact = log(xmax)-log(xmin);
  err = area-exact;
  printf("xmin=%g, xmax=%g, area=%g, exact=%g, err=%g \n\n",
         xmin, xmax, area, exact, err);

  xmin = 0.001;
  area = aquad3(f, xmin, xmax, eps);
  exact = log(xmax)-log(xmin);
  err = area-exact;
  printf("xmin=%g, xmax=%g, area=%g, exact=%g, err=%g \n\n",
         xmin, xmax, area, exact, err);

  xmin = 0.0001;
  area = aquad3(f, xmin, xmax, eps);
  exact = log(xmax)-log(xmin);
  err = area-exact;
  printf("xmin=%g, xmax=%g, area=%g, exact=%g, err=%g \n\n",
         xmin, xmax, area, exact, err);

  xmin = 0.00001;
  area = aquad3(f, xmin, xmax, eps);
  exact = log(xmax)-log(xmin);
  err = area-exact;
  printf("xmin=%g, xmax=%g, area=%g, exact=%g, err=%g \n\n",
         xmin, xmax, area, exact, err);

  printf("1.0/((x-0.3)*(x-0.3)+0.01) + 1.0/((x-0.9)*(x-0.9)+0.04) -6 \n");
  xmin = 0.0;
  xmax = 1.0;
  eps = 0.001;
  area = aquad3(f1, xmin, xmax, eps);
  exact = 29.85832540;
  err = area - exact;  
  printf("xmin=%g, xmax=%g, area=%g, exact=%g, err=%g \n\n",
         xmin, xmax, area, exact, err);
   
  printf("test_aquad3.c finished \n");
  return 0;
}