自适应辛普森（Simpson）积分

最新推荐文章于 2025-02-04 10:13:42 发布

徒手装机甲

最新推荐文章于 2025-02-04 10:13:42 发布

阅读量1.4k

点赞数

分类专栏：数论文章标签：算法

本文链接：https://blog.csdn.net/weixin_45757507/article/details/127349665

版权

数论专栏收录该内容

11 篇文章

订阅专栏

本文介绍了一种数值积分方法——自适应辛普森积分法。该方法通过将待积分区间细分为多个小段并利用二次多项式逼近，进而计算积分的近似值。为提高精度，采用递归分治策略，确保误差在允许范围内。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

自适应辛普森（Simpson）积分

很多时候，我们会面临一些求积分的问题，无论是直接给你函数，让你想办法求解积分，还是对于一些计算几何问题，无法直接推导积分，我们都可以用这种方法来求一段区域的积分，积分的相关基础概念这里就不再赘述，今天主要就是说明他的大致原理，和他的用法。

他的本质就是把函数看作一个二次函数，如果区间够小，那么函数就足够近似，我们可以直接按照二次函数来计算每个小段的积分。

推导其实就是最基本的积分过程

$\int_{a}^{b} f(x) dx=\int_{a}^{b}(Ax^2+Bx+C)dx$

$=\frac{A}{3}(b^3-a^3)+\frac{B}{2}(b^2-a^2)+C(b-a)$

$=\frac{2A(b^3-a^3)+3B(b^2-a^2)+6C(b-a)}{6}$

$=\frac{(b-a)(Aa^2+Ba+C+Ab^2+Bb+C+Aa^2+2Aab+Ab^2+2Bb+2Ba+4C)}{6}$

$=\frac{(b-a)(f(a)+f(b)+A(a+b)^2+2B(a+b)+4C)}{6}$

$=\frac{(b-a)(f(a)+f(b)+4(A(\frac{(a+b)}{2})^2+B(\frac{a+b}{2})+C)}{6}$

$=\frac{(b-a)(f(a)+f(b)+4f(\frac{a+b}{2})}{6}$

于是我们就得到了一个近似的辛普森公式。

b-a越小，他就越接近正确答案。

一个普通的思路就是，把要求的积分区间分解成一个个小的区间，然后累加在一起，但是这样做的缺点就是，如果为了保证精度，每个区间都会很小，所以这样就很容易会TLE。那么如何控制这个区间范围比较好呢？

其实就是把问题不断分解，比如求一个区间（l，r）的积分。

先取一个mid=（l+r）/2，在对两部分分别计算辛普森积分，并且对整个区间计算辛普森积分，然后判断他们的误差是否在精度要求内，如果在就直接返回，否则继续递归下去，很朴素的一个分治思想。代码也很短。

#include<bits/stdc++.h>
using namespace std;
using ll=long long;
double a,b,c,d,l,r;
const double eps=1e-8;
double f(double x)//原函数
{  
	return  ;       
}
double simpson(double l,double r) //Simpson公式
{ 
	double mid=(l+r)/2;
	return (f(l)+4*f(mid)+f(r))*(r-l)/6;
}
double solve(double l,double r,double ans) 
{
	double mid=(l+r)/2;
	double ll=simpson(l,mid),rr=simpson(mid,r);
	if(fabs(ll+rr-ans)<=eps) return ll+rr;     //确认精度
	return solve(l,mid,ll)+solve(mid,r,rr);     //精度不够则递归调用
}
int main()
{
    cin.tie(nullptr)->sync_with_stdio(false);
    //cin>>a>>b>>c>>d>>l>>r;
    cout<<fixed<<setprecision(6)<<solve(l,r,simpson(l,r))<<"\n";
    return 0;
}

不过需要注意的是，在实际应用是，要灵活调整精度，否则可能在wa和tle之间反复徘徊。

这里在介绍一个优化技巧。

double solve(double l,double r,double eps,double ans) 
{
	double mid=(l+r)/2;
	double ll=simpson(l,mid),rr=simpson(mid,r);
	if(fabs(ll+rr-ans)<=15*eps) return ll+rr+(ll+rr-ans)/15.0;     //确认精度
	return solve(l,mid,eps/2,ll)+solve(mid,r,eps/2,rr);     //精度不够则递归调用
}

这样对于每个区间，我都可以保证，整个区间的精度可以继续保持。

至于为啥要弄15，这里有个论文给了解释，感兴趣可以了解，不感兴趣，存个板子直接用就好了。

AdaptiveQuadProof.pdf