1、close(10)2 理论简述本文采用了三种方式进行matlab仿真gauss消元法、Jacob迭代、SOR超松弛迭代,经过实验可知道Jacob迭代是失败的,并不收敛。2.1 gauss消元法对于线性方程组Ax=b,经过n-1次行初等变换,将B=(A,b)变成上三角阵。对上式进行回带,从而化简成为单位I阵,从而得到解形式如下所示。其中,。该方法对线性方法是最原始的最准确的代数解法。2.2 Jacobi迭代法如果假设如下两系数矩阵,那么可以得到Jacobi迭代公式式子中A为方程系数矩阵,b为右端常数项。2.3 SOR迭代公式为可选择的松弛因子,那么得到如下迭代公式。在松弛迭代公式里有向前和向后的
2、超松弛迭代方法,这里就不详细阐述。3 仿真实验3.1 数据提取程序通过matlab程序提取txt文档里将其保存在A,b中,从而提取Ax=b的方程,程序如下所示。% 初始化clear all;clc;% 数据提取A=textread();% A=textread(test.txtN=A(1,1);for i=1: N for j=1:N+1 B(i,j)=A(i-1)*(N+1)+j+1,1); endend% 清除不必要数据节约内存C=B(:,1:N);f=B(:,N+1);clear A B i j;得到的矩阵C和f3.2 Gauss消元法对于一个3x3的矩阵Ax=b进行验证计算,首先我们已
3、知该方程的矩阵A,b,以及解x;理论解: 用gauss解出结果如所示图 1 方程解分布图通过计算得到解如图 2所示,程序如gaus_program.m所示。图 2 线性Ax=b对应解3.3 Jacobi迭代法取eps0.001 X=Sw*Y+kw; eps_temp=sum(abs(X-Y) eps(kk)=eps_temp; kk=kk+1; Y=X;% 画出解的图,便于点数多后查看figure(1);plot(X);grid on;xlabel(X(n)中的nfigure(2);plot(eps);ylabel(X的数值解w=1.5;P=inv(D-w*L);Gw=P*(1-w)*D+w*
4、U;fw=w*P*b;% X的初值 X=Gw*Y+fw; eps(kk)=sum(abs(X-Y)% % 清除数据% clear all;% clc;% % 读取数据% % A=textread(% N=A(1,1);% for i=1:% for j=1:% B(i,j)=A(i-1)*(N+1)+j+1,1);% end% %提取系数矩阵和常数矩阵Cx=f% C=B(:% f=B(:% clear A B i j;% A=C;% B=A f;% b=f;% 系数矩阵Bj=inv(D)*(L+U);fj=inv(D)*b;Y=0.01*ones(N,1); X=Bj*Y+fj; eps=su
5、m(abs(X-Y)% 高斯解法计算clear C f;n=length(b);RA=rank(A);RB=rank(B);judge=RB-RA;if judge0, disp(因为AB秩不同,所以此方程组无解. returnif RA=RB if RA=n因为AB秩相同且为N,所以此方程组有唯一解. % RA=n X=zeros(n,1); C=zeros(1,n+1); for p= 1:n-1 for k=p+1:n m= B(k,p)/ B(p,p); B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1); b=B(1:n,n+1); A=B(1:n,1:n); X(n)=b(n)/A(n,n); for q=n-1:-1:1 X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)/A(q,q); else因为R(A)=R(B)n,所以此方程组有无穷多解.
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1