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

实验题目：线性代数求解方程组一、实验目的1. 理解线性代数中方程组的求解方法。

2. 掌握利用计算机求解线性方程组的算法。

3. 熟悉数学软件（如MATLAB、Python等）在数学问题中的应用。

二、实验内容本次实验主要利用数学软件求解线性方程组。

线性方程组是线性代数中的一个基本问题，其求解方法有很多种，如高斯消元法、矩阵求逆法等。

本实验以高斯消元法为例，利用MATLAB软件求解线性方程组。

三、实验步骤1. 编写高斯消元法算法程序。

2. 输入方程组的系数矩阵和常数项。

3. 调用程序求解方程组。

4. 输出解向量。

四、实验代码及分析1. 高斯消元法算法程序```matlabfunction x = gaussElimination(A, b)[n, m] = size(A);assert(n == m, 'The matrix A must be square.');assert(n == length(b), 'The length of b must be equal to the number of rows in A.');% 初始化解向量x = zeros(n, 1);% 高斯消元for i = 1:n-1% 寻找最大元素[~, maxIdx] = max(abs(A(i:n, i)));maxIdx = maxIdx + i - 1;% 交换行A([i, maxIdx], :) = A([maxIdx, i], :);b([i, maxIdx]) = b([maxIdx, i]);% 消元for j = i+1:nfactor = A(j, i) / A(i, i);A(j, i:n) = A(j, i:n) - factor A(i, i:n); b(j) = b(j) - factor b(i);endend% 回代求解for i = n:-1:1x(i) = (b(i) - A(i, i+1:n) x(i+1:n)) / A(i, i); endend```2. 输入方程组的系数矩阵和常数项```matlabA = [2, 1, -1; 1, 2, 1; -1, 1, 2];b = [8; 5; 2];```3. 调用程序求解方程组```matlabx = gaussElimination(A, b);```4. 输出解向量```matlabdisp('解向量为：');disp(x);```五、实验结果与分析实验结果：```解向量为：2-13```实验分析：通过高斯消元法，我们成功求解了给定的线性方程组。

Gauss 列主消元法与Gauss-Seidel 迭代法一、实验目的：通过本次实验，学会用Gauss 列主消元法对线性方程组的解的情况进行判断；并且学会利用Gauss-Seidel 迭代法求线性方程组的解。

通过本次实验，加强用Matlab 编程的能力和上机调试能力。

二、实验基本原理:Gauss 列主消元法:根据线性方程组任意交换两个方程的次序，方程组的同解性不变，且解的分量次序也不变。

于是，第k 步在顺序消元法进行之前，从()k A 的第k 列元素()k kk a ，()k k k a 1+，…()k nk a 中选取绝对值最大者，并记录所在行，即()()k ik n i k k k i a a k ≤≤=max ，记k i l =,如果k l ≠，则交换矩阵()()[]k k b A 的第k 行与第l 列所有对应的元素。

然后，再进行第k 步顺序消元算法。

Gauss-Seidel 迭代法：在迭代的每一步计算过程中用()k x 的全部分量来计算()1+k x 的所有分量，进而在计算第i 个分量()1+k i x 时，已经计算出最新的分量()11+k x …()11++k i x 没有被利用。

因而由线性方程组的矩阵形式f x G x J +=（1）,任取初始解向量()()()()()Tn x x x x 002010,,, =代入(1)式右端的即可得到迭代格式： ()()()ii i j n i j k j ij k jij i k i a x a x a b x /11111⎪⎪⎭⎫ ⎝⎛--=∑∑-=+=++()n i ,,2,1 =。

三、实验步骤：（一）Gauss 列主消元法：1、编写Gauss 列主消元法以及Gauss 回代法的Matlab 程序。

具体步骤如下：(1)输入A,b(2)计算B=[A b]; n=length(b); RA=rank(A);RB=rank(B);zhica=RB-RA;(3)如果zhica>0,方程组无解如果if RA==RB ，if RA==n ，方程有唯一解且 X=backsub(B(1:n,1:n), B(1:n,n+1));否则RA=RB<n ,方程有无穷多解。

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

