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<bits/stdc++.h>
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
#define db double
db 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<bits/stdc++.h>
#define db double
db 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)<eps*15)return (left+right+(left+right-ans)/15.0);
    else return asr(l,mid,eps/2,left)+asr(mid,r,eps/2,right);
}
int main()
{
    scanf("%lf",&a);
    if(a<0)puts("orz");
    else printf("%.5lf\n",asr(1e-8,20,1e-8,simpson(1e-8,20)));
    return 0;
} 
 

Blog來自PaperCloud，未經允許，請勿轉載，TKS！

自適應辛普森法