赞
踩
用Jacobi迭代计算如下方程组
在求解线性方程组时, 程序中输入的是线性方程组的增广矩阵,并设定送
代结果的精度要求, 在开始迭代之前则可以进行以下判断:
function函数:
function solution = Jacobi(Ab,epsilon) %% Ab为输入的增广矩阵 %% epsilon为输入的精度要求 %% 输入参数检查 if nargin ==1 disp('请输入精度要求epsilon') return end row = size(Ab,1); col = size(Ab,2); if ndims(Ab)~= 2 | col - row ~= 1 disp('矩阵的大小有误,不能采用Jacobi迭代法') return end A = Ab(:,1:row); b = Ab(:,col); ddet = abs(det(A)); ddiag = abs(det(diag(diag(A)))); if ddet < eps | ddiag < eps disp('该方程的系数矩阵行列式不为0,方程组无解或有无穷解,或系数矩阵的对角线有零元,不能采用Jacobi迭代法') return end %%提取上下三角矩阵及对角矩阵 U = -triu(A,1); %提取对角元素为0的上三角矩阵 L = -tril(A,-1); %提取对角元素为0的下三角矩阵 Dinv = diag(diag(A).^(-1)); B = Dinv * (L + U); f = Dinv * b; if max(abs(eig(B))) > 1 disp('迭代法不收敛!!!') return end %% 迭代过程 error = 10; n = 0; start = zeros(row,1); xk = start; xknext = B * xk +f; error = norm(xk - xknext); while error > epsilon xk = xknext; xknext = B * xk + f; error = norm(xk - xknext); n = n+1; end fprintf('Jacobi迭代次数:%d',n) solution = xknext;
主函数:
Ab = [10 -2 -1 3;-2 10 -1 15;-1 -2 5 10]
Jacobi(Ab,1e-8)
**注意:**function函数和主函数是在同一个文件夹路径下的两个文件!
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。