赞
踩
雅克比(Jacobi)迭代法求解线性方程组
简单的讲其实就是我们平时求解的方法(最常用的方法)
以下是Jacobi的迭代过程:
线性方程组为:
将A分裂:
从而得到迭代公式为:
由于D为对角元素,从而有:
即编程语言可以写为:
推到过程为:
完整程序代码为:
function [x,error,iter]=GJJacobi_solve(A,b,epsilon,max_iter) % x , error , iter 为输出变量 A ,b epsilon,max_iter为输入变量 n=length(b); x=zeros(n,1); % 首先编写第一步代码目的 x =-(L+U)*x/D+b/D %第二部明确符号的意义 D表示对角线的元素 L表示下三角的元素 U表示上三角的元素 b表示列向量的长度 %确定误差error = norm(x-y) %输入矩阵A 以及 列向量b %表示L U D L=tril(A); %tril(A)表示取矩阵A的下三角元素, 包含了对角线元素 U=triu(A); %triu(A) 表示取矩阵A的上三角元素, 包含了对角线元素 D=diag(diag(A)); %表示取矩阵A的主对角元素 %想办法使得L,U对角线元素都为0 % L(logical(eye(size(L))))=0; %由于维度不对,从而只是等号右边为零不对 % U(logical(eye(size(U))))=0; %改正 %B = zeros(1,n); %[0,0,0] L(logical(eye(size(L))))=0; %[8;11;12] U(logical(eye(size(U))))=0; %需要看看是否输出的是对角线元素为0 disp(L) disp(U) % 先写一个循环(发现运用一个变量时需要先初始化error=1) % 设置最小误差为10e-6 error = 1 ;%初始化误差变量 iter = 0; %初始化迭代步数变量 while error>epsilon && iter<max_iter % 给出循环条件,以及函数输入变量epsilon,max_iter for i =1:max_iter y=x;%迭代前的矩阵,用来计算迭代误差 x=-D\(L+U)*x+D\b ; %从这里知道需要把L,U,D,b表示出来 x=-(L+U)*x/D+b/D 这个维度错误 error=norm(x-y); %知道迭代前的数据y,以及迭代后的数据x,设定误差为(x-y)的范数 if error>epsilon %给出判断误差是否减小 iter =iter+1; break end %if 需要 end 结束 end %for 需要 end 结束 end %while 需要 end 结束测试代码为:
A=[8 -3 2 4 11 -1 6 3 12]; b=[20;33;36]; epsilon=10e-6; max_iter=15; [x,error,iter]=GJJacobi_solve(A,b,epsilon,max_iter); disp('程序计算的精确解为:'); disp(x); disp('最大迭代次数下的误差:'); disp(error); disp('最小迭代次数:'); disp(iter+1);运行结果为:
0 0 0 4 0 0 6 3 0 0 -3 2 0 0 -1 0 0 0 程序计算的精确解为: 3.0000 2.0000 1.0000 最大迭代次数下的误差: 4.1060e-12 最小迭代次数: 14
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。