高斯消元法MATLAB实现(可编辑修改word版)

合集下载

MATLAB中高斯消去法

MATLAB中高斯消去法

clear allclcdata_fname='E:\学习\计算方法\fun005.dat';%读入二进制数据a=5;%输入头文件%打开方式为二进制打开file_id=fopen(data_fname,'rb');%读取二进制文件head_data=fread(file_id,a,'uint');%头文件第一个元素为该文件的标识。

id=dec2hex(head_data(1));%头文件第二个元素为数据的版本号,十六进制ver=dec2hex(head_data(2));%其中101表示稀疏矩阵;102表示条状对角阵;201表示压缩格式下的条状对角矩阵%头文件第三个元素表示为该方程的介绍n=head_data(3);p=head_data(4);%条状矩阵下带宽q=head_data(5);%条状矩阵上带宽jump_distance=20;%跳过字节数fseek(file_id,jump_distance,'bof');[raw_data,count]=fread(file_id,inf,'float');%读入所有元素fclose(file_id);%关闭二进制文件%形成b矩阵b=raw_data(count-n+1:count); A1=raw_data(1:count-n);A=(reshape(A1,p+q+1,n))';%GAUSS算法¨%消去过程for k=1:n-1t=A(k+1,p)/A(k,p+1);for j=1:p+qA(k+1,j)=A(k+1,j)-t*A(k,j+1);endb(k+1)=b(k+1)-t*b(k);end%回待过程x(n)=b(n)/A(n,p+1);x(n-1)=(b(n-1)-A(n-1,p+2)*x(n))/A(n-1,p+1);x(n-2)=(b(n-2)-A(n-2,p+3)*x(n)-A(n-2,p+2)*x(n-1))/A(n-2,p+1);for i=n-3:-1:1s=b(i);s=s-A(i,3)*x(i+1)-A(i,4)*x(i+2)-A(i,5)*x(i+3);x(i)=s/A(i,p+1);end%输出结果fprintf('The solution of the equation set is:\n');for k=1:10fprintf('x[%d]=%f',k,x(k));fprintf('x[%d]=%f',k+10,x(k+10));fprintf('x[%d]=%f',k+n-20,x(k+n-20)); fprintf('x[%d]=%f\n',k+n-10,x(k+n-10)); end。

“线性方程组高斯消去法”实验报告(内含matlab程序)

“线性方程组高斯消去法”实验报告(内含matlab程序)

实验一实验报告一、实验名称:线性方程组高斯消去法。

二、实验目的:进一步熟悉理解Guass 消元法解法思路,提高matlab 编程能力。

三、实验要求:已知线性方程矩阵,利用软件求解线性方程组的解。

