龍貝格積分

**：

龍貝格積分對於程式設計實現來說,一開始不太好懂.

對於不易直接用積分公式計算的原函式，通常用

復合梯形求積公式或復合拋物線求積公式等方法，但這些方法精度不高，收斂的速度緩慢。為了提高收斂速度，減少計算量，人們尋求其他方法.

romberg方法也稱為逐次分半加速法。它是在梯形公式、辛卜生公式和柯特斯公式之間的關係的基礎上，構造出一種加速計算積分的方法。

作為一種外推算法, 它在不增加計算量的前提下提高了誤差的精度.

在等距基點的情況下，用計算機計算積分值通常都採用把區間逐次分半的方法進行。這樣，前一次分割得到的函式值在分半以後仍可被利用,且易於程式設計。

推導和證明資料：

以下是一些對理解romberg演算法很關鍵的資訊.

對區間[a， b]，令h=b-a構造梯形值序列。

t1=h[f(a)+f(b)]/2

把區間二等分，每個小區間長度為 h/2=(b-a)/2，於是

t2 =t1/2+[h/2]f(a+h/2）

把區間四(2

2)等分，每個小區間長度為h/2

2 =（b-a）/4，於是

t4 =t2/2+[h/2^2][f（a+h/4）+f（a+3h/4).

把[a，b]

2k等分，分點xi=a+（b-a）/ 2

k ·i （i =0，1，2 · · · 2k）每個小區間長度為（b-a）/ 2

k ，由歸納法可得下圖的第乙個公式.

整個程式就是循著這四個公式進行計算的.

sn,cn, rn

分別代表特例梯形積分,拋物線積分,龍貝格積分.當然,程式設計的時候統一處理即可.

應用中,算到rn級別即可(7階的精度), 精度再高則會與累計誤差相衝突.

下面舉乙個例,**執行一下：取e=0.0001,用龍貝格方法計算積分

i = ∫0

1( 4/(1+x2)) dx

解 按上述五步計算，此處

f(x)=4/(1+x2)

a=0 b=1

f(0)=4

f(1)=2

由梯形公式得

t1=1/2[f(0)+f(1)]=3

計算f(1/2)=16/5

用變步長梯形公式得

t2=1/2[t1+f(1/2)]=3.1

由加速公式得

s1=1/3(4t2-t1)=3.133333333

求出f(1/4)

f(3/4) 進而求得

t4=1/2

=3.131176471

s2=1/3(4t4-t2)=3.141568628

c1=1/15(16s2-s1)=3.142117648

計算f(1/8)

f(3/8)

f(5/8)

f(7/8)進而求得

t8=1/2

=3.138988495

s4=1/3(4t8-t4)=3.141592503

c2=1/15(16s4-s2)=3.141594095

r1=1/63(64c2-c1)=3.141585784

把區間再二分，重複上述步驟算得

t16=3.140941613

s8=3.141592652

c4=3.141592662

r2=3.141592640

由於 |r1-r2|<=0.00001,計算可停止，取r2=3.14159

題目：poj2369

數值積分 龍貝格 Romberg 積分

數值積分在工程上是個比較有用的數學工具。在工程上有很多數學問題，看似簡單，計算所用的數學公式不算複雜，但是求解起來卻很困難，很難獲得解析解的公式，這個時候就需要用到數值求解的辦法來獲取滿足工程需要的近似數值結果。舉個簡單的例子，求乙個變截面圓柱 或者叫圓台 的等效截面慣性矩，計算公式為ie h 3 ...

龍貝格演算法

龍貝格演算法是以復化梯形公式，復化希普森，復化牛頓 科特斯為基礎採用的是逐漸平衡誤差的做法，如下 include 龍貝格 include define maxlen 100 double a double b double f double x else double t int n else el...

用C 寫的龍貝格（Romberg）積分法

很久之前了，群上有人問有沒有用c 寫的積分函式，因為我自己以前做過乙個積分計算器，就跟他說搜尋一下romberg積分法吧，但是後來他說網上寫的都是用c 寫的方法，沒有c 的，問我有沒有，但是我之前是在家裡做的，學校裡沒有副本，可惜。今天翻一下郵箱，發現居然是我忘記貼標籤，導致當時沒有找到，感覺有點對...