例

台形公式と同じ をシンプソンの公式で, 計算したのが次です. 少し結果が良くなります.

/* File name 9-2.c */

#include <stdio.h>

double simpson(double a, double b, int n);
double fun(double x);

main()
{
     printf("Pi is about %1.15f.\n", simpson(0.0, 1.0, 65536));
}

double simpson(double a, double b, int n)
{
     double sum=fun(a);           
     double interval=(b-a)/n;   /* The length of the small interval */
     int k;
    
     for (k=1; k < n ; k++){ 
         sum = sum + 2.0*fun(a+k*interval)+4.0*fun(a+k*interval-0.5*interval); 
     }
     return interval*(sum + fun(b)+4.0*fun(b-0.5*interval))/6.0;
}

double fun(double x)
{
      return 4.0/(1.0 + x*x);
}