Gauss-Seidel迭代法
求解线性方程组 (MATLAB)
(数值分析作业自编)
代码如下:
*将代码分块,很容易理解。
左端是输出项,右端输入自定义项
function [R,k,x] = Gauss_Seidel(A,b,error)
% Gauss-Seidel迭代
% R为迭代矩阵B的谱半径;k为迭代步数;x为解向量
% A为系数矩阵;b为常数项;error为指定相邻两次迭代精度
%% 初始化
n=length(b);
x=zeros(n,1);x0=zeros(n,1);
k=0;
%%循环
while (1)
% 用谱半径判断收敛性 不满足敛散性退出程序并报错
B=(diag(diag(A))-tril(-A,-1))\triu(-A,1);
R=max(sqrt(eig(B'*B)));
if R>=1
disp('错误,迭代不收敛');
break;
end
%x0上一次迭代的解向量 算法主体
for i=1:n
x(i)=(b(i)-A(i,1:i-1)*x(1:i-1)-A(i,i+1:n)*x0(i+1:n))/A(i,i);
end
if norm(x-x0,Inf)<error
break;
end
x0=x;
k=k+1;
end
end
输入
A=[5 2 1; -1 4 2; 2 -3 10];
b=[-12 20 3]';
error=1e-4;
调用
[x,k,R]= Gauss_Seidel(A,b,error)
运行结果
>> test_Gauss_Seidel
x =
0.6512
k =
7
R =
-4.0000
3.0000
2.0000
>>