摘要
本文叙述了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）进行回代求解，即得方程组的解。手算的过程是非常清楚的，现在需回答的是计算机求解，如何实现上述计算过程。
设线性方程组为
1111221331121122223322
112233n n n n n n n nn n n
a x a x a x a x
b a x a x a x a x b a x a x a x a x b +++???+=??+++???+=??
?????????????
?+++???+=? 写成矩阵形式为
111211121222222122
2m m m n n a a a x b a
a a x
b 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 =。
第一次消元。设(1)
110a ≠，
为将第二个方程至第n 个方程的1x 系数(1)1i a 消成零，构造乘数
(1)1
1(1)11
i i a l a = (2,3,,i n =
用1i l -乘以矩阵(1)[A b ]（1）
的第一行加到第i 行上(2,3,
,)i n =得
(1)
(1)
(1)(1)(1)11121311(2)(2)(2)(2)(2)22
23
22
(2)
(2)
(2)
(2)230
[A
b ]0n
n
n n nn
n a a a a b a a a b a a a b ?????
?=???????
?
（2）
其中，
(2)(1)
(1)11(2)(1)(1)11
a a (,2,3,,)(,2,3,
,)
ij ij i j i i i l a i j n b b l b i j n ?=-=?=-=?
假设经过1k -次消元得同解方程组()()k k A x b =，此时
(1)(1)
(1)(1)(1)11121311(2)(2)(2)(2)k ()22
23
22
()()()
()()n
k n k k kn k k k k nk
nn
n a a a a b a a a b A
b a b a a b ????????=???
??????
?
（）
(1.1.3 假若()
0k kk a ≠，则第k 次消元如下：构造乘数
()
()
(1,,)k ik
ik k kk
a l i k n a ==+
用ik l -乘以矩阵()()k k A b ???
?第k 行加到第i 行上(1,
,)i k n =+得同解方程
组(1)(1)k k A x b ++=，其中，(1)(1)[]k k A b ++中的元素计算公式为
(1)()()(1)
()()a a (,2,3,,)(,2,3,
,)k k k ij ij ik kj k k k i
i ik k l a i j n b b l b i j n ++?=-=?=-=?
如此进行1n -次消元，即1,2,,1k n =-，最后得同解方程组
()()n n A x b =
其中，
(1)
(1)
(1)(1)(1)11121311(2)(2)(2)(2)()()
k ()22
23
22
()()n
n n k n
n n nn
n a a a a b a a a b A x b
A
b a b ????????==?????????
?
（）
假若()
0n nn a ≠，对式（3.1.3）进行回代计算得方程组的解
()
()()()1/()/(1,2,1)
n n n nn n i i i i i ij j ii
j i x b a x b a x a i n =+?=??=-=-??
∑
上述全过程称为Gauss 顺序消元法。
1.2Gauss 主元素消元法
Gauss 顺序消元法每一步总有一个要求()
0k kk a ≠，否则算法就失败。然而从方
程组的理论，只要det()0A ≠，则方程组解存在且唯一。由此可知，Gauss 顺序
消元法可执行条件苛刻。不仅如此，从数值计算角度，当()0k kk a ≠，但()k kk
a 很小时。用()k kk a 做分母运算，会引起误差界的放大，数值不稳定现象产生，严重时导
致解失真。为克服这一缺点，常采用选主元的消元法。
Gauss 列主元消元法主要依据线性方程组任意交换两个方程的次序，方程组的同解性不变，且解的分量次序也不变。于是，第k 步在顺序消元法进行之前，
从()k A 的第k 列元素()k kk a ，()
1k k k a +，
，()
k nk a 中选取绝对值最大者，并记录所在行，
即
()()
max k k k i k ik k i n
a a ≤≤=，记k l i =
如果l k ≠，则交换矩阵()()k k A b ???
?的第k 行与第l 行所有对应的元素，然后，再
进行第k 步顺序消元法。
1.3Gauss 全主元消元法
Gauss 全主元消元法是在Gauss 列主元消元法的基础上进一步改进，目的是更好的提高数值稳定性和解的精度。对那些矩阵性态不太理想的方程组能给出最
大可能性的满意解。
全主元算法的第k 步是从()
(,1,
,)k ij a i k k n =+中选取绝对值最大者作为主元，
并记录主元所在列与所在行，即
()()
,max k k k k i j ij k i j n
a a ≤≤=，记,k k l i t j ==
如果l k ≠，则交换矩阵()()k k A b ????的第k 行与第l 行所有对应的元素；
如果t k ≠，则交换矩阵()
()k k A
b ????的第k 行与第t 列所有对应的元素。值得注意的是，进行
行元素的交换，解的分量次序不变，但是进行列元素交换时，解的分量次序将有所改变。
1.4 Gauss Seidel -迭代法
将方程组Ax b =改写成等价方程组
112
11111111
12221222222221200
0n n n n n n n nn nn
nn
a a
b a a a x x a b a x x a a a x x b a
a a a a ????
-
-
??
??????????????????--????
????=+????????????????????
??????
??--
???????
? (1.4.1)
矩阵形式为
J x G x f =+ (1.4.2)
其中
112
111122122221200
0n n J n n nn
nn
a a a a a a a a G a a a a ?
?-
-
??
????--??
=????????
--
????，111222n nn b a b a f b a ??
????????=??????
?????? (1.4.3)
任取初始解向量(0)(0)(0)
(0)12(,,
,),T
n x x x x =带入式(1.4.2)右端得迭代格式
(1)()(0,1,2,)k k J x G x f
k +=+=
其分量形式为
1
()1()/(1,2,,)(0,1,2,)n
k k i
i ij j ii
j j i
x
b a x a i n k +=≠=-==∑ (1.4.4)
此式称为Jacobi 迭代法，简称J 法。 将上式修改为
1
1(1)()1
1
()/(1,2,,)i n
k k k i
i ij j
ij
j
ii
j j i x
b a x
a x
a i n -++==+=--
=∑∑
则为Gauss Seidel -迭代法。
二、计算机的基本配置
版本：7win 旗舰版 安装内存()RAM :2GB
处理器：()Intel R ()Core TM 3i CPU M 350 2.27GHz 系统类型：32位操作系统
三、参数设置及计算结果
在本文中选取一个六元方程组
1234561
234561234561
23456123456123456234564423457873
45469792
123454653444122445894523433415235677435645
x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x +++++=??+++++=??+++++=??
+++++=??+++++=?+++++=?? 这是本人随机写的一个方程组，本文将以之为例，用摘要里说的方法来解方程组。
3.1Gauss 全主元消元法 3.1.1算法流程图：
1) 输入(,1,2,,),(1,2,
,)ij i a i j n b i n ==；
2) ()(1,2,,)z i i i n ==；
3) 对1,2,
,1k n =-做
a) 选主元：,max k k i j ij k i j n
a a ≤≤=，并记,k k l i t j ==，即
max ;;;kk a a l k t k === 对1,,i k n =+做 对1,
,j k n =+做
if max ij a a > then max ,,ij a a l i t j === b) if max 0a =，输出奇异标志，停机； c) if l k ≠ then (,1,
,);kj lj k l a a j k k n b b ?=+?；
if t k ≠ then (1,2,,);()()ik it a a i n z k z t ?=?；
d) 对1,
,i k n =+做
i. /ik ik ik kk a l a a ←=；
ii. i i ik k b b l b =- iii.
对1,
,j k n =+做ij ij ik kj a a l a =-；
4) if 0nn a = then 输出奇异标志，并停机 else 做
a) /n n nn b b a =； b) 1
()/(1,
,1);n
i i ij
j
ii
j i b b a b a i n =+=-
=-∑
c) (())(1,2,,)i x z i b i n ==；
5) 输出解(1,2,
,)i x i n =
3.1.2计算结果及分析（程序见附录）：
123 x =18.5115 x 7.2791 x = -12.0967= 456 x = -4.6715 x -5.9736 x = 15.9624=
系数矩阵为：
1 2 3 4 5 62 3 4 5 7 84 5 4 6 9 7a=12 34 54 65 34 4424 45 89 45 23 4315 23 56 77 43 56?????????????????? 447392b=123445??
????????????????
构成增广矩阵为
1 2 3 4 5 6
442 3 4 5 7 8734 5 4 6 9 79212 34 54 65 34 441224 45 89 45 23 433415 23 56 77 43 56
45A ?????
???=?
??????
???
经过Gauss 全主元消元法消元后矩阵a 变为
89.00 45.0000 45.0000 23.0000 43.0000 24.000 0 48.6854 -5.3146 28.5281 28.9438 -0.1011 0 0 10.8117 -2.0441 -4.5008 -2.4835 0 0 0 6.28a =06 4.1230 3.7133 0 0 0 0 1.9065 -1.4148 0 0 0 0 0 -0.2568??????????????????????
消元后的b 变为
34.000 23.6067 -26.9077 97.0342 4.2430 -4.7529b ????????=??????????
最终将解得的x 带入原方程组，计算出误差矩阵 error AX b =-
下面给出误差矩阵
120.01420.02840.0284100.11370.11370.1137error -????????
=???????????
由误差矩阵可知误差已接近于0，可知：这种算法的具有数值稳定，精度高
的优点，是比较好的算法，但是缺点是在选主元的过程中，消耗了较多的时间。 该实例的运行时间为0.048152 seconds. 3.2Gauss Seidel -迭代法
3.2.1Gauss Seidel -迭代法算法流程图：
1) 输入：(,1,2,
,),(1,2,
,),;ij i a i j n b i n ε===
2) 输入：初始解向量(1,2,,)i x i n =；
3) 对1,2,
,i n =做
a) 1()/n
i i ij j ii j j i
x b a x a =≠=-∑；
b) i i i e y x =-； c) i i x y =； 4)
max{}i i i n
f
e ε≤≤< then 输出：i x (1,2,
,)i n =，停机else 返回第三步。
3.2.2 参数的设置：
本题初始解：(0)(0,0,0,0,0,0)T x =; 误差精度0.0001ε=;
3.2.3计算结果及分析（程序见附录）： 针对这个线性方程组，用迭代法无法求解，因为这个系数矩阵的谱半径不小于1，在本题中谱半径为161.1222远大于1。这个问题的出现从侧面反映了迭代
法的适应面会很小，因为现实生活中的很多矩阵都不是会符合谱半径小于1这个条件的，虽然这种方法迭代速度回很快，但是这个谱半径要求小于1弊端是其致命的弱点。
四、计算结果分析及遇到的问题
4.1计算结果的分析
根据误差矩阵来看，此方法解决此类问题具有很好的效果，误差数量级在12
10 上面，误差是非常小的，几乎可以说是没有，而且运行时间也很快，而且当矩阵变得很大时，运行时间并没有发生明显的变化，下满给出运行时间与矩阵阶数的变化关系的图表（注：由于编程的方便，下表给出的运行时间的程序中使用的a，b均是用随机数来产生的）。
运行时间（s）矩阵阶数
0.042626 2
0.046684 4
0.047888 8
0.044901 16
0.046488 32
0.050995 64
0.084790 128
0.475078 256
6.029018 512
55.884309 1024
通过上表可知:当矩阵阶数太大时，此方法不是一个很好的方法，当矩阵的
阶数小于256阶时，采用此方法比较好，因为上述数据是用服从正态分布的
随机数的出来的，结果固然具有一定的波动性，但是这种趋势是正确的，我
们画出其变化曲线，结果会更加一目了然。
从这幅图像中，我们知道运行时间并不与阶数成正比，而是一种类似二次曲线的效果，所以阶数越大，此种方法越不可行。 4.2遇到的问题
在编这个算法的程序时有一个语句1:1for
j n =- 这语句一直没有发现毛
病，运行是没有报错，但是因为得出了错误的结果，便开始查程序。查程序时发现没有进入到循环里面去，突然意识到应该写成1:1:1for j n =--，
犯了一个很低级的错误，害的我找了半天的错误。
五、算法的改进
针对算法的改进，本人没有什么好的见解，一般的Gauss 消元法（包括本文叙述的这三种），都要采用回代的方法解出x ，在这里推荐使用Gauss Jordan -消元法，该算法是不需要回代过程的。Gauss Jordan -消元法在形式上要比Gauss 消元法简单一些，没有回代过程，消元完毕就获得解答，
但从工作量角度来看，它大约要31()2O n ，比有回代的Gauss 消元法31()3O n 要
多出31
()6
O n 。
附录
%%%%%%%%%%%%%%%%%%Gauss 全主元消元法%%%%%%%%%%%%%%%%% tic clc; clear all ;
a=[1 2 3 4 5 6;2 3 4 5 7 8;4 5 4 6 9 7;12 34 54 65 34 44;24 45 89 45 23 43;15 23 56 77 43 56]; b=[44;73;92;12;34;45]; A=a; B=b;
n=size(a,1); %A=[a b]; z=1:n; for i=1:n-1
[I J]=find(abs(a)==max(max(abs(a(i:end,i:end)))));
l=I(1);t=J(1);
number=a(l,t);
if number(1)==0
disp('method failed');
break;
end
k=i;
if i~=l
temp=a(l,:);temp1=b(l);
a(l,:)=a(i,:);b(l)=b(i);
a(i,:)=temp;b(i)=temp1;
end
if i~=t
temp=a(:,t);
a(:,t)=a(:,i);
a(:,i)=temp;
temp1=z(i);
z(i)=z(t);
z(t)=temp1;
end
k=i;
ll=a(k+1:end,k)./a(k,k);
a(k+1:end,:)=a(k+1:end,:)-ll*a(k,:);
b(k+1:end)=b(k+1:end)-b(k).*ll;
end
if a(n,n)==0
disp('sigular');
else
b(n)=b(n)./a(n,n);
for j=n-1:-1:1
b(j)=(b(j)-a(j,j+1:end)*b(j+1:end))./a(j,j);
end
end
x(z)=b;
disp(x);
error=A*x'-B;
disp(error);
toc
%%%%%%%%%%%%%%%%%%%% Gauss-seidel法%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clc;
clear all;
a=[1 2 3 4 5 6;2 3 4 5 7 8;4 5 4 6 9 7;12 34 54 65 34 44;24 45 89 45 23 43;15 23 56 77 43 56];
b=[44;73;92;12;34;45];
lamada=0.00001%?ó2????è
x=[0 0 0 0 0 0]';%3?ê??a
A=a;
n=size(a,1);
e=ones(n,1);
[U V]=eig(A);
ll=diag(V);
G=max(abs(ll));
if G<1
while max(e)>lamada
for i=1:n
s=A(i,:)*x-A(i,i)*x(i); y=(b(i)-s)/A(i,i);
e(i)=y-x(i);
x(i)=y;
end
disp(s);
end
disp(x);
else
disp('method failed');
end
解线性方程组的迭代法 1.rs里查森迭代法求线性方程组Ax=b的解 function[x,n]=rs(A,b,x0,eps,M) if(nargin==3) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值elseif(nargin==4) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-A)*x0+b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 2.crs里查森参数迭代法求线性方程组Ax=b的解 function[x,n]=crs(A,b,x0,w,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-w*A)*x0+w*b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x;
if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 3.grs里查森迭代法求线性方程组Ax=b的解 function[x,n]=grs(A,b,x0,W,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1;%前后两次迭代结果误差 %迭代过程 while(tol>eps) x=(I-W*A)*x0+W*b;%迭代公式 n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 4.jacobi雅可比迭代法求线性方程组Ax=b的解 function[x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200; elseif nargin<3 error return elseif nargin==5 M=varargin{1}; end D=diag(diag(A));%求A的对角矩阵 L=-tril(A,-1);%求A的下三角阵
高斯列主元消元法解线性方程组 一、题目：用Gauss 列主元消去法解线性方程组Ax b =，其中， A=17.031 -0.615 -2.991 1.007 -1.006 0.000-1.000 34.211 -1.000 -2.100 0.300 -1.7000.000 0.500 13.000 -0.500 1.000 -1.5004.501 3.110 -3.907 -61.705 12.170 8.9990.101 -8.012 -0.017 -0.910 4.918 0.1001.000 2.000 3.000 4.500 5.000 21.803?? ? ? ? ? ? ? ? ??? 0.230 -52.322 54.000 240.236 29.304 -117.818b ?? ? ? ?= ? ? ? ? ??? T X=(0.907099 -1.961798 3.293738 -4.500708 3.029344 -5.255068) 二、原理及步骤分析 设 n n ij R a A ?∈=][)1(，n n R b b b b ∈=],,,[)1()2(2)1(1 。若约化主元素 ),,2,1(0)(n k a k kk =≠，则通过高斯消元法将方程b AX =约化为三角形方程组求解。 如果在消元过程中发现某个约化主元0) (=k kk a , 则第K 次消元就无法进行。此外，即 使所有约化主元全不为零，虽然可以完成方程组的求解，但也无法保证结果的可靠性，因为计算过程中存在舍入误差。 为减少计算过程中的舍入误差对解的影响，在每次消元前，应先选择绝对值尽可能大的元作为约元的主元，如果在子块的第一列中选取主元，则相应方法称为列主元消元法。相应过程为： （1）选主元：在子块的第一列中选择一个元) (k k i k a 使) (max k ik n i k k k i a a k ≤≤= 并将第k 行元与第k i 行元互换。 （2）消元计算：对k=1,2,……n-1依次计算 ()()()?? ?? ?????++=-=++=-=++==++n k k i b m b b n k k j i a m a a n k k i a a m k k ik k i k i k kj ik k ij k ij k kk k ik k ik ,,2,1,,2,1,,,2,1) ()()1() ()()1()() ()( （3）回代求解
用高斯消元法求解线性代数方程组 1234111 5 -413-2823113-2104151 3-21719x x x x ??????????????????=?????? ?????? ?????? 1111X *??????=?????? （X*是方程组的精确解） 1 高斯消去法 1.1 基本思想及计算过程 高斯（Gauss ）消去法是解线性方程组最常用的方法之一，它的基本思想是通过逐步消元，把方程组化为系数矩阵为三角形矩阵的同解方程组，然后用回代法解此三角形方程组得原方程组的解。 为便于叙述，先以一个三阶线性方程组为例来说明高斯消去法的基本思想。 ??? ??=++II =++I =++III) (323034)(5 253)(6432321 321321x x x x x x x x x 把方程（I ）乘（2 3 - ）后加到方程（II ）上去，把方程（I ）乘（2 4- ）后加到方程（III ）上 去，即可消去方程（II ）、（III ）中的x 1，得同解方程组 ?? ? ??=+-II -=-I =++III) (20 223)(445.0)(6 4323232321x x x x x x x 将方程（II ）乘（ 5 .03 ）后加于方程（III ），得同解方程组： ?? ? ??-=-II -=-I =++III) (42)(445.0)(6432332321x x x x x x 由回代公式（3.5）得x 3 = 2，x 2 = 8，x 1 = -13。 下面考察一般形式的线性方程组的解法，为叙述问题方便，将b i 写成a i , n +1，i = 1, 2,…,n 。
作业六：分别编写用Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组Ax=B的标准程序，并求下列方程组的解。 可取初始向量 X(0) =(0，0，0)’； 迭代终止条件||x(k+1)-x(k)||<=10e-6 （1） = （2） = Jacobi迭代法： 流程图 开 始 判断b中的最大值 有没有比误差大 给x赋初值 进行迭代 求出x，弱到100次还没到，警告不收 结束
程序 clear;clc; A=[8,-1,1;2,10,01;1,1,-5]; b=[1;4;3]; e=1e-6; x0=[0;0;0]'; n=length(A); x=zeros(n,1); k=0; r=max(abs(b)); while r>e for i=1:n d=A(i,i); if abs(d)
程序结果（1）
（2）
Gauss-Seidel迭代法： 程序 clear;clc; %A=[8,-1,1;2,10,01;1,1,-5]; %b=[1;4;3]; A=[5,2,1;-1,4,2;2,-3,10]; b=[-12;20;3]; m=size(A); if m(1)~=m(2) error('矩阵A不是方阵'); end n=length(b); %初始化 N=0;%迭代次数 L=zeros(n);%分解A=D+L+U,D是对角阵，L是下三角阵，U是上三角阵U=zeros(n); D=zeros(n); G=zeros(n);%G=-inv(D+L)*U d=zeros(n,1);%d=inv(D+L)*b x=zeros(n,1); for i=1:n%初始化L和U for j=1:n if i
求解线性方程组——超松弛迭代法 #include
cout<<"输入松弛因子w (1 Gauss-Seidel迭代法求解线性方程组 一． 问题描述 用Gauss-Seidel 迭代法求解线性方程组 由Jacobi 迭代法中，每一次的迭代只用到前一次的迭代值。使用了两倍的存储空间，浪费了存储空间。若每一次迭代充分利用当前最新的迭代值，即在计算第i 个分量 ) 1(+k i x 时，用最新分量 ) 1(1 +k x ， ???+) 1(2 k x ) 1(1 -+k i x 代替旧分量 ) (1 k x ， ???) (2 k x ) (1 -k i x ，可以起 到节省存储空间的作用。这样就得到所谓解方程组的Gauss-Seidel 迭代法。 二． 算法设计 将A 分解成U D L A --=，则b x =A 等价于b x =--U)D (L 则Gauss-Seidel 迭代过程 ) ()1()1(k k k Ux Lx b Dx ++=++ 故 ) ()1()(k k Ux b x L D +=-+ 若设1 )(--L D 存在，则 b L D Ux L D x k k 1)(1)1()()(--+-+-= 令 b L D f U L D G 11)()(---=-=， 则Gauss-Seidel 迭代公式的矩阵形式为 f Gx x k k +=+) () 1( 其迭代格式为 T n x x x x ) ()0()0(2)0(1)0(，，，???= (初始向量)， ) (1 1 1 1 1 )()1()1(∑∑-=-+=++--=i j i i j k j ij k j ij i ii i i x a x a b a x )210i 210(n k ???=???=，，，；，，， 或者 ?? ???--=???=???==?+=∑∑-=-+=+++) (1)210i 210(111 1)()1()1()()1(i j i i j k j ij k j ij i ii i i i k i k i x a x a b a x n k k x x x ，，，；，，， 三． 程序框图 迭代法解线性方程组作业 沈欢00986096 北京大学工学院，北京100871 2011年10月12日 摘要 由所给矩阵生成系数矩阵A和右端项b，分析系数矩阵A，并用Jacobi迭代法、GS迭代法、SOR(逐步松弛迭代法)解方程组Ax=b 1生成系数矩阵A、右端项b,并分析矩阵A 由文件”gr900900c rg.mm”得到了以.mm格式描述的系数矩阵A。A矩阵是900?900的大型稀 疏对称矩阵。于是，在matlaB中，使用”A=zeros(900,900)”语句生成900?900的零矩阵。再 按照.mm文件中的描述，分别对第i行、第j列的元素赋对应的值，就生成了系数矩阵A，并 将A存为.mat文件以便之后应用。 由于右端项是全为1的列向量，所以由语句”b=ones(900,1)”生成。 得到了矩阵A后，求其行列式，使用函数”det(A)”,求得结果为”Inf”,证明行列式太大，matlaB无法显示。由此证明，矩阵A可逆，线性方程组 Ax=b 有唯一解。 接着，判断A矩阵是否是对称矩阵(其实，这步是没有必要的，因为A矩阵本身是对称矩阵，是.mm格式中的矩阵按对称阵生成的)。如果A是对称矩阵，那么 A?A T=0 。于是，令B=A?A T,并对B求∞范数。结果显示： B ∞=0,所以，B是零矩阵，也就是：A是对称矩阵。 然后，求A的三个条件数： Cond(A)= A ? A?1 所求结果是，对应于1范数的条件数为：377.2334;对应于2范数的条件数为：194.5739;对应 于3范数的条件数为：377.2334; 1 从以上结果我们看出，A是可逆矩阵，但是A的条件数很大，所以，Ax=b有唯一解并且矩阵A相对不稳定。所以，我们可以用迭代方法来求解该线性方程组，但是由于A的条件数太大迭代次数一般而言会比较多。 2Jacobi迭代法 Jacobi迭代方法的程序流程图如图所示： 图1:Jacobi迭代方法程序流程图 在上述流程中，取x0=[1,1,...,1]T将精度设为accuracy=10?3，需要误差满足： error= x k+1?x k x k+1 消元法解线性方程组 学校：青海师范大学 院系：数学系 专业：数学与应用数学 班级：10B 指导教师：邓红梅 学号：20101611218 姓名：梅增旺 摘要：线性方程组在数学的各个分支，在自然科学，工程技术，生产实际中经常遇到，而且未知元的个数及方程的个数可达成百上千，因此它的理论是很重要的，其应用也很广泛。本篇将就解线性方程组在此做一浅谈，以消元法为主要方法。消元法是解一般线性方程组行之有效的方法，早在中学大家都已经有接触，消元法的基本思想是通消元变形把方程组化成容易求解的同解方程组进行求解。 关键字：线性方程组消元法求解 Abstract: linear equations in various branches of mathematics, natural science,engineering technology, often encountered in actual production, and the unknown element number and the number of equations can be hundreds, so itis important in the theory, its application is very extensive. This article on thesolution of linear equations based on a discussion, mainly by means ofelimination method. Elimination method is the general linear equations ofeffective early in high school, everyone has a contact, the basic idea ofelimination method is through the elimination of the equations of deformationinto easy to solve with the solution of equations. Keywords:elimination method for solving linear equations 西京学院数学软件实验任务书 课程名称数学软件实验班级数学0901 学号0912020112 姓名*** 实验课题 线性方程组高斯消去法，高斯列主元消去法，高斯全 主元消去法 实验目的熟悉线性代数方程组高斯消去法，高斯列主元消去法，高斯全主元消去法 实验要求运用Matlab/C/C++/Java/Maple/Mathematica等其中一种语言完成 实验内容线性方程组高斯消去法 线性方程组高斯列主元消去法线性方程组高斯全主元消去法 成绩教师 实 验 报 告 实验名称：Guass 消元法编程求解线性方程 实验目的：进一步熟悉理解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步消元。 对线性方程组中有第2，3，.。。。N 个方程组成的n —1元线性方程组做同样的处理，消去其除第一个方程组之外的所有变元2x ，可得到 ???? ?? ? ??????=++=++=++=++)3()3(3)3(3)2(3)2(33)2(33)1(2)1(22)1(22)0(1)0(11)0(11... . . . ... ... ...n n nn n n n n n n n b x a x a b x a x a b x a x a b x a x a 依次类推，当做到n-1步消元后，就完成了Guass 消元过程，得到上三角方程组 实验内容：利用Guass 消元操作的原理，求解线性方程组 ?? ?? ? ????==++=++--) 1()1()1(2)1(22)1(22) 0(1)0(11)0(11 . . ... ...n n n n nn n n n n b x a b x a x a b x a x a 回代过程： 在最后的一方程中解出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 b x 流程图如下： 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; 解线性方程组的消元法及其应用 （朱立平 曲小刚） ● 教学目标与要求 通过本节的学习，使学生熟练掌握一种求解方程组的比较简便且实用的方法—高斯消元法，并能够熟练应用消元法将矩阵化为阶梯形矩阵和求矩阵的逆矩阵. ● 教学重点与难点 教学重点：解线性方程组的高斯消元法，利用消元法求逆矩阵. 教学难点：高斯消元法，利用消元法求逆矩阵. ● 教学方法与建议 先向学生说明由于运算量的庞大，克莱姆法则在实际应用中是很麻烦的，然后通过解具体的方程组，让学生自己归纳出在解方程组的时候需要做的三种变换，从而引出解高阶方程组比较简便的一种方法—高斯消元法，其三种变换的实质就是对增广矩阵的初等行变换，最后介绍利用消元法可以将矩阵化为阶梯形矩阵以及求矩阵的逆。 ● 教学过程设计 1.问题的提出 由前面第二章的知识，我们知道当方程组的解唯一的时候，可以利用克莱姆法则求出方程组的解，但随着方程组阶数的增高，需要计算的行列式的阶数和个数也增多，从而运算量也越来越大，因此在实际求解中该方法是很麻烦的. 引例 解线性方程组 ??? ??=+-=+=++132724524321 21321x x x x x x x x )3()2()1( 解 （1）???→??)2()1(?????=+-=++=+13245247 232132121x x x x x x x x )3()2()1(????→?+-?+-?) 3()2()1()2()4()1(?????-=+-=+=+133524567232 3221x x x x x x )3()2()1( ????→?+-?)3()65 ()2(??????? =--=+=+76 724567233221x x x x x )3()2()1( 用回代的方法求出解即可. 问题：观察解此方程组的过程，我们总共作了三种变换：（1）交换方程次序，（2）以不等于零的数乘某个方程，（3）一个方程加上另一个方程的k 倍.那么对于高阶方程组来说，是否也可以考虑用此方法. 2.矩阵的初等变换 定义1 阶梯形矩阵是指每一非零行第一个非零元素前的零元素个数随行序数的增加而增加的矩阵. 定义2 下面的三种变换统称为矩阵的初等行变换： i. 互换矩阵的两行（例如第i 行与第j 行，记作j i r r ?）， ii. 用数0≠k 乘矩阵的某行的所有元素（例如第i 行乘k ，记作i kr ）， iii. 把矩阵某行的所有元素的k 倍加到另一行的对应元素上去（例如第j 行的k 倍加到第i 行上，记作j i kr r +）. 同理可以定义矩阵的初等列变换. 定义 3 如果矩阵A 经过有限次初等变换变为矩阵B ，则称矩阵A 与B 等价，记作 A ~ B . 注：任意一个矩阵总可以经过初等变换化为阶梯形矩阵. 3. 高斯消元法 对于一般的n 阶线性方程组 ?????? ?=++=+++=+++n n nn n n n n n n b x a x a x a b x a x a x a b x a x a x a ΛΛΛΛΛΛΛΛΛΛΛΛΛΛΛ22112 22221211 1212111 )()2()1(n （3.1） 若系数行列式0det ≠A ，即方程组有唯一解，则其消元过程如下： 第一步，设方程（1）中1x 的系数01≠l a 将方程)(l 与（1）对调，使对调后的第一个方程1x 的系数不为零.作)1(11 1 a a i i - ),3,2(n i Λ=，得到同解方程组 ?? ? ????=++=++=+++)1()1(2)1(2) 1(2 )1(22)1(22)0(1)0(12)0(121)0(11n n nn n n n n n b x a x a b x a x a b x a x a x a ΛΛΛΛΛΛΛΛΛΛΛΛΛ （3.2） 第二步，设0) 1(22≠a ，保留第二个方程，消去它以下方程中的含2x 的项，得 数值分析实验报告 0 a 12 K a 1,n 1 K a 2,n 1 U O M 则有： 第一步： Jacobi 迭代法 a 1n a 2n M , 则有： A D L U a n 1,n Ax b A A x D b L U (D L U)x b Dx (L U)x b x D (L U)x D b 令 J D (L U) 则称 J 为雅克比迭代矩阵 f D b 由此可得雅克比迭代的迭代格式如下： x (0) , 初始向量 x (k 1) Jx (k) f ,k 0,1,2,L 第二步 Gauss-Seidel 迭代法 Ax b (D L U )x b (D L)x Ux b x (D L) Ux (D L) b A D L U a 11 a 12 L a 1n a 11 A a 21 a 22 L a 2n a 22 M MM MO a n1 a n2 L a nn a 11 得到 D a 22 O a nn 由 a 21 0 M M O a n 1,1 a n 1,2 L 0 a nn a n1 a n2 L a n,n a 21 L M M O a n 1,1 a n 1,2 L a n1 a n2 L a n,n 1 a 12 K a 1,n 1 a 1n 0 K a 2,n 1 a 2n O M M a n 1,n 10 令 G (D L) U ，则称G 为Gauss-Seidel 迭代矩阵 f (D L) b 由此可得 Gauss-Seidel 迭代的迭代格式如下： x (0) , 初始向量 第三步 SOR 迭代法 w0 AD L U 1 ( D 1 wL ((1 w)D wU )) (D 1 wL) ((1 w)D wU ) w w w 令M w 1 (D wL), N 1 ((1 w)D wU )则有：A MN w w Ax b AM L W N M (M N )x b Mx Nx b x M Nx M b N M, 令W f Mb 带入 N 的值可有 L W ((1 w)D wU) (D wL) 1((1 w)D wU) (D wL) f 1 b w 1(D wL) 1b 1 (D wL) w 称 L W 为 SOR 迭代矩阵，由此可得 SOR 迭代的迭代格式如下： x (0) ,初始向量 二、算法程序 Jacobi 迭代法的 M 文件： function [y,n]=Jacobi(A,b,x0,eps) %************************************************* %函数名称 Jacobi 雅克比迭代函数 %参数解释 A 系数矩阵 % b 常数项 % x0 估计解向量 x (k 1) Gx (k) f ,k 0,1,2,L (k 1) f,k 0,1,2,L 一． 问题描述 用Gauss-Seidel 迭代法求解线性方程组 由Jacobi 迭代法中，每一次的迭代只用到前一次的迭代值。使用了两倍的存储空间，浪 费了存储空间。若每一次迭代充分利用当前最新的迭代值，即在计算第i 个分量) 1(+k i x 时， 用最新分量) 1(1 +k x ，???+) 1(2 k x ) 1(1 -+k i x 代替旧分量)(1k x ，???) (2 k x ) (1-k i x ，可以起到节省存储 空间的作用。这样就得到所谓解方程组的Gauss-Seidel 迭代法。 二． 算法设计 将A 分解成U D L A --=，则b x =A 等价于b x =--U)D (L 则Gauss-Seidel 迭代过程 ) ()1()1(k k k Ux Lx b Dx ++=++ 故 )()1()(k k Ux b x L D +=-+ 若设1 )(--L D 存在，则 b L D Ux L D x k k 1)(1)1()()(--+-+-= 令 b L D f U L D G 11)()(---=-=， 则Gauss-Seidel 迭代公式的矩阵形式为 f Gx x k k +=+)()1( 其迭代格式为 T n x x x x )()0()0(2)0(1)0(，，，???= (初始向量)， )(1111 1 )() 1()1(∑∑-=-+=++--=i j i i j k j ij k j ij i ii i i x a x a b a x )210i 210(n k ???=???=，，，；，，， 或者 ?? ???--=???=???==?+=∑∑-=-+=+++) (1)210i 210(111 1)() 1()1()()1(i j i i j k j ij k j ij i ii i i i k i k i x a x a b a x n k k x x x ，，，；，，， 三． 程序框图 河北工业大学计算机软件技术基础（VC）课程设计报告 学院管理班级管理104班姓名杨立宝 __ 学号 101707____ 成绩 __ ____ 一、题目： 求线性代数方程组的解（高斯消去法）（C13） 二、设计思路 1、总体设计 1）分析程序的功能 第一：编写输入程序，通过键盘先输入对应的已知量及函数的大小n和系数a[i]和得数b[i]。 第二：编写中间程序，通过函数的调用先定义线性代数方程，然后通过程序求出方程的梯形矩阵系数，并最终得出结果。 第三编写输出程序，输出最终结果。 2）系统总体结构：设计程序的组成模块，简述各模块功能。 模块一：各函数的具体内容 A：三个输入函数，分别输入n，一维数组，二维数组。即输入已知量。 B：中间运算函数，计算是使得方程系数所成的矩阵成梯形矩阵，未知数的结果。即计算中间变量及结果。 C：最后输出函数，输出最后计算结果。 模块二：各函数原型的声明 a写头文件。 b变量声明:存放输入数据的数组的声明，存放中间变量的数组的声明，存放运算结果的数组的声明。分别存放对应数据。 c输入有关操作的文字 d函数调用，在运算中自动调用对应的函数解决对应问题。 模块三:主函数 2、各功能模块的设计：说明各功能模块的实现方法 模块一：各个函数的声明，直接声明。 模块二：各函数都通过for循环来实现各个数组之间的基本运算。 3、设计中的主要困难及解决方案 在这部分论述设计中遇到的主要困难及解决方案。 1）困难1 函数调用是怎么用？ 解决方案： 仔细阅读课本，以及同学之间的讨论，和老师的帮助。 4、你所设计的程序最终完成的功能 1）说明你编制的程序能完成的功能 输入线性代数的系数后，运行程序即可得到梯形矩阵和结果。 2）准备的测试数据及运行结果 三、程序清单 如果是使用一个文件完成的程序，只需列出程序代码。 如果是使用多文件完成的程序，首先说明程序中的代码存放在哪些文件中，说明文件名（例如：本程序包含first.cpp、second.cpp、third.cpp和all.h四个文件）；然后依次给出每个文件名及该文件清单，例如： #include 用高斯消元法和列主元消去法求解线性代数方程组 （X*是方程组的精确解） 1 高斯消去法 1.1 基本思想及计算过程 高斯（Gauss ）消去法是解线性方程组最常用的方法之一，它的基本思想是通过逐步消元，把方程组化为系数矩阵为三角形矩阵的同解方程组，然后用回代法解此三角形方程组得原方程组的解。 为便于叙述，先以一个三阶线性方程组为例来说明高斯消去法的基本思想。 ??? ??=++II =++I =++III) (323034)(5 253)(6 432321 321321x x x x x x x x x 把方程（I ）乘（2 3 - ）后加到方程（II ）上去，把方程（I ）乘（2 4- ）后加到方程（III ）上 去，即可消去方程（II ）、（III ）中的x 1，得同解方程组 ?? ? ??=+-II -=-I =++III) (20 223)(445.0)(6 4323232321x x x x x x x 将方程（II ）乘（ 5 .03 ）后加于方程（III ），得同解方程组： ?? ? ??-=-II -=-I =++III) (42)(445.0)(6432332321x x x x x x 由回代公式（3.5）得x 3 = 2，x 2 = 8，x 1 = -13。 下面考察一般形式的线性方程组的解法，为叙述问题方便，将b i 写成a i , n +1，i = 1, 2,…,n 。 ??? ?? ??=++++=++++=+++++++1,3322111 ,223232221211,11313212111n n n nn n n n n n n n n n a x a x a x a x a a x a x a x a x a a x a x a x a x a （1-1） 如果a 11 ≠ 0，将第一个方程中x 1的系数化为1，得 ) 1(1,1)1(12)1(121+=+++n n n a x a x a x 其中)0(11 ) 0()1(1a a a ij j = ， j = 1, …, n + 1（记ij ij a a =) 0(，i = 1, 2, …, n ; j = 1, 2, …, n + 1） 从其它n –1个方程中消x 1，使它变成如下形式 ?? ? ????=++=++=++++++)1(1,)1(2)1(2) 1(1 ,2)1(22)1(22) 1(1,1)1(12)1(121n n n nn n n n n n n n a x a x a a x a x a a x a x a x （1-2） 其中n i a m a a ij i ij ij ,,2)1(1) 1( =?-=，1,,3,211 )1(1 1+== n j a a m i i 由方程（1-1）到（1-2）的过程中，元素11a 起着重要的作用，特别地，把11a 称为主元素。 如果（1-2）中0) 1(22≠a ，则以) 1(22a 为主元素，又可以把方程组（1-2）化为： ?? ? ??????=++=++=+++=+++++++)2(1 ,)2(3)2(3) 3(1,3)2(33)2(33) 2(1 ,2)2(23)2(232) 1(1,1)1(12)1(121 n n n nn n n n n n n n n n n a x a x a a x a x a a x a x a x a x a x a x （1-3） 针对（1-3） 继续消元，重复同样的手段，第k 步所要加工的方程组是： ?? ?? ?? ?? ? ????=++=++=+++=+++=++++-+---+---+-----++) 1(1,)1()1() 1(1,)1()1() 1(1,1)1()1(11) 2(1 ,2)2(23)2(232) 1(1,1)1(13)1(132)1(121 k n n n k nn k k nk k n k n k nn k k kk k n k n k kn k k k k n n n n n n a x a x a a x a x a a x a x a x a x a x a x a x a x a x a x 用Gauss-Seidel迭代法解线性方程组的C语言源代码：#include ptr2->link=ptr1; ptr2=ptr1; } } j++; } if(ptr2!=NULL){ ptr2->link=NULL; ptr2=NULL; } if(Rhead!=NULL){ p1=(struct Line*)malloc(sizeof(Line)); if(p1==NULL){ printf("内存分配错误!\n"); exit(1); } p1->L=i; p1->head=Rhead; if(p2==NULL){ Lhead=p1; p2=p1; } else{ p2->next=p1; p2=p1; } } i++; Rhead=NULL; j=1; } if(p2!=NULL) p2->next=NULL; return Lhead; } struct Line *Change(struct Line*Lhead,int n){ struct Line*p1,*p2,*p3,*p; struct Row*ptr; int i=1,k,j; float max,t; if(Lhead==NULL){ printf("链表为空!\n");Gauss-Seidel迭代法求解线性方程组
迭代法解线性方程组
消元法解线性方程组
线性方程组的Guass消元法求解
MATLAB之GAUSS消元法解线性方程组
(完整版)解线性方程组的消元法及其应用
实验解线性方程组的基本迭代法实验
Gauss-Seidel迭代法求解线性方程组
C++课程设计高斯消元法求线性代数方程组的解
作业一 高斯消元法和列主元消元法
迭代法解线性方程组(C语言描述)
相关主题
文本预览