数值分析算法 MATLAB 实践 线性方程组迭代法
Jacobi迭代法
雅可比迭代法保证收敛的条件是矩阵A(Ax=b)为严格的行对角占优矩阵，对于每一行，对角线上的元素之绝对值大于其余元素绝对值的和。需要说明的是：即使不满足此条件，雅可比法有时仍可以收敛。
%% 雅可比迭代法 [x,k,index] = Jacobimethod(A,b,ep)
% A为方程组的系数矩阵；
% b为方程组的右端项;
% ep为精度要求，缺省值为1e-5;
% it_max为最大迭代次数，缺省值为100;
% x为方程组的解;
% k为迭代次数;
% index为指标变量，index=0表示迭代失败，index=1表示收敛到指定要求
A = [10 3 1; 2 -10 3 ;1 3 10];
b = [14 -5 14 ]'; %b = [14; -5; 14 ];
eps = 0.005;
[x_0,k0_cnt,index] = Jacobimethod(A,b,eps);
disp('迭代次数:k0_cnt=')
disp(k0_cnt)
disp(['方程组的解：x_0 = '])
disp(x_0)
%% 求线性方程组的Jacobi迭代法，调用格式为[x, k] = JacobiFunc(A,b,x0,eps,it_max)
% 其中, A 为线性方程组的系数矩阵，b 为常数项，eps 为精度要求，默认为1e-6,x0迭代初始值
% it_max 为最大迭代次数，默认为1000
% x 为线性方程组的解，k迭代次数
x0=[0,0,0]';%[x1;x2;x3]列向量
it_max = 1000;eps=1e-6;
[x1, k1_cnt] = JacobiFunc(A,b,x0,eps,it_max);
disp('迭代次数:k1_cnt=');
disp(k1_cnt)
disp(['方程组的解：x1 = ']);
disp(x1)
%% [x,k]=JacobiFunmethod(A,b,x0,N,emg)
% A:线性方程组左端矩阵，b:线性方程组右端向量，x0:迭代初值
% N:迭代次数上界，若迭代次数大于n，则迭代失败， emg:精度指标
% k:迭代次数，
% x:用迭代法求得的线性方程组的近似解
x0=[0,0,0]';%[x1;x2;x3]列向量
it_max = 1000;eps=1e-6;
[x2,k2_cnt]=JacobiFunmethod(A,b,x0,it_max,eps);
disp('迭代次数:k2_cnt=');
disp(k2_cnt)
disp(['方程组的解：x2 = ']);
disp(x2)
function [x,k] = JacobiFunc(A,b,x0,eps,it_max)
% 求线性方程组的Jacobi迭代法，调用格式为[x, k] = JacobiFunc(A,b,x0,eps,it_max)
% 其中, A 为线性方程组的系数矩阵，b 为常数项，eps 为精度要求，默认为1e-6,x0迭代初始值
% it_max 为最大迭代次数，默认为200
% x 为线性方程组的解，k迭代次数
if nargin == 3
eps = 1.0e-6;
M = 200;
elseif nargin<3
disp('输入参数数目不足3个');
return
elseif nargin ==5
M = it_max;
end
D = diag(diag(A));%求A的对角矩阵
L = -tril(A,-1);%求A的下三角矩阵
U = -triu(A,1);%求A的上三角矩阵
B = D\(L+U);
f = D\b;
x = B*x0+f;
k = 1;%迭代次数
while norm(x-x0)>=eps
x0 = x;
x = B*x0+f;
k = k+1;
if(k>=M)
disp('Warning:迭代次数太多,可能不收敛!');
return;
end
end
end
function [ x,k,index]=Jacobimethod(A,b,ep,it_max)
% 求线性方程组的雅可比迭代法，其中，
% A为方程组的系数矩阵；
% b为方程组的右端项;
% ep为精度要求，缺省值为1e-5;
% it_max为最大迭代次数，缺省值为100;
% x为方程组的解;
% k为迭代次数;
% index为指标变量，index=0表示迭代失败，index=1表示收敛到指定要求,
[n,m] = size(A);nb = length(b);
%当方程组行与列的维数不相等时，停止计算，并输出出错信息。
if n ~=m
error('The rows and columns of matrix A must be equal! ');
return;
end
% 当方程组与右端项的维数不匹配时，停止计算，并输出出错信息。
if m~=nb
error ('The columns of A must be equal the length of b! ');
return;
end
if nargin<4
it_max =100;
end
if nargin<3
ep = 1e-5;
end
k=0;x = zeros (n,1);y=zeros (n,1);index=1;
while 1
for i=1 :n
y(i) =b(i) ;
for j=1:n
if j~=i
y(i) =y(i) -A(i,j)*x(j);
end
end
if abs(A(i,i))<1e-10 &&k==it_max % abs绝对值函数
index =0 ;return;
end
y(i) =y(i)/A(i,i);
end
k = k +1;
if norm(y-x,inf) <ep
break;
end
x = y;
end
function [x,k]=JacobiFunmethod(A,b,x0,N,emg)
% A:线性方程组左端矩阵，b:线性方程组右端向量，x0:迭代初值
% N:迭代次数上界，若迭代次数大于n，则迭代失败， emg:精度指标
% k:迭代次数，
% x:用迭代法求得的线性方程组的近似解
n=length(A);
x=zeros(n,1); %设置变量
X=zeros(n,1); % X
x=x0; k=0;
r=max(abs(b-A*x));
while (r>emg) % 迭代循环过程
for i=1:n
sum=0;
for j=1:n
if i~=j
sum=sum+A(i,j)*x(j);
end
end
X(i)=(b(i)-sum)/A(i,i);
end
r=max(abs(X-x));
x=X;
k=k+1;
if k>N
disp('迭代失败，返回');
return;
end
end