第三讲 Matlab 求解代数方程组理论介绍：直接法+迭代法，简单介绍相关知识和应用条件及注意事项 软件求解：各种求解程序讨论如下表示含有n 个未知数、由n 个方程构成的线性方程组：11112211211222221122n n n n n n nn n na x a x a xb a x a x a x b a x a x a x b +++=⎧⎪+++=⎪⎨⎪⎪+++=⎩ （1）一、直接法 1.高斯消元法：高斯消元法的基本原理： 在（1）中设110,a ≠将第一行乘以111,k a a -加到第(2,3,,),k k n = 得： (1)(1)(1)(1)11112211(2)(1)(2)22112(2)(2)(2)22n n n n n nn n n a x a x a x b a x a x b a x a x b ⎧+++=⎪++=⎪⎨⎪⎪++=⎩（2）其中(1)(1)1111,.k k aa b b ==再设(2)220,a ≠将（2）式的第二行乘以(2)2(2)22,(3,,)k a k n a -= 加到第k 行，如此进行下去最终得到：(1)(1)(1)(1)11112211(2)(1)(2)22112(1)(1)(1)1,111,1()()n n n n n n n n n n n n n n n n nn n n a x a x a x b a x a x b a x a x b a x b --------⎧+++=⎪++=⎪⎪⎨⎪+=⎪⎪=⎩（3） 从（3）式最后一个方程解出n x ，代入它上面的一个方程解出1n x -，并如此进行下去，即可依次将121,,,,n n x x x x - 全部解出，这样在()0(1,2,,)k kk a k n ≠= 的假设下，由上而下的消元由下而上的回代，构成了方程组的高斯消元法. 高斯消元法的矩阵表示：若记11(),(,,),(,,)T T ij n n n n A a x x x b b b ⨯=== ，则（1）式可表为.Ax b =于是高斯消元法的过程可用矩阵表示为：121121.n n M M M Ax M M M b --=其中：(1)21(1)111(1)1(1)11111n a a M a a ⎛⎫ ⎪ ⎪- ⎪=⎪ ⎪ ⎪ ⎪- ⎪⎝⎭ (2)32(2)222(2)2(2)221111n a a M a a ⎛⎫⎪⎪ ⎪-⎪=⎪ ⎪ ⎪⎪- ⎪⎝⎭高斯消元法的Matlab 程序： %顺序gauss 消去法,gauss 函数 function[A,u]=gauss(a,n) for k=1:n-1%消去过程 for i=k+1:n for j=k+1:n+1%如果a(k,k)=0,则不能削去 if abs(a(k,k))>1e-6 %计算第k 步的增广矩阵 a(i,j)=a(i,j)-a(i,k)/a(k,k)*a(k,j); else%a(k,k)=0,顺序gauss 消去失败 disp (‘顺序gauss 消去失败‘)； pause; exit; end end end end%回代过程 x(n)=a(n,n+1)/a(n,n); for i=n-1:-1:1 s=0; for j=i+1:n s=s+a(i,j)*x(j); endx(i)=(a(i,n+1)-s)/a(i,i); end%返回gauss 消去后的增广矩阵 A=triu(a); %返回方程组的解 u=x ;练习和分析与思考： 用高斯消元法解方程组：12345124512345124512452471523814476192536x x x x x x x x x x x x x x x x x x x x x x ++++=⎧⎪+++=⎪⎪++++=⎨⎪+++=⎪+++=⎪⎩2.列主元素消元法在高斯消元法中进行到第k 步时，不论()k ik a 是否为0，都按列选择()||(,,)k ik a i k n = 中最大的一个，称为列主元，将列主元所在行与第k 行交换再按高斯消元法进行下去称为列主元素消元法。

陆韶琦 3110000441说明：本程序需要将方程组的增广矩阵保存在程序.m文件相同的目录下，命名为mytxt.txt.运行程序时输入系数矩阵的行数即可进行方程组的求解。

1.列主元高斯消元法流程图：程序代码%列主元高斯消元法¨%控制输入增广矩阵行数N=input('what is the rank of the matrix:');%从mytxt.txt文件中读入增广矩阵f=fopen('mytxt.txt','r');A=fscanf(f,'%g',[N+1 N]);fclose(f);A=A';disp('A|b=');disp(A)%列主元的寻找for i=1:Nmax=abs(A(i,i));for j=i+1:Nif A(j,i)>maxflag=j;max=A(j,i);endend%判断矩阵是否满秩if A(flag,i)==0disp('the martrix is not full rank');flag=0;break;end;%实现两行的交换for k=i:N+1B=A(flag,k);A(flag,k)=A(i,k);A(i,k)=B;end%实现消元for kh=i+1:Nm=-A(kh,i)/A(i,i);A(kh,i)=0;for kl=i+1:N+1A(kh,kl)=A(kh,kl)+m*A(i,kl);endendendif flag~=0disp('U|b=');%输出此时得到的上三角矩阵disp(A);X=zeros(N,1);%回代得到解for i=N:-1:1for j=i-1:-1:1m=-A(j,i)/A(i,i);A(j,N+1)=A(j,N+1)+m*A(i,N+1);endX(i,1)=A(i,N+1)/A(i,i);end disp('X=');disp(X);end运行equation.m文件，结果如下：2.程序代码：%列主元三角分解法N=input('what is the rank of the matrix:');f=fopen('mytxt.txt','r');A=fscanf(f,'%g',[N+1 N]);fclose(f);A=A';disp('A|b=');disp(A);%形成Li1,U1jmax=A(1,1);flag1=1;for i=2:Nif A(i,1)>max;max=A(i,1);flag=i;endend%判断矩阵是否满秩if A(flag,1)==0flag1=0;disp('the matrix is not full rank');elsefor i=1:N+1m=A(1,i);A(1,i)=A(flag,i);A(flag,i)=m;endfor i=2:NA(i,1)=A(i,1)/A(1,1);end% 计算aij-∑Lij*Ujkfor k=2:Nfor i=k:Nfor j=1:k-1 A(i,k)=A(i,k)-A(i,j)*A(j,k);endend%寻找主元max=A(k,k);for i=k:Nif A(i,k)>max;max=A(i,k);flag=i;endend%判断矩阵是否满秩if A(flag,i)==0flag1=0;disp('the matrix is not full rank');break;else%换行for i=1:N+1m=A(k,i);A(k,i)=A(flag,i);A(flag,i)=m;end%形成Uijfor i=k+1:N+1for j=1:k-1A(k,i)=A(k,i)-A(k,j)*A(j,i);endend%形成Lijfor i=k+1:NA(i,k)=A(i,k)/A(k,k);endendendif flag1~=0disp('LU=');disp(A);X=zeros(N,1);%回代得到解for i=N:-1:1for j=i-1:-1:1m=-A(j,i)/A(i,i);A(j,N+1)=A(j,N+1)+m*A(i,N+1);endX(i,1)=A(i,N+1)/A(i,i);enddisp('X=');disp(X);endend运行equation2.m文件，结果如下：编程感想：本程序利用matlab平台，实现了线性方程组的求解，但是这是很底层的方法，如果将其中的语句稍加修改，就可以用c编译器进行计算。

