复化辛普森公式

/*
 * 复化梯形公式
 * 复化辛普森公式
 */
#include <iostream>
#include <iomanip>
#include <cmath>
using namespace std;
/*
 * 复化梯形公式
 * a表示区间起点
 * b表示区间终点
 * n表示区间长度
 */
double Trapezoid(double(*fun)(double x), double a, double b, int n)
{
	double result = 0.0;
	if(a > b)
		swap(a, b);
	result = (fun(a) + fun(b))/2;
	double h = (1.0/n) * (b-a);		//步长
	cout << setw(5) << a << setw(10) << fun(a) << endl;
	for(int k=1; k<n; k++)
	{
		cout << setw(5) << a+k*h << setw(10)<< fun(a+k*h) <<endl;
		result += fun(a + k*h);
	}
	cout << setw(5) << b << setw(10) << fun(b) << endl;
	result *= h;
	return result;
}

double Simpson(double(*fun)(double x), double a, double b, int n)
{
	double result = 0.0;
	if(a > b)
		swap(a, b);
	double h = (1.0/n) * (b-a);	//步长
    result += (fun(a) - fun(b));
    double x = a;
    for(int k=1; k<=n; k++)
	{
		x += h/2;
		result += 4*fun(x);
		x += h/2;
		result += 2*fun(x);
	}
	return (result*h)/6;
}

//待积分函数
double function(double x)
{
	if(x == 0)
		return 1;
	else
	return sin(x)/x;
}

int main()
{
	cout << "---复化梯形公式---" << endl;
	cout << "Result = " << Trapezoid(function, 0, 1, 8) << endl;
	cout << endl;
	cout << "---复化辛普森公式---" << endl;
	cout << "Result = " << Simpson(function, 0, 1, 4) << endl;
	return 0;
}