/*
  Simple trapezoidal integration in 1D.
  Dag Lindbo, dag@kth.se, 2009-08-26
*/

#include <stdio.h>
#include <math.h>
double trapez(double (*func)(double), double a, double b, int N)
{
    double h=(b-a)/N, sum, x;
    int i;
    for(i=1,sum=0,x=a+h; i<N; i++,x+=h)
	sum += func(x);
    return h*(sum+0.5*(func(a)+func(b)));
}

int main(void)
{
    const int N = 102400000;
    const double a = 1;
    const double b = 2;
    
    double quad = trapez(sin, a, b, N);
    printf("  %.16f\n", quad);
    return 0;
}
