实验五 解线性方程组的直接方法

数值分析

实验五 解线性方程组的直接方法

实验5.1 （主元的选取与算法的稳定性）

问题提出：Gauss消去法是我们在线性代数中已经熟悉的。但由于计算机的数值运算是在一个有限的浮点数集合上进行的，如何才能确保Gauss消去法作为数值算法的稳定性呢？Gauss消去法从理论算法到数值算法，其关键是主元的选择。主元的选择从数学理论上看起来平凡，它却是数值分析中十分典型的问题。

实验内容：考虑线性方程组

Ax b,A Rn n,b Rn

编制一个能自动选取主元，又能手动选取主元的求解线性方程组的Gauss消去过程。

实验要求：

61 7 861 15 （1）取矩阵A ,b ，则方程有解x* (1,1, ,1)T。 861 15 86 14

取n=10计算矩阵的条件数。让程序自动选取主元，结果如何？

（2）现选择程序中手动选取主元的功能。每步消去过程总选取按模最小或按模尽可能小的元素作为主元，观察并记录计算结果。若每步消去过程总选取按模最大的元素作为主元，结果又如何？分析实验的结果。

（3）取矩阵阶数n=20或者更大，重复上述实验过程，观察记录并分析不同的问题及消去过程中选择不同的主元时计算结果的差异，说明主元素的选取在消去过程中的作用。

（4）将上述矩阵A中的主元改为0.00006再重新作一次数值实验看看。

（5）选取其他你感兴趣的问题或者随机生成矩阵，计算其条件数。重复上述实验，观察记录并分析实验结果。

实验5.2（线性代数方程组的性态与条件数的估计）

问题提出：理论上，线性代数方程组Ax b的摄动满足

xcon(dA) x1 A 1 A A b A b

矩阵的条件数确实是对矩阵病态性的刻画，但在实际应用中直接计算它显然不现实，因为计算A 1通常要比求解方程Ax b还困难。

实验内容：MATLAB中提供有函数“condest”可以用来估计矩阵的条件数，它给出的是按1-范数的条件数。首先构造非奇异矩阵A和右端，使得方程是可以精确求解的。再人为地引进系数矩阵和右端的摄动 A和 b，使得 A b充

数值分析

分小。

实验要求：

b b，以1-范数，（1）假设方程Ax=b的解为x，求解方程(A A)x

给出 x

x xx

x的计算结果。

（2）选择一系列维数递增的矩阵（可以是随机生成的），比较函数“condest”所需机器时间的差别.考虑若干逆是已知的矩阵，借助函数“eig”很容易给出cond2(A)的数值。将它与函数“cond(A,2)”所得到的结果进行比较。

（3）利用“condest”给出矩阵A条件数的估计，针对（1）中的结果给出 x

x

的理论估计，并将它与（1）给出的计算结果进行比较，分析所得结果。注意，如果给出了cond(A)和A的估计，马上就可以给出A 1的估计。

（4）估计著名的Hilbert矩阵的条件数。

H (hi,j)n n,hi,j 1,i j 1i,j 1,2, ,n

思考题一：（Vadermonde矩阵）设

1 1

A 1

1 x0x1x2 xn2x0 x12 2x2 2xn ni x0 0n i x0n xi n 1x1 i 0n x2，b ni ， x2 i 0

n xn ni xn i 0

其中，xk 1 0.1k,k 0,1, ,n，

（1）对n=2,5,8，计算A的条件数；随n增大，矩阵性态如何变化？

（2）对n=5，解方程组Ax=b；设A的最后一个元素有扰动10-4，再求解Ax=b

（3）计算（2）扰动相对误差与解的相对偏差，分析它们与条件数的关系。

（4）你能由此解释为什么不用插值函数存在定理直接求插值函数而要用拉格朗日或牛顿插值法的原因吗？

数值分析

相关 MATLAB 函数提示： zeros(m,n) 生成 m 行，n 列的零矩阵 ones(m,n) 生成 m 行，n 列的元素全为 1 的矩阵 生成 n 阶单位矩阵 eye(n) rand(m,n) 生成 m 行,n 列(0,1)上均匀分布的随机矩阵 diag(x) 返回由向量 x 的元素构成的对角矩阵 tril(A) 提取矩阵 A 的下三角部分生成下三角矩阵 triu(A) 提取矩阵 A 的上三角部分生成上三角矩阵 返回矩阵 A 的秩 rank(A) det(A) 返回方阵 A 的行列式 inv(A) 返回可逆方阵 A 的逆矩阵 [V,D]=eig(A) 返回方阵 A 的特征值和特征向量 norm(A,p) 矩阵或向量的 p 范数 cond(A,p) 矩阵的条件数 [L,U,P]=lu(A) 选列主元 LU 分解 R=chol(X) 平方根分解 Hi=hilb(n) 生成 n 阶 Hilbert 矩阵