数值积分

数值积分方法,对复化simpson和Romberg哪个更好的问题,我也不知道怎么比较
Romberg更耗内存,但是它可以对精度控制的相当好…
这个优点简直就是让我没有办法不选它…
复化simpson相当的直接,就相当于把被积区间给分成n份,然后对每份都用simpson公式…
如果是分成n份的话,就需要调用2n + 2次func…所以对算法的复杂度可以控制的很好…
不过,浙大的算法库上用的是Romberg…


//这只是一个函数的示例,使用时按照需求自行修改(增加系数个数,修改函数内容等等)
double func(double x){
return log(x + cos(x));
}

//最简单的simpson积分,积分端点a, b
//只是随便写写,基本没有什么很好的精度

double Simpson(double a, double b){
return (b - a) * (func(a) + func(b) + 4 * func((a+b)/2)) / 6.0;
}

//复化simpson积分,积分端点a, b,以及你想把区间分成的份数n
//这个在于积分区间不变或者变数不大的情况下特别好用
//如果一定要用,则要随着积分区间的的增大,按比例调整n的大小

double C_simpson(double a, double b, int n){
double h = (b - a) / (2 * n), s1 = 0, s2 = 0, x;
for(int i = 1; i < n; i++){
x = a + 2 * i * h;
s1 += func(x);
}
for(int i = 0; i < n; i++){
x = a + (2 * i + 1) * h;
s2 += func(x);
}
return h * (func(a) + 4 * s2 + 2 * s1 + func(b))/3;
}
//Romberg数值积分,输入积分端点a, b,精度要求deta
//计算次数相对较多,效率并不很高,占用内存也相对较大
//优点是可以控制精度
//当数组t的大小不够算出所需的精度的时候(此时函数会返回100000),可以把t开大

double Romberg(double a, double b, double deta){
double h = b - a, t[10][10]; //数组t的大小可以按照所需的精度自行控制
memset(t, 0, sizeof(t));
unsigned int m = 1;;
t[0][0] = h * (func(a) + func(b))/2;
for(int k = 1; k < 10;k++){
h /= 2;
for(int i = 0; i < m; i++)
t[k][0] += func(a + (2 * i + 1) * h);
t[k][0] = 0.5 * (t[k - 1][0] + t[k][0] * h * 2);
for(int i = 1; i <= k; i++){
t[k][i] = t[k][i - 1] + (t[k][i - 1] - t[k - 1][i - 1]) /(pow(2.0, 2.0 * i) - 1);
if(fabs(t[k][i] - t[k - 1][i - 1]) < deta) return t[k][i];
}
m <<= 1;
}
return 100000;
}