XX大
學
綜
合設
計
報
告
綜合設計五
多方法求解數值積分
學生姓名:
學
號:
年級專業:
指導老師:
學
院:
評閱成績:
評閱意見:
成績評定教師簽名:
時間:
提交日期:2014年X月
多方法求解數值積分
具體題目要求:用不同數值方法計算積分
(1)
取不同的步長,分別用復合梯形及復合辛普森公式計算積分,給出誤差中關于的函數,并與積分精確值比較兩個公式的精度,是否存在一個最小的,使得精度不能再被改善?
(2)
用龍貝格求積計算完成問題(1);
(3)
用自適應辛普森積分,使其精度達到。
1設計目的、要求
由積分學基本理論,定積分可由公式計算,但是對于一些無法找到原函數的函數(如等)不能通過牛頓—萊布尼茲公式計算,就必須得另尋它法。因此需要我們能夠熟練地應用常用的數值積分計算方法(如機械求積、公式等)并掌握結合數值計算軟件(、等)及計算機高級語言進行對應算法實現的技能。
熟練數學軟件求解數學問題,掌握各種數學問題的求解方法。本設計主要是通過多種復合求積公式求解積分,主要包括復化梯度法、復化辛普森法、龍貝格以及自適應辛普森法等求解方法,利用軟件編寫相對應的算法進行求解,大大地提高了解題的速度。
2設計原理
由積分中值定理我們可以知道在積分區間內存在一點,使得式子
成立。這個式子在于對于點的具體位置一般是不知道的,因此難以準確算出的值。也就是不同算法求得平均高度,對應的就是一種不同的數值求積方法。更一般地,我們可以在區間上適當選取某些節點,然后用的加權平均得到平均高度的近似值,這樣構造出的求積公式具有下列形式:
稱為機械求積公式。
復合梯形公式、復合辛普森公式、龍貝格求積公式以及自適應辛普森公式都以此公式的基礎,對積分區間進行變步長的劃分求得近似的平均高度值,得到積分函數的近似值。也由于牛頓—柯特斯公式在時不具有穩定性,所以不可能通過提高階的方法來提高求積精度。為了我提高精度通常可把積分區間等分成為若干個子區間,再在每個子區間上用低價求積公式,這就是復合求積方法。但是這樣的積分求解方法也是存在不容忽視的誤差。因此需要在設計算法時考慮到算法存在的誤差(舍入誤差、截斷誤差等),并對誤差作出分析。
3采用軟件、設備
軟件
4設計內容
第一步:復合梯形公式、復合辛普森公式算法
(一)、復合梯形公式計算積分
復化梯形公式的主要思想是利用若干小梯形的面積代替原方程的積分,利用微元法,可以求出坐標面上由函數與坐標軸圍城的圖像的面積的近似值,符合了計算機計算存儲的思想。下面,我們在探討復化梯形公式的計算規律:
設將求積區間分成等份,則一共有個分點,按梯形公式
計算積分值,需要提供個函數值。
這里代表步長,分點為,其中
(二)、復合辛普森公式計算積分
算法的基本思想是:把積分區間等分成若干個子區間,而在每一個子區間上用辛普森
求積公式:
得到復合辛普森求積公式:
并且用軟件來求解。
第二步:龍貝格算法
考慮積分,欲求其近似值,通常有復化的梯形公式、公式和公式。但是給定一個精度,這些公式達到要求的速度很緩慢。如何提高收斂速度,自然是人們極為關心的課題。為此,記,為將區間進行等分的復化的梯形公式計算結果,記,為將區間進行等分的復化的公式計算結果,記,為將區間進行等分的復化的公式計算結果。根據外推加速方法,可以得到收斂速度較快的積分法。其具體的計算公式為:
1、準備初值,計算
2、按梯形公式的遞推關系,計算
3、按Romberg積分公式計算加速值,4、精度控制。對給定的精度,若
則終止計算,并取為所求結果;否則返回2重復計算,直至滿足要求的精度為止。
第三步:自適應辛普森算法
復合求積方法通常適用于被積函數變化不太大的積分,如果在積分區間被積函數變化很大有的部分函數值變化劇烈而有的部分則是變化平緩,如果此時還是將積分區間等分后用復合求積公式的話計算量很大。而采用針對被積函數在區間上的不同情形而采用不同的步長,使得在滿足精度的前提下積分計算量減少,這就是自適應積分方法,能自動地在被積函數變化劇烈的區域增多節點,而在被積函數變化平緩的地方減少節點。因此它是一種不均勻區間的積分方法。題目要求使相鄰兩個區間的誤差達到一定的要求,即用自適應辛普森公式來求積分,先算出積分區間的左右端點函數值,求出區間中點函數值與左右端點的函數差值,再與所要求的精度比較,不滿足的對所在區間二等分,接著算出每個子區間端點的函數值判斷時否符合精度要求,直到積分每個子區間內都滿足精度要求,最后所得各個區間端點的函數值之和即為積分的近似值。
第四步:誤差余項以及精度分析:
由插值型的求積公式我們得到求積公式誤差余項的表達式:
其中表示求積公式的代數精度,為不依賴于的待定系數,這個結果表明當是次數小于等于的多項式時由于,故此時,即前面的求積公式精確成立。而當時,故可求得:
因為梯形公式的代數精度為1,可以得到的值為
于是得到梯形公式的余項為:
又因為復合梯形公式要滿足
綜上所述,就得到了復合梯形公式的余項表達式:
同理可得復合辛普森公式的余項表達式:
結果分析:從以上余項的表達式可以看出復合辛普森公式的代數精度為3,而復合梯形公式的代數精度為1,所以復合辛普森比復合梯形精確度更高。對于算法的精度,是通過對設計所得值與準確值之間的誤差值來評判。將變步長的復合求積方法每次求得計算結果與準確值進行比較求出誤差值,通過畫出誤差值的變化趨勢圖比較復合梯形公式與復合辛普森公式這兩種算法的精度。經過實驗中驗證,也表明自己的初步推理是正確的,無論是復合梯形公式還是復合辛普森公式它們最終結果都會隨著步長
值的減小而更加精確。復合梯形公式和復合辛普森公式計算出的結果進行比較,發現復合辛普森公式計算出的結果更加的精確。
5原始程序、數據
文件:f.m
function
y=f(x);
y=sqrt(x)*log(x);
1、復合梯形公式求解算法:
文件:trapezoid.m
clc
a=0;
%積分下限
b=1;
%積分上限
T=[];
%用來裝不同n值所計算出的結果
R=[];
G=[];
m=120;
%等分數
true=-(4/9);
for
n=2:m;
h=(b-a)/n;
%步長
x=zeros(1,n+1);
%給節點定初值
y=zeros(1,n+1);
for
i=1:n+1
x(i)=a+(i-1)*h;
end
x(1,1)=0.000000001;
for
i=1:n+1
y(i)=x(i).^(1/2)*log(x(i))
%
g=(-(b-a)/12*h^
2)*(-log(x(i))/(4*x(i)*x(i)^(1/2)))
%準確的積分余項(計算誤差)
end
%
G=[G,g];
t=0;
r=0;
for
i=1:n
format
long
t=t+h/2*(y(i)+y(i+1))
;
%利用復化梯形公式求值
err=t-floor(t);
digits(7);
%
此處為需要的小數位+1
t=floor(t)+vpa(err,6)
;
%
此處控制顯示的小數點位數,更改顯示的小數位數
r=t-true;
%計算的值與真實值之差(實際誤差)
end
T=[T,t]
;
%把不同n值所計算出的結果裝入
T中
R=[R,r];
end
x=linspace(0,1,m-1);
plot(x,R,'*')
%將計算誤差與實際誤差用圖像畫出來
2、復合辛普森積分求解算法:
simpon.m
clc
clear
a=0;
%積分下限
b=1;
%積分上限
T=[];
%用來裝不同n值所計算出的結果
R=[];
true=-(4/9);
m=20;
%等分數
for
n=2:m
h=(b-a)/(2*n);
%步長
x=zeros(1,2*n+1);
%給節點定初值
y=zeros(1,2*n+1);
for
i=1:2*n+1
x(i)=a+(i-1)*h;
%給節點賦值
end
x(1,1)=0.000000001;
for
i=1:2*n+1
y(i)=x(i).^(1/2)*log(x(i));
%給相應節點處的函數值賦值
end
t=0;
r=0;
for
i=1:n
format
long
t=t+h/3*(y(2*i-1)+4*y(2*i)+y(2*i+1));
%利用復化simpson公式求值
err=t-floor(t)
digits(7);
%
此處為需要的小數位+1
t=floor(t)+vpa(err,6);
r=t-true;
end
T=[T,t]
%把不同n值所計算出的結果裝入
T中
R=[T,r];
end
%
R=(-(b-a)/180*((b-a)/2).^4*24)
%積分余項(計算誤差)
%
true=quad(@fx1,0,1);
%積分的真實值
x=linspace(0,1,2*m-1);
plot(x,R,'*')
3、龍貝格算法
rebeg.m
%%龍貝格
clear
clc
a=0;
b=1;
%確定積分上下限
eps=10^(-4);
err=1;
k=1;
a=0.0000001;
T(1,1)=(b-a)/2*(f(a)+f(b));
while(err>eps)
h=(b-a)/2^(k-1);
S=0;
for
x=a:h:b-h
S=S+f(x+h/2);
end
T(k+1,1)=1/2*T(k,1)+h/2*S;
k=k+1;
for
i=2:k
T(k,i)=(4^(k-1)*T(k,i-1)+T(k-1,i-1))/(4^(k-1)-1);
end
err=abs(T(i,i)+4/9);
end
fprintf('龍貝格求積算法積分值為%.10f\n',T(k,k));
disp(T)
T
龍貝格求積算法積分值為-0.44438207534、自適應辛普森算法:
%自適應辛普森算法
Self_Adaptive_integral.m
function
s=Self_Adaptive_integral(a,b,tol)
k=0;
w=0;
x=a;
y=b;
t=0;
h=(b-a)/2;
s=0;
i=0;
to=abs(simpson_integral(x,y,2)-simpson_integral(x,y,1));
while
to>=tol
i=i+1;
while
to>=tol
t=x;
if
k==0
x=t;
y=t+h;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
k=1;
w=0;
end
if
w==0
x=t+h;
y=t+2*h;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
k=0;
w=1;
end
end
s=s+simpson_integral(x,y,2);
if
k==0
x=t;
y=t+h;
h=h/2;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
end
if
w==0
x=t+h;
y=t+2*h;
h=h/2;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
end
if
to s=s+simpson_integral(x,y,2); to=tol-10; end end %自適應辛普森算法 simpson_integral.m function s=simpson_integral(a,b,m) h=(b-a)/(2*m); s1=0; s2=0; s=0; if m>1 for i=1:(m-1) x=a+2*i*h; s1=s1+f(x); end for i=1:m x=a+(2*i-1)*h; s2=s2+f(x); end s=h/3*(f(a)+f(b)+2*s1+4*s2); else s=s+h/3*(f(a)+f(b)++4*f((a+b)/2)); end 6結果分析與設計總結 結果分析: 初步分析:通過對步長h的值的改變,只要h值越小(等分數n的值越大),即等分的區間越小,結果應該更加精確,精確度越高。實驗結果分析: 1、復合梯度算法: 通過算法的運行結果可得:當等分數n從2開始變化到50時實驗計算結果以及與準確值之間的誤差可以達到-0.4410968和0.003347665; 當等分數更改為80時實驗計算結果以及與準確值之間的誤差可以達到-0.4426581和0.001786344; 結果分析:復合梯形求積公式隨著區間數不斷增加,積分的誤差不斷減小。運算所得結果如下圖所示: 2、復合辛普森算法: 當等分區間數目達到50時,實驗計算出來的結果以及與準確值(-4/9)之間的誤差值分別為:-0.443796和0.0006484617; 而當n為80時實驗結果分別為-0.4441055和0.0003389765; 結果分析:復合辛普森求積公式也是隨著等分數不斷增加,積分的誤差不斷減小。算法計算結果如圖所示: 將這兩種插值形的積分算法所得結果與準確值比較,得出兩種方法產生的誤差趨勢圖如下: 經過算法設計驗證,也表明自己的初步推理是正確的,無論是復合梯形公式還是復合辛普森公式計算的最終結果都會隨著步長值的減小而更加精確,更加趨近于準確值。復合梯形公式和復合辛普森公式計算出的結果進行比較,發現復合辛普森公式計算出的結果更加的精確。 3、龍貝格算法: 通過龍貝格算法,當要求積分計算結果達到精度為時:實驗所得結果為-0.***,而當要求的計算精度達到時,計算所得結果為-0.***。 計算結果為: k h …… 0 …… 0 0 …… …… 0 …… …… -0.444***4 0 …… …… -0.*** -0.*** 0 -0.*** -0.*** -0.*** 4、自適應辛普森算法: 輸入s1=Self_Adaptive_integral(0.0000001,1,0.01)時 運算所得結果為-0.***; 輸入s2=Self_Adaptive_integral(0.0000001,1,0.0001)時,運算結果也是-0.***; 設計總結: 雖然此次課程設計時間不是很長,但是還是讓我學會了不少。不僅是運籌學知識的應用還是對于數值分析中多種數值計算方法的回顧,都讓我對于專業知識得到進一步地加深理解。本課程設計也讓我更加熟練地掌握了應用MATLAB編寫相應的算法求解相應的數學問題,將理論知識與實際應用想結合,提高了自身的算法設計能力以及編程程序的技能。這一過程也得利于到老師、同學以及組員的幫助,才能如期完成自己的任務。這一次的課程設計更讓學會對問題的分析以及思考,以及查找算法中的不足并作出改進。 參考文獻 [1] 李慶揚,王能超,易大義.數值分析第四版[M].北京:清華大學出版社.2001.[2] 董霖.MATLAB使用詳解第一版[M].科學出版社.2008 [3] 龔純,王正林.MATLAB語言常用算法程序集[M].電子工業出版社.2008 [4] 李慶楊.科學計算方法基礎[M].北京:清華大學出版社,2006 [5] 白峰杉.數值計算引論[M].北京:高等教育出版社,2004