计算方法实验报告-线性方程组的数值解法

重庆大学
学生实验报告实验课程名称计算方法
开课实验室DS1421
学院年级专业
学生姓名学号
开课时间至学年第学期
1．实验目的
（1）高斯列主元消去法求解线性方程组的过程
（2）熟悉用迭代法求解线性方程组的过程
（3）设计出相应的算法，编制相应的函数子程序
2．实验内容
分别用高斯列主元消去法 ,Jacobi 迭代法,Gauss--Saidel 迭代法,超松弛迭代法求解线性方程组
⎥⎥⎥⎥⎦
⎤⎢⎢⎢⎢⎣⎡-=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡------725101391444321131243301024321x x x x 3.实验过程
解：（1）高斯列主元消去法
编制高斯列主元消去法的M 文件程序如下：
%高斯列主元消元法求解线性方程组Ax=b
%A 为输入矩阵系数,b 为方程组右端系数
%方程组的解保存在x 变量中
format long;%设置为长格式显示,显示15位小数
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
b=[10,5,-2,7]'
[m,n]=size(A);
%先检查系数正确性
if m~=n
error('矩阵A 的行数和列数必须相同');
return;
end
if m~=size(b)
error('b 的大小必须和A 的行数或A 的列数相同');
return;
end
%再检查方程是否存在唯一解
if rank(A)~=rank([A,b])
error('A 矩阵的秩和增广矩阵的秩不相同,方程不存在唯一解');
return;
end
c=n+1;
A(:,c)=b; %(增广)
for k=1:n-1
[r,m]=max(abs(A(k:n,k))); %选主元
m=m+k-1; %修正操作行的值
if(A(m,k)~=0)
if(m~=k)
A([k m],:)=A([m k],:); %换行
end
A(k+1:n, k:c)=A(k+1:n, k:c)-(A(k+1:n,k)/ A(k,k))*A(k, k:c); %消去 end
end
x=zeros(length(b),1); %回代求解
x(n)=A(n,c)/A(n,n);
for k=n-1:-1:1
x(k)=(A(k,c)-A(k,k+1:n)*x(k+1:n))/A(k,k);
end
disp('X=');
disp(x);
format short;%设置为默认格式显示,显示5位
运行，结果如下所示：
（2）Jacobi迭代法
编制迭代计算的M文件程序如下:
%Jacobi迭代法求解
% A为方程组的增广矩阵
clc;
A=[2,10,0,-3,10;-3,-4,-12,13,5;1,2,3,-4,-2;4,14,9,-13,7] MAXTIME=50;%最多进行50次迭代
eps=1e-5;%迭代误差
[n,m]=size(A);
x=zeros(n,1);%迭代初值
y=zeros(n,1);
k=0;
%进入迭代计算
disp('迭代过程X的值情况如下:')
disp('X=');
while 1
disp(x');
for i=1:1:n
s=0.0;
for j=1:1:n
if j~=i
s=s+A(i,j)*x(j);
end
y(i)=(A(i,n+1)-s)/A(i,i);
end
end
for i=1:1:n
maxeps=max(0,abs(x(i)-y(i))); %检查是否满足迭代精度要求 end
if maxeps<=eps%小于迭代精度退出迭代
for i=1:1:n
x(i)=y(i);%将结果赋给x
end
return;
end
for i=1:1:n%若不满足迭代精度要求继续进行迭代
x(i)=y(i);
y(i)=0.0;
end
k=k+1;
if k>MAXTIME%超过最大迭代次数退出
error('超过最大迭代次数,退出');
return;
end
end
运行该程序结果如下:
（3）Gauss--Saidel迭代法
编制求解程序Gauss_Seidel.m如下:
%Gauss_Seidel.m
% A为方程组的增广矩阵
clc;
format long;
A=[2,10,0,-3,10;-3,-4,-12,13,5;1,2,3,-4,-2;4,14,9,-13,7]
[n,m]=size(A);
%最多进行50次迭代
Maxtime=50;
%控制误差
Eps=10E-5;
%初始迭代值
x=zeros(1,n);
disp('x=');
%迭代次数小于最大迭代次数,进入迭代
for k=1:Maxtime
disp(x);
for i=1:n
s=0.0;
for j=1:n
if i~=j
s=s+A(i,j)*x(j);%计算和
end
end
x(i)=(A(i,n+1)-s)/A(i,i);%求出此时迭代的值
end
%因为方程的精确解为整数,所以这里将迭代结果向整数靠近的误差作为判断迭代是否停止的条件
if sum((x-floor(x)).^2)<Eps
break;
end;
end;
X=x;
disp('迭代结果:');
X
format short;
运行结果如下所示：
（4）超松弛迭代法
编写函数M文件如下：%SOR法求解
%w=1.45
%方程组系数矩阵
clc;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
%方程组右端系数
b=[10,5,-2,7]'
w=1.45;
%最大迭代次数
Maxtime=100;
%精度要求
Eps=1E-5;
%以15位小数显示
format long;
n=length(A);
k=0;
%初始迭代值
x=ones(n,1);
y=x;
disp('迭代过程:');
disp('x=');
while 1
y=x;
disp(x');
%计算过程
for i=1:n
s=b(i);
for j=1:n
if j~=i
s=s-A(i,j)*x(j);
end
end
if abs(A(i,i))<1E-10 | k>=Maxtime
error('已达最大迭代次数或矩阵系数近似为0,无法进行迭代'); return;
end
s=s/A(i,i);
x(i)=(1-w)*x(i)+w*s;
end
if norm(y-x,inf)<Eps%达到精度要求退出计算
break;
end
k=k+1;
end
disp('最后迭代结果:');
%最后的结果
X=x'
%设为默认显示格式
format short; 结果如下：
4．实验环境及实验文件存档名
实验环境：Matlab7.0
文件存档名：Gauss.m，Jacobi.m，Gauss_Seidel.m，SOR.m
5.实验结果及分析
=1.0000
=2.0000
=3.0000
=4.0000
经过验证，高斯列主元消结果正确。

对2，3，4方法不收敛。

错误分析
迭代法收敛的充要条件ρ（Bj）<1
ρ（Bj）= max(abs(eig(Bj)))=3.9412>1
因此迭代法不收敛。

6.总结与体会
1.通过该实验弄懂了高斯列主元消去法求解线性方程组的过程，练习并掌握了高斯列主元消去法求解线性方程组的方法。

2. 熟悉了用迭代法求解线性方程组的过程，学会了Jacobi迭代法,Gauss--Saidel迭代法,超松弛迭代法求解线性方程组的方法。

教师签名
年月日。