數(shù)值分析作業(yè)MATLAB.doc

1.用二分法解方程 x-lnx=2 在區(qū)間【2 ，4】內(nèi)的根 方法: 二分法 算法：f=inline(x-2-log(x);a=2;b=4;er=b-a; ya=f(a);er0=.00001;while erer0 x0=.5*(a+b); y0=f(x0); if ya*y0 answer1 3 4 3.0000 3.5000 3.0000 3.2500 3.1250 3.2500 3.1250 3.1875 3.1250 3.1563 3.1406 3.1563 3.1406 3.1484 3.1445 3.1484 3.1445 3.1465 3.1455 3.1465 3.1460 3.1465 3.1460 3.1462 3.1461 3.1462 3.1462 3.1462 3.1462 3.1462 3.1462 3.1462 3.1462 3.1462最終結果為: 3.14622.試編寫MATLAB函數(shù)實現(xiàn)Newton插值，要求能輸出插值多項式。對函數(shù) 在區(qū)間-5，5上實現(xiàn)10次多項式插值。Matlab程序代碼如下：%此函數(shù)實現(xiàn)y=1/(1+4*x2)的n次Newton插值，n由調用函數(shù)時指定%函數(shù)輸出為插值結果的系數(shù)向量（行向量）和插值多項式 算法：function t y=func5(n)x0=linspace(-5,5,n+1);y0=1./(1.+4.*x0.2);b=zeros(1,n+1);for i=1:n+1 s=0; for j=1:i t=1; for k=1:i if k=j t=(x0(j)-x0(k)*t; end; end; s=s+y0(j)/t; end; b(i)=s;end; t=linspace(0,0,n+1);for i=1:n s=linspace(0,0,n+1); s(n+1-i:n+1)=b(i+1).*poly(x0(1:i); t=t+s;end;t(n+1)=t(n+1)+b(1);y=poly2sym(t);10次插值運行結果： b Y=func5(10)b = Columns 1 through 4 -0.0000 0.0000 0.0027 -0.0000 Columns 5 through 8 -0.0514 -0.0000 0.3920 -0.0000 Columns 9 through 11 -1.1433 0.0000 1.0000Y = - (7319042784910035*x10)/147573952589676412928 + x9/18446744073709551616 + (256*x8)/93425 - x7/1152921504606846976 - (28947735013693*x6)/562949953421312 - (3*x5)/72057594037927936 + (36624*x4)/93425 - (5*x3)/36028797018963968 - (5148893614132311*x2)/4503599627370496 + (7*x)/36028797018963968 + 1b為插值多項式系數(shù)向量，Y為插值多項式。插值近似值： x1=linspace(-5,5,101); x=x1(2:100); y=polyval(b,x)y = Columns 1 through 12 2.7003 3.9994 4.3515 4.0974 3.4926 2.7237 1.9211 1.1715 0.5274 0.0154 -0.3571 -0.5960 Columns 13 through 24 -0.7159 -0.7368 -0.6810 -0.5709 -0.4278 -0.2704 -0.1147 0.0270 0.1458 0.2360 0.2949 0.3227 Columns 25 through 36 0.3217 0.2958 0.2504 0.1915 0.1255 0.0588 -0.0027 -0.0537 -0.0900 -0.1082 -0.1062 -0.0830 Columns 37 through 48 -0.0390 0.0245 0.1052 0.2000 0.3050 0.4158 0.5280 0.6369 0.7379 0.8269 0.9002 0.9549 Columns 49 through 60 0.9886 1.0000 0.9886 0.9549 0.9002 0.8269 0.7379 0.6369 0.5280 0.4158 0.3050 0.2000 Columns 61 through 72 0.1052 0.0245 -0.0390 -0.0830 -0.1062 -0.1082 -0.0900 -0.0537 -0.0027 0.0588 0.1255 0.1915 Columns 73 through 84 0.2504 0.2958 0.3217 0.3227 0.2949 0.2360 0.1458 0.0270 -0.1147 -0.2704 -0.4278 -0.5709 Columns 85 through 96 -0.6810 -0.7368 -0.7159 -0.5960 -0.3571 0.0154 0.5274 1.1715 1.9211 2.7237 3.4926 4.0974 Columns 97 through 994.3515 3.9994 2.7003繪制原函數(shù)和擬合多項式的圖形代碼：plot(x,1./(1+4.*x.2)hold allplot(x,y,r)xlabel(X)ylabel(Y)title(Runge現(xiàn)象)gtext(原函數(shù))gtext(十次牛頓插值多項式)繪制結果：誤差計數(shù)并繪制誤差圖： hold off ey=1./(1+4.*x.2)-yey = Columns 1 through 12 -2.6900 -3.9887 -4.3403 -4.0857 -3.4804 -2.7109 -1.9077 -1.1575 -0.5128 -0.0000 0.3733 0.6130 Columns 13 through 24 0.7339 0.7558 0.7010 0.5921 0.4502 0.2943 0.1401 0.0000 -0.1169 -0.2051 -0.2617 -0.2870 Columns 25 through 36 -0.2832 -0.2542 -0.2053 -0.1424 -0.0719 -0.0000 0.0674 0.1254 0.1696 0.1971 0.2062 0.1962 Columns 37 through 48 0.1679 0.1234 0.0660 0.0000 -0.0691 -0.1349 -0.1902 -0.2270 -0.2379 -0.2171 -0.1649 -0.0928 Columns 49 through 60 -0.0271 0 -0.0271 -0.0928 -0.1649 -0.2171 -0.2379 -0.2270 -0.1902 -0.1349 -0.0691 0.0000 Columns 61 through 72 0.0660 0.1234 0.1679 0.1962 0.2062 0.1971 0.1696 0.1254 0.0674 0.0000 -0.0719 -0.1424 Columns 73 through 84 -0.2053 -0.2542 -0.2832 -0.2870 -0.2617 -0.2051 -0.1169 0.0000 0.1401 0.2943 0.4502 0.5921 Columns 85 through 96 0.7010 0.7558 0.7339 0.6130 0.3733 0.0000 -0.5128 -1.1575 -1.9077 -2.7109 -3.4804 -4.0857 Columns 97 through 99 -4.3403 -3.9887 -2.6900 plot(x,ey) xlabel(X) ylabel(ey) title(Runge現(xiàn)象誤差圖)3.應用牛頓迭代法于方程f(x)-1=0,導出平方根的迭代公式，用此公式計算. 算法：f=inline(1-115/x2);f1=inline(230/x3);x0=10;er=1;k=0;while er0.00001 x=x0-f(x0)/f1(x0); er=abs(x-x0) x0=x; disp(x);end求解結果： answer9er =0.652210.6522er =0.070910.7231er =7.1604e-0410.7238er =7.1729e-0810.7238最終結果：10.72384.實驗數(shù)據(jù)使用次數(shù)x容積y使用次數(shù)x容積y2106.4211110.593108.2612110.605109.5814110.726109.5016110.907109.8617110.769110.0019111.1010109.9320111.30選用雙曲線對數(shù)據(jù)進行擬合，使用最小二乘法求出擬合函數(shù)，做出擬合曲線圖?！窘狻縞lear,clc;%題目條件x=2 3 5 6 7 9 10 11 12 14 16 17 19 20;y=106.42,108.26,109.58,109.50,109.86,110.00,109.93. 110.59,110.60,110.72,110.90,110.76,111.10,111.30;%使用最小二乘法求出1次多項式擬合系數(shù)a=polyfit(1./x,1./y,1);%繪制擬合圖像xx=0.04:0.01:0.5;yy=a(1)*xx + a(2);plot(1./xx,1./yy,x,y,*);hold on;xx=-0.5:0.01:-0.04;yy=a(1)*xx + a(2);plot(1./xx,1./yy);使用最小二乘法擬合的曲線方程為下圖為繪制出的擬合曲線，并同時將一直點用“*”表示到圖中。5、利用LU分解法解方程組首先，編輯一個LU分解函數(shù)如下functionL,U=Lu(A)% 求解線性方程組的三角分解法% A為方程組的系數(shù)矩陣%L和U為分解后的下三角和上三角矩陣n,m=size(A);if n=m error(The rows and columns of matrix A must be equal!); return;end%判斷矩陣能否LU分解for ii=1:n for i=1:ii for j=1:ii AA(i,j)=A(i,j); end endif (det(AA)=0) error(The matrix can not be divided by LU!) return;endend%開始計算，先賦初值L=eye(n);U=zeros(n,n);%計算U的第一行，L的第一列for i=1:n U(1,i)=A(1,i); L(i,1)=A(i,1)/U(1,1);end%計算U的第r行，L的第r列for i=2:n for j=i:n for k=1:i-1 M(k)=L(i,k)*U(k,j); end U(i,j)=A(i,j)-sum(M); end for j=i+1:n for k=1:i-1 M(k)=L(j,k)*U(k,i); end L(j,i)=(A(j,i)-sum(M)/U(i,i); endend然后，編輯一個通過LU分解法解線性方程組的函數(shù)如下function L,U,x=Lu_x(A,d)%三角分解法求解線性方程組，LU法解線性方程組Ax=LUx=d%A為方程組的系數(shù)矩陣%d為方程組的右端項%L和U為分解后的下三角和上三角矩陣%x為線性方程組的解n,m=size(A);if n=m error(The rows and columns of matrix A must be equal!); return;end%判斷矩陣能否LU分解for ii=1:n for i=1:ii for j=1:ii AA(i,j)=A(i,j); end endif (det(AA)=0) error(The matrix can not be divided by LU!) return;endendL,U=Lu(A); %直接調用自定義函數(shù)，首先將矩陣分解，A=LU%設Ly=d由于L是下三角矩陣,所以可求y(i)y(1)=d(1);for i=2:n for j=1:i-1 d(i)=d(i)-L(i,j)*y(j); end y(i)=d(i);end%設Ux=y,由于U是上三角矩陣，所以可求x(i)x(n)=y(n)/U(n,n);for i=(n-1):-1:1 for j=n:-1:i+1