您好,欢迎访问三七文档
当前位置:首页 > 电子/通信 > 综合/其它 > Jacobi迭代和gs迭代的程序代码
1234512123523443253283103381050102515015450530010JacobiGaussSeideliiiiiiiiiiiiiiiiiiJacobiGaussSeidel用和方法解线性方程组。有一给定的电路图,其中电流,,,,满足方程组试用和迭代法解之,使误差小于。Jacobi迭代的程序代码:%Jacobi迭代法解方程组function[s,i,y]=jacobi(a,b,r)s=1;i=0;y=0;[m,n]=size(a);%矩阵a的维数ifm~=n|rank(a)~=n%矩阵不满秩,不可解,输出错误信息s=0;elsed=diag(a);%d为矩阵a对角线上的元素%若矩阵a的对角线元素存在0元,则将该元素所在列的任一不为0的元素所在的行加到该行z=find(d==0);%搜索0元l=size(z,1);ifl%存在0元fori=1:lj=find(abs(a(:,z(i))),1,'first');a(z(i),:)=a(z(i),:)+a(j,:);%两行相加b(z(i))=b(z(i))+b(j);endendd=1./diag(a);g=eye(n)-diag(d)*a;%迭代矩阵y=zeros(n,1);%初始向量b=diag(d)*b;fori=1:10^5%进行迭代x=y;y=g*x+b;ifmax(abs(y-x))rbreak%满足精度要求,退出迭代endendendGauss-Seidel迭代的程序代码:%Gauss-Seidel迭代法解方程组function[s,i,y]=gaussseidel(a,b,r)s=1;i=0;y=0;[m,n]=size(a);%矩阵a的维数ifm~=n|rank(a)~=n%矩阵不满秩,不可解,输出错误信息s=0;elsed=diag(a);%d为矩阵a对角线上的元素%若矩阵a的对角线元素存在0元,则将该元素所在列的任一不为0的元素所在的行加到该行z=find(d==0);%搜索0元l=size(z,1);ifl%存在0元fori=1:lj=find(abs(a(:,z(i))),1,'first');a(z(i),:)=a(z(i),:)+a(j,:);%两行相加b(z(i))=b(z(i))+b(j);endendd=diag(a);tl=tril(a);f=eye(n);e=[tl,f];tu=diag(d)-triu(a);forj=1:n-1e(j,:)=e(j,:)/e(j,j);e(j+1:n,j:end)=e(j+1:n,j:end)-e(j+1:n,j)*e(j,j:end);%进行消去ende(n,:)=e(n,:)/e(n,n);e=e(:,n+1:2*n);g=e*tu;%迭代矩阵b=e*b;%迭代向量y=zeros(n,1);%初始向量fori=1:10^5%进行迭代x=y;y=g*x+b;ifmax(abs(y-x))rbreak%满足精度要求,退出迭代endendend主程序为:%计算方法上机第四题clear;clc;a=[28-3000;-338-100-5;0-1025-150;00-15450;0-50030];%输入矩阵ab=[100000]';%输入矩阵br=10^-3;%指定迭代误差%采用jacobi迭代进行计算[s,i,x]=jacobi(a,b,r);if~sdisp('Error!Jacobiiteratingmethodcannotgo!');elseifi==10^5disp('Error!ItisdivergentwhenusingJacobiiteratingmethod!');else%输出结果disp(['After',num2str(i),'timesJacobiiteration,thesolutionconvergeon:']);disp(x);end%采用Gauss-Seidel迭代进行计算[s,i,x]=gaussseidel(a,b,r);if~sdisp('Error!Gauss-seideliteratingmethodcannotgo!');elseifi==10^5disp('Error!ItisdivergentwhenusingGauss-seideliteratingmethod!');else%输出结果disp(['After',num2str(i),'timesGauss-seideliteration,thesolutionconvergeon:']);disp(x);end
本文标题:Jacobi迭代和gs迭代的程序代码
链接地址:https://www.777doc.com/doc-5481713 .html