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

maleesky

<>积分的源代码,效率很好的
' U3 x5 F* k) b- R! K  d

$ l2 H$ O9 C1 {; K. o# x1 \//////////////////////////////////////////////////////////////////////
// 龙贝格求积法
//
0 j: m4 |6 W. k2 T/ p, D" |// 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x)
3 i! ~) H6 V( H5 c& H" Y, B//
5 o) k8 S5 Q9 `1 X% r// 参数：
$ k! o% o3 M, G// 1. a - Double型变量，积分下限
6 l3 {+ C. W6 ~$ B: p3 B// 2. b - Double型变量，积分上限，要求b&gt;a
// 3. eps - Double型变量，积分精度要求
//
( M" p2 t" {3 L! R  D// 返回值：double 型，积分值
//////////////////////////////////////////////////////////////////////
$ Y3 r9 H0 ?' f8 \0 Tdouble Integral_Romberg(double a, double b, double eps /*= 0.000001*/)
5 f7 w" K4 Y0 M' u{
    int m,n,i,k;
! }( X( s' N5 K' Y    double y[10],h,ep,p,x,s,q;</P>
9 \& V* A* \! I2 S, W" R9 m<> // 迭代初值
+ E0 R, L/ P2 w* g6 E    h=b-a;
    y[0]=h*(Func(a)+Func(b))/2.0;
; f3 o( _  v  o    m=1;
n=1;
: N0 n0 y+ U7 L6 H8 M7 m1 [ ep=eps+1.0;
   
// 迭代计算
, X; x9 u0 m: z$ f7 P; L" _) Q while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
. b2 r  t  n, c) X    {
+ f# @; `# N2 J& P& o4 z  p=0.0;
1 z; O+ D- g0 x3 ]/ {1 o9 ~+ _- ~        for (i=0;i&lt;=n-1;i++)
' }/ [- z/ }0 Z( {3 ~# v; s        {
   x=a+(i+0.5)*h;
            p=p+Func(x);
6 W# d1 ?! h; ?# v        }
4 @, d/ x0 v" f% ?) @        
  p=(y[0]+h*p)/2.0;
5 [9 ~' E! O. n+ |/ v  C6 v3 r) T        s=1.0;
        for (k=1;k&lt;=m;k++)
        {
. m9 m2 L4 L: ]$ g' }   s=4.0*s;
9 t* B2 _% t6 k/ Y0 c6 h            q=(s*p-y[k-1])/(s-1.0);
            y[k-1]=p; p=q;
        }</P>
<>        ep=fabs(q-y[m-1]);
. w( `6 @, ?* @8 h! J4 q        m=m+1;
; g' r8 s9 m" v) }1 k1 V  y[m-1]=q;
  n=n+n;
5 O0 F2 x+ d6 P2 M% A$ k  h=h/2.0;
    }
* E8 V# L% J- [5 H/ E3 u$ j; m3 Q1 }   
3 t( H* Q, i( x) I1 {9 p3 R, T return(q);
}</P>