摘要本文叙述了Gauss 顺序消元法解线性方程的算法思想以及其求解过程，同时简要叙述了Gauss 主元素消元法以及Gauss 全主元消元法。

紧接着给出了Gauss Seidel -迭代法的算法思想，本文给出了这三个消元方法以及一个迭代法的算法流程图，由于全主元消元法是前两个算法的基础上改进而来，故本文采用第三种方法进行编程计算，前两种方法不再重复编程，然后给出一个实例的计算结果，运行时间，在文章最后分析该实例的计算结果，针对同一实例，又采用Gauss Seidel -方法编程实现，然后对结果进行分析和对比。

最后给出了本人在编程时遇到的一些问题和解决办法。

关键词：Gauss 顺序消元法 Gauss 主元素消元法 Gauss 全主元消元法一、算法的简要描述1.1Gauss 顺序消元法Gauss 消元法在中学里已经学习过，其方法实质，就是运用初等变换，将线性方程组Ax b =转化为同解的上三角矩阵方程组1Ux L b -=（1.1.1）其中，U 为上三角矩阵，L 为下三角矩阵。

然后对式（1.1.1）进行回代求解，即得方程组的解。

手算的过程是非常清楚的，现在需回答的是计算机求解，如何实现上述计算过程。

设线性方程组为1111221331121122223322112233n n n n n n n nn n na x a x a x a xb a x a x a x a x b a x a x a x a x b +++⋅⋅⋅+=⎧⎪+++⋅⋅⋅+=⎪⎨⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⎪⎪+++⋅⋅⋅+=⎩ 写成矩阵形式为1112111212222221222m m m n n a a a x b aa a xb a a a x b ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦（1.1.2）设线性方程组如上式所示，记(1)A A =，(1)b b =，与是增广矩阵具有形式(1)(1)[][]A b A b =，此时方程组为(1)(1)A x b =。

