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);endx

2、2.運行結果：高斯選列主元消元法1. 程序:clearformat ratA=input('輸入增廣矩陣A=')m,n=size(A);for i=1:(m-1) numb=int2str(i);disp('第',numb,'次選列主元后的增廣矩陣')temp=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

3、 j=(i+1):m A(j,:)=A(j,:)-A(i,:)*A(j,i)/A(i,i); end Aend%回代過程 disp('回代求解') x(m)=A(m,n)/A(m,m);for i=(m-1):-1:1 x(i)=(A(i,n)-A(i,i+1:m)*x(i+1:m)')/A(i,i); end x2.運行結果：追趕法1. 程序：function x,L,U=zhuiganfa(a,b,c,f) a=input('輸入矩陣-1對角元素a=');b=input('輸入矩陣對角元素b=');c=input('輸入矩陣+

4、1對角元素c=');f=input('輸入增廣矩陣最后一列元素f=');n=length(b);% 對A進行分解 u(1)=b(1); for i=2:n if(u(i-1)=0) l(i-1)=a(i-1)/u(i-1); u(i)=b(i)-l(i-1)*c(i-1); else break; end end L=eye(n)+diag(l,-1); U=diag(u)+diag(c,1); x=zeros(n,1); y=x; % 求解Ly=b y(1)=f(1); for i=2:n y(i)=f(i)-l(i-1)*y(i-1); end % 求解Ux=y i

5、f(u(n)=0) x(n)=y(n)/u(n); end for i=n-1:-1:1 x(i)=(y(i)-c(i)*x(i+1)/u(i); end 2.運行結果：高斯-塞德爾迭代格式1.程序：function x=Gauss_Seidel(a,b)a=input('輸入系數矩陣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 sum=0; E=0; t(1:n)=x(1:n); for i=1:n x(i)=(b(i)-a(i

6、,1:(i-1)*x(1:(i-1)-a(i,(i+1):n)*t(i+1):n)/a(i,i); end if norm(x-t)<e k break; endend2. 運行結果：雅戈比迭代格式1.程序：function x=Jocabi(a,b)a=input('輸入系數矩陣a=');b=input('輸入增廣矩陣最后一列b=');e=0.5e-7;n=length(b);N=100;x=zeros(n,1);y=zeros(n,1);for k=1:N sum=0; for i=1:n y(i)=(b(i)-a(i,1:n)*x(1:n)+a(i,

7、i)*x(i)/a(i,i); end for i=1:n sum=sum+(y(i)-x(i)2; end if sqrt(sum)<e k break; else for i=1:n x(i)=y(i); end endendif k=N warning('未能找到近似解');end2.運行結果：逐次超松弛法（SOR）1.程序：function n,x=sor22(A,b,X,nm,w,ww)%用超松弛迭代法求解方程組Ax=b%輸入：A為方程組的系數矩陣，b為方程組右端的列向量，X為迭代初值構成的列向量，nm為最大迭代次數，w為誤差精度，ww為松弛因子%輸出：x為求得

8、的方程組的解構成的列向量，n為迭代次數A=input('輸入系數矩陣A=');b=input('輸入方程組右端的列向量b=');X=input('輸入迭代初值構成的列向量X=');nm=input('輸入最大迭代次數nm=');w=input('輸入誤差精度w=');ww=input('輸入松弛因子ww=');n=1;m=length(A);D=diag(diag(A); %令A=D-L-U,計算矩陣DL=tril(-A)+D; %令A=D-L-U,計算矩陣LU=triu(-A)+D; %令A=D-

9、L-U,計算矩陣UM=inv(D-ww*L)*(1-ww)*D+ww*U); %計算迭代矩陣g=ww*inv(D-ww*L)*b; %計算迭代格式中的常數項%下面是迭代過程while n<=nm x=M*X+g; %用迭代格式進行迭代 if norm(x-X,'inf')<w disp('迭代次數為');n disp('方程組的解為');x return; %上面：達到精度要求就結束程序，輸出迭代次數和方程組的解 end X=x;n=n+1;end%下面：如果達到最大迭代次數仍不收斂，輸出警告語句及迭代的最終結果（并不是方程組的解）d

10、isp('在最大迭代次數內不收斂!');disp('最大迭代次數后的結果為');2.運行結果：二分法求解方程的根1.程序：%其中a，b表示查找根存在的范圍，M表示要求解函數的值%f(x)表示要求解根的方程%eps表示所允許的誤差大小function y=er_fen_fa(a,b,M)k=0;eps=0.05while b-a>eps x=(a+b)/2; %檢查是否大于值 if (x3)-3*x-1)>M b=x else a=x end k=k+1end2.運行結果：Newton 迭代法（切線法）1.程序：function x=nanewton(

11、fname,dfname,x0,e,N)%newton迭代法解方程組%fname和dfname分別表示F(x)及其導函數的M函數句柄或內嵌函數，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)>e&k<N, k=k+1; x0=x;x=x0-feval(fname,x0)/feval(dfname,x0); disp(x)endif k=N,warning('已達迭代次數上限');end2.運行結果：割線方式迭代法1

12、.程序：function x=ge_xian_fa(fname,dfname,x0,x1,e,N)%割線方式迭代法解方程組%fname和dfname分別表示F(x)及其導函數的M函數句柄或內嵌函數，x0,x1分別為迭代初值，e為精度要求k=0;a=x1;b=x0;if nargin<5,N=500;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)

13、*feval(fname,x0)>0, x0=x;b=x0; else x1=x;a=x1; end x=x1-(x1-x0)/(feval(fname,x1)-feval(fname,x0)*feval(fname,x1); numb=int2str(k); disp('第',numb,'次計算后x=') fprintf('%fnn',x); endif k=N,warning('已達迭代次數上限');end2.運行結果：Newton插值1.程序：%保存文件名為New_Int.m%Newton基本插值公式%x為向量，全部的

14、插值節點%y為向量，差值節點處的函數值%xi為標量，是自變量%yi為xi出的函數估計值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!');return;end%計算均差表YY=zeros(n);Y(:,1)=y'for k=1:n-1for i=1:n-kif abs(x(i+k)-x(i)<epserror('the DATA is error!');return;endY(i,k

15、+1)=(Y(i+1,k)-Y(i,k)/(x(i+k)-x(i);endend%計算牛頓插值公式yi=0;for i=1:nz=1;for k=1:i-1z=z*(xi-x(k);endyi=yi+Y(1,i)*z;end2.運行結果：Lagrange插值1.程序：function y0 = Language(x,y,x0)syms t l;if length(x)=length(y) n = length(x);else disp('x和y的維數不相等！'); return; %檢錯endh=sym(0);for i=1:n l=sym(y(i); for 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; h=h+l;endsimplify(h);if nargin = 3 y0 = subs (h,'t',x0); %計算插值點的函數值else y0=collect(h); y0 = vpa(y0,6); %將插值多項式的系數化成6位精度的小數end2.運行結果：最小二乘法1.程序：function p=nafit(x,y,m)%多項式擬合% x, y為已知數據點向量, 分別表示