龙贝格求积法源代码,效率很好的

maleesky

<>积分的源代码,效率很好的

9 N+ a/ R& A0 H8 S1 V
8 \% \& Q1 j/ U3 z7 x8 n//////////////////////////////////////////////////////////////////////
// 龙贝格求积法
//
// 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x)
7 A! I3 I6 b' _( q//
9 ]! F+ @& T) }) Z: j// 参数：
3 t0 _' K+ H0 Q% `( i// 1. a - Double型变量，积分下限
; I- e2 y9 Q, c3 F1 t// 2. b - Double型变量，积分上限，要求b&gt;a
! }( ^2 |: Q/ D// 3. eps - Double型变量，积分精度要求
//
$ x$ [5 ~, [- X2 @+ V# Y// 返回值：double 型，积分值
//////////////////////////////////////////////////////////////////////
" T+ o/ U0 z* H# v9 p/ u1 idouble Integral_Romberg(double a, double b, double eps /*= 0.000001*/)
! _. d+ j" y9 M" n{
    int m,n,i,k;
  H! r' D9 j9 F; S/ L    double y[10],h,ep,p,x,s,q;</P>
5 W7 B/ r  X1 t4 Z<> // 迭代初值
  m5 [: C0 [4 b6 E" p% g8 @3 o* j    h=b-a;
0 H0 F" K0 |/ G6 r+ g3 W  O* ~; E    y[0]=h*(Func(a)+Func(b))/2.0;
9 V( B* |2 S$ C" i$ V6 H# I    m=1;
* }/ R" S7 P4 _/ P5 [% q n=1;
% v3 _, ]( `' u( F ep=eps+1.0;
6 v& C+ a7 x2 W2 f+ C   
// 迭代计算
) H+ ^1 U9 t2 e while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
) m2 k( k* L6 N2 E    {
  p=0.0;
        for (i=0;i&lt;=n-1;i++)
        {
   x=a+(i+0.5)*h;
            p=p+Func(x);
        }
/ X& _7 U7 Z9 m/ _$ A9 }, L5 {        
, o  Z1 ^# g+ Q! j( w, M) x2 C5 h" @) u  p=(y[0]+h*p)/2.0;
        s=1.0;
5 U4 ~7 D4 {( O        for (k=1;k&lt;=m;k++)
+ h- n1 E; {8 C" }. W        {
   s=4.0*s;
            q=(s*p-y[k-1])/(s-1.0);
+ A* Z/ s, W, V( B2 N& Q  j; R            y[k-1]=p; p=q;
        }</P>
<>        ep=fabs(q-y[m-1]);
2 G  J9 r* A* j6 V5 K- ]6 s5 T        m=m+1;
  y[m-1]=q;
  n=n+n;
1 R' q: @0 _/ G% C& D6 G. ^) u5 |7 t* s  h=h/2.0;
( L2 |4 h% W/ {" @% l) c# O: M    }
   
return(q);
5 F) m, O8 n0 R$ u$ [+ s+ [}</P>