在这章中我们要学习线性方程组的直接法，特别是适合用数学软件在计算机上求解 的方法•3.1 方程组的逆矩阵解法及其MATLA 程序3.1.3 线性方程组有解的判定条件及其 MATLAB?序判定线性方程组 A m n X 二b 是否有解的MATLA 醉序fun ction [RA,RB, n]=jiepb(A,b) B=[A b]; n=le ngth(b); RA=ra nk(A); RB=ra nk(B);zhica=RB-RA;if zhica>0,disp （'请注意：因为RA~=RB ，所以此方程组无解.'） return endif RA==RB if RA==ndisp （'请注意：因为RA=RB=n ，所以此方程组有唯一解 .'） elsedisp （'请注意：因为RA=RB<n ，所以此方程组有无穷多解 .'）end解 在MATLABT 作窗口输入程序>> A=[2 3 -1 5;3 1 2 -7;4 1 -3 6;1 -2 4 -7]; b=[ 0; 0; 0; 0]; [RA,RB, n]=jiepb(A,b)运行后输出结果为请注意：因为 RA=RB=n ，所以此方程组有唯一解 RA = 4,RB =4, n =4在MATLABT 作窗口输入>>X=A\b,运行后输出结果为 X =(0 0 0 0)'.(2) 在MATLABT 作窗口输入程序>> A=[3 4 -5 7;2 -3 3 -2;4 11 -13 16;7 -2 1 3];b=[ 0; 0; 0; 0];[RA,RB, n]=jiepb(A,b)运行后输出结果请注意：因为RA=RB<n ，所以此方程组有无穷多解• RA =2,RB =2, n =4(3)在MATLABT 作窗口输入程序例 3.1.4判断下列线性 方程组解 军的情况”2X 1+ 3X 2 - X 3+ 5X 4 =0, (1)丿3x1十x 2十2X 3—7 X 4 ",(2)4X 1 + x ? - 3 X 3+ 6x 4 =0,J1 --2X 2 + 4X 3—7 X 4 =0;<4x 1+ 2X 2 -X 3 =2,⑶ 丿 3X 1 —X 2 + 2 X 3 = 10,(4)11x 1 +3x 2 =8.如果有唯一解，则用表 3-2方法求解3x 1 + 4X 2 - 5X3 + 7X 4 = 0,2x 1 -3x 2 +3x 3 - 2x 4 =0, 4x 111x 2「13X 316x 4 = 0,7x 1 - 2X 2 X 3 3X 4 = 0; 2x y -z w=1,« 4x + 2y - 2z + w = 2, 2x + y — z — w=1.end>> A=[4 2 -1;3 -1 2;11 3 0]; b=[2;10;8]; [RA,RB, n]=jiepb(A,B)运行后输出结果请注意：因为RA~=RB，所以此方程组无解•RA =2,RB =3 ,n =3(4) 在MATLABT作窗口输入程序>> A=[2 1-1 1;4 2 -2 1;2 1-1-1];b=[1; 2; 1]; [RA,RB, n]=jiepb(A,b)运行后输出结果请注意：因为RA=RB <n,所以此方程组有无穷多解.RA =2,RB =2, n =33.2 三角形方程组的解法及其MATLA程序3.2.2 解三角形方程组的MATLAB?序解上三角形线性方程组AX二b的MATLAB^序fun ction [RA,RB ,n, X]=sha ngsa n(A,b)B=[A b]; n=le ngth(b); RA=ra nk(A);RB=ra nk(B);zhica=RB-RA;if zhica>0,disp('请注意：因为RA~=RB，所以此方程组无解.') returnendif RA==RBif RA==ndisp('请注意：因为RA=RB=n，所以此方程组有唯一解.') X=zeros(n,1);X(n)=b(n)/A(n,n);for k=n-1:-1:1X(k)=(b(k)-sum(A(k,k+1: n)*X(k+1: n)))/A(k,k);endelsedisp('请注意：因为RA=RB<n，所以此方程组有无穷多解.') endend例3.2.2 用解上三角形线性方程组的MATLAB程序解方程组馭- x2+ 2 x3+ 3x4 = 20,* -2X2 +7X3 -4X4 = -7,| 6x3 +5X4 =4,J 3X4 = 6.解在MATLABT作窗口输入程序>>A=[5 -1 2 3;0 -2 7 -4;0 0 6 5;0 0 0 3];b=[20; -7; 4;6];[RA,RB, n,X]=sha ngsa n(A,b)运行后输出结果请注意：因为RA=RB=n，所以此方程组有唯一解.RA = RB =4, 4,3.3高斯(GausS 消元法和列主元消元法及其MATLA 程序3.3.1高斯消元法及其MATLAB?序用高斯消元法解线性方程组AX 二b 的MATLA 龍序fun ctio n [RA,RB ,n, X]=gaus(A,b) B=[A b]; n=le ngth(b); RA=ra nk(A); RB=ra nk(B);zhica=RB-RA; if zhica>0,disp('请注意：因为RA~=RB ，所以此方程组无解.') return endif RA==RB if RA==ndisp('请注意：因为RA=RB=n ，所以此方程组有唯一解.') X=zeros (n ,1); C=zeros(1, n+1); for p= 1:n-1for 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 endb=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); end elsedisp('请注意：因为RA=RB<n ，所以此方程组有无穷多解 .')end end 例3.3.2用高斯消元法和 MATLAB^序求解下面的非齐次线性方程组，并且用逆矩阵解方程组的方法验证.乙一x 2x 3 — 3X 4 = 1,-X 2 - X3 + X 4 = 0, 2X 1 - 2x ? - 4X 3 ' 6X 4 = — 1,I x 1 - 2X 2 - 4X 3 x 4 - -1. 解在MATLABT 作窗口输入程序>> A=[1 -1 1 -3; 0 -1 -1 1;2 -2 -4 6;1 -2 -4 1]; b=[1;0; -1;-1]; [RA,RB, n,X] =gaus (A,b)运行后输出结果请注意：因为RA=RB=n ，所以此方程组有唯一解. X =0 -0.5000 0.5000 3.3.2 列主元消元法及其MATLAB?序用列主元消元法解线性方程组 AX二b 的MATLAB^序fun ction [RA,RB ,n, X]=liezhu(A,b) B=[A b]; n=le ngth(b); RA=ra nk(A); RB=ra nk(B);zhica=RB-RA; if zhica>0,disp('请注意：因为RA~=RB ，所以此方程组无解.') return endif RA==RB if RA==ndisp('请注意：因为RA=RB=n ，所以此方程组有唯一解.') X=zeros (n ,1);4,X =[2.4 -4.0 -1.02.0]RA =RB =4n =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例333 用列主元消元法解线性方程组的MATLAB程序解方程组—X 2 —X3 + X4 = 0,X1 - X2 + X3 -3X4 = 1,_2X2_4X36X4 = -1,x1 _ 2x2 _ 4x3 x4 = _1.解在MATLABT作窗口输入程序>> A=[0 -1 -1 1;1 -1 1 -3;2 -2 -4 6;1 -2 -4 1];b=[0;1;-1;-1]; [RA,RB, n,X]=liezhu(A,b)运行后输出结果请注意：因为RA=RB=n，所以此方程组有唯一解•RA = 4 , RB = 4 , n = 4 , X =[0 - 0.5 0.5 0]'3.4 LU分解法及其MATLA程序3.4.1判断矩阵LU分解的充要条件及其MATLAB?序判断矩阵A能否进行LU分解的MATLAB^序fun ctio n hl=pdLUfj(A)[n n] =size(A); RA=ra nk(A); if RA~=ndisp('请注意：因为A的n阶行列式hl等于零，所以A不能进行LU分解.A的秩RA如下：'),RA,hl=det(A); returnend if RA==nfor p=1: n,h(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;RA, returnend end if h(1,i)~=0disp('请注意：因为A 的各阶主子式都不等于零，所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值 hl 依次如下：')hl;RA end end>> A=[1 2 3;1 12 7;4 5 6];hl=pdLUfj(A)运行后输出结果为请注意：因为A 的各阶主子式都不等于零，所以 各阶顺序主子式值hl 依次如下：RA = 3 ， hl = 1 10 -48(2 )在MATLA 工作窗口输入程序>> A=[1 2 3;1 2 7;4 5 6];hl=pdLUfj(A) 运行后输出结果为请注意：因为A 的r 阶主子式等于零，所以 A 不能进行LU 分解.A 的秩RA 和各阶 顺序主子式值hl 依次如下：RA = 3 ， hl =1 0 12(3 )在MATLA 工作窗口输入程序>> A=[1 2 3;1 2 3;4 5 6];hl=pdLUfj(A) 运行后输出结果为请注意：因为A 的n 阶行列式hl 等于零，所以A 不能进行LU 分解.A 的秩RA 如下RA = 2 ， hl = 03.4.2 直接LU 分解法及其MATLA 程序将矩阵A 进行直接LU 分解的MATLAB^序fun ction hl=zhjLU(A)[n n] =size(A); RA=ra nk(A); if RA~=ndisp('请注意：因为A 的n 阶行列式hl 等于零，所以A 不能进行LU 分解.A的秩 RA 如下:’),RA,hl=det(A);returnend if RA==n for 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;RAreturn end endif h(1,i)~=0disp('请注意：因为A 的各阶主子式都不等于零，所以 A 能进行LU 分解.A123、1 2 3 'q 2 3、 1 12 7 ；(2) 1 2 7 ；(3) 12 3 I 45 6I 4 5 6>1<5 6判断下列矩阵能否进行 LU 分解，并求矩阵的秩(1) 解 (1)在MATLA 工作窗口输入程序例341A 能进行LU 分解.A 的秩RA 和的秩RA 和各阶顺序主子式值hl 依次如下：') for j=1: nU(1,j)=A(1,j); endfor k=2: n for i=2: nfor j=2: n L(1,1)=1;L(i,i)=1; if i>j L(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); else U(k,j)=A(k,j)-L(k,1:k-1)*U(1:k-1,j); end end end endhl;RA,U ,L end end例343 用矩阵进行直接 LU 分解的MATLABq 0 20 1 0 A =1 2 4 (0 1 0 解 在MATLA 工作窗口输入程序>> A=[1 0 2 0;0 1 0 1;1 2 4 3;0 1 0 3];hl=zhjLU(A)运行后输出结果请注意：因为A 的各阶主子式都不等于零， 所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值 hl 依次如下：RA = =4L = 1 0 0 0 U = 1 0 2 0 0 1 0 0 0 1 0 1 1 2 1 0 0 0 2 1 0 1 0 1 0 0 0 2 hl = 1 1 2 43.4.4 判断正定对称矩阵的方法及其 MATLAB?序 判断矩阵A 是否是正定对称矩阵的 MATLAB^序function hl=zddc(A)[n n] =size(A); for p=1: nh(p)=det(A(1:p, 1:p)); endhl=h(1: n);zA=A'; for i=1: nif h(1,i)<=0disp('请注意：因为A 的各阶顺序主子式hl 不全大于零，所以A 不是正 定的.A的转置矩阵zA 和各阶顺序主子式值 hl 依次如下：'),hl;zA, returnend endif h(1,i)>0 disp('请注意：因为A 的各阶顺序主子式hl 都大于零，所以A 是正定的.A 的 转置矩阵zA 和各阶顺序主子式值hl 依次如下：')hl;zA end程序分解矩阵 01 . 33丿>> A=[0.1 2 3 4;-1 2 -3 4;11 21 13 41;5 7 8 9];hl=运行后输出结果请注A 不是对称矩阵请注意：因为A 的各阶顺序主子式 hl 不全大于零，所以A 不是正定的.A 的转 置矩阵zA 和各阶顺序主子式值 hl 依次如下：zA =1/10 -1 11 52 2 2173 -3 13 84 4 41 9hl = :1/10 11/5 -1601/10 3696/5因此，A 即不是正定矩阵，也不是对称矩阵 •（2）在MATLA 工作窗口输入程序>> A=[1 -1 2 1;-1 3 0 -3;2 0 9 -6;1 -3 -6 19],hl=zddc （A ）运行后输出结果A = 1 -1 2 1 -1 3 0-3 2 09 -6 1-3 -6 19 请注意： A 是对称矩阵请注意：因为A 的各阶顺序主子式hl 都大于零，所以A 是正定的.A 的转置矩阵zA和各阶顺序主子式值 hl 依次如下： zA = 1 -1 2 1 -1 3 0-3 2 09 -6 1-3 -6 19 hl = 1 2 6 24(3) 在MATLA 工作窗口输入程序>> A=[1/sqrt(2) -1/sqrt(2) 0 0; -1/sqrt(2) 1/sqrt (2) 0 0; 00 1/sqrt(2) -1/sqrt(2); 0 0 -1/sqrt(2) 1/sqrt(2)], hl= zddc (A)运行后输出结果A= 985/1393 -985/1393 0 0 -985/1393 985/1393 0 0 0 0 985/1393 -985/1393 0 0 -985/1393 985/1393 请注意：A 是对称矩阵请注意：因为A 的各阶顺序主子式 hl 不全大于零，所以A 不是正定的.A 的转置 矩阵zA和各阶顺序主子式值 hl 依次如下：zA = 985/1393 -985/1393 0 0 -985/1393 985/1393 0 0 0 0 985/1393 0 0 -985/1393hl = 985/1393 0 0可见，A 不是正定矩阵，是半正定矩阵;（4）在MATLA 工作窗口输入程序 >> A=[-2 1 1;1 -6 0;1 0 -4];hl=运行后输出结果A = -2110.1 2 3 4、Z1 -12 1、(1)二 2 d 4 ；(2) 4 3 0；；11 21 13 412 0 9 -6(5 78 9』J ； -6 19』1 一2 120 0<-2 1 「 1 -6 0 J 0 Tzddc (A)-985/1393985/1393 0 因为 A = A T因此, A 是对称矩阵zddc (A)例345 判断下列矩阵是否是正定对称矩阵:解 （1）在MATLA 工作窗口输入程序<1 1 ■2 0 01 -6 0 1 0 -4请注意：A 是对称矩阵请注意：因为A 的各阶顺序主子式 hl 不全大于零，所以A 不是正定的.A 的转置 矩阵zA和各阶顺序主子式值 hl 依次如下：zA = -2 1 1 hl = -2 11 -381 -6 0 1 0 -4可见A 不是正定矩阵，是负定矩阵；因为 A = A T 因此，A 是对称矩阵.3.5 求解线性方程组的LU 方法及其MATLA 程序3.5.1 解线性方程组的楚列斯基(Cholesky )分解法及其MATLAB?序例3.5.1先将矩阵A 进行楚列斯基分解,然后解矩阵方程 AX =b ，并用其他方法验证.-1 2 1、-1 3 0 -33 A =2 0 9-6,b = 5J-3—6 19」<7解在工作窗口输入>>A=[1 -1 2 1;-1 3 0 -3; 2 0 9 -6;1 -3 -6 19];b1=1:2:7; b=b1： R=chol(A);C=A-R'*R,R 仁in v(R);R2=R1：x=R1*R2*b,Rx=A\b-x运行后输出方程组的解和验证结果x = Rx =1.0e-014 * C =1.0e-015 *-8.0000 -0.7105 0 0 0 0 0.3333 -0.0833 0 -0.4441 0 0 3.6667 0.2220 0 0 00 2.0000 0.1332 0 0 0 0A 分解为一个单位下三角形矩阵 L 和一个上三角形矩阵 U 的积 A = LU . (2)在工作窗口输入>> U=[1 0 2 0;0 1 0 1;0 0 2 1;0 0 0 2]; L=[1 0 0 0;0 1 0 0;13.5.2 解线性方程组的直接 LU 分解法及其MATLAB?序例3.5.2 首先将矩阵A 直接进行分解，然后解矩阵方程 AX =b解(1)首先将矩阵 >> A=[1 0 A 直接进行2 0;0 10 1;1 2 4 3;0 1 0 3];b=[1;2;-1;5]; hl=zhjLU(A),A-L*U 运行后输出LU 分解 请注意：因为A 的各阶主子式都不等于零, 所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值 hl 依次如下:RA = 4U = 1 0 2 00 1 0 1 0 0 2 1 0 0 0 2 L = 1 0 0 0 0 1 0 0 1 2 1 0 0 1 0 1 hl = 1 1 2 4LU A =2 1 0;0 1 0 1];b=[1;2;-1;5];U 仁i nv(U); L 仁 in v(L); X=U1*L1*b,x=A\b 运行后输出方程组的解 X = 8.50000000000000 0.50000000000000 -3.75000000000000 1.500000000000003.5.3 解线性方程组的选主元的 LU 方法及其MATLAB?序先将矩阵A 进行LU 分解，然后解矩阵方程 AX =b 其中根据（3.28 ）式编写MATLA 程序，然后在工作窗口输入>> A=[0.1 2 3 4;-1 2 -3 4;11 21 13 41;5 7 8 9];b=[1;2;-1;5]; [L U P]=LU(A), U1=inv(U);L1=inv(L);X=U1*L1*P*b运行后输出结果 L = 1.0000 0 0 0 -0.0909 1.0000 0 0 0.0091 0.4628 1.0000 0 0.4545 -0.6512 0.2436 1.0000 U =11.0000 21.0000 13.0000 41.00000 3.9091 -1.8182 7.7273 0 0 3.7233 0.05120 0 0 -4.6171方法2根据(3.29 )式编写MATLA 程序，然后在工作窗口输入>> A=[0.1 2 3 4;-1 2 -3 4;11 21 13 41;5 7 8 9];b=[1;2;-1;5]; [F U]=LU(A), U1=inv(U); F1=inv(F); X=U1*F1*b 运行后输出结果用LU 分解法解线性方程组 A m n X = b 的MATLAB^序fun ction [RA,RB, n,X,Y]=LUjfcz(A,b)[n n] =size(A);B=[A b]; RA=ra nk(A); RB=ra nk(B); for 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;RA return end end if h （1,i ）~=0disp （'请注意：因为A 的各阶主子式都不等于零，所以A 能进行LU 分 解.A 的秩RA 和各阶顺序主子式值hl 依次如下：'）3 4、 「1、 -3 4 ，b = 213 41-189<5丿巾.1 2 -1 2 A =11 21 <57P =0 0 1 0 0 1 0 0 1 0 0 00 0 0 1 X =[-1 .2013 3.36770.0536例 3.5.3解方法1F=0.0091 0.4628 1.0000 0 -0.0909 1.0000 0 01.0000 0 0 00.4545 -0.6512 0.2436 1.0000X =[-1.2013 3.3677 0.0536 -U=11.0000 21.0000 13.0000 41.00000 3.9091 -1.8182 7.7273 0 0 3.7233 0.0512 0 0 0 -4.6171 1.4440]' -1.4440]X=zeros (n ,1); Y=zeros( n,1); C=zeros(1, n);r=1:1; for p=1: n-1 [max1,j]=max(abs(A(p: n,p))); C=A(p,:);A(p,：)= A(j+P1,:); C= A(j+P1,:); g=r(p); r(p)= r(j+P1); r(j+P1)=g; for k=p+1: nH= A(k,p)/A(p,p); A(k,p) = H; A(k,p+1: n)=A(k,p+1: n)- H* A(p,p+1: n); end endY(1)=B(r(1)); for k=2: nY(k)= B(r(k))- A(k,1:k-1)* Y(1:k-1); endX( n)= Y( n)/ A( n,n); for i=n-1:-1:1X(i)= (Y(i)- A(i, i+1:n) * X (i+1:n))/ A(i,i);end end[RA,RB,n,X,Y]'；解有何变化•其中以用下列命令计算n 阶希尔伯特矩阵>>h=hilb(n) % 输出 h 为n 阶 Hilbert 矩阵 在MATLA 工作窗口输入程序>> A=hilb(7);b=[1;2;2;2;2;2;2];X=A\b运行后输出 AX=b 的解为 X =( -35 504 -1260 -4200 20790 -27720 12012 )T .(2 )在MATLAB ：作窗口输入程序>> B =[0.001,zeros(1,6);zeros(6,1),zeros(6,6)]; A=(B+hilb(7)); b=[1;2;2;2;2;2;2];X=A\b运行后输出方程的解为 X= (-33 465 -966 -5181 22409 -29015 12413 )T .在MATLA 工作窗口输入程序>> X =[-33, 465,-966,-5181,22409,-29015,12413]';3.6 误差分析及其两种MATLA 程序3.6.1用MATLA 软件作误差分析例362解下列矩阵方程 AX 二b ，并比较方程（1）和（2）有何区别，它们的(1) (2)『1 /2 1 / 31 /4 1 / 21/ 31 / 4 1 /5 1 / 3 1 / 41 / 5 1 /6 1 / 41 / 51 / 6 1 /7 1 / 51 / 61 / 7 1 /8 1 / 61 / 71 / 81 /9/ 71 /81 / 9 1/10 ■1. 001 1 /2 1 /3 11 / 21 / 3 1 / 41 / 31 / 4 1/5 1 = 1 / 41 / 5 1/ 6 11 / 51 / 6 1/ 71 / 5 1 / 6 1 /7 1 / 61 /71 / 81 / 7 1 /8 1 /9 1 /8 1 /9 1 /10 1 / 9 1/101/11 1 /101 / 11 1 /12 1/111 /12 1 / 131/91 /101/11 1 / 121 / 8 1 / 9 1 / 10 1/11 1 /12 1 /13解（1）AX= b 的系数矩阵A 为7阶希尔伯特（Hilbert）矩阵, 广1 \2 2 2 . 2 2我们可A A 1 / 61 /71 /8 W 7 矩阵方程 /4 1 /5 1 / 6 1/7 / 5 1 /6 1 / 7 1 /8 / 6 1 / 71 / 8 1 / 9 / 7 1 / 8 1 / 91 / 10/ 8 1 /9 1/10 1 / 11 b =2 2 2 2X1 =[-35,504,-1260,-4200,20790,-27720,12012]'; wu=X1'- X'运行后输出方程(1 )和(2)的解的误差为X =(-2 39 -294 981 -1619 1295 -401 )1©001 0仅、方程(1)和(2)的系数矩阵的差为5A= 減,常数向量相同，则Ax = b1。

数值计算实验——解线性方程组西南交通大学2012级茅7班20123257 陈鼎摘要本报告主要介绍了基于求解线性方程组的高斯消元法和列主消元法两种数值分析方法的算法原理及实现方法。

运用matlab数学软件辅助求解。

实验内容1．编写用高斯消元法解线性方程组的MATLAB程序，并求解下面的线性方程组，然后用逆矩阵解方程组的方法验证。

2．编写用列主消元法解线性方程组的MATLAB程序，并求解下面的线性方程组，然后用逆矩阵解方程组的方法验证。

给定方程组如下：①0.325x1+2.564x2+3.888x3+5x4=1.521②-1.548x1+3.648x2+4.214x3-4.214x4=2.614③-2.154x1+1.647x2+5.364x3+x4=3.978④0x1+2.141x2-2.354x3-2x4=4.214A.高斯消元法一、算法介绍高斯消元法是一种规则化的加减消元法。

基本思想是通过逐次消元计算把需要求解的线性方程组转化成为上三角方程组，即把现形方程组的系数矩阵转化为上三角矩阵，从而使一般线性方程组的求解转化为等价的上三角方程组的求解。

二、matlab程序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，所以此方程组无解.')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和常系数向量b；输出的量：系数矩阵A和增广矩阵B的秩RA、RB，方程中未知量的个数n和有关方程组解X及其解的信息。

高斯消去法是一种用于求解线性方程组的经典方法，它可以通过矩阵的初等变换将方程组化为上三角形式，然后通过回代的方式求解方程组。

在Matlab中，我们可以利用高斯消去法求解方程组，这样可以更加高效地进行数值计算。

下面我们将简要介绍高斯消去法的原理，并通过Matlab代码演示如何使用高斯消去法求解方程组。

一、高斯消去法原理及步骤高斯消去法是一种通过矩阵的初等变换将线性方程组化为上三角形式的方法，其求解过程主要包括以下几个步骤：1. 将系数矩阵增广为增广矩阵；2. 首先通过初等行变换将增广矩阵化为上三角矩阵；3. 然后通过回代的方式求解方程组。

通过这样的步骤，我们可以将原始的线性方程组化简为上三角形式，从而更容易求解方程组。

二、Matlab代码演示在Matlab中，我们可以通过编写代码实现高斯消去法来求解线性方程组。

下面是一个简单的例子代码，用来演示如何在Matlab中使用高斯消去法求解方程组：```matlabfunction x = gauss_elimination(A, b)[n, m] = size(A);if n ~= merror('A must be a square matrix');endAb = [A, b];for k = 1 : n - 1for i = k + 1 : nfactor = Ab(i, k) / Ab(k, k);Ab(i, k : n + 1) = Ab(i, k : n + 1) - factor * Ab(k, k : n + 1); endendx = zeros(n, 1);x(n) = Ab(n, n + 1) / Ab(n, n);for i = n - 1 : -1 : 1x(i) = (Ab(i, n + 1) - Ab(i, i + 1 : n) * x(i + 1 : n)) / Ab(i, i); endend```通过以上的Matlab代码，我们可以实现高斯消去法的求解过程，并得到方程组的解。

Matlab方程组解1. 引言方程组是数学中一个重要的概念，它描述了多个未知数之间的关系。

解方程组的过程在科学、工程和计算机科学等领域中有着广泛的应用。

Matlab作为一种高级数值计算环境，提供了丰富的工具和函数来解决方程组的求解问题。

本文将介绍如何使用Matlab解方程组，包括线性方程组和非线性方程组的求解方法。

2. 线性方程组的求解2.1 利用矩阵求解线性方程组可以表示为矩阵形式，例如：Ax = b，其中A是系数矩阵，x是未知数向量，b是常数向量。

在Matlab中，可以使用线性代数工具箱中的函数来求解线性方程组。

2.1.1 使用inv函数求解如果系数矩阵A是可逆的，可以使用inv函数求解线性方程组。

具体步骤如下： 1. 计算A的逆矩阵：A_inv = inv(A) 2. 计算解向量：x = A_inv * b2.1.2 使用linsolve函数求解linsolve函数可以直接求解线性方程组，无需计算逆矩阵。

具体步骤如下： 1. 调用linsolve函数：x = linsolve(A, b)2.2 利用高斯消元法求解高斯消元法是一种常用的求解线性方程组的方法，它通过矩阵的行变换将方程组转化为上三角矩阵，然后通过回代得到解。

在Matlab中，可以使用lu函数来进行高斯消元法求解。

2.2.1 使用lu函数求解lu函数可以将方程组的系数矩阵A分解为一个下三角矩阵L和一个上三角矩阵U，即A = LU。

具体步骤如下： 1. 调用lu函数：[L, U] = lu(A) 2. 解得方程组：x = U \ (L \ b)3. 非线性方程组的求解非线性方程组是指未知数与其函数之间存在非线性关系的方程组。

与线性方程组不同，非线性方程组的求解通常需要借助数值方法。

Matlab提供了多种函数和工具箱来解决非线性方程组的求解问题。

3.1 利用fsolve函数求解fsolve函数是Matlab中用于求解非线性方程组的函数，它通过迭代的方式逼近方程组的解。

matlab求线性方程组的解求解线性方程分为两种方法–直接法和迭代法常见的方法一共有8种直接法Gauss消去法Cholesky分解法迭代法Jacobi迭代法Gauss-Seidel迭代法超松弛迭代法共轭梯度法Bicg迭代法Bicgstab迭代法这里我从计算代码的角度来解释一下，代码按以下顺序给出。

把方程组直接带入已知条件，就可以得到答案。

适用条件Gauss消去法：求解中小规模线性方程（阶数不过1000），一般用于求系数矩阵稠密而且没有任何特殊结构的线性方程组Cholesky分解法：对称正定方程优先使用，系数矩阵A是n 阶对称正定矩阵Jacobi迭代法非奇异线性方程组，分量的计算顺序没有关系Gauss-Seidel迭代法与Jacobi迭代法相似，但计算的分量不能改变超松弛迭代法Jacobi迭代法和Gauss-Seidel迭代法的加速版，由Gauss-Seidel迭代法改进而来，速度较快共轭梯度法需要确定松弛参数w,只有系数矩阵具有较好的性质时才可以找到最佳松弛因子。

但好处是不用确定任何参数，他是对称正定线性方程组的方法也是求解大型稀疏线性方程组最热门的方法Bicg迭代法本质是用双共轭梯度求解线性方程组的方法，对求解的方程没有正定性要求Bicgstab迭代法本质是用稳定双共轭梯度求解线性方程组的方法，对求解的方程没有正定性要求Gauss消去法第一、二个函数ltri、utri是一定要掌握的，后面的几乎每个函数都要用到ltri简单来说，当Ly=bb,L(非奇异下三角矩阵）已知求yfunction y =ltri(L,b)n=size(b,1);y=zeros(n,1);for j =1:n-1y(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-y(j)*L(j+1:n,j); endy(n)=b(n)/L(n,n);utri简单来说，当Ux=yy,U(非奇异上三角矩阵）已知求xfunction x =utri(U,y)n=size(y,1);x=zeros(n,1);for j = n:-1:2x(j)=y(j)/U(j,j);y(1:j-1)=y(1:j-1)-x(j)*U(1:j-1,j);endx(1)=y(1)/U(1,1);gauss算法，计算时粘贴过去就好function[L,U]=gauss(A)n=size(A,1);for k =1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k +1:n)-A(k+1:n,k)*A(k,k+1:n);endL=tril(A,-1)+eye(n);U=triu(A);使用例子已经知道一个线性方程组，这里我就不写出数学形式了，A是系数矩阵，直接把上面写好的函数复制过来在运算就可以。

§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： 用高斯消元法求解下面的非齐次线性方程组。