自适应辛普森

辛普森定理:

对于一段函数f在[l,r]范围内的面积∫(r,l)f(x)dx,我们可以由辛普森定理得到下列式子(∫是积分的意思)。

即:∫(r,l)f(x)dx=(r-l)/6*(f(l)+4f((l+r)/2)+f(r))。

但问题是这样做是把f(x)在[l,r]范围内看做一个二次函数来算的,所以精度并不高。因此,我们需要一种可以在较复杂函数中正确算出面积的自适应希普森。

自适应辛普森其实就是二分递归用辛普森定理来求面积。

S(l,r)=S(l,m)+S(m,r)

首先我们先算出S(l,r),S(l,m),S(m,r)

如果|S(l-r)-S(l,m)-S(m,r)|<你需要的精度,就返回S(l,r),否则返回F(l,m)+F(m,r)。

还蛮好理解的,下面是代码:

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define ll long long
#define il inline
#define db double
using namespace std;
il db f(db x)
{
    return (f(x));//f函数 
}
il db simpson(db l,db r)
{
    db m=(l+r)/2;
    return (f(l)+f(r)+4*f(m))*(r-l)/6;//辛普森定理 
}
db asr(db l,db r,db eps,db sum)//eps是精度,sum是simpson(l,r) 
{
    db m=(l+r)/2;
    db L=simpson(l,m),R=simpson(m,r);//二分的左边和右边 
    if(fabs(sum-L-R)<=15*eps)//在精度范围内 
    return L+R+(L+R-sum)/15.0;
    else
    return asr(l,m,eps/2,L)+asr(m,r,eps/2,R);
}
int main()
{
    return 0;
}

 

posted @ 2017-07-31 15:28  GSHDYJZ  阅读(284)  评论(0编辑  收藏  举报