数值分析报告线性方程组
0 0 0 0 0 0 1 3 2 0
0 0 0 0 0 0 0 1 3 2
0 0 0 0 0 0 0 0 1 3];
[n,n]=size(A);
L=zeros(n,n);
U=zeros(n,n);
fori=1:n
L(i,i)=1;
end
fork=1:n
forj=k:n
U(k,j)=A(k,j)-sum(L(k,1:k-1).*U(1:k-1,j)');
0 0 0 1 1 0 0 0 0 0
0 0 0 0 1 1 0 0 0 0
0 0 0 0 0 1 1 0 0 0
0 0 0 0 0 0 1 1 0 0
0 0 0 0 0 0 0 1 1 0
0 0 0 0 0 0 0 0 1 1
2、求解下面代数方程组的解
x =
1
2
3
4
5
6
7
8
9
10
3、jacobi迭代求下面方程组的解
数值分析实验报告
班级
姓名
学号
课程名称
数值分析
时间
地点
一、实验名称:
实验四、线性代数方程组直接解法、迭代法
二、实验目的:
1.高斯消去法的原理和计算步骤,理解迭代法基本原理及步骤;
2.了解顺序消去法能够实现的条件;
3.掌握矩阵的三角分解(即LU分解),平方分解的原理与步骤;
4. 利用LU分解、平方分解求解线性代数方程组.
5.利用Jacobi,Gauss-seidel迭代求解方程组
三、实验内容及要求:
1.编写LU分解matlab通用程序,并用该程序求下面矩阵的LU分解
试写出 , 阵。
2.编写平方分解matlab通用程序,并用该程序求解下面代数方程组的解:
其中
3.编写Jacobi、Gauss-seidel迭代法求解线性代数方程组通用程序,并利用Jacobi、Gauss-seidel迭代求下面方程组的解,精度控制为 (误差最大范数).
End
以Fjacobi.m保存
a=[28 1 3;2 8 1; 2 -3 15];
b=[32 11 14]';
x0=[0 0 0]';
[x,k]=Fjacobi(a,b,x0,1e-3)
以jacobi保存
再运行程序jacobi
max1=300;%默认迭代300次
D=diag(diag(A));%求a的对角矩阵
end
fori=k+1:n
L(i,k)=(A(i,k)-sum(L(i,1:k-1).*U(1:k-1,k)'))/U(k,k);
end
end
2、a=[1 2 0 0 0 0 0 0 0 0
2 5 2 0 0 0 0 0 0 0
0 2 5 2 0 0 0 0 0 0
0 0 2 5 2 0 0 0 0 0
x =
1.0e+005 *
-0.3211
-0.3203
1.1929
k =
300
Gauss-seidel迭代求下面方程组的解
x =
14.7758
-19.5433
-5.8788
k =
300
六、讨论
直接法可以求得精确解是指就计算公式而言保证得到精确解,但计算机计算过程中的舍入误差是不可避免的,这种误差对解的精度影响会不会太大,也就是计算的稳定性,是要考虑的问题。直接法普遍适用,但要求计算机有较大的存储量,迭代法要求的存储量较小,但必须在收敛性得以保证的情况下才能使用。对于迭代法,其收敛性则是要考虑的问题。
U=-tril(A,1);%求a的上角
B=D\(L+U);%用左除比求矩阵的逆
f=D\b;
x=B*x0+f;
k=1;%迭代次数
while norm (x-x0)>=tol
x0=x;
x=B*x0+f;
k=k+1;
if(k>=maxl)
disp('迭代超过300次,方程可能不收敛');
return;
End%[k x']如果想显示每一步的迭代结果
return;
end
%[k x']如果想显示每一步的迭代结果
End
以Fgseid.m保存
a=[28 1 3;2 8 1; 2 -3 15];
b=[32 11 14]';
x0=[0 0 0]';
[x,k]=Fgseid(a,b,x0,1e-3)
以gseid.m保存
再运行程序gseid
八、成绩评定
教师签字:年月日
L=-tril(A,-1);%求a的下三角
U=-tril(A,1);%求a的上三角
G=(D-L)\U;%用左除比求矩阵的逆
f=(D-L)\b;
x=G*x0+f;
k=1;%迭代次数
while norm (x-x0)>=tol
x0=x;
x=G*x0+f;
k=k+1;
if(k>=max1)
disp('迭代超过300次,方程可能不收敛');
七、程序附录
1、function[L,U]=myLU(A)
A=[1 2 0 0 0 0 0 0 0 0
1 3 2 0 0 0 0 0 0 0
0 1 3 2 0 0 0 0 0 0
0 0 1 3 2 0 0 0 0 0
0 0 0 1 3 2 0 0 0 0
0 0 0 0 1 3 2 0 0 0
0 0 0 0 0 1 3 2 0 0
[l,u]=lu(a);
x=u\(l\b)
以LUfj.m保存文件
在执行LUfj
3、function [x,k]=Fjacobi(A,b,x0,tol)
%tol是输入误差容限,x0是初值
maxl=300;%默认迭代300次
D=diag(diag(A));%求a的对角矩阵
L=-tril(A,-1);%求a的下三角
四、直接法的实验步骤:
步骤1:打开matlab,存储行列式值,并选主元
步骤2:对矩阵进行消元,再回代
步骤3:编写函数代码,保存并运行
…
五、直接法的实验结果:
1、矩阵分解结果
ans =
1 0 0 0 0 0 0 0 0 0
1 1 0 0 0 0 0 0 0 0
0 1 1 0 0 0 0 0 0 0
0 0 1 1 0 0 0 0 0 0
0 0 0 2 5 2 0 0 0 0
0 0 0 0 2 5 2 0 0 0
0 0 0 0 0 2 5 2 0 0
0 0 0 0 0 0 2 5 2 0
0 0 0 0 0 0 0 2 5 2
0 0 0 0 0 0 0 0 2 5];
b=[5 18 27 36 45 54 63 72 81 68]';