1、西南交通大學數值分析2015上機實習報告2015年11月 11 數值分析2015上機實習報告目 錄第2題11. 程序12. 結果分析2第3題51. 程序62. 結果分析7第5題81. 程序82. 結果分析9第2題2. 某過程測涉及兩變量x 和y, 擬分別用插值多項式和多項式擬合給出其對應規律的近似多項式，已知xi與yi之間的對應數據如下，xi=1,2,10yi = 34.6588 40.3719 14.6448 -14.2721 -13.3570 24.8234 75.2795 103.5743 97.4847 78.2392（1）請用次數分別為3，4，5，6的多項式擬合并給出最好近似結果f(

2、x)。（2）請用插值多項式給出最好近似結果下列數據為另外的對照記錄，它們可以作為近似函數的評價參考數據。xi = Columns 1 through 7 1.5000 1.9000 2.3000 2.7000 3.1000 3.5000 3.9000 Columns 8 through 14 4.3000 4.7000 5.1000 5.5000 5.9000 6.3000 6.7000 Columns 15 through 17 7.1000 7.5000 7.9000yi = Columns 1 through 7 42.1498 41.4620 35.1182 24.3852 11.27

3、32 -1.7813 -12.3006 Columns 8 through 14 -18.1566 -17.9069 -11.0226 2.0284 19.8549 40.3626 61.0840 Columns 15 through 17 79.5688 93.7700 102.36771. 程序（1）多項式擬合程序n=input('輸入所要擬合的階數n=');hold on;x=1:10;y=34.6588 40.3719 14.6448 -14.2721 -13.3570 24.8234 75.2795 103.5743 97.4847 78.2392;P=polyfit

4、(x,y,n)xi=1:.2:10; yi=polyval(P,xi); plot(xi,yi,x,y,'r*');legend('擬合曲線','原始數據')（2）拉格朗日插值多項式擬合程序clc;x=1:10;y=34.6588 40.3719 14.6448 -14.2721 -13.3570 24.8234 75.2795 103.5743 97.4847 78.2392;plot(x,y,'r*');hold on; syms t; n=length(x); f = 0.0; for i =1:n l = y(i); fo

