自适应辛普森积分 学习笔记

# 1. Simpson 公式 [百度百科上的解释](https://baike.baidu.com/item/%E8%BE%9B%E6%99%AE%E6%A3%AE%E7%A7%AF%E5%88%86%E6%B3%95/23337870) 原理是通过函数 $f(x)$ 在 $l$ , $r$ 和 $mid$ 上的值将 $f(x)$ 简化为一个二次函数进行处理。 于是就有了： $$\int_a^bf(x)dx\approx\frac{(b-a)(f(a)+4f(\frac{a+b}{2})+f(b))}{6}$$ 推导过程可以看[此博客](https://www.luogu.com.cn/blog/beili233/solution-p4525) # 2.Adaptive 主过程 是一个分治的过程，采用递归实现 大致思想是： - `左区间积分` + `右区间积分` 与 `总积分` 之差不大于eps，则直接返回 - 否则分治处理 具体看代码 ```cpp db adaptive(db l,db r,db pre) { db mid=(l+r)/2.0; db ansl=simpson(l,mid),ansr=simpson(mid,r); db ans2=ansl+ansr; if(abs(pre-ans2)<=eps)return ans2; return adaptive(l,mid,ansl)+adaptive(mid,r,ansr); } ``` # 3.模板 因为时间相隔久远，码风会有些不一样qwq * [P4525 【模板】自适应辛普森法1](https://www.luogu.com.cn/problem/P4526) ```cpp #include<bits/stdc++.h> using namespace std; #define fo(v,a,b) for(int v=a;v<=b;v++) #define fo2(v,a,b) for(v=a;v<=b;v++) #define fr(v,a,b) for(int v=a;v>=b;v--) #define fe(v,a) for(int v=head[a];v;v=Next[v]) #define rg register #define il inline typedef double db; const db eps=1e-8; db a,b,c,d,L,R; il db f(db x){ return (c*x+d)/(a*x+b); } il db simpson(db l,db r) { db mid=(l+r)*0.5; return (f(l)+4*f(mid)+f(r))*(r-l)/6; } db solve(db l,db r,db All) { db mid=(l+r)*0.5; db l_num=simpson(l,mid),r_num=simpson(mid,r); if(fabs(l_num+r_num-All)<eps)return l_num+r_num; return solve(l,mid,l_num)+solve(mid,r,r_num); } il db Int(){ return solve(L,R,simpson(L,R)); } int main() { cin>>a>>b>>c>>d>>L>>R; printf("%.6lf\n",Int()); return 0; } ``` * [P4526 【模板】自适应辛普森法2](https://www.luogu.com.cn/problem/P4526) ```cpp #include<bits/stdc++.h> using namespace std; #define fo(v,a,b) for(int v=a;v<=b;v++) #define fo2(v,a,b) for(v=a;v<=b;v++) #define fr(v,a,b) for(int v=a;v>=b;v--) #define rg register #define il inline typedef double db; const db eps=1e-9; db a; il db f(db x) { return pow(x,(a/x)-x); } db simpson(db l,db r) { db mid=(l+r)/2.0; return (f(l)+4.0*f(mid)+f(r))*(r-l)/6.0; } db adaptive(db l,db r,db pre) { db mid=(l+r)/2.0; db ansl=simpson(l,mid),ansr=simpson(mid,r); db ans2=ansl+ansr; if(abs(pre-ans2)<=eps)return ans2; return adaptive(l,mid,ansl)+adaptive(mid,r,ansr); } int main() { cin>>a; if(a<0) { puts("orz"); exit(0); } printf("%.5lf\n",adaptive(eps,15.0,simpson(eps,15.0))); return 0; } ```