博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
自适应辛普森法
阅读量:5860 次
发布时间:2019-06-19

本文共 1639 字,大约阅读时间需要 5 分钟。

Simpson公式

\(f(x)\)为原函数,$g(x)=Ax2+Bx+C $ 为拟合后的函数,则有:

\[ \int_{a}^{b}f(x)dx \approx \int_{a}^{b}Ax^2+Bx+C = \frac{A}{3}(b^3-a^3)+\frac{B}{2}(b^2-a^2)+C(a-b) =\frac{(b-a)}{6}(f(a)+f(b)+4f(\frac{a+b}{2})) \]
然后就得到了Simpson公式
\[ \int_{a}^{b}f(x)dx \approx \frac{(b-a)}{6}(f(a)+f(b)+4f(\frac{a+b}{2})) \]
所以

inline double simpson(double l,double r){    double mid=(l+r)/2.;    return (r-l)*(f(l)+f(r)+4.*f(mid))/6.;}

洛谷P4525 【模板】自适应辛普森法1

Solution

考虑二分区间,当精度满足时,终止递归

Code 

#include
#define max(a,b) ((a)>(b)?(a):(b))#define min(a,b) ((a)<(b)?(a):(b))#define db doubledb a,b,c,d;inline db f(db x){return (c*x+d)/(a*x+b);}inline db simpson(db l,db r){return (f(l)+4.*f((l+r)/2.)+f(r))*(r-l)/6.;}inline db asr(db l,db r,db ans,db eps){ db mid=(l+r)/2,L=simpson(l,mid),R=simpson(mid,r); if(std::fabs(L+R-ans)<=14*eps) return L+R+(L+R-ans)/14; return asr(l,mid,L,eps/2)+asr(mid,r,R,eps/2);}int main(){ db L,R; scanf("%lf%lf%lf%lf%lf%lf",&a,&b,&c,&d,&L,&R); printf("%.6lf",asr(L,R,simpson(L,R),1e-6)); return 0;}

洛谷P4526 【模板】自适应辛普森法2

Solution

  • \(a<0\)时原积分发散

    \[ f(x)=x^{\frac{a-x^2}{x}}=\frac{1}{x^{\frac{x^2-a}{x}}} \]
    \(x\)趋近于0时,\(\frac{x^2-a}{x}\) 趋近于\(\infty\),\(f(x)\)趋近于\(\infty\)

  • \(a>0\)时原积分收敛

不再赘述。。。

然后我们假装是\(( eps , 20 )\)的定积分

Code 

#include
#define db doubledb a;inline db f(db x){return pow(x,a/x-x);}inline db simpson(db l,db r){return (f(l)+f(r)+4.0*f((l+r)/2.))*(r-l)/6.0;}inline db asr(db l,db r,db eps,db ans){ db mid=(l+r)/2.0,left=simpson(l,mid),right=simpson(mid,r); if(std::fabs(left+right-ans)


Blog来自PaperCloud,未经允许,请勿转载,TKS!

转载于:https://www.cnblogs.com/PaperCloud/p/10201156.html

你可能感兴趣的文章
[na]office 2010 2013卸载工具
查看>>
Dynamic Performance Tables not accessible Automatic Statistics Disabled for this session
查看>>
Linux中使用vim乱码
查看>>
MR1和MR2的工作原理
查看>>
Eclipse中修改代码格式
查看>>
关于Keytool创建服务器自签名证书
查看>>
GRUB Legacy
查看>>
iOS开发之常用的那些工具类和方法
查看>>
关于 error: LINK1123: failure during conversion to COFF: file invalid or corrupt 错误的解决方案...
查看>>
linix下用keepalived搭建高可用myqsl-ha
查看>>
我的友情链接
查看>>
hexo博客解决不蒜子统计无法显示问题
查看>>
python实现链表
查看>>
java查找string1和string2是不是含有相同的字母种类和数量(string1是否是string2的重新组合)...
查看>>
Android TabActivity使用方法
查看>>
java ShutdownHook介绍与使用
查看>>
在Deferred框架中轻松实现Decal
查看>>
Eclipse的 window-->preferences里面没有Android选项
查看>>
《麦田里的守望者》--[美]杰罗姆·大卫·塞林格
查看>>
[置顶] 深入探析Java线程锁机制
查看>>