计算方法大作业2 高斯消元解线性方程组
- 格式：pdf
- 大小：467.07 KB
- 文档页数：7
种列主元的消去法的主要步骤如下： （1）消元过程 对 k =1,2,…,n-1,进行如下步骤： 1）选主元，记
| a rk | max a ik
ik
若 | a rk | 很小，这说明方程的系数矩阵严重病态，给出警告，提示结果可能不对。
2
《计算方法》大作业（二）
Gauss 消元法解线性方程组
2） 交换增广阵 A 的 r，k 两行的元素。
arj akj ( j k ,..., n 1)
3） 计算消元
aij aij aik akj / akk (i k 1,..., n; j k 1,..., n 1)
（2）回代过程 对 k = n, n-1,…,1,进行如下计算
x k ( a k , n 1
4
《计算方法》大作业（二）
Gauss 消元法解线性方程组
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 X end
（5）
k n 1, n 2,..., 2,1
这种通过消元再回代的求解方法即称为 Gauss 消元法。 由于一般线性方程在使用 Gauss 消去法求解时，从求解的过程中可以看到，若
( k 1) ( k 1) akk 0 ，则必须进行行交换，才能使消去过程进行下去。有的时候即使 akk 0 ，但
1020 时，不选主元时，程序报错：Warning: Divide by zero.这是因为机器计算的最小
精度为 10-15 ，所以此时的 1020 就认为是 0，故出现了错误现象。而选主元时则没有 这种现象，将求得的解代入原方程，可以发现符合得较好，说明了列主元消去法计算的 有效性。在此基础上，如果考虑使用按比例列主元消去法，效果会更优，在算法上也相 对较易实现。
6
采用 Gauss 消元法求解。其中 为一小数，当 105 ,1010 ,1015 ,1020 时，分别采用 列主元和不列主元的 Gauss 消去法求解，并比较结果。 二、 数学原理 对于线性方程组，作如下运算不会改变方程组的解： (1) 交换任意两个方程的顺序； (2) 方程组中任意一个方程乘上一个非零数； (3) 方程组中任一方程加上另一个方程的倍数。 所以，可以对方程组的系数矩阵进行变换，由此有以下的计算方法： 首先，对于方程 AX=b 有
(i 1, 2,..., n; j 1, 2,..., n) ，则原方程组改写为
(0) a1(0) n x n a1n 1 (0) (0) a2 n x n a 2 n 1
（1）
（0） （0 ） a 11 x1 a12 x2 （0） （0 ） a21 x1 a22 x2 （0） （0 ） a n1 x1 an2 x2
是其绝对值非常小，由于机器舍入误差的影响，消去过程也会出现不稳定的现象，导致 结果不正确。因此有必要进行列主元技术，以最大可能的消除这种现象。这一技术要寻 找行 r，使得
( k 1) ( k 1) | a rk | max a ik ik
( k 1) 并将第 r 行和第 k 行的元素进行交换，以使得当前的 akk 的数值比 0 要大的多。这
（2）
(0) (0) ann xn ann 1
(0) 0 ，在增广矩阵中每一行减去第一行的 不妨设 a11
ai(10 )
(0) a11
倍，则方程组变为
（0 ） （0 ） (0) a a1(0) 11 x1 a12 x2 n x n a1 n 1 (1) (1) (1) a22 x2 a 2 n xn a 2 n 1 (1) (1) (1) an ann xn ann 2 x2 1
（3）
1
《计算方法》大作业（二）
Gij aij
ai(0) 1
(0) a11
a1(0) j
(i 2, 3,..., n; j 2, 3,..., n 1)
(1) (1) 若式（3）中 a22 为主元素，消去第 3 至 n 个方程中的第 2 个未知量。 0 ，则以 a22
五、 计算结果 经过编程求解，得到计算结果如下表：
x x1
不选主元
选主元
0.999993478281169 2.00002173901126 2.99997608710494 1.00000363545405 1.99999502380626 3.00000165856458
0.999993478303404 2.00002173898866 2.99997608711248 0.999999999934784 2.00000000021739 2.99999999976087
（4）
将方程组（1）化至（4）的过程为消元过程，最后由回代过程可以求得原方程组的解为
( n 1) ann 1 x n ( n 1) ann n ( k 1) ( k 1) a xj akj kn 1 j k 1 xk ( k 1) akk
j k 1
a
n
kj
x j / a kk )
至此，完成整个方程组的求解。 三、 算法框图
开始
输入增广矩阵 A
Y
ark=0? N
k =1 交换 A 中的 r、k 两行 输出奇异标志
k<n-1?
结束
3
《计算方法》大作业（二）
Gauss 消元法解线性方程组
四、 程序实现 编写 Gauss 消元法函数如下：
1015
x2 x3 x1
1020
x2 x3
六、 结果分析 从以上结果可以看到， 当 105 时， 此时， 由于 不是很小， 机器误差就不是很大， 可以看出不选主元的计算结果精度还可以，因此此时可以考虑不选主元以减少计算量。 当 1010 时， 此时可以看出不选主元的计算精度不太好， 误差开始增大。 当 1015 时， 可以看出不选主元的结果误差比较大，是不正确的结果，这是由机器误差引起的。当
编写列主元Gauss消元法函数如下：
function [RA,RB,n,X]=lgauss(A,b) B=[A b]; n=length(b); RA=rank(A); RB=rank(B); zhica=RB-RA; if zhica>0 disp('请注意：因为RA~=RB，所以此方程组无解.') return
function [RA,RB,n,X]=gauss(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 X end
以此类推，重复上述过程 n-1 次后，我们得到与原方程组等价的系数矩阵为三角阵的方 程组：
（0 ） （0 ） （0 ） (0) a a1(0) 11 x1 a12 x2 a13 x3 n x n a1 n 1 (1) (1) (1) (1) a22 x2 a23 x3 a2 n xn a 2 n 1 (2) (2) (2) a33 x3 a3 n xn a3 n 1 ( n 1) ( n 1) ann xn ann 1
a11 x1 a12 x2 a21 x1 a22 x2 an1 x1 an 2 x2
(0) (0) aij , ain 为了符号统一，记 aij 1 bi
a1n xn b1 a2 n xn b2 ann xn bn
《计算方法》
上机大作业（二）
题目： Gauss 消元法解线性方程组
姓名： 班级：
李楠 B1603193
学号： 116031910075
2016 年 11 月
《计算方法》大作业（二）
Gauss 消元法解线性方程组
一、 题目描述 本实验对于线性方程组