5、r j = 1:i-1 l = l*(t-x(j)/(x(i)-x(j); end; for j = i+1:n l = l*(t-x(j)/(x(i)-x(j); end; f = f + l; simplify(f); f = collect(f); f = vpa(f,6); end ti=1.0:0.4:10; f=subs(f,'t',ti); plot(ti,f)legend('擬合曲線','原始數據')title '拉格朗日插值多項式擬合' 2. 結果分析（1）請用次數分別為3，4，5，6的多項式擬合并給出最好近似結

6、果f(x)。圖2-1 3次多項式擬合結果 圖2-2 4次多項式擬合結果圖2-3 5次多項式擬合結果圖2-4 6次多項式擬合結果從繪制的圖形來看，當采用6次多項式擬合的時候，擬合的曲線已經與所給出的點非常逼近了。6次多項式擬合曲線為：f(x)=0.0194x6-0.5408x5+5.1137x4-16.8973x3-0.8670x2+66.3750x-18.6991（2）拉格朗日插值多項式給出最好近似結果圖1-4 拉格朗日插值多項式擬合第3題3用雅格比法與高斯賽德爾迭代法解下列方程組Ax=b，研究其收斂性，上機驗證理論分析是否正確，比較它們的收斂速度，觀察右端項對迭代收斂有無影響。(1)A行分別

7、為A1=6,2,-1,A2=1,4,-2,A3=-3,1,4； b1=-3,2,4T, b2=100,-200,345T,(2) A行分別為A1=1,0,8,0.8,A2=0.8,1,0.8,A3=0.8,0.8,1；b1=3,2,1 T, b2=5,0,-10T,(3)A行分別為A1=1,3,A2=-7,1；b=4,6T,1. 程序（1）雅格比法程序function jacobi( )clc;clear;A=6 2 -1;1 4 -2;-3 1 4;B=-3;2;4;Err_user=0.01;N=500;m,n=size(A);X=zeros(n,1);k=1;while k<=N;

8、 Xk=X; for i=1:n for j=1:n if i=j AX(j)=A(i,j)*Xk(j); end end Sum_AX=sum(AX); AX=0; X(i)=(B(i)-Sum_AX)/A(i,i); end E=max(abs(Xk-X); if E<Err_user break; end k=k+1;enddisp(X); %顯示迭代結果disp(k); %顯示迭代次數（2）高斯賽德爾迭代法程序function GS( )clc;clear;A=6 2 -1;1 4 -2;-3 1 4;B=-3;2;4;A=1 0.8 0.8;0.8 1 0.8;0.8 0.8

9、1;B1=3;2;1;B2=5;0;-10;Err_user=0.0001;N=50;m,n=size(A);X=zeros(n,1);k=1;while k<=NXk=X; for i=1:n for j=1:n if i=j AX(j)=A(i,j)*X(j); %Óë Jacobi ·¨Ö÷ÒªÇø±ð end endsum_AX=sum(AX);AX=0;X(i)=(B1(i)-sum_AX)/A(i,i); endEr=max(abs(Xk-X); if Er

10、<Err_user break; endk=k+1;enddisp(X);disp(k); end2. 結果分析（1）1）A行分別為A1=6,2,-1,A2=1,4,-2,A3=-3,1,4； b1=-3,2,4T雅克比迭代：計算結果 x=-0.6363 0.5960 0.3737T,.迭代次數16次。高斯賽德迭代：計算結果 x=-0.6363 0.5959 0.3738T,.迭代次數10次。2）A行分別為A1=1,0,8,0.8,A2=0.8,1,0.8,A3=0.8,0.8,1；b1=3,2,1 T 雅克比迭代：計算不收斂高斯賽德迭代：計算結果 x=5.7690 0.7694 -4.

11、2307T,.迭代次數30次。（2）A行分別為A1=1,3,A2=-7,1；b1=4,6T。雅克比迭代：計算結果不收斂高斯賽德迭代：計算結果不收斂通過對雅克比法和高斯法的上機編程實習，分析對比實驗結果可得：在方程組Ax=b中，右端項對迭代收斂是有影響的，即當b增大時，迭代次數增加，收斂速度降低。并且通過對比可知，在相同條件下，高斯賽德爾迭代法比雅克比迭代法收斂速度快；方法的選擇也很重要，比如第二問，用雅克比迭代法是發散的，而用高斯迭代法則是收斂的。通過上機驗證，我們也得出結論：理論分析是正確的，即當迭代矩陣的譜半徑小于1時，迭代法是收斂的，而當迭代矩陣譜半徑大于1的時，迭代法都是發散的。第5題

12、5. 用Runge-Kutta 4階算法對初值問題y/=-20*y，y(0)=1按不同步長求解，用于觀察穩定區間的作用，推薦兩種步長h=0.1,0.2。注：此方程的精確解為：y=e-20x1. 程序%Runge-Kutta4階算法%f=-20y%y（0）=1%clc;clear;N=10;%設定節點個數h=0.05;%設定步長x(1)=0;%x0=0y(1)=1;%y（0）=1yy(1)=exp(-20*x(1);%y（0）的精確解%開始用runge-kutta法計算for i = 2:N K1 = -20*(y(i-1);%(xi,yi)點的導數為f=-20*y K2 = -20*(y(i-

13、1)+K1*h/2); K3 = -20*(y(i-1)+K2*h/2); K4 = -20*(y(i-1)+K3*h); delta = h*(K1 + 2*K2 + 2*K3 + K4)/6; y(i) = y(i-1) + delta;%計算y（i）值 x(i) = x(i-1) + h;%計算下個節點的x（i)值 yy(i) = exp(-20*x(i);%計算y(i)的精確值end;%存儲計算結果dlmwrite('f:計算結果.xls', x, 'delimiter', 't', 'precision', 8);dlm

14、write('f:計算結果.xls', y,'-append', 'delimiter', 't', 'precision', 8);dlmwrite('f:計算結果.xls', yy,'-append', 'delimiter', 't', 'precision', 8);%畫圖功能 plot(x,y,'o',x,yy,'*'); xlabel('x'); ylabel('y'); legend('四階Runge-Kutta法','精確解'); pause; close;%2. 結果分析（1）步長h = 0.1 時結果當步長取為0.1時，計算結果如圖1.2-2：圖5-1 h=0.05時