自适应辛普森积分

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
const double eps=1e-10;
double a,b,c,d,l,r;

double f(double x){ //积分函数
  return (c*x+d)/(a*x+b);
}
double simpson(double l,double r){//辛普森公式
  return (r-l)*(f(l)+f(r)+4*f((l+r)/2))/6;
}//二次函数特性

double asr(double l,double r,double ans){//自适应
//分段simpson，如果划分足够小，低于误差就可以
  auto m=(l+r)/2,a=simpson(l,m),b=simpson(m,r);
  if(fabs(a+b-ans)<eps) return ans;
  return asr(l,m,a)+asr(m,r,b);
}
int main(){
  scanf("%lf%lf%lf%lf%lf%lf",&a,&b,&c,&d,&l,&r);
  printf("%.6lf",asr(l,r,simpson(l,r)));
  return 0;
}

时间复杂度：O(log(n/eps))O(log(n/eps))

tips:要注意保证给的初始区间的积分是收敛的并且不要出现无定义点

∫0∞xax−xdx\displaystyle{\int_0^\infty x^{\frac{a}{x}-x}\mathrm{d}x}

保留至小数点后55位。若积分发散，请输出orz\text{orz}。

∣a∣≤50|a|\le50。

Solution:首先把反常积分的发散部分特判，也就是a<0.

对于初始区间，显然不能直接赋值0和无穷大，左端点复制成eps。

考虑右端点，根据a的取值，当x=20的时候代入发现已经远小于eps了故右端点设计为20.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
const double eps = 1e-10;
int n, m;
double a;
double f(double x){
	  return pow(x,(a/x)-x);
}
double sim(double l,double r){
	return (r-l)*(f(l)+f(r)+4*f((l+r)/2))/6;
}
double asr(double l,double r,double ans){
	auto m=(l+r)/2,a=sim(l,m),b=sim(m,r);
	if(fabs(a+b-ans)<eps)return ans;
	return asr(l,m,a)+asr(m,r,b);
}

void solve(){
	cin>>a;
	if(a<0){
		cout<<"orz";return ;
	}
	baoliu(asr(eps,20,sim(eps,20)),5);
}