數(shù)值分析報告編程及運(yùn)行結(jié)果高斯順序消元法

1、高斯消元法1 .程序:clearformat ratA=input('輸入增廣矩陣A=')m,n=size(A);for i=1:(m-1)numb=int2str(i);disp('第',numb,'次消元后的增廣矩陣')for j=(i+1):mA(j,:)=A( j,:)-A(i,:)*A( j,i)/A(i,i);endAend%同代過程disp('同代求解,)x(m)=A(m,n)/A(m,m);for i=(m-1):-1:1x(i)=(A(i,n)-A(i,i+1:m)*x(i+1:m)')/A(i,i);end實用

2、標(biāo)準(zhǔn)文案2 .運(yùn)行結(jié)果:Command Window-'Fils Edit Debug Desktop '' W indow Help和蕨增廣晦加二i11 T-1 3 L 42 -6 1 -54-54-82-1第L次消元后的增廣拒腥第2:次消元后的母廣指陣回代求解精彩文檔實用標(biāo)準(zhǔn)文案高斯選列主元消元法1 .程序：clearformat ratA=input('輸入增廣矩陣A=')m,n=size(A);for i=1:(m-1)numb=int2str(i);disp('第',numb,'次選列主元后的增廣矩陣')temp

3、=max(abs(A(i:m,i);a,b=find(abs(A(i:m,i)=temp);tempo=A(a(1)+i-1,:);A(a(1)+i-1,:)=A(i,:);A(i,:)=tempodisp('第',numb,'次消元后的增廣矩陣')for j=(i+1):mA(j,:)=A( j,:)-A(i,:)*A( j,i)/A(i,i);end精彩文檔實用標(biāo)準(zhǔn)文案Aend%同代過程disp('同代求解,)x(m)=A(m,n)/A(m,m);for i=(m-1):-1:1x(i)=(A(i,n)-A(i,i+1:m)*x(i+1:m)'

4、;)/A(i,i);endx2 .運(yùn)行結(jié)果：精彩文檔實用標(biāo)準(zhǔn)文案小 Command WindowEIb Edit Debug desktop Window H,elp輸入增L矩陌/二1116-13 142 6 L 5A =111-1312-51第1次選則主元后的增廣指院A =2-ai-131111第1次消元后的憎廣拒降A(chǔ) =2-S1003/204L/2第2次選刊主元后的塔廣植應(yīng)3/21772A =精彩文檔精彩文檔-50X3/21/2- D3/217/2第2次選到主元后的增廣矩降策2次消元后的增廣柜隨回代求解1/23/21/23/217/23/2-q17/23/2OVR .：實用標(biāo)準(zhǔn)文案, Co

5、mmand Window XFile Edit D&bucj Desktop Vindow Help第L次消元后的增廣拒院實用標(biāo)準(zhǔn)文案追趕法1 .程序：function x,L,U=zhuiganfa(a,b,c,f)a=input('輸入矩陣-1對角元素a=');b=input('輸入矩陣對角元素b=');c=input('輸入矩陣+1對角元素c=');f=input('輸入增廣矩陣最后一列元素f=');n=length(b);%對A進(jìn)行分解u(1)=b(1);for i=2:nif(u(i-1)=0)l(i-1)=a(

6、i-1)/u(i-1);u(i)=b(i)-l(i-1)*c(i-1);elsebreak;endend精彩文檔實用標(biāo)準(zhǔn)文案L=eye(n)+diag(l,-1);U=diag(u)+diag(c,1);x=zeros(n,1);y=x;%求解Ly=by(1)=f(1);for i=2:ny(i尸f-l(i-1)*y(i-1);end%求解Ux=yif(u(n)=0)x(n)=y(n)/u(n);endfor i=n-1:-1:1x(i)=(y(i)-c(i)*x(i+1)/u(i);end2 .運(yùn)行結(jié)果：精彩文檔實用標(biāo)準(zhǔn)文案小 Comrnand 'WindowFile Edit De

7、bu9 Desktop1 Vindovy-XHa Ip，輸入矩陣T對角元素a=二乙乙11 輸入指B主肅晶元熟式豆2,3.3 輸入矩陣十1對角元素5 n, 1,1輸入增廣矩陣最后一列元素f= 2, L 2,期 ans =15/19-7/1910/1922/19»高斯-塞德爾迭代格式1 .程序：function x=Gauss_Seidel(a,b)a=input('輸入系數(shù)矩陣a=')b=input('輸入增廣矩陣最后一列 b=');e=0.5e-7;n=length(b);N=50;x=zeros(n,1);t=zeros(n,1);for k=1:N

8、sum=0;精彩文檔實用標(biāo)準(zhǔn)文案E=0;t(1:n)=x(1:n);for i=1:nx(i)=(b(i)-a(i,1:(i-1)*x(1:(i-1)-a(i,(i+1):n)*t(i+1):n)/a(i,i);endif norm(x-t)<ekbreak;endend2 .運(yùn)行結(jié)果：精彩文檔實用標(biāo)準(zhǔn)文案,Comrnand Window XFile Edit Debus Desktop ' Vindow Help輸入系數(shù)矩陣口 =X, 3i 0 ;3 4H,Q, 443034-1-104輸入憎廣指牌最后一列b=24, 30,-24：k =24ans -96/31120/31-1

9、62731OVR .；精彩文檔實用標(biāo)準(zhǔn)文案雅戈比迭代格式1 .程序：function x=Jocabi(a,b)a=input('輸入系數(shù)矩陣a=');b=input('輸入增廣矩陣最后一列b=');e=0.5e-7;n=length(b);N=100;x=zeros(n,1);y=zeros(n,1);for k=1:Nsum=0;for i=1:ny(i)=(b(i)-a(i,1:n)*x(1:n)+a(i,i)*x(i)/a(i,i);endfor i=1:nsum=sum+(y(i)-x(i)八2;endif sqrt(sum)<ekbreak;e

10、lse精彩文檔實用標(biāo)準(zhǔn)文案for i=1:nx(i)=y(i);endendendif k=N warning(, 未能找到近似解,);end2 .運(yùn)行結(jié)果： *l Command Window XFile Edit Debua D 噌 5kt0P Vjndow Help ans =J3617/116S42S9/U0S-19163/3667系數(shù)一蔣工。317 ;-L CL i 輸入揩廣矩陣最后一子AH24,翔,-Mk -62一L ans -96/31 120/31 -162/3111 » |OW精彩文檔實用標(biāo)準(zhǔn)文案逐次超松弛法（SOR）1.程序：function n,x=sor22(

11、A,b,X,nm,w,ww)%用超松弛迭代法求解方程組 Ax=b %輸入：A為方程組的系數(shù)矩陣，b為方程組右端的列向量，X為 迭代初值構(gòu)成的列向量，nm為最大迭代次數(shù)，w為誤差精度，ww 為松弛因子%輸出：x為求得的方程組的解構(gòu)成的列向量， n為迭代次數(shù)A=input（'輸入系數(shù)矩陣A='）;b=input（'輸入方程組右端的列向量b='）;X=input（'輸入迭代初值構(gòu)成的列向量X='）;nm=input（' 輸入最大迭代次數(shù) nm='）;w=input（'輸入誤差精度w='）;ww=input（' 輸

12、入松弛因子ww='）;n=1;m=length(A);D=diag(diag(A);%令人=口1-5 計算矩陣 DL=tril(-A)+D;%令 A=D-L-U,計算矩陣 LU=triu(-A)+D;%令人=口1-5 計算矩陣 UM=inv(D-ww*L)*(1-ww)*D+ww*U);% 計算迭代矩陣g=ww*inv(D-ww*L)*b;%計算迭代格式中的常精彩文檔實用標(biāo)準(zhǔn)文案數(shù)項%下面是迭代過程while n<=nmx=M*X+g;%用迭代格式進(jìn)行迭代if norm（x-X,'inf）<wdisp（'迭代次數(shù)為,）;ndisp（'方程組的解為,）

13、;xreturn;%上面：達(dá)到精度要求就結(jié)束程序，輸出迭代次數(shù)和方程組 的解endX=x;n=n+1;end%下面：如果達(dá)到最大迭代次數(shù)仍不收斂，輸出警告語句及迭代的最終結(jié)果（弁不是方程組的解）disp（'在最大迭代次數(shù)內(nèi)不收斂!,）;disp（'最大迭代次數(shù)后的結(jié)果為,）;精彩文檔實用標(biāo)準(zhǔn)文案2.運(yùn)行結(jié)果:4 Command Window XFile Edit D 金 buq Desktop Vindow HelpA輸入系數(shù)柜除3二船3,0:3, 4,-l;-l,0,0輸入方程蛆右端的列向堂M24: 30 -241輸入迭代初值構(gòu)成的列向裝; 1 ;口輸入最大誄代次效nm=l叩

14、輸入誤差蓿度w=Ck 5?-7輸入松弛因子w-L 25迭代次數(shù)為 n =方程蛆的解為x 3, oses工 871C一5. 225SOVR精彩文檔實用標(biāo)準(zhǔn)文案二分法求解方程的根1 .程序：%其中a, b表示查找根存在的范圍，M表示要求解函數(shù)的值%f(x)表示要求解根的方程%eps表示所允許的誤差大小function y=er_fen_fa(a,b,M)k=0;eps=0.05while b-a>epsx=(a+b)/2;%檢查是否大于值if (xA3)-3*x-1)>Mb=xelsea=xendk=k+1end精彩文檔實用標(biāo)準(zhǔn)文案2 .運(yùn)行結(jié)果:，Command WindowFie

15、Edit Debu q Desktop Window Helpfaval1 fnaae, si ：fb= feral (fnaasa b > ：? Input argument -fnameM is undefined.E1TDT Liler f en f a a± 6fa=f*val (fnaae, a) :fb«f evalifnaj!)f = b):» funFinliiM fi'S-S*!-1* )：» nabiscct (funn 1，2, Q. 05>? Undefined function ci method *nabi

16、sect for input arguments of type>> =nabiEect (futi3 1, 2, 0. 05)1.9C63OVR精彩文檔實用標(biāo)準(zhǔn)文案Newton 迭代法(切線法)1 .程序：function x=nanewton(fname,dfname,x0,e,N)%newton迭代法解方程組%fname和dfname分別表示F(x)及其導(dǎo)函數(shù)的M函數(shù)句柄或內(nèi)嵌函數(shù)，x0為迭代初值，e為精度要求x=x0;x0=x+2*e;k=0;if nargin<5,N=500;endif nargin<4 e=1e-4;endwhile abs(x0-x)&g

17、t;e&k<N,k=k+1;x0=x;x=x0-feval(fname,x0)/feval(dfname,x0);disp(x)endif k=N,warning('已達(dá)迭代次數(shù)上限');end精彩文檔實用標(biāo)準(zhǔn)文案2 .運(yùn)行結(jié)果: Command WindowFie Edit Debu q Desktop Window Help» fun-inlLae Ca-costs)1) ;dfun=inlinaC l+sin(K)!);>> s=natiewt on. (ftm, dfim, 13 le-5, 1OCJ? Undefined funct

18、ion ci method ,nantwion, fcr input ar gun ent s of type» x-nanewr an (fun, dfim, L le-5, 10CD0. 7504 1391Q, 73910,73910.7391精彩文檔實用標(biāo)準(zhǔn)文案割線方式迭代法1 .程序：function x=ge_xian_fa(fname,dfname,x0,x1,e,N)%割線方式迭代法解方程組%fname和dfname分別表示F(x)及其導(dǎo)函數(shù)的M函數(shù)句柄或內(nèi)嵌函數(shù)，x0,x1分別為迭代初值，e為精度要求k=0;a=x1;b=x0;if nargin<5,N=500

19、;endif nargin<4 e=1e-4;endwhile abs(a-b)>e&k<N,k=k+1;x=x1-(x1-x0)/(feval(fname,x1)-feval(fname,x0)*feval(fname, x1);if feval(fname,x)*feval(fname,x0)>0,x0=x;b=x0;elsex1=x;a=x1;endx=x1-(x1-x0)/(feval(fname,x1)-feval(fname,x0)*feval(fname,x1);精彩文檔實用標(biāo)準(zhǔn)文案numb=int2str(k);disp('第',

20、numb,'次計算后 x=') fprintf('%fnn',x);endif k=N,warning('已達(dá)迭代次數(shù)上限,);end2 .運(yùn)行結(jié)果：精彩文檔實用標(biāo)準(zhǔn)文案* Command Window XFie Edit Debug Desktop '.Viridow Help» fuii=inline (' x-cos (x)3 ) ;dfun= inline (? Hsin(x)?);» ge_xian_f a(funa dfun,13 2, Le-5, 100)第1次計算后*0. 742293第七次計茸扁=0.

21、 739491第之次計算后行。 73»13«第1次計算后宜工0. 739092第壬次計算后工=0. 739086第6次計算后=0. 739085第7次計算后x=0,739035第4次計算房0. 739085第9次計宜后所0. 739085第L。次計算后K=0. 739096精彩文檔實用標(biāo)準(zhǔn)文案Command WindowFie Edit Debug Desktop Window Help第L。次計算后收5 739085第"次計算后瀘0. 729095第12次計算后。0. 7390S5第國次計算后M=0. 739085竟乜次計算后x0. 739086第15次計算后爐

22、0. 739085第L6次計箕后行 0.739085第1：次計算后會0. 739085ans0. 7391OVR .：精彩文檔實用標(biāo)準(zhǔn)文案精彩文檔實用標(biāo)準(zhǔn)文案Newton 插值1 .程序：%保存文件名為 New_Int.m%Newton 基本插值公式%x為向量，全部的插值節(jié)點%y為向量，差值節(jié)點處的函數(shù)值%xi為標(biāo)量，是自變量%yi為xi出的函數(shù)估計值function yi=newton_chazhi(x,y,xi)n=length(x);m=length(y);if n=merror('The lengths of X ang Y must be equal!');retur

23、n;end%計算均差表YY=zeros(n);Y(:,1)=y'for k=1:n-1for i=1:n-k if abs(x(i+k)-x(i)<epserror('the DATA is error!');精彩文檔實用標(biāo)準(zhǔn)文案return;endY(i,k+1)=(Y(i+1,k)-Y(i,k)/(x(i+k)-x(i);endend%計算牛頓插值公式y(tǒng)i=0;for i=1:nz=1;for k=1:i-1z=z*(xi-x(k);endyi=yi+Y(1,i)*z;end精彩文檔實用標(biāo)準(zhǔn)文案2 .運(yùn)行結(jié)果:一L Command WindowFie Edit

24、 Debuq Desktop Window Help» cleax» yi=nTtoii_eha2hL (pi/6 :pt/4 T 0, 5:0. 7071, 2*pi/9)? Undefined function ex mo%hod 1 newxon_chashin for inpul ar&uaents of iyyi-newtchazhi (pi/6 ;pi/4 9 0« 5; D« 7071 i 2*pi/9)0.6381y i =ne wt on_ chaahi ( p i/ 6 ;pi/4;pi/3'T ' 0. 5

25、; 0. 7071; 0_ 8650 T 2 MpM9)0. 6434；OVR精彩文檔實用標(biāo)準(zhǔn)文案Lagrange 插值1 .程序：function y0 = Language(x,y,x0)syms t l;if length(x)=length(y)n = length(x);elsedisp('x和y的維數(shù)不相等！');return; % 檢錯endh=sym(0);for i=1:nl=sym(y(i);for j=1:i-1l=l*(t-x( j)/(x(i)-x( j);end;for j=i+1:nl=l*(t-x( j)/(x(i)-x( j);end;h=h+l;end精彩文檔實用標(biāo)準(zhǔn)文案simplify(h);if narg