四、实验原理:消元过程:设0)0(11≠a ,令乘数)0(11)0(11/a a m i i -=,做(消去第i 个方程组的i x )操作1i m ×第1个方程+第i 个方程(i=2,3,.....n )则第i 个方程变为1)1(2)1(2...i n in i b x a x a =++ 这样消去第2,3,。

,n 个方程的变元i x 后。

原线性方程组变为:⎪⎪⎪⎩⎪⎪⎪⎨⎧=++=++=++)1()1(2)1(2)1(2)1(22)1(22)0(1)0(11)0(11... . .... ...n n nn n n n n n b x a x a b x a x a b x a x a 这样就完成了第1步消元。

回代过程:在最后的一方程中解出n x ,得:)1()1(/--=n nn n n n a b x再将n x 的值代入倒数第二个方程,解出1-n x ,依次往上反推,即可求出方程组的解:其通项为3,...1-n 2,-n k /)()1(1)1()1(=-=-+=--∑k kk n k j j k kj k k k a x a bx五、实验内容:A=[1 1 1;0 4 -1;2 -2 1];%ϵÊý¾ØÕó b=[6 5 1]'%³£ÊýÏînum=length(b)for k=1:num-1for i=k+1:numif A(k,k)~=0l=A(i,k)/A(k,k); A(i,:)=A(i,:)-A(k,:).*l; b(i)=b(i)-b(k)*l; endendendAb%»Ø´úÇóxx(num)=b(num)/A(num,num); for i=num-1:-1:1sum=0;for j=i+1:numsum=sum+A(i,j)*x(j); endx(i)=(b(i)-sum)/A(i,i); endx六、实验结果:A =1 1 1 0 4 -10 0 -2b =65-6x =1 2 3。

(完整word版)Gauss消去法Matlab

(完整word版)Gauss消去法Matlab

实验一 用列主元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。

gauss列主元素消去法matlab

gauss列主元素消去法matlab

高斯列主元素消去法是一种解线性方程组的常用方法,特别在数值分析和线性代数中应用广泛。

在Matlab中,我们可以使用该方法来解决大规模的线性方程组,包括矩阵的求解和矩阵的反转。

一、高斯列主元素消去法的基本原理高斯列主元素消去法是一种基于矩阵消元的方法,它通过一系列的矩阵变换将原始的线性方程组转化为上三角形式,然后再进行回代求解。

这个方法的核心就是通过矩阵的变换来简化原始的线性方程组,使得求解过程更加简单高效。

在Matlab中,我们可以利用矩阵运算和函数来实现高斯列主元素消去法,如`lu`分解函数和`\"`运算符等。

通过这些工具,我们能够快速地求解各种规模的线性方程组并得到准确的结果。

二、高斯列主元素消去法在Matlab中的实现在Matlab中,我们可以通过调用`lu`函数来实现高斯列主元素消去法。

该函数返回一个上三角矩阵U和一个置换矩阵P,使得PA=LU。

通过对U进行回代求解,我们可以得到线性方程组的解。

除了`lu`函数之外,Matlab还提供了一些其他的函数和工具来帮助我们实现高斯列主元素消去法,比如`\"`运算符和`inv`函数等。

通过这些工具的组合使用,我们能够更加灵活地进行线性方程组的求解,并且可以方便地处理特殊情况和边界条件。

三、高斯列主元素消去法的应用与局限性高斯列主元素消去法在实际应用中具有广泛的适用性,特别是对于大规模的线性方程组或者稀疏矩阵的求解。

通过Matlab中的工具和函数,我们可以快速地求解各种规模的线性方程组,并得到高精度的数值解。

然而,高斯列主元素消去法也存在一些局限性,比如对于奇异矩阵或者接近奇异矩阵的情况时,该方法的求解精度可能会下降。

在实际应用中,我们需要结合具体的问题和矩阵特性来选择合适的求解方法,以确保得到准确的结果。

四、个人观点和总结作为一种经典的线性方程组求解方法,高斯列主元素消去法在Matlab 中具有较好的实现和应用效果。

通过对其原理和实现细节的深入理解,我们能够更加灵活地应用该方法,并且能够更好地理解其适用性和局限性。

(完整)高斯消去法解线性方程的Matlab程序

(完整)高斯消去法解线性方程的Matlab程序

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« .。

高斯消元法,列主元素消元法及LU分解法的matlab程序

高斯消元法,列主元素消元法及LU分解法的matlab程序

§2.2.1高斯消元法的MATLAB程序f unction [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,所以此方程组无解.')returnendif RA==RBif RA==ndisp('请注意:因为RA=RB=n,所以此方程组有唯一解.')X=zeros(n,1); C=zeros(1,n+1);for p= 1:n-1for k=p+1:nm= B(k,p)/ B(p,p);B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1);endendb=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n);for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);endelsedisp('请注意:因为RA=RB<n,所以此方程组有无穷多解.') endend运行命令及结果a=[2.51 1.48 4.53;1.48 0.93 -1.30 ;2.68 3.04 -1.48];b=[0.05;1.03;-0.53];[RA,RB,n,X] =gaus (A,b)请注意:因为RA=RB=n,所以此方程组有唯一解.RA =3RB =3n =3X =1.4531-1.5892-0.2749§2.2.2 列主元素消元法的MATLAB程序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,所以此方程组无解.')returnendif RA==RBif RA==ndisp('请注意:因为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:nm= B(k,p)/ B(p,p);B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1);endendb=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n);for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);endelsedisp('请注意:因为RA=RB<n,所以此方程组有无穷多解.') endend运行命令及结果a=[2.51 1.48 4.53;1.48 0.93 -1.30 ;2.68 3.04 -1.48];b=[0.05;1.03;-0.53];[RA,RB,n,X]=liezhu(A,b)请注意:因为RA=RB=n,所以此方程组有唯一解.RA =3RB =3n =3X =1.4531-1.5892-0.2749§2.2.3 LU分解法的MATLAB程序function hl=zhjLU(A)[n n] =size(A); RA=rank(A);if RA~=ndisp('请注意:因为A的n阶行列式hl等于零,所以A不能进行LU分解.A的秩RA如下:'), RA,hl=det(A);returnendif RA==nfor p=1:nh(p)=det(A(1:p, 1:p));endhl=h(1:n);for i=1:nif h(1,i)==0disp('请注意:因为A的r阶主子式等于零,所以A不能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:'), hl;RAreturnendendif h(1,i)~=0disp('请注意:因为A的各阶主子式都不等于零,所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:')for j=1:nU(1,j)=A(1,j);endfor k=2:nfor i=2:nfor j=2:nL(1,1)=1;L(i,i)=1;if i>jL(1,1)=1;L(2,1)=A(2,1)/U(1,1);L(i,1)=A(i,1)/U(1,1);L(i,k)=(A(i,k)- L(i,1:k-1)*U(1:k-1,k))/U(k,k);elseU(k,j)=A(k,j)-L(k,1:k-1)*U(1:k-1,j);endendendendhl;RA,U,Lendend运行命令及结果a=[2.51 1.48 4.53;1.48 0.93 -1.30 ;2.68 3.04 -1.48];hl=zhjLU(A)请注意:因为A的各阶主子式都不等于零,所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:RA =3U =2.5100 1.4800 4.53000 0.9300 -3.97110 0 -0.0837L =1.0000 0 00.5896 1.0000 01.0677 1.5696 1.0000hl =2.5100 0.1439 13.6410>> U=[2.5100 1.4800 4.53000 0.9300 -3.97110 0 -0.0837];>>L= [1.0000 0 00.5896 1.0000 01.0677 1.5696 1.0000];>> b=[0.05;1.03;-0.53];U1=inv(U); L1=inv(L); X=U1*L1*b,x=A\bX =-111.8440110.953125.7324x =1.4531-1.5892-0.2749例2.1: 用高斯消元法求解下面的非齐次线性方程组。

gauss在matlab的程序编写

gauss在matlab的程序编写

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经验证是正确的,太长了只节了部分。

MATLAB之GAUSS消元法解线性方程组

MATLAB之GAUSS消元法解线性方程组
matlab之gauss消元法解线性方程组 matlab之gauss消元法解线性方程组 1.gauss消元法 function x=delgauss(a,b) %gauss消去法 [n,m]=size(a); nb=length(b); det=1;%存储行列式值 x=zeros(n,1); for k=1:n-1 for i=k+1:n if a(k,k)==0 return end m=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-m*a(k,j); end b(i)=b(i)-m*b(k); end det=det*a(k,k);%计算行列式 end det=det*a(n,n); for k=n:-1:1%回代求解 for j=k+1:n b(k)=b(k)-a(k,j)*x(j); end x(k)=b(k)/a(k,k); end example: >>a=[1.0170-0.00920.0095;-0.00920.99030.0136;0.00950.0136 0.9898]; >>b=[101]'; >>x=delgauss(a,b) x= 0.9739 -0.0047 1.0010 2.列主元gauss消去法: function x=detgauss(a,b) %gauss列主元消去法 [n,m]=size(a); nb=length(b); det=1;%存储行列式值 x=zeros(n,1); for k=1:n-1 amax=0;%选主元 for i=k:n if abs(a(i,k))>amax amax=abs(a(i,k));r=i; end end if amax<1e-10 return; end if r>k%交换两行 for j=k:n z=a(k,j);a(k,j)=a(r,j);a(r,j)=z; end z=b(k);b(k)=b(r);b(r)=z;det=-det; end for i=k+1:n%进行消元 m=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-m*a(k,j); end b(i)=b(i)-m*b(k); end det=det*a(k,k); end det=det*a(n,n); for k=n:-1:1%回代求解 for j=k
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

《数值分析》实验报告
一、实验目的与要求
1.掌握高斯消去法的基本思路和迭代步骤;
2.培养编程与上机调试能力。

二、实验内容
1.编写用高斯消元法解线性方程组的MATLAB程序,并求解下面的线性方程组,然后用逆矩阵解方程组的方法验证.
⎧0.101x1+ 2.304x2+ 3.555x3=1.183 ⎧5x1+ 2x2+x3=8(1) ⎪-1.347x + 3.712x + 4.623x = 2.137 (2)⎪2x + 8x - 3x = 21
⎨ 1 2 3 ⎨ 1 2 3
⎪-2.835x +1.072x + 5.643x = 3.035 ⎪x - 3x - 6x =1
⎩ 1 2 3 ⎩ 1 2 3
2.编写用列主元高斯消元法解线性方程组的MATLAB程序,并求解下面的线性方程组,然后用逆矩阵解方程组的方法验证.
⎧0.101x1+ 2.304x2+ 3.555x3=1.183 ⎧5x1+ 2x2+x3=8(1) ⎪-1.347x + 3.712x + 4.623x = 2.137 (2)⎪2x + 8x - 3x = 21
⎨ 1 2 3 ⎨ 1 2 3
⎪-2.835x +1.072x + 5.643x = 3.035 ⎪x - 3x - 6x =1
⎩ 1 2 3 ⎩ 1 2 3 三.MATLAB 计算源程序
1.用高斯消元法解线性方程组AX =b 的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 b 的MATLAB程序
输入的量:系数矩阵 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 =
RB = n = X = 3
3
3
-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;], 致使程序不能运行,无法得到预期的实验结果,后经多番仔细查找,最终发现分号应为英文输入法,以后在编程时,一定更加认真仔细,不犯同样的错误!。

相关文档
最新文档