Gauss-Seidel迭代法

设线性方程组为AX=B,则Gauss-Seidel迭代法的迭代公式为

\[{x}_{m}^{k+1}=\frac{1}{{a}_{mm}}(-\sum_{i=1}^{m-1} {{a}_{mi}{x}_{i}^{k+1}-\sum_{i=m+1}^{n} {{a}_{mi}{x}_{i}^{k}}+{b}_{m}})\]

对应的矩阵表达式为

\[{x}^{\{k+1\}}={B}_{0}{x}^{k}+f\]

其中,B0为Gauss-Seidel迭代矩阵,B0=(D+L) 1U,f=(D+L)1b,D为对角矩阵,L和U分别为严格下三角矩阵和严格上三角矩阵。据此,可以编写实现Gauss-Seidel迭代法的M文件Gauss.m,如下所示。

code.matlab
function s=Gauss(a,b,x0,eps)
  % 用Gauss-Seidel迭代法解线性方程组
  % a为系数矩阵,b为方程组ax=b的右端项,x0为初值
  if nargin==3
      eps=1.0e-6;
  elseif nargin<3
      error
      return
  end
  D=diag(diag(a));	%求对角矩阵
  L=tril(a,-1);		%求严格下三角矩阵
  U=triu(a,1); 		%求严格上三角矩阵
  C=inv(D+L);
  B=-C*U;
  f=C*b;
  s=B*x0+f;
  while norm(s-x0)>=eps
      x0=s;
      s=B*x0+f;
  end
  return

【例4】用上面编写的Gauss函数求解下列方程组

\[\left\{ \begin{matrix} 10{x}_{1}-2{x}_{2}-{x}_{3}=6 \\ -2{x}_{1}+2{x}_{2}-{x}_{3}=10 \\ -{x}_{1}-2{x}_{2}+5{x}_{3}=10 \end{matrix} \right.\]

下面用MATLAB实现,在命令窗口中输入

code.matlab
>> a=[10 -2 -1;-2 2 -1;-1 -2 5];
>> b=[6 10 10]';
>> x0=[0 0 0]';
>> gauss(a,b,x0)
ans =
    4.0000
   13.0000
    8.0000