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

maleesky

<>积分的源代码,效率很好的
4 c* U. r6 M$ z- Y* v
  U. ]- y# q, y4 R& l/ Q
3 C3 B/ J# k- p+ F//////////////////////////////////////////////////////////////////////
+ R& m7 Z/ P& k' D, E// 龙贝格求积法
//
+ k: i" s  e) Z$ }) i2 m// 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x)
//
  ?) J& ]5 A- L1 \8 f0 b5 A// 参数：
  t7 v1 b9 O' s) J! M% l8 s0 H6 E3 M// 1. a - Double型变量，积分下限
// 2. b - Double型变量，积分上限，要求b&gt;a
// 3. eps - Double型变量，积分精度要求
//
: p9 g/ H( W( K- G// 返回值：double 型，积分值
//////////////////////////////////////////////////////////////////////
double Integral_Romberg(double a, double b, double eps /*= 0.000001*/)
/ H! m& S; R6 V- f{
    int m,n,i,k;
    double y[10],h,ep,p,x,s,q;</P>
# M/ ?! M8 x) e$ Y$ R6 }( U/ N9 q% W<> // 迭代初值
5 ]" V4 ?) d- I5 I& V5 q+ q    h=b-a;
    y[0]=h*(Func(a)+Func(b))/2.0;
    m=1;
- a% ]$ @6 [( c  p: ? n=1;
* G9 ]6 @) i  K3 B. O ep=eps+1.0;
( W, F) P/ f$ [7 w/ c4 z   
// 迭代计算
- m" ?( Y* n2 z6 F while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
& w3 ^* ]$ X3 W    {
  p=0.0;
3 n  u8 J( T3 N6 _1 |& |        for (i=0;i&lt;=n-1;i++)
        {
   x=a+(i+0.5)*h;
3 t7 l4 }' u& L            p=p+Func(x);
        }
; m9 D6 E4 [! p3 J3 k        
  p=(y[0]+h*p)/2.0;
3 ^9 p) W  `) _5 c7 k. H1 {! R        s=1.0;
* `! U/ _8 V' w        for (k=1;k&lt;=m;k++)
' w, e) v$ x4 b3 g$ q9 H        {
8 }' B& I8 n5 r# F# z3 s/ S   s=4.0*s;
            q=(s*p-y[k-1])/(s-1.0);
            y[k-1]=p; p=q;
        }</P>
<>        ep=fabs(q-y[m-1]);
$ W) D4 S" R& A9 ~9 Z( i% P: ?        m=m+1;
) B7 [9 J0 E6 |6 L  y[m-1]=q;
  n=n+n;
  h=h/2.0;
    }
; k* G# g2 f! j" E   
  T/ L9 V5 m' P; p8 S8 x* @" \ return(q);
% T: M; L9 j# a1 r' o}</P>