// Approximating the area under y = f(x) on [a, b] using the Trapezoidal Rule.
// Version 4:  Approximating area under sqrt(4-x^2) on [0, 2].
//    Here sqrt is computed via Newton's Method, rather than through math.h.

#include <iostream.h>
#include <assert.h>
const double a = 0.0;
const double b = 2.0;

double sqrt(double r) 
/* function to compute square root of r */
{ if (r == 0.0)
     return 0.0;         /* the square root of 0.0 is a special case */
  assert (r > 0.0);      /* negative square roots are not defined */
  double change = 1.0;    /* square roots are found using Newton's method */
  double x = r;
  while ((change > 0.00005) || (change < -0.0005))
    { change = (x*x - r) / (2*x);
      x -= change;
    }
  return x;
}

double f(double x) 
/* function to be used in the area approximation */
{  return sqrt(4.0 - x*x);
}

int main (void)
{  int  n;
   double width, sum, xvalue, area;

   cout << "Program approximates the area under a function using the "
        << "Trapezoidal Rule." << endl;
   cout << "Enter number of subintervals to be used: ";
   cin  >> n;

   width = (b - a) / n; 

   /* compute the sum in the area approximation */
   sum = (f(a) + f(b)) / 2.0;            /* first and last terms in sum */
   for (xvalue = a + width; xvalue < b; xvalue += width)
      sum += f(xvalue);

   area = sum * width;
   cout << "The approximate area is " << area << endl;
   return 0;
}