MATLAB解线性方程组的直接方法
MATLAB计算方法3解线性方程组计算解法
直到(n-1) 原方程组化为
a11 x1 a12 x2 a1n xn a1,n1 a22 x2 a2 n xn a2 ,n1
ann xn an ,n1
(上三角方程组) （3.2) 以上为消元过程。
(n) 回代求解公式
a n ,n1 xn a nn n x k 1 [a k ,n1 a kj x j ] a kk j k 1 ( k n 1, n 2,...,1)
由矩阵乘法 (1) 1) l11 a11 l11
umj 1 ukj a kj ukj a kj l km umj
m 1
k 1
2 求L的第k列：用L的第i行 u的第k列
(i k 1, , n)，即 ( l i 1 , , l ik , l kk , 0 0) ( u1k , u2 k , , ukk , 0 0)' a ik
( 2) 1)求u的第2行：用L的第2行 u的第j列 （j 2, , n) l 21 u1 j 1 u2 j a 2 j u2 j a 2 j l 21u1 j 2)求L的第2列：用L的第i行 u的第2列 （i 3,4, , n) l i 1 u12 l i 2 u22 a i 2 l i 2 (a i 2 l i 1 u12 ) / u22
m 1
l
k 1
im
umk l ik ukk a ik
k 1
l ik a ik l im umk ukk m 1
LU分解式： u1 j a1 j ( j 1,2, n) l i 1 a i 1 u11 ( i 2,3, , n) k 1 ukj a kj l km umj a kj m 1 ( j k , k 1, , n) k 1 l ik a ik l im umk ukk a ik m 1 ( i k 1, , n) ( k 2, 3, , n )
利用Matlab求线性方程组的通解
一、 引一 线性代数是数学中的一个重要分支．很多理论问题和实际问题都需要借 助于线性代数的理论工具来分析解决。学习线性代数有两大难点：一是概 念、理论抽象，二是计算量大。由于这两个难点，初学者往往很难掌握好线 性代 数的 知识理 论。若 能掌 握好用 于解 数学问 题的 l I at l ab软件 ，则能 轻松 快 捷的解决很多线 性代数问题。 ( 一 ) 线性方程组有关定理
X x=
0 ．7530
- 0 ．4167 _o ．3043 -t ) ．3 04 3
0 ．0176
- 0 ．7464
-t }．0 00 0
基础解系含有 刀一，个解向量，( 1) 若舌 ，参，…，靠，为导出组的 基础解系。
- 0 ．0000
- 0．7071 0 ．70 7l -- 0．0 00 0
则毛毒十如最+…+ko就是导出组朋=0的全部解，也称为通解：
ab方法
200 6．
参考文献： [ 1】胡良剑、孙晓君。 ( 1l at l ab数学实验， [ 1I 】，北京：高等教育出版社．
[ 2] 赵 静、但 琦． 出版社， 2003 ． [ 3] 同济大学应用数学系， 出版社． 2004 ．
‘教学 建模与 数学实验 ’( 第2版 ) [ 岫。北 京：高等 教育 ‘高等教学' ( 第5版) [ M ] ．北京t 高等教育
1
1 毛+ 毛+而 +毛+ 毛 =7
例 3求线性方程组
3而+2 毛+气+‘一 3而=—2的 通解
毛+勰+2 ‘+6毛=23
5而 +4b +3而 +3丘 一 毛=1 2
MATLAB程序设计第六讲
MATLAB程序设计杨凯2010 . 11主要内容自学））*MATLAB解方程与函数极值解方程与函数极值（（自学（自学）自学）线性方程组求解（一、线性方程组求解二、非线性方程组求解三、函数极值四、常微分方程初值问题的数值解法*MATLAB符号计算一、符号计算基础二、微积分三、简化方程表达式四、解方程一、线性方程组求解（自学）1.1 直接解法1．利用左除运算符的直接解法对于线性方程组Ax=b，可以利用左除运算符“\”求可以利用左除运算符“解：x=A\b例*：用直接解法求解下列线性方程组用直接解法求解下列线性方程组。
命令如下命令如下：：A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];b=[13,-9,6,0]';x=A\b2．利用矩阵的分解求解线性方程组矩阵分解是指根据一定的原理用某种算法将一个矩阵分解成若干个矩阵的乘积算法将一个矩阵分解成若干个矩阵的乘积。
常见的矩阵分解有LU 分解分解、、QR 分解分解、、Cholesky 分解分解，，以及Schur 分解分解、、Hessenberg 分解分解、、奇异分解等奇异分解等。
(1) LU 分解矩阵的LU 分解就是将一个矩阵表示为一个交换下三角矩阵和一个上三角矩阵的乘积形式角矩阵和一个上三角矩阵的乘积形式。
线性代数中已经证明，只要方阵A 是非奇异的是非奇异的，，LU 分解总是可以进行的分解总是可以进行的。
MATLAB 提供的lu 函数用于对矩阵进行LU 分解分解，，其调用格式为式为：：[L,U]=lu(X)：产生一个上三角阵U 和一个变换形式的下三角阵L(行交换)，使之满足X=LU 。
注意注意，，这里的矩阵X 必须是方阵是方阵。
[L,U,P]=lu(X)：产生一个上三角阵U 和一个下三角阵L 以及一个置换矩阵P ，使之满足PX=LU 。
当然矩阵X 同样必须是方阵方阵。
实现LU 分解后分解后，，线性方程组Ax=b 的解x=U\(L\b)或x=U\(L\Pb)，这样可以大大提高运算速度这样可以大大提高运算速度。
matlab最小二乘解方程
matlab最小二乘解方程最小二乘法是求解线性方程组的一种有效方法，可以通过最小化误差平方和来得到最优解。
在MATLAB中，我们可以使用“\”操作符或者使用“pinv”函数来求解一个线性方程组的最小二乘解。
以下是关于如何在MATLAB中使用最小二乘法来求解线性方程组的详细内容:1. 使用“\”操作符使用“\”操作符可以很方便地求解一个线性方程组的最小二乘解。
例如，假设我们有一个由n个方程组成的线性方程组：Ax = b其中，A是一个m ×n的矩阵，x是一个n维向量，b是一个m维向量。
则它的最小二乘解为：x = (A' A)^(-1) A' b在MATLAB中，我们可以通过以下代码实现最小二乘解：A = [1 1 1; 2 3 4; 4 5 7; 5 6 8];b = [1; 2; 3; 4];x = A \ b;其中，反斜杠符号“\”表示求解线性方程组的最小二乘解。
2. 使用“pinv”函数除了使用“\”操作符，我们也可以使用MATLAB中的“pinv”函数来求解一个线性方程组的最小二乘解。
例如，我们可以通过以下代码实现最小二乘解：A = [1 1 1; 2 3 4; 4 5 7; 5 6 8];b = [1; 2; 3; 4];x = pinv(A) * b;其中，pinv函数表示求矩阵A的伪逆矩阵。
使用“pinv”函数来求解线性方程组的最小二乘解与使用“\”操作符的结果是等价的。
需要注意的是，在使用最小二乘法来求解线性方程组时，矩阵A的列应该是线性无关的，否则可能会出现唯一最小二乘解不存在的情况。
综上所述，MATLAB中使用最小二乘法来求解线性方程组非常简单。
我们可以通过“\”操作符或者“pinv”函数来求解一个线性方程组的最小二乘解。
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 行交换再按高斯消元法进行下去称为列主元素消元法。
matlab中快速求解xa=b的方法
matlab中快速求解xa=b的方法在Matlab中，要快速求解线性方程组xa=b，可以使用以下几种方法：1. 直接求解法（\）：直接使用斜杠操作符（\）可以求解线性方程组。
例如，对于方程组xa=b，可以直接使用x = A\b来解决，其中A是系数矩阵，b是常数向量。
这种方法使用了高效的LU分解算法，并且能够自动适应方程组的类型（如稀疏矩阵或密集矩阵），因此是一种快速求解线性方程组的常用方法。
2. QR分解法：QR分解是一种将矩阵分解为正交矩阵和上三角矩阵的方法。
在Matlab中，可以使用qr函数对系数矩阵进行QR分解，然后使用这个分解求解线性方程组。
具体而言，可以使用[q,r] = qr(A)将系数矩阵A分解为正交矩阵q和上三角矩阵r，然后使用x = r\(q'*b)求解方程组。
这种方法通常适用于方程组的系数矩阵具有较大的条件数或者方程组数目较多的情况。
3. Cholesky分解法：如果线性方程组的系数矩阵是对称正定的，那么可以使用Cholesky分解来求解方程组。
在Matlab中，可以使用chol函数对系数矩阵进行Cholesky分解，然后使用这个分解求解线性方程组。
具体而言，可以使用R = chol(A)将系数矩阵A分解为上三角矩阵R，然后使用x = R'\(R\b)求解方程组。
Cholesky分解法通常适用于系数矩阵具有良好的性质（如对称正定）的情况。
4. 迭代法：如果线性方程组的系数矩阵是稀疏的，那么可以使用迭代法来求解方程组。
迭代法的基本思想是通过迭代改进解的逼近值。
在Matlab中，可以使用pcg函数（预处理共轭梯度法）或者bicg函数（双共轭梯度法）来求解稀疏线性方程组。
这些函数需要提供一个预处理矩阵，用于加速迭代过程。
预处理矩阵可以根据具体问题进行选择，常见的预处理方法包括不完全LU分解（ilu）和代数多重网格（amg）等。
通过使用上述方法，可以在Matlab中快速求解线性方程组xa=b。
MATLAB-平方根法和改进平方根法求解线性方程组例题与程序演示教学
M A T L A B-平方根法和改进平方根法求解线性方程组例题与程序（2）设对称正定阵系数阵线方程组12345678424024000221213206411418356200216143323218122410394334411142202531011421500633421945x x x x x x x x -⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥---⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥----⎢⎥⎢⎥⎢⎥----⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥----⎢⎥⎢⎥⎢⎥----⎢⎥⎢⎥⎢⎢⎥⎢⎥⎢---⎢⎥⎢⎥⎢--⎢⎥⎢⎢⎥⎣⎦⎣⎦⎣⎦⎥⎥⎥⎥ (1,1,0,2,1,1,0,2)T x *=--二、数学原理 1、平方根法解n 阶线性方程组Ax=b 的choleskly 方法也叫做平方根法，这里对系数矩阵A 是有要求的，需要A 是对称正定矩阵，根据数值分析的相关理论，如果A 对称正定，那么系数矩阵就可以被分解为的T A=L L •形式，其中L 是下三角矩阵，将其代入Ax=b 中，可得：T LL x=b 进行如下分解：T L xL by y ⎧=⎨=⎩ 那么就可先计算y,再计算x ，由于L 是下三角矩阵，是T L 上三角矩阵，这样的计算比直接使用A 计算简便，同时你应该也发现了工作量就转移到了矩阵的分解上面，那么对于对称正定矩阵A 进行Cholesky 分解，我再描述一下过程吧： 如果你对原理很清楚那么这一段可以直接跳过的。
设T A=L L •，即1112111112112122221222221212....................................n n n n n n nn n n nn nn a a a l l l l aa a l l l l a a a l l l l ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦其中,,1,2,...,ij ji a a i j n ==第1步，由矩阵乘法，211111111,i i a l a l l ==g 故求得111111,2,3,...i i a l l i n a === 一般的，设矩阵L 的前k-1列元素已经求出 第k 步，由矩阵乘法得112211k k kk kmkkik im km ik kkm m a l l a l l l l --===+=+∑∑, 于是11(2,3,...,n)1(),1,2,...kk k ik ik im km m kk l k l a l l i k k n l -=⎧=⎪⎪=⎨⎪=-=++⎪⎩∑ 2、改进平方根法在平方根的基础上，为了避免开方运算,所以用TLDL A =计算；其中，11122.........n d D D D d ⎤⎤⎡⎤⎥⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥===⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎢⎢⎥⎣⎦⎣⎣；得1121212212111111n n n n n d l l l d l A l l d ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦L L M MO O O M L按行计算的L 元素及D 对元素公式 对于n i ,,2,1Λ=11(1,21)j ij ij ik jk k t a t l j i -==-=-∑…，./(1,2,)ij ij j l t d j ==…，i-1.11i i ii ik ikk d a t l -==-∑计算出LD T =的第i 行元素(1,2,i-1)ij t j =…，后，存放在A 的第i 行相置，然后再计算L 的第i 行元素，存放在A 的第i 行.D 的对角元素存放在A 的相应位置.对称正定矩阵A 按T LDL 分解和按T LL 分解计算量差不多，但T LDL 分解不需要开放计算。
利用matlab解线性方程组
数值计算实验——解线性方程组西南交通大学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及其解的信息。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
在这章中我们要学习线性方程组的直接法，特别是适合用数学软件在计算机上求解的方法.3.1 方程组的逆矩阵解法及其MATLAB 程序3.1.3 线性方程组有解的判定条件及其MATLAB 程序 判定线性方程组A n m ⨯b X =是否有解的MATLAB 程序function [RA,RB,n]=jiepb(A,b)B=[A b];n=length(b); RA=rank(A); RB=rank(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 end例3.1.4 判断下列线性方程组解的情况.如果有唯一解，则用表 3-2方法求解.(1) ⎪⎪⎩⎪⎪⎨⎧=-+-=+-+=-++=+-+;0742,0634,0723,05324321432143214321x x x x x x x x x x x x x x x x (2) ⎪⎪⎩⎪⎪⎨⎧=++-=+-+=-+-=+-+;0327,01613114,02332,075434321432143214321x x x x x x x x x x x x x x x x (3) ⎪⎩⎪⎨⎧=+=+-=-+;8311,1023,22421321321x x x x x x x x (4) ⎪⎩⎪⎨⎧=--+=+-+=+-+.12,2224,12w z y x w z y x w z y x解 在MATLAB 工作窗口输入程序>> 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 在MATLAB 工作窗口输入>>X=A\b,运行后输出结果为 X =(0 0 0 0)’.(2) 在MATLAB 工作窗口输入程序>> 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) 在MATLAB 工作窗口输入程序>> 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)在MATLAB 工作窗口输入程序>> 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 三角形方程组的解法及其MATLAB 程序3.2.2 解三角形方程组的MATLAB 程序 解上三角形线性方程组b AX =的MATLAB 程序function [RA,RB,n,X]=shangsan(A,b)B=[A b]; n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0,disp('请注意：因为RA~=RB ，所以此方程组无解.') return endif RA==RB if 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);end elsedisp('请注意：因为RA=RB<n ，所以此方程组有无穷多解.')end end例3.2.2 用解上三角形线性方程组的MATLAB 程序解方程组⎪⎪⎩⎪⎪⎨⎧==+-=-+-=++-.63,456,7472,203254434324321x x x x x x x x x x . 解 在MATLAB 工作窗口输入程序>>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]=shangsan(A,b)运行后输出结果请注意：因为RA=RB=n ，所以此方程组有唯一解. RA = RB =4, 4, n =4,X =[2.4 -4.0 -1.0 2.0]’3.3 高斯（Gauss ）消元法和列主元消元法及其MATLAB 程序3.3.1 高斯消元法及其MATLAB 程序用高斯消元法解线性方程组b AX =的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 ，所以此方程组无解.') 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:nm= 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); endelsedisp('请注意：因为RA=RB<n ，所以此方程组有无穷多解.') end end例3.3.2 用高斯消元法和MATLAB 程序求解下面的非齐次线性方程组，并且用逆矩阵解方程组的方法验证.⎪⎪⎩⎪⎪⎨⎧-=+---=+--=+--=-+-.142,16422,0,13432143214324321x x x x x x x x x x x x x x x 解 在MATLAB 工作窗口输入程序>> 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 ，所以此方程组有唯一解. RA =4RB =4n =43.3.2 列主元消元法及其MATLAB 程序用列主元消元法解线性方程组b AX =的MATLAB 程序function [RA,RB,n,X]=liezhu(A,b)B=[A b]; n=length(b); RA=rank(A);X = 0 -0.5000 0.5000 0RB=rank(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-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); 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); endelsedisp('请注意：因为RA=RB<n ，所以此方程组有无穷多解.') end end例3.3.3 用列主元消元法解线性方程组的MATLAB 程序解方程组⎪⎪⎩⎪⎪⎨⎧-=+---=+--=-+-=+--.142,16422,13,0432143214321432x x x x x x x x x x x x x x x . 解 在MATLAB 工作窗口输入程序>> 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 分解法及其MATLAB 程序3.4.1判断矩阵LU 分解的充要条件及其MATLAB 程序 判断矩阵A 能否进行LU 分解的MATLAB 程序function hl=pdLUfj(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: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 endif h(1,i)~=0disp('请注意：因为A 的各阶主子式都不等于零，所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值hl 依次如下：')hl;RA end end例3.4.1 判断下列矩阵能否进行LU 分解，并求矩阵的秩.（1）⎪⎪⎪⎭⎫ ⎝⎛6547121321；（2）⎪⎪⎪⎭⎫ ⎝⎛654721321；（3）⎪⎪⎪⎭⎫ ⎝⎛654321321.解 （1）在MATLAB 工作窗口输入程序>> A=[1 2 3;1 12 7;4 5 6];hl=pdLUfj(A)运行后输出结果为请注意：因为A 的各阶主子式都不等于零，所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值hl 依次如下：RA = 3， hl = 1 10 -48（2）在MATLAB 工作窗口输入程序>> 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）在MATLAB 工作窗口输入程序>> 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 分解法及其MATLAB 程序 将矩阵A 进行直接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);return endif 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 分解.A的秩RA 和各阶顺序主子式值hl 依次如下：')for j=1:nU(1,j)=A(1,j); endfor k=2:n for i=2:n for 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); end end end endhl;RA,U,L end end例3.4.3 用矩阵进行直接LU 分解的MA TLAB 程序分解矩阵⎪⎪⎪⎪⎪⎭⎫⎝⎛=3010342110100201A . 解 在MATLAB 工作窗口输入程序>> 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 = 4U = 1 0 2 00 1 0 10 0 2 10 0 0 2 3.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)>0disp('请注意：因为A 的各阶顺序主子式hl 都大于零，所以A 是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下：')hl;zA endL = 1 0 0 0 0 1 0 0 1 2 1 0 0 1 0 1 hl = 1 1 2 4例3.4.5 判断下列矩阵是否是正定对称矩阵：（1）⎪⎪⎪⎪⎪⎭⎫ ⎝⎛--98754113211143214321.0；(2) ⎪⎪⎪⎪⎪⎭⎫ ⎝⎛------19631690230311211； (3) ⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛----212100212100002121002121；（4）⎪⎪⎪⎭⎫⎝⎛---401061112. 解 （1）在MATLAB 工作窗口输入程序>> A=[0.1 2 3 4;-1 2 -3 4;11 21 13 41;5 7 8 9];hl=zddc (A)运行后输出结果请注意： A 不是对称矩阵请注意：因为A 的各阶顺序主子式hl 不全大于零，所以A 不是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下：zA = 1/10 -1 11 5 2 2 21 7 3 -3 13 8 4 4 41 9 hl = 1/10 11/5 -1601/10 3696/5因此,A 即不是正定矩阵，也不是对称矩阵.（2）在MATLAB 工作窗口输入程序>> 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 0 9 -6 1 -3 -6 19 请注意： A 是对称矩阵请注意：因为A 的各阶顺序主子式hl 都大于零，所以A 是正定的.A 的转置矩阵zA和各阶顺序主子式值hl 依次如下：zA = 1 -1 2 1 -1 3 0 -3 2 0 9 -6 1 -3 -6 19 hl = 1 2 6 24 （3）在MATLAB 工作窗口输入程序>> 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 -985/1393 0 0 -985/1393 985/1393 hl = 985/1393 0 0 0可见，A 不是正定矩阵，是半正定矩阵；因为A = A T因此，A 是对称矩阵.（4）在MATLAB 工作窗口输入程序>> A=[-2 1 1;1 -6 0;1 0 -4];hl=zddc (A)运行后输出结果A = -2 1 11 -6 0 1 0 -4 请注意： A 是对称矩阵请注意：因为A 的各阶顺序主子式hl 不全大于零，所以A 不是正定的.A 的转置矩阵zA 和各阶顺序主子式值hl 依次如下：zA = -2 1 1 hl = -2 11 -38 1 -6 0 1 0 -4可见A 不是正定矩阵，是负定矩阵；因为A = A T因此，A 是对称矩阵.3.5 求解线性方程组的LU 方法及其MATLAB 程序3.5.1 解线性方程组的楚列斯基（Cholesky ）分解法及其MATLAB 程序例3.5.1 先将矩阵A 进行楚列斯基分解,然后解矩阵方程b AX =，并用其他方法验证.⎪⎪⎪⎪⎪⎭⎫⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛------=7531,19631690230311211b A . 解 在工作窗口输入>>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,R1=inv(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 0 0 2.0000 0.1332 0 0 0 03.5.2 解线性方程组的直接LU 分解法及其MATLAB 程序例3.5.2 首先将矩阵A 直接进行LU 分解，然后解矩阵方程b AX =⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=3010342110100201A ，⎪⎪⎪⎪⎪⎭⎫⎝⎛-=5121b . 解 (1) 首先将矩阵A 直接进行LU 分解.在MATLAB 工作窗口输入程序>> A=[1 0 2 0;0 1 0 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 10 0 2 10 0 0 2A 分解为一个单位下三角形矩阵L 和一个上三角形矩阵U 的积 LU A =.(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;1L = 1 0 0 0 0 1 0 0 1 2 1 0 0 1 0 1 hl = 1 1 2 42 1 0;0 1 0 1];b=[1;2;-1;5];U1=inv(U); L1=inv(L); X=U1*L1*b,x=A\b运行后输出方程组的解X = 8.50000000000000 0.50000000000000 -3.75000000000000 1.500000000000003.5.3 解线性方程组的选主元的LU 方法及其MATLAB 程序例3.5.3 先将矩阵A 进行LU 分解，然后解矩阵方程b AX = 其中⎪⎪⎪⎪⎪⎭⎫ ⎝⎛--=98754113211143214321.0A ，⎪⎪⎪⎪⎪⎭⎫⎝⎛-=5121b . 解 方法1 根据（3.28）式编写MATLAB 程序，然后在工作窗口输入>> 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.0000 0 3.9091 -1.8182 7.72730 0 3.7233 0.05120 0 0 -4.6171方法2 根据（3.29）式编写MATLAB 程序，然后在工作窗口输入>> 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运行后输出结果F=0.0091 0.4628 1.0000 0 -0.0909 1.0000 0 0 1.0000 0 0 0 0.4545 -0.6512 0.2436 1.0000 X =[-1.2013 3.3677 0.0536 -1.4440]’ 用LU 分解法解线性方程组A n m ⨯b X =的MATLAB 程序function [RA,RB,n,X,Y]=LUjfcz(A,b)[n n] =size(A);B=[A b]; RA=rank(A); RB=rank(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 endif h(1,i)~=0disp('请注意：因为A 的各阶主子式都不等于零，所以A 能进行LU 分解.A 的秩RA 和各阶顺序主子式值hl 依次如下：')P = 0 0 1 0 0 1 0 0 1 0 0 0 0 0 0 1 X =[-1.2013 3.3677 0.0536 -1.4440]’ U=11.0000 21.0000 13.0000 41.0000 0 3.9091 -1.8182 7.7273 0 0 3.7233 0.0512 0 0 0 -4.6171X=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]’;3.6 误差分析及其两种MATLAB 程序3.6.1 用MATLAB 软件作误差分析例3.6.2 解下列矩阵方程b AX =，并比较方程（1）和（2）有何区别，它们的解有何变化.其中,13/112/111/110/19/18/17/112/111/110/19/18/17/16/111/110/19/18/17/16/15/110/19/18/17/16/15/14/19/18/17/16/15/14/13/18/17/16/15/14/13/12/17/16/15/14/13/12/11)1(⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=A ;2222221⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=b ,13/112/111/110/19/18/17/112/111/110/19/18/17/16/111/110/19/18/17/16/15/110/19/18/17/16/15/14/19/18/17/16/15/14/13/18/17/16/15/14/13/12/17/16/15/14/13/12/1001.1)2(⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=A .2222221⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=b 解 (1) 矩阵方程b AX =的系数矩阵A 为7阶希尔伯特（Hilbert ）矩阵，我们可以用下列命令计算n 阶希尔伯特矩阵>>h=hilb(n) % 输出h 为n 阶Hilbert 矩阵 在MATLAB 工作窗口输入程序>> A=hilb(7);b=[1;2;2;2;2;2;2];X=A\b运行后输出b AX =的解为 X =（-35 504 -1260 -4200 20790 -27720 12012T).（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 12413T).在MATLAB 工作窗口输入程序>> X =[-33, 465,-966,-5181,22409,-29015,12413]';X1 =[-35,504,-1260,-4200,20790,-27720,12012]'; wu=X1'- X' 运行后输出方程（1）和（2）的解的误差为=δX 401- 1295 1619- 981 294- 39 -2(T ).方程（1）和（2）的系数矩阵的差为⎪⎪⎭⎫⎝⎛=δ⨯⨯⨯661661001.0O O O A ，常数向量相同，则b Ax =的解的差为=δX 40112951619981294392(----T ).A 的微小变化，引起X 的很大变化，即X 对A 的扰动是敏感的.3.6.2 求P 条件数和讨论b AX =解的性态的MATLAB 程序求P 条件数和讨论b AX =解的性态的MATLAB 程序function Acp=zpjxpb(A)Acw = cond (A, inf);Ac1= cond (A,1);Ac2= cond (A,2);Acf = cond (A,'fro');dA=det(A);if (Acw>50)&(dA<0.1)disp('请注意：AX=b 是病态的，A 的∞条件数,1条件数,2条件数, 弗罗贝尼乌斯条件数和A 的行列式的值依次如下：')Acp=[Acw Ac1 Ac2 Acf dA]';elsedisp(' AX=b 是良态的，A 的∞条件数,1条件数,2条件数,弗罗贝尼乌斯条件数和A 的行列式的值依次如下：')Acp=[Acw Ac1 Ac2 Acf dA]';end例3.6.3 根据定理3.10，讨论线性方程组b AX =解的性态，并且求出A 的4种条件数.其中（1）A 为7阶希尔伯特矩阵；（2）⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-----=7421631472135132A . 解 （1）首先将求P 条件数和讨论b AX =解的性态的MATLAB 程序保存名为zpjxpb.m 的M 文件，然后在MATLAB 工作窗口输入程序>> Acp =zpjxpb(hilb(7)); Acp',det(hilb(7))运行后输出结果请注意：AX=b 是病态的，A 的∞条件数,1条件数,2条件数, 弗罗贝尼乌斯条件数和A 的行列式的值依次如下：ans = 1.0e+008 *9.8519 9.8519 4.7537 4.8175 0.0000ans = 4.8358e-025（2）在MATLAB 工作窗口输入程序>> A=[2 3 -1 5;3 1 2 -7;4 1 -3 6;1 -2 4 -7];Acp=zpjxpb(A); Acp' 运行后输出结果 AX=b 是良态的，A 的∞条件数,1条件数,2条件数, 弗罗贝尼乌斯条件数和A 的行列式的值依次如下：ans =14.1713 19.4954 8.2085 11.4203 327.00003.6.3 用P 范数讨论b AX =解和A 的性态的MATLAB 程序用P 范数讨论b AX =解和A 的性态的MATLAB 程序function Acp=zpjwc(A,jA,b,jb,P)Acp = cond (A,P);dA=det(A); X=A\b;dertaA=A-jA;PndA=norm(dertaA, P);dertab=b-jb;Pndb=norm(dertab, P);if Pndb>0jX=A\jb; Pnb= norm(b, P);PnjX = norm(jX,P); dertaX=X-jX;PnjdX= norm(dertaX, P);jxX= PnjdX/PnjX; PnjX =norm(jX,P);PnX = norm(X,P); jxX= PnjdX/PnjX; xX= PnjdX/PnX;Pndb=norm(dertab,P);xAb=Pndb/Pnb;Pnbj=norm(jb,P); xAbj=Pndb/Pnbj;Xgxx= Acp*xAb;endif PndA>0jX=jA\b; dertaX=X-jX;PnX = norm(X,P);PnjdX= norm(dertaX, P);PnjX = norm(jX,P); jxX= PnjdX/PnjX;xX= PnjdX/PnX;PnjA=norm(jA,P); PnA=norm(A,P);PndA=norm(dertaA,P);xAbj= PndA/PnjA;xAb= PndA/PnA;Xgxx= Acp*xAb;endif (Acp >50)&(dA<0.1)disp('请注意：AX=b 是病态的，A 的P 条件数Acp,A 的行列式值dA ，解X ，近似解jX ，解的相对误差jxX ，解的相对误差估计值Xgxx ，b 或A 的相对误差xAb 依次如下：')Acp,dA,X',jX',xX',jxX',Xgxx',xAb',xAbj'elsedisp('请注意： AX=b 是良态的，A 的P 条件数Acp,A 的行列式值dA ，解X ，近似解jX ，解的相对误差jxX ，解的相对误差估计值Xgxx ，b 或A 的相对误差xAb 依次如下：')Acp,dA,X',jX',xX',jxX',Xgxx',xAb',xAbj'end例3.6.4 根据定理3.10，讨论线性方程组b AX =解的性态，并利用（3.32）式讨论当A 的每个元都取4位有效数字时，其解的相对误差.其中A 为7阶希尔伯特矩阵，()22224311=b T .解 （1）取∞范数和∞条件数，线性方程组b AX =的b 不变时，取∞范数和∞条件数，系数矩阵A 为7阶希尔伯特矩阵，A 中的每个元素取4位有效数字.用P 范数讨论b AX =解和A 的性态的MATLAB 程序保存名为zpjwc.m 的文件,然后在工作窗口输入MATLAB 程序>> jA =[1.0000 0.5000 0.3333 0.2500 0.2000 0.1667 0.14290.5000 0.3333 0.2500 0.2000 0.1667 0.1429 0.12500.3333 0.2500 0.2000 0.1667 0.1429 0.1250 0.11110.2500 0.2000 0.1667 0.1429 0.1250 0.1111 0.10000.2000 0.1667 0.1429 0.1250 0.1111 0.1000 0.09090.1667 0.1429 0.1250 0.1111 0.1000 0.0909 0.08330.1429 0.1250 0.1111 0.1000 0.0909 0.0833 0.0769]; A=hilb(7);b=[1;1/3;4;2;2;2;2];jb=[1;0.3333;4;2;2;2;2]; Acp=zpjwc(A,jA,b,jb,inf)运行后输出结果请注意：AX=b 是病态的，A 的P 条件数Acp,A 的行列式值dA ，解X ，近似解jX ，解的相对误差jxX ，解的相对误差估计值Xgxx ，b 或A 的相对误差xAb 依次如下：Acp = dA =9.8519e+008 4.8358e-025ans =1.0e+007 *0.0020 -0.0697 0.6243 -2.3054 4.0677 -3.41231.0943ans =1.0e+004 *0.0349 -0.4807 2.1126 -5.1087 7.6557 -6.32392.1112xX = jxX = Xgxx =0.9981 530.3248 4.9291e+004 xAb = xAbj =5.0032e-005 5.0031e-005由此可见，因为∞条件数(Cond ≈∞)A 985 194 889.719 848 31>>，所以此方程组为病态的b AX =的解T 932)94210320,12334-790,676380,402.305-300,2436 256,697- 565,19( =X ， b X jA =)(的解为T )11221239,63-,55767087,51-126,21 079,48- ,493( =jX 解的相对误差120.998≈∞∞X Xδ,291.2749≤+∞∞X X X δδ , 530.32≈+∞∞X X X δδ,,1025.0035-∞∞⨯≈A Aδ即相对误差放大了约985 194 889.72倍.(2) 如果取2范数和2条件数计算，在MATLAB 工作窗口输入程序>> Acp =zpjwc(A,jA,b,jb,2)运行后输出结果请注意：AX=b 是病态的，A 的P 条件数Acp,A 的行列式值dA ，解X ，近似解jX ，解的相对误差jxX ，解的相对误差估计值Xgxx ，b 或A 的相对误差xAb 依次如下：Acp = dA =4.7537e+008 4.8358e-025ans =1.0e+007 *0.0020 -0.0697 0.6243 -2.3054 4.0677 -3.41231.0943ans =1.0e+004 *0.0349 -0.4807 2.1126 -5.1087 7.6557 -6.32392.1112xX = jxX = Xgxx=0.9981 511.0640 2.9951e+004xAb = xAbj =6.3006e-005 6.3005e-005因为2条件数Cond ≈2)(A 475 367 356.591>>，所以此方程组为病态的.解的相对误差10.998 22≈X X δ,,105300.6522-⨯≈A A δ,06.51122≈+X X X δδ.85.9502922≤+A A A δδ 即相对误差放大了约475 367 356.59倍.。