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

double Quadrature11Points(double(*func)(double),double a,double b,int n)
			{
				double h,x,sum; int i;
				if(func == NULL || n == 0) return 0.0;

				h = (b-a)/(10*n);
				sum = 16067.0*(-func(a) + func(b));

				for(i = 0; i < n; i++){
					x = a +10*h*i;
					sum += 32134.0*func(x)      + 
                     106300.0*(func(x+h)   + func(x+9*h)) -
                      48525.0*(func(x+2*h) + func(x+8*h)) + 
                     272400.0*(func(x+3*h) + func(x+7*h)) - 
                     260550.0*(func(x+4*h) + func(x+6*h)) + 
                     427368.0*func(x+5*h);
				}
				return 5*h*sum/299376;
			}
double PI;
         
double funC(double x){return 1.0;}
         
double funNP(double x){ return exp(-x*x*0.5)/sqrt(2*PI); }

int main(){

double a=0.0, b=0.2;
double x[1],y[1];

PI=4*atan(1.0);
printf("F(x) = 1;          Exact result = %e; Numerical result = %e\n",
   (b-a),Quadrature11Points(funC,a,b,1));
x[0]=b/sqrt(2.0); vdErf(1,x,y);   
printf("F(x) = Norm.p.d.f; Exact result = %e; Numerical result = %e\n",
   (1+y[0])*0.5,0.5+Quadrature11Points(funNP,a,b,1));

}