gauss.m

《MATLAB实用指南》系列丛书DE源代码

function s=Gauss()
  % Gauss-Seidel迭代法解线性方程组
  % a为系数矩阵，b为方程组ax=b的右端项，x0为初值
  global a b x0 eps;
  b=b';
  x0=x0';
  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
  s=s';
  return