高斯消元法MATLAB实现
- 格式:doc
- 大小:77.50 KB
- 文档页数:5
matlab高斯-约旦消去法
高斯-约旦消去法是一种线性代数中的消元法,常用于求解线性方程组。
该方法通过矩阵的初等变换将方程组转化为阶梯型矩阵,从而求解出未知数的值。
具体步骤如下:
设有n个未知数,m个方程,方程组的系数矩阵为A,右端常数为b。
1. 将系数矩阵A和右端常数b组合成增广矩阵Ab。
2. 从第一行开始,将该行的第一个非零元素(称为主元)作为消元元素,用该元素将下面所有行的对应列元素消为零。
3. 重复以上步骤,依次将每一行的主元素作为消元元素,直到将整个矩阵消成阶梯型矩阵。
4. 倒序回代,求出每个未知数的值。
以上就是高斯-约旦消去法的主要步骤。
在实际应用中,需要注意判断矩阵是否可逆,以及主元素是否为零等情况,以保证求解的正确性。
实验一 用列主元Gauss 消去法求解线性方程组实验目的会使用Matlab 语言编程使用列主元Gauss 消去法求解线性方程组.实验原理1、 列主元Gauss 消去法记线性方程组1112111212222212n n n n nn n n a a a x b a a a x b a a a x b ⎛⎫⎛⎫⎛⎫⎪⎪ ⎪ ⎪⎪ ⎪=⎪⎪ ⎪⎪⎪ ⎪⎝⎭⎝⎭⎝⎭ 为Ax=b, 其中A =111212122212n n n n nn a a a a a a a a a ⎛⎫ ⎪ ⎪ ⎪⎪⎝⎭,x=12n x x x ⎛⎫⎪ ⎪ ⎪⎪⎝⎭, b=12n b b b ⎛⎫ ⎪ ⎪ ⎪ ⎪⎝⎭, 记其增广矩阵为()(1)(1)(1)1111(1)(1)(1)(1)(1)2122(1)(1)(1)1n nn nnn a a b aa b Ab a a b ⎛⎫ ⎪ ⎪= ⎪ ⎪ ⎪⎝⎭。
设主元(1)11a 0≠,记(1)11(1)11(2,3,,)i i a l i n a =-=,用1i l 乘增广矩阵()(1)(1)A b 的第1行,再分别与第i 行相加,得()(1)(1)(1)(1)111211(1)(1)(2)(2)(2)2222(2)(2)(2)2b 00n nn nnn a a a a a b Ab a a b ⎛⎫ ⎪⎪= ⎪ ⎪ ⎪⎝⎭, 其中(2)(1)(1)1,ij ij i ij a a l a =+ i ,j=2,3,,n(2)(1)(1)11,i i i b b l b =+ i=2,3,,n又设主元(2)(2)i222i2(2)22a 0,l =-a a≠用乘矩阵()(2)(2)A b 的第二行,再与第i 行相加(i=3,4,,n ),得()(1)(1)(1)(1)(1)1112131n 1(2)(2)(2)(2)22232n 2(3)(3)(3)(3)(3)333n3(3)(3)(3)n3nnn b 0b Ab =0b 00b a a a a a a a a a a a ⎛⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎝⎭. 经过n-1步消去后,增广矩阵最终变为()=)()( n n b A实验程序function x=gaussc(A ,b ) [n,m ]=size(A); A=[A ,b]; for k=1:n —1for p=k+1:nif abs(A (p ,k))>abs(A(k ,k )) for j=k :n+1 t=A(k,j );A (k,j)=A(p,j ); A(p,j)=t; end endend %搜索主元并交换 for i=k+1:nl=-A(i ,k )/A(k,k); for j=k+1:n+1A (i,j )=A(i ,j )+l*A(k,j); end endend %消去过程结束 x(n)=A (n,n+1)/A (n ,n); for i=n —1:-1:1 s=0;for j=i+1:ns=s+A (i,j)*x(j ); endx (i)=(A (i,n+1)-s )/A (i,i); end实验结果设A=[2,5;4,6],b=[3;4],求解线形方程组Ax=b.实验步骤:1) 先在matlab 里输入上面的程序;2) 然后输入A=[2,5;4,6] b=[3;4]3)再输入x=gaussc(A,b)命令即得出结果.由以上程序可求解得到x=( 0.2500 0。
高斯列主元素消去法是一种解线性方程组的常用方法,特别在数值分析和线性代数中应用广泛。
在Matlab中,我们可以使用该方法来解决大规模的线性方程组,包括矩阵的求解和矩阵的反转。
一、高斯列主元素消去法的基本原理高斯列主元素消去法是一种基于矩阵消元的方法,它通过一系列的矩阵变换将原始的线性方程组转化为上三角形式,然后再进行回代求解。
这个方法的核心就是通过矩阵的变换来简化原始的线性方程组,使得求解过程更加简单高效。
在Matlab中,我们可以利用矩阵运算和函数来实现高斯列主元素消去法,如`lu`分解函数和`\"`运算符等。
通过这些工具,我们能够快速地求解各种规模的线性方程组并得到准确的结果。
二、高斯列主元素消去法在Matlab中的实现在Matlab中,我们可以通过调用`lu`函数来实现高斯列主元素消去法。
该函数返回一个上三角矩阵U和一个置换矩阵P,使得PA=LU。
通过对U进行回代求解,我们可以得到线性方程组的解。
除了`lu`函数之外,Matlab还提供了一些其他的函数和工具来帮助我们实现高斯列主元素消去法,比如`\"`运算符和`inv`函数等。
通过这些工具的组合使用,我们能够更加灵活地进行线性方程组的求解,并且可以方便地处理特殊情况和边界条件。
三、高斯列主元素消去法的应用与局限性高斯列主元素消去法在实际应用中具有广泛的适用性,特别是对于大规模的线性方程组或者稀疏矩阵的求解。
通过Matlab中的工具和函数,我们可以快速地求解各种规模的线性方程组,并得到高精度的数值解。
然而,高斯列主元素消去法也存在一些局限性,比如对于奇异矩阵或者接近奇异矩阵的情况时,该方法的求解精度可能会下降。
在实际应用中,我们需要结合具体的问题和矩阵特性来选择合适的求解方法,以确保得到准确的结果。
四、个人观点和总结作为一种经典的线性方程组求解方法,高斯列主元素消去法在Matlab 中具有较好的实现和应用效果。
通过对其原理和实现细节的深入理解,我们能够更加灵活地应用该方法,并且能够更好地理解其适用性和局限性。
gauss消去法matlabGauss消去法是一种常用的线性方程组求解方法,它可以通过消元和回代的方式,将一个复杂的线性方程组转化为一个简化的三角形方程组,从而得到方程组的解。
在MATLAB中,我们可以使用高斯消去法函数来求解线性方程组。
我们需要明确线性方程组的形式。
一个典型的线性方程组可以表示为:Ax = b其中,A是一个n×n的系数矩阵,x是一个n×1的未知向量,b是一个n×1的常数向量。
接下来,我们可以使用MATLAB中的高斯消去法函数来求解线性方程组。
在MATLAB中,我们可以使用“[L,U,P] = lu(A)”函数来进行高斯消去法的分解,其中L是单位下三角矩阵,U是上三角矩阵,P 是置换矩阵。
通过高斯消元法的分解,我们可以得到三角形方程组:L(Ux) = b然后,我们可以使用“y = L\b”函数来求解下三角方程Ly = b,再使用“x = U\y”函数来求解上三角方程Ux = y。
最终,我们可以得到线性方程组的解x。
除了使用MATLAB中的高斯消去法函数,我们还可以手动实现高斯消去法。
首先,我们可以通过消元操作将系数矩阵A转化为上三角矩阵U。
消元操作的基本步骤如下:1.选择主元:选择第一列中绝对值最大的元素作为主元,并将其所在的行交换到第一行。
2.消元操作:对于第一行以下的每一行,将其第一列元素消为0。
具体操作是,将第一行乘以一个适当的倍数,然后从当前行中减去第一行的倍数。
3.重复以上步骤,直到所有的主元都不为0或者所有的行都消元结束。
接下来,我们可以使用回代操作将上三角矩阵U转化为解向量x。
回代操作的基本步骤如下:1.确定最后一个未知量:将最后一行的最后一个元素设为1。
2.回代计算:从最后一行开始,依次计算每个未知量的值。
具体操作是,将当前行的右侧元素减去已知的未知量的倍数,然后除以当前行对角线上的系数。
通过手动实现高斯消去法,我们可以更好地理解高斯消去法的原理和过程。
matlab⾼斯消gf,⾼斯消去、追赶法matlab(⽰例代码)1. 分别⽤Gauss消去法、列主元Gauss消去法、三⾓分解⽅法求解⽅程组程序:(1)Guess消去法:function x=GaussXQByOrder(A,b)%Gauss消去法N = size(A);n = N(1);x = zeros(n,1);for i=1:(n-1)for j=(i+1):nif(A(i,i)==0)disp(‘对⾓元不能为0‘);return;endm = A(j,i)/A(i,i);A(j,i:n)=A(j,i:n)-m*A(i,i:n);b(j)=b(j)-m*b(i);endendx(n)=b(n)/A(n,n);for i=n-1:-1:1x(i)=(b(i)-sum(A(i,i+1:n)*x(i+1:n)))/A(i,i);end命令⾏输⼊:A=[1 -1 2 1;-1 3 0 -3 ;2 0 9 -6;1 -3 -6 19];b=[1 3 5 7];b=b‘;x=GaussXQByOrder(A,b)-8.0000000000000000.3333333333333333.6666666666666672.000000000000000(2)列主元Gauss消去法程序:function x=GaussXQLineMain(A,b)%列主元Gauss消去法N = size(A);n = N(1);x = zeros(n,1);zz=zeros(1,n);for i=1:(n-1)[~,p]=max(abs(A(i:n,i)));zz=A(i,:);A(i,:)=A(p+i-1,:);A(p+i-1,:)=zz;temp=b(i);b(i)=b(i+p-1);b(i+p-1)=temp;for j=(i+1):nm = A(j,i)/A(i,i);A(j,i:n)=A(j,i:n)-m*A(i,i:n);b(j)=b(j)-m*b(i);endendx(n)=b(n)/A(n,n);for i=n-1:-1:1x(i)=(b(i)-sum(A(i,i+1:n)*x(i+1:n)))/A(i,i); end命令⾏:A=[1 -1 2 1;-1 3 0 -3 ;2 0 9 -6;1 -3 -6 19];x=GaussXQLineMain(A,b)运⾏结果:x =-8.0000000000000050.3333333333333323.6666666666666682.000000000000000(3)三⾓分解⽅法程序:function x = LU(A,b)%三⾓分解N = size(A);n = N(1);L = eye(n,n);U = zeros(n,n);x = zeros(n,1);y = zeros(n,1);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i);endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k); endendy(1)=b(1)/L(1,1);for i=2:ny(i)=b(i)-sum(L(i,1:i-1)*y(1:i-1));endx(n)=y(n)/U(n,n);x(i)=(y(i)-sum(U(i,i+1:n)*x(i+1:n)))/U(i,i);end命令⾏:A=[1 -1 2 1;-1 3 0 -3 ;2 0 9 -6;1 -3 -6 19];b=[1 3 5 7];b=b‘;x=LU(A,b)运⾏结果:x =-8.0000000000000000.3333333333333333.6666666666666672.000000000000000程序:function [times,wucha]=zhuiganfa(a,b,c,f)%追赶法:x为所求解,times为所有乘除运算次数(即时间),wucha为误差的2-范数。
1151091 杨晨辉高斯消去法解线性方程的Matlab 程序方法一:function x = gauss(A,b) n = length(b);for k = 1 : n-1if A(k,k)==0fprintf( 'Error: the %dth pivot element equal to zero!\n' return ;endindex = [k+1:n];m = -A(index,k)/A(k,k);A(index,index) = A(index,index) + m*A(k,index); b(index) = b(index) + m*b(k);endx = zeros(n,1);x(n) = b(n)/A(n,n);for i = n-1:-1:1x(i) = ( b(i) - A(i,[i+1:n])*x([i+1:n]) )/A(i,i); end运行结果:>> A=[1 1.355 1.4 2; 3 3.5 0.22 1; 0.5 2 2.1 3;>> b=[2.00,1.00,0.55,3.00]' b =2.00001.00000.55003.0000 >> gauss(A,b)ans =2.5225-2.23130.01771.2381 方法二:矩阵求逆:function [ B ] = qiuni( A )%UNTITLED Summary of this function goes here% Detailed explanation goes here n=numel(A);r=rank(A);B=eye(r);if n==L2for k=1:rfor i=1:r ,k);0.3 0.1 -0.55 2];for j=1:rii=i-1;jj=j-1;if ii==0 && jj==0;ii=r;jj=r;B(ii,jj)=1/A(1,1);elseif ii==0 && jj~=0;ii=r;B(ii,jj)=-A(1,j)/A(1,1);elseif jj==0 && ii~=0;jj=r;B(ii,jj)=A(i,1)/A(1,1);elseB(ii,jj)=A(i,j)-A(i,1)*A(1,j)/A(1,1);endendendA=B;B=eye(r);endB=A;elsemsgbox( ' 矩阵不可逆' , 'message' , 'warn' ); end endGr^l Jr«t jskiALM13LII4I3 S ]J. IXI f 8 4»LEF3f* ll.Tfil™'.*?a. nr-ii I. >3 DOC-1- E^:I.MQP n rml.»HGil・a血«.4WtOMHi超if iMr MHfH 曰ITfar !■ i: rIlliUliNJTL2in kLIU#•十*minilu-dzr齡 4 J 2.4 7^3.9 il-l-Z U • (14 M,1 1*1,1 i Ih.i 他[S 4 J 2 4 2 it 3.3 I ? 5.4 3 3 lh”gg飞『可选1M^4T<E AfrrE喩耳.庇賈■斗『.*玛11册V ttii< .M ■帥?脅 Q* M It 当方程不可逆时:第二种:fun ction [ B ] = lyxq inv( A )%UNTITLED Summary of this fun ctio n goes here % Detailed expla nati on goes here n=n umel(A);r=ra nk(A);B=eye(r);if n==L2E=eye(r);A=[A E];for k=1:(r-1)for i=(k+1):rfor j=(k+1):2*rA(i,j) = A(i,j)-A(i,k)* A(k,j) / A(k,k);endj=k;A(i,j) = A(i,j)-A(i,k)* A(k,j) / A(k,k);endendfor k=2:rpuspus !( .UJBM , £ .sBesseiu, £ ,廡址k 捌目# . )xoq6siu9S|9 :(j^:(i,+j )i:)a=apue pue©!)▼/(「!)皆(「!)日J^:L=[」oj」j=!」ojpue pue:(>l 1>l )V/(r>l )V.(>l 1!)V-(r!)V=(r!)V■>l=f pue:(>i i >i )v/(r>i )v.(>i i !)v-(r!)v=(r!)v j^:(i /+>i )=r 」oj(宀)J=!」OJ!■« ^'1 1 I t'C 0 £ S't t fe- S]-*-fiN1=5 3 If• ◎ ID * ・E -*w J£ 口1M^!PIS^'D-EMJ'DgfWP Sfft F IMF CGiK S^£j£K 加屜iz "UP EWtT>- UKH MX®MJFV- HQFU岡Z 伽巾w:MMsmn a t i 11 ft g i t z t if 町■ V —1 J£*frdrj_ E-C-^'P ■fp XklMC I «f*E|«ul3 5 ¥ 1 i I r'C 9 £ C I ・ ** epiiJbcMH L ^btLJ=Ul] Ci 1 J't 1 I L'L 9 ■£ R* L I 3]^j ",|Q 鼻 H-R r ft ・ 4 ■卩 * C- [- T'P- f ■门碑 :\|T- i TI :-q r R!- 3 I 旷£- « 1- 4 1- f E- F]-» 14肯1肚丫 F1 41 J. I 时4T ■ = FIT IF P 町贰F 3»rl-qll 5 I 3 -5 A I E :f B 5 P C E > 5]<*V] *■*iia t i i i E C g i P S I t ・町詐(q Fri "T<e fc-xsr[i —we*Fip - . jr r ■'w■I« !J►5 ft O EWK陆羽 LfEMLD?K :fl.«ll r miEI f lfiMlbiM.T -[| g 盅 E fi J E £'L 5> E »'Z € t £]=9 «IH -MnhiQui'u^iJEi 1I in-M ka« tfun cti on X=jie(A,b); [m,n ]=size(A);B=[A,b];RA=ra nk(A);RB=ra nk(B);formatratif RA==RB & RA==n %判断方程有唯一解 X=A\b;else if RA==RB & RA<n %判断方程有无穷解X=A\b %求特解 C=null(A,'R' );%求AX=0勺基础解系else X='方程无解’;%判断方程无解 end end上图的方程有两个,第一个有解,显示结果;第二个无解,显示 方程无解当方程有无穷解时,显示其特解■ -1^ u* sf Bi fclw JtJf• 6 •CiJ □ £J 4■谍 1 3 7 fl I 7 >7 13 D<L. rlsipc-hn- >Ji>; «»EV i > 1.4 > H >J < i T >■ t iipiM -caim 勒诵 < 3 j,4 g ].] I 7 5.? 3 fi« j a m 半・ jj.^(]!4M r A ).i i Fv.? > t 盯・・口©・ *-C5 4 3 a pfUkb)A F IB 4 ft 1J ji«L4ib)X|t 4 1 1,4 •EE i 3 2;a Jt Ct < 5 1.4 ■£L -2 J -L J -i & -J,2 I 2 -2P r «^(b 2 ]]',1 i 1J l !*«/■■ 1]X14N3 g 3.3 I T 5.S 3 5 l.|.b-[14Mi 1 6 3.3 I r 9.2 3 1 d|il>±(14rt if山n 4 I 2 4 I 4 3,9 3 r t,3 ) 6 lixim KJ Hip.EDM 酉52IDiU|i -a v -i ,s -i b o s ? ar■AT“n4 a M H 4 3 IQ:\・B I 各聲a&14+3 3l :gi52図邸H - g ><k 申r 5.71 %门.・口0医2 d ).] i. ?3 S i|,b-{^94 I■IWr=iIMIS» M 【l I 3 ・l"・l •> 4.1 5 V -«),tell < •】•.”心”」UTIMntrtE ex J 9 Q, " Si・ V12M)l»724?»n« Ut9・ i/i2«e )in^4rmi« m|g 门“t—SJ. *Mb 4 ) 1.4 2 ・).3 I r e.4 ) 1MH 4 >:.< 2 • X> I r «.? ) t ibMOHM “ 4 > 2.4 2 e X )I r 5.2 ) 3ibUOUMA^li 4 ) 2.4 2 « >.> I I C.2 ) < l)#b-(l«M IJWU.WMIS 4 > 2.4 7 • J.5 I f 5.2 S 1 IJUU.b>A F 5 < > >.4 < 3.S I > 6.J > « l).^(l«M 1 Mb <>:.<? « XI I r o.2 ) 9 lbW (l«M I "4 2 6 3.3 I 7 6.J 1 5l),W(l<M I4 > 2.4 2 • >.> t r B.2 > ■ IA=!l -2 3 -l.> -I 5 -3.2 I 2 -2!.bs(> 2 «* J A*|l 1 J -I i ・(・1 4.1 i -• -l|.b-(l 4 •)•*nu *" >«w “5*a□ Jk fi W 2f » -ute«twec»>x - €>»»<» M S^*9V«9<aa >4X<t * * J (J)Mg ■ r.Tp■"daStn m 一「一 . M4te~ 3144-18S44 3B • r *H» " MIL X14 4-189 15 M SmM*lte 3314^21719X3M4・ X144-3215152® .dbr«4M»U RMUMXH4^259O4C N 1 Qf* n M 低 an 心 2( aooTXT Mr201H-2 Id 1)518)0 &1Q @6 ■.U, _ B Ulw« Mi . : MUI 4iwAM« .。
3线性代数方程组数值解法39.(上机题)列主元高斯消去法对于某电路的分析,归结为求线性方程组RI=V ,其中⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−−−−−−−−−−−−−−−−−−−=2920009000227005000000413000000003047700000507573000090003079100000000103190000011093513000100001331R ()15,27,23,0,20,12,7,7,10TT V =−−−−(1)编制解n 阶线性方程组Ax b =的列主元高斯消去法的通用程序;(2)用所编程序解线性方程组RI V =,并打印出解向量,保留5位有效数;(3)本题编程中,你提高了那些编程能力?本程序用matlab 编写(1)通用程序如下function [x,det,flag ]=Gauss(A,b)%A 为方程组的系数矩阵%b 为方程组的右端项%x 为方程组的解%det 为系数矩阵A 的行列式的值%flag 为指标向量,flag=‘failure ’表示失败,flag=‘OK ’表示成功[n,m]=size(A);nb=length(b);if n~=merror('A 不是方阵')return;endif m~=nberror('b 的长度不等于A 的阶数')return;endflag='OK';det=1;x=zeros(n,1);for k=1:n-1max1=0;for i=k:nif abs(A(i,k))>max1max1=abs(A(i,k));r=i;endendif max1<1e-10flag='failure';return;endif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfor i=k+1:nm=A(i,k)/A(k,k);for j=k+1:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);enddet=det*A(k,k);enddet=det*A(n,n)if abs(A(n,n))<1e-10flag='failure';return;endfor k=n:-1:1for j=k+1:nb(k)=b(k)-A(k,j)*x(j);endx(k)=b(k)/A(k,k);endx(k)=b(k)/A(k,k);endvpa(x)digits(5)(2)在命令栏输入矩阵,并执行guass命令如下>>A=[31-13000-10000-1335-90-1100000-931-100000000-1079-30000-9000-3057-70-500000-747-300000000-3041000000-50027-2000-9000-229];>>b=[-15;27;-23;0;-20;12;-7;7;10];>>x=Gauss(A,b);fprintf('%.5g\n',x)运行结果如下:det=6.1817e+013X=-0.289230.34544-0.71281-0.22061-0.43040.15431-0.0578230.201050.29023(3)通过此次编程,让我更加熟悉了matlab程序的编辑,对选择语句,循环语句,循环嵌套等的运用更加熟练,也对列主元Guass消去法解方程组的思想理解更加深刻。
1.高斯消元法function [x]=mgauss(A,b,flag)ticif nargin<3,flag=0;endn=length(b);for k=1:(n-1)m=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-m*A(k,k+1:n);b(k+1:n)=b(k+1:n)-m*b(k);A(k+1:n,k)=zeros(n-k,1);if flag~=0,Ab=[A,b],endendx=zeros(n,1);x(n)=b(n)/A(n,n);for k=n-1:-1:1x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k);endtoc2.程序验证n=6;>> r=rand(n);a=r'*r+n*eye(n);b=a*ones(n,1);>> x=mgauss(a,b)Elapsed time is 0.000000 seconds.x =1.00001.00001.00001.00001.00001.0000n=100;r=rand(n);a=r'*r+n*eye(n);b=a*ones(n,1);x=mgauss(a,b) Elapsed time is 0.015000 seconds.x =1.00001.00001.00001.00001.00001.00001.0000function [RA,RB,x]=gaus(A,b,flag)ticif nargin<3,flag=0;endB=[A b]; n=length(b); RA=rank(A);RB=rank(B);if RB-RA>0disp('请注意:因为RA~=RB,所以此方程组无解.')returnendif RA==RBif RA==ndisp('请注意:因为RA=RB=n,所以此方程组有唯一解.')for k=1:(n-1)m=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-m*A(k,k+1:n);b(k+1:n)=b(k+1:n)-m*b(k);A(k+1:n,k)=zeros(n-k,1);if flag~=0,Ab=[A,b],endendx=zeros(n,1);x(n)=b(n)/A(n,n);for k=n-1:-1:1x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k);endelsedisp('请注意:因为RA=RB<n,所以此方程组有无穷多解.')endendtoc程序验证n=100; r=rand(n);a=r'*r+n*eye(n);b=a*ones(n,1);[ra,rb,x]=gaus(a,b) 请注意:因为RA=RB=n,所以此方程组有唯一解.Elapsed time is 0.032000 seconds.ra =100rb =100x =1.00001.00001.00001.00001.00001.0000图形2.LU分解程序function[l,u,x]=mlu(a,b)ticn=length(b);l=zeros(n);u=zeros(n);x=zeros(n,1);y=zeros(n,1); for k=1:nu(k,k:n)=a(k,k:n)-l(k,1:k-1)*u(1:k-1,k:n);l(k:n,k)=(a(k:n,k)-l(k:n,1:k-1)*u(1:k-1,k))/u(k,k);l(k,k)=1;endfor k=1:ny(k)=b(k)-l(k,1:k-1)*y(1:k-1);endfor k=n:-1:1x(k)=(y(k)-u(k,k+1:n)*x(k+1:n))/u(k,k);endtoc验证n=6,r=rand(n),a=r'*r+n*eye(n);b=a*ones(n,1);[l,u,x]=mlu(a,b)n =6r =0.8462 0.6813 0.3046 0.1509 0.4966 0.34200.5252 0.3795 0.1897 0.6979 0.8998 0.28970.2026 0.8318 0.1934 0.3784 0.8216 0.34120.6721 0.5028 0.6822 0.8600 0.6449 0.53410.8381 0.7095 0.3028 0.8537 0.8180 0.72710.0196 0.4289 0.5417 0.5936 0.6602 0.3093 Elapsed time is 0.000000 seconds.l =1.0000 0 0 0 0 00.2303 1.0000 0 0 0 00.1367 0.1246 1.0000 0 0 00.2291 0.1977 0.1438 1.0000 0 00.2676 0.2622 0.1441 0.2122 1.0000 00.1814 0.1540 0.0926 0.1288 0.1104 1.0000u =8.1875 1.8854 1.1195 1.8760 2.1912 1.48510 7.8060 0.9728 1.5430 2.0464 1.20180 0 6.7424 0.9694 0.9715 0.62430 0 0 7.5994 1.6123 0.97890 0 0 0 7.6472 0.84410 0 0 0 0 6.4954 x =1.00001.00001.00001.00001.00001.0000n=90,r=rand(n),a=r'*r+n*eye(n);b=a*ones(n,1);[l,u,x]=mlu(a,b)n =90r =Columns 1 through 90.8385 0.0164 0.4199 0.2731 0.8518 0.6859 0.8686 0.5152 0.19300.5681 0.1901 0.7537 0.6262 0.7595 0.6773 0.6264 0.6059 0.90960.3704 0.5869 0.7939 0.536x =1.00001.00001.00001.00001.0000经验证是正确的,太长了只节了部分。
《数值分析》实验报告
一、实验目的与要求
1.掌握高斯消去法的基本思路和迭代步骤;
2.培养编程与上机调试能力。
二、实验内容
1.编写用高斯消元法解线性方程组的MATLAB程序,并求解下面的线性方程组,然后用逆矩阵解方程组的方法验证.
(1)
123
123
123
0.101 2.304 3.555 1.183
1.347 3.712 4.623
2.137
2.835 1.072 5.643
3.035
x x x
x x x
x x x
++=
⎧
⎪
-++=
⎨
⎪-++=
⎩
(2)
123
123
123
528
28321
361
x x x
x x x
x x x
++=
⎧
⎪
+-=
⎨
⎪--=
⎩
2.编写用列主元高斯消元法解线性方程组的MATLAB程序,并求解下面的线性方程组,然后用逆矩阵解方程组的方法验证.
(1)
123
123
123
0.101 2.304 3.555 1.183
1.347 3.712 4.623
2.137
2.835 1.072 5.643
3.035
x x x
x x x
x x x
++=
⎧
⎪
-++=
⎨
⎪-++=
⎩
(2)
123
123
123
528
28321
361
x x x
x x x
x x x
++=
⎧
⎪
+-=
⎨
⎪--=
⎩
三.MATLAB计算源程序
1. 用高斯消元法解线性方程组b
AX=的MATLAB程序
输入的量:系数矩阵A和常系数向量b;
输出的量:系数矩阵A和增广矩阵B的秩RA,RB, 方程组中未知量的个数n 和有关方程组解X及其解的信息.
function [RA,RB,n,X]=gaus(A,b)
B=[A b]; n=length(b); RA=rank(A);
RB=rank(B);zhica=RB-RA;
if zhica>0,
disp('请注意:因为RA~=RB,所以此方程组无解.')
return
end
if RA==RB
if RA==n
disp('请注意:因为RA=RB=n,所以此方程组有唯一解.')
X=zeros(n,1); C=zeros(1,n+1);
for p= 1:n-1
for k=p+1:n
m= B(k,p)/ B(p,p); B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1);
end
end
b=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n);
for q=n-1:-1:1
X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);
end
else
disp('请注意:因为RA=RB<n,所以此方程组有无穷多解.') End
End
2.列主元消元法及其MATLAB程序
AX 的MATLAB程序
用列主元消元法解线性方程组b
输入的量:系数矩阵A和常系数向量b;
输出的量:系数矩阵A和增广矩阵B的秩RA,RB, 方程组中未知量的个数n和有关方程组解X及其解的信息.
function [RA,RB,n,X]=liezhu(A,b)
B=[A b]; n=length(b); RA=rank(A);
RB=rank(B);zhica=RB-RA;
if zhica>0,
disp('请注意:因为RA~=RB,所以此方程组无解.')
return
end
if RA==RB
if RA==n
disp('请注意:因为RA=RB=n,所以此方程组有唯一解.')
X=zeros(n,1); C=zeros(1,n+1);
for p= 1:n-1
[Y,j]=max(abs(B(p:n,p))); C=B(p,:);
B(p,:)= B(j+p-1,:); B(j+p-1,:)=C;
for k=p+1:n
m= B(k,p)/ B(p,p); B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1);
end
end
b=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n);
for q=n-1:-1:1
X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);
end
else
disp('请注意:因为RA=RB<n,所以此方程组有无穷多解.') end
end
三.实验过程:
1(1)编写高斯消元法的MATLAB文件如下:
clear;
A=[0.101 2.304 3.555;-1.347 3.712 4.623;-2.835 1.072 5.643];
b=[1.183;2.137;3.035];
[RA,RB,n,X] =gaus (A,b)
运行结果为:
请注意:因为RA=RB=n,所以此方程组有唯一解.
RA =
3
RB =
3
n =
3
X =
-0.3982
0.0138
0.3351
(2)编写高斯消元法MATLAB文件如下:
clear;
A=[5 2 1;2 8 -3;1 -3 -6];
b=[8;21;1;];
[RA,RB,n,X] =gaus (A,b)
运行结果为:
请注意:因为RA=RB=n,所以此方程组有唯一解.
RA =
3
RB =
3
n =
3
X =
1
2
-1
在MATLAB中利用逆矩阵法检验结果:
(1) 在command windows中直接运行命令:
A=[0.101 2.304 3.555;-1.347 3.712 4.623;-2.835 1.072 5.643];
b=[1.183;2.137;3.035];X=A\b
运行结果为:
X =
-0.3982
0.0138
0.3351
(2) 在command windows中直接运行命令:
A=[5 2 1;2 8 -3;1 -3 -6];
b=[8;21;1;];X=A\b
运行结果为:
X =
1
2
-1
两小题所得结果相同,检验通过
2(1)编写列组高斯消元法MATLAB文件如下:
clear;
A=[0.101 2.304 3.555;-1.347 3.712 4.623;-2.835 1.072 5.643]; b=[1.183;2.137;3.035];
[RA,RB,n,X] =liezhu(A,b)
运行结果:
请注意:因为RA=RB=n,所以此方程组有唯一解.
RA =
3
RB =
3
n =
3
X =
-0.3982
0.0138
0.3351
(2)编写列组高斯消元法的MATLAB文件如下:
clear;
A=[5 2 1;2 8 -3;1 -3 -6];
b=[8;21;1;]
[RA,RB,n,X] =liezhu(A,b)
运行结果为:
请注意:因为RA=RB=n,所以此方程组有唯一解.
RA =
3
RB =
3
n =
3
X =
1
2
-1
与题 1 中逆矩阵计算所得结果相同,检验通过
四.实验体会:
通过实验我掌握了消元法解方程的一些基本算法以及用matlab实现矩阵的几种基本计算。
对MATLAB软件有了更深的了解。
同时,在实验中,在输入向量b=[8;21;1;]时错误的输成b=[8;21;1;],致使程序不能运行,无法得到预期的实验结果,后经多番仔细查找,最终发现分号应为英文输入法,以后在编程时,一定更加认真仔细,不犯同样的错误!
如有侵权请联系告知删除,感谢你们的配合!。