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

1.用二分法解方程 x-lnx=2 在區(qū)間【2 ，4】?jī)?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最終結(jié)果為: 3.14622.試編寫MATLAB函數(shù)實(shí)現(xiàn)Newton插值，要求能輸出插值多項(xiàng)式。對(duì)函數(shù) 在區(qū)間-5，5上實(shí)現(xiàn)10次多項(xiàng)式插值。Matlab程序代碼如下：%此函數(shù)實(shí)現(xiàn)y=1/(1+4*x2)的n次Newton插值，n由調(diào)用函數(shù)時(shí)指定%函數(shù)輸出為插值結(jié)果的系數(shù)向量（行向量）和插值多項(xiàng)式 算法：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次插值運(yùn)行結(jié)果： 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為插值多項(xiàng)式系數(shù)向量，Y為插值多項(xiàng)式。插值近似值： 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ù)和擬合多項(xiàng)式的圖形代碼：plot(x,1./(1+4.*x.2)hold allplot(x,y,r)xlabel(X)ylabel(Y)title(Runge現(xiàn)象)gtext(原函數(shù))gtext(十次牛頓插值多項(xiàng)式)繪制結(jié)果：誤差計(jì)數(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.應(yīng)用牛頓迭代法于方程f(x)-1=0,導(dǎo)出平方根的迭代公式，用此公式計(jì)算. 算法：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求解結(jié)果： answer9er =0.652210.6522er =0.070910.7231er =7.1604e-0410.7238er =7.1729e-0810.7238最終結(jié)果：10.72384.實(shí)驗(yàn)數(shù)據(jù)使用次數(shù)x容積y使用次數(shù)x容積y2106.4211110.593108.2612110.605109.5814110.726109.5016110.907109.8617110.769110.0019111.1010109.9320111.30選用雙曲線對(duì)數(shù)據(jù)進(jìn)行擬合，使用最小二乘法求出擬合函數(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次多項(xiàng)式擬合系數(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);使用最小二乘法擬合的曲線方程為下圖為繪制出的擬合曲線，并同時(shí)將一直點(diǎn)用“*”表示到圖中。5、利用LU分解法解方程組首先，編輯一個(gè)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%開(kāi)始計(jì)算，先賦初值L=eye(n);U=zeros(n,n);%計(jì)算U的第一行，L的第一列for i=1:n U(1,i)=A(1,i); L(i,1)=A(i,1)/U(1,1);end%計(jì)算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然后，編輯一個(gè)通過(guò)LU分解法解線性方程組的函數(shù)如下function L,U,x=Lu_x(A,d)%三角分解法求解線性方程組，LU法解線性方程組Ax=LUx=d%A為方程組的系數(shù)矩陣%d為方程組的右端項(xiàng)%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); %直接調(diào)用自定義函數(shù)，首先將矩陣分解，A=LU%設(shè)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%設(shè)Ux=y,由于U是上三角矩陣，所以可求x(i)x(n)=y(n)/U(n,n);for i=(n-1):-1:1 for j=n:-1:i+1