/**
復化Simpson積分
m-1 m-1
復化Simpson公式: Tn = (h/3) * { f(a)/2 + 4 * ∑ f ( x2i + 2 * ∑ f ( x2i+1 ) + f(b)/2 )}
i=0 i=1
步長: h=b-a/n, n=2*m
屬性: 數值積分法
《數值計算方法與算法》-2 Editon -科學出版社 P55
代碼維護:2007.04.20 pengkuny
**/
#include<iostream>
#include<cmath>
using namespace std;
#define f(x) (sin(x)) //舉例函數
#define epsilon 0.00001 //精度
//復化Simpson公式
float Simpson(float aa, float bb, long n)

{//aa,bb:端點 n,分區數
if (!(n%2)) //保證n為偶數
n = n+1;
long m = n/2;
float sum = 0;
float h = (bb-aa)/n; //步長
float sum1 = 0, sum2 = 0;
for (long i=0; i<m; i++)
{
sum1 += f(aa + (2*i+1)*h);
}
for (long i=1; i<m; i++)
{
sum2 += f(aa + (2*i)*h);
}
sum = f(aa) + f(bb) + 2*sum2 + 4*sum1;
return (h*sum/3);
}

int main()

{
float a,b;
long n;
cout<<"復化Simpson積分,請輸入積分范圍a,b:"<<endl;
cin>>a>>b;
cout<<"請輸入分割區間數n:"<<endl;
cin>>n;
cout<<"積分結果:"<<Simpson(a, b, n)<<endl;
system("pause");
return 0;
}

