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