luoguP4526 【模板】自适应辛普森法2

Description

试计算积分${\int }_0^{\infty }{x}^{\frac{a}{x}-x}dx(|a|\le 50)$，保留至小数点后$5$位。若积分发散，请输出orz。

Solution

如何分析题目中给定积分

设题目给定的无界（广义）积分中被积函数为$f(x)$；首先，上述积分中被积函数自变量$x$取不到$0$且无上界；其次，由于$f(x)$是“幂指函数”，故先“指对互化”：

$f(x)=x^{\frac{a}{x}-x}=e^{(\frac{a}{x}-x)lnx}$ 设$g(x)=(\frac{a}{x}-x)lnx$，则$f(x)=e^{g(x)}$

通常情况下，对于此题，我们判断$x\to 0$这种情况下积分是否收敛可试着求${\lim }_{\epsilon \to 0^{+}}{\int }_{0+\epsilon }^{\infty }{x}^{\frac{a}{x}-x}dx$看其是否存在（即是否等于一常数）；对于$x\to +\infty$的情况，则看$f(x)$原函数在$x\to +\infty$时极限是否存在；但此积分着实没法用常规方法求原函数，故考虑审敛法。

（注意：被积函数在某处趋近于正无穷（常数）并不一定代表其积分发散（收敛）！即：不能简单根据被积函数的取值判断收敛发散。

Eg：1、${\int }_0^{1}lnxdx$：此积分可用上面设$\epsilon$的方法分部积分后得到其$=-1$，收敛；2、虽然$x\to +\infty$时$\frac{1}{x}$趋近于$0$，但${\int }_1^{+\infty }\frac{1}{x}dx$显然发散。）

关于审敛法

一、无穷界积分的审敛法

设$f(x)$在积分区间$[a,+\infty )$上恒有$f(x)\ge 0$：

（1）、若存在常数$p>1$使得${\lim }_{x\to +\infty }{x}^{p}f(x)=c$，那么${\int }_a^{+\infty }f(x)dx$收敛。

（2）、若${\lim }_{x\to +\infty }xf(x)=d>0$，那么${\int }_a^{+\infty }f(x)dx$发散。

将积分区间换为$(-\infty ,a]$时定理仍成立。

二、无穷间断点的审敛法

设$a$是$f(x)$的无穷间断点，$f(x)$在$(a,b]$上有$f(x)\ge 0$：

1、若存在常数$0<p<1$使得${\lim }_{x\to a^{+}}(x-a{)}^{p}f(x)=c$，那么${\int }_a^{b}f(x)dx$收敛。

2、若${\lim }_{x\to a^{+}}(x-a)f(x)=d>0$，那么${\int }_a^{b}f(x)dx$发散。

类似的，设$b$是$f(x)$的无穷间断点，$f(x)$在$[a,b)$上有$f(x)\ge 0$：

1、若存在常数$0<p<1$使得${\lim }_{x\to b^{-}}(b-x{)}^{p}f(x)=c$，那么${\int }_a^{b}f(x)dx$收敛。

2、若${\lim }_{x\to b^{-}}(b-x)f(x)=d>0$，那么${\int }_a^{b}f(x)dx$发散。

（注意：以上$c$均为常数，$d$可为常数也可取无穷）

上述两种审敛法的证明回头拿到同济版高数课本再补。

对于此题：$x\ge 0$，故$f(x)\ge 0$恒成立。

1、$x\to 0^{+}$：${\lim }_{x\to 0^{+}}f(x)=e^{{\lim }_{x\to 0^{+}}\frac{(a-x^2)lnx}{x}}=e^{{\lim }_{x\to 0^{+}}(a-x^2){\lim }_{x\to 0^{+}}\frac{lnx}{x}}$

$=e^{a{\lim }_{x\to 0^{+}}\frac{lnx}{x}}=$(1)$e^{+\infty }=+\infty ,a<0$

​ (2)$e^0=1,a=0$

​ (3)$e^{-\infty }=0,a>0$

所以只有$a<0$时${\int }_0^{\infty }f(x)dx$可能在$x=0$处发散。

又$\because$${\lim }_{x\to 0^{+}}(x-0)f(x)=e^{{\lim }_{x\to 0^{+}}\frac{(a-x^2+x)lnx}{x}}=e^{{\lim }_{x\to 0^{+}}lnx(1-x+\frac{a}{x})}=e^{+\infty }=+\infty$，故${\int }_0^{\infty }f(x)dx$在$x=0$处发散。

2、$x\to +\infty$：可利用洛必达法则求得：

${\lim }_{x\to +\infty }{x}^{2}f(x)=e^{{\lim }_{x\to +\infty }\frac{(a-x^2+2x)lnx}{x}}=e^{{\lim }_{x\to +\infty }[(-2x+2)lnx+\frac{a}{x}-x+2]}=e^{-\infty }=0$，故${\int }_0^{\infty }f(x)dx$在$x\to +\infty$处收敛。

$a\ge 0$情况代码实现

由于$a>0$时原积分在$x\to +\infty$处收敛且$f(10)=10^{\frac{a}{10}-10},f(20)=20^{\frac{a}{20}-20}$。即便$a$取到题目给定上限$50$，$f(10)=10^{-5}=0.0001,f(20)=20^{-15}$小的可以在保留$5$位小数情况下忽略不记，故自适应$Simpson$法的左右端点可取$0,20$（注意：右端点取值过大会导致初始状态用$Simpson$公式近似定积分值时得到$f(mid)$与$f(r)$等于$0$从而导致拟合误差较大）。同时，为了提高准确性，初始精度$eps$可适当调高。

总结

1、若写的自适应Simpson法实际运行的时候出了误差，可从修改初次调用Simpson函数时积分区间的端点$l$和$r$、长度$eps$的值入手；

2、求极限时：能求提则求提！！（再次重申：求提原理：若$f(a)、f(b)$在$x\to a$时极限均存在，则$li{m}_{x\to a}f(a)f(b)=li{m}_{x\to a}f(a)li{m}_{x\to a}f(b)$！）；

3、公式推导务求过程严谨。
Code

#include<cstdio>
#include<cmath>
#include<iostream>
using namespace std;

double A;

inline double dcmp(double x,double y) { return fabs(x-y)<=1e-6; }
inline double f(double x) 
{ 
	if(dcmp(x,0.0))
	{
		if(dcmp(A,0.0))	return 1.0;
		return 0.0;		
	}
	return exp((A-x*x)*log(x)/x); 
}
inline double Simpson(double a,double b) { return (f(a)+4*f((a+b)/2.0)+f(b))*(b-a)/6.0; }

inline double queryans(double l,double r,double eps,double ans)
{
	register double mid=(l+r)/2.0;
	double ansl=Simpson(l,mid),ansr=Simpson(mid,r);
	if(fabs(ansl+ansr-ans)<=15*eps)	return ansl+ansr+(ansl+ansr-ans)/15.0;
	return queryans(l,mid,eps/2.0,ansl)+queryans(mid,r,eps/2.0,ansr);
}

int main()
{
	scanf("%lf",&A);
	if(A<0.0)	{ cout<<"orz"; return 0; }
	printf("%.5lf",queryans(0,20,1e-7,Simpson(0,20)));
	return 0;
}