#include <iostream.h>
#include <math.h>
// the function we want to integrate
float f (float x)
{ 
	return (1.0/x);               
}
// Numerical integration using Simpson's rule 
float simpson (int n, float min, float max)
{  
   int i;				 
   float interval, sum, x;
   interval = ((max -min) /n);
   sum=0;
   for (i=1; i<n; i=i+2)           //loop for odd points
   {
       x = min + interval * i;
       sum = sum+ 4 * f(x);
   }
   for (i=2; i<n; i=i+2)           // loop for even points 
   {
      x = min+interval * i;
      sum = sum + 2 * f(x);
   }   
   sum = sum + f(min) + f(max);	 	// add first and last value          
   sum = sum * interval/3.;        	// then multilpy by interval
   return sum;
} 
void main()
{
	cout<<simpson (4, 1,2)<<endl;
}