2013230073
逸群
昆明理工大学工科研究生
《数值分析》上机实验报告
姓 名:王逸群
学 号:2013230073
专 业:材料加工工程
学 院:材料科学与工程学院
授课教师:李玉兰
1
一、课题名称 矩阵求逆与矩阵行列式
二、班级、姓名、学号
材料加工工程 王逸群 2013230073
三、目的和意义
方法的理论意义和实用价值。
1、通过该课题的实验,主要掌握Gauss-Jordan消去法求非奇异矩阵的逆阵的程序设计方法;
2、体会Gauss-Jordan消去法求解线性方程组的关键步骤;
3、提高科学计算和编程的能力。
问题提出
应用列主元Gauss-Jordan 消去法求满秩矩阵的逆矩阵,并计算的行列式的值,如下列矩阵之一:
(1)
(2)
(3)
(4)
要求
1、分析列主元Gauss-Jordan消去法的计算公式;
2、确定选主元,换行,计算行(约化非主元行)和交换列序等四个子程序;
3、*应用结构程序设计编出计算n阶非奇异方阵的通用程序*;
4、*计算考核题
*
四、计算公式
列主元素法就是在待消元的所在列中选取主元,经方程的行交换,置主元素于对角线位置后进行消元的方法。即初等行变换方法。
如第一题
(1)
所以
流程应该是:选主元→换行→消去,将主元置于对角线→左边化为单位矩阵
用高斯消去法解方程组时,小主元可能导致计算失败,因为用绝对值很小的数作除数,乘数很大,引起约化中间结果数量级严重增长,再舍入就使得计算结果不可靠了,故避免采用绝对值很小的主元素。以便减少计算过程中舍入误差对计算解的影响。
本题目计算中因为1和3的绝对值相差不大,因此主元变换不受影响。但是下面程序设计中设计了首先选择大主元。
五、结构程序设计
Matlab程序设计:
%矩阵赋值%
>> A=[1 2 0 0;3 4 0 0;0 0 4 1;0 0 3 2]
A =
1 2 0 0
3 4 0 0
0 0 4 1
0 0 3 2
>> E=eye(4)
E =
1 0 0 0
0 1 0 0
0 0 1 0
0 0 0 1
>> A1=[A E]
A1 =
1 2 0 0 1 0 0 0
3 4 0 0 0 1 0 0
0 0 4 1 0 0 1 0
0 0 3 2 0 0 0 1
%比较第一列,选最大值%
>> v=max(A(1,1),A(2,1))
v =
3
%交换第一行和最大值的行%
>> t1=A1(1,:)
t1 =
1 2 0 0 1 0 0 0
>> A1(1,:)=A1(2,:)
A1 =
3 4 0 0 0 1 0 0
3 4 0 0 0 1 0 0
0 0 4 1 0 0 1 0
0 0 3 2 0 0 0 1
>> A1(2,:)=t1
A1 =
3 4 0 0 0 1 0 0
1 2 0 0 1 0 0 0
0 0 4 1 0 0 1 0
0 0 3 2 0 0 0 1
%通过行变换消去非主元%
>> t=A1(2,:)-A1(1,:)/v
t =
0 0.6667 0 0 1.0000 -0.3333 0 0
>> A1(2,:)=t
A1 =
3.0000 4.0000 0 0 0 1.0000 0 0
0 0.6667 0 0 1.0000 -0.3333 0 0
0 0 4.0000 1.0000 0 0 1.0000 0
0 0 3.0000 2.0000 0 0 0 1.0000
%继续选择主元,消去非主元%
>> v=A1(3,3)/A1(4,3)
v =
1.3333
>> t=A1(4,:)-A1(3,:)/v
t =
0 0 0 1.2500 0 0 -0.7500 1.0000
>> A1(4,:)=t
A1 =
3.0000 4.0000 0 0 0 1.0000 0 0
0 0.6667 0 0 1.0000 -0.3333 0 0
0 0 4.0000 1.0000 0 0 1.0000 0
0 0 0 1.2500 0 0 -0.7500 1.0000
%此时可计算行列式的值%
>> AA=-A1(1,1)*A1(2,2)*A1(3,3)*A1(4,4)
AA =
-10
%与matlab自带程序结果计算相同%
>> det(A)
ans =
-10
%继续消去上三角内非主元元素%
>> v=A1(3,4)/A1(4,4)
v =
0.8000
>> t=A1(3,:)-A1(4,:)*v
t =
0 0 4.0000 0 0 0 1.6000 -0.8000
>> A1(3,:)=t
A1 =
3.0000 4.0000 0 0 0 1.0000 0 0
0 0.6667 0 0 1.0000 -0.3333 0 0
0 0 4.0000 0 0 0 1.6000 -0.8000
0 0 0 1.2500 0 0 -0.7500 1.0000
>> v=A1(1,2)/A1(2,2)
v =
6.0000
>> t=A1(1,:)-A1(2,:)*v
t =
3.0000 0 0 0 -6.0000 3.0000 0 0
>> A1(1,:)=t
A1 =
3.0000 0 0 0 -6.0000 3.0000 0 0
0 0.6667 0 0 1.0000 -0.3333 0 0
0 0 4.0000 0 0 0 1.6000 -0.8000
0 0 0 1.2500 0 0 -0.7500 1.0000
%将对角矩阵化为单位矩阵%
>> A1(1,:)=A1(1,:)/A1(1,1)
A1 =
1.0000 0 0 0 -2.0000 1.0000 0 0
0 0.6667 0 0 1.0000 -0.3333 0 0
0 0 4.0000 0 0 0 1.6000 -0.8000
0 0 0 1.2500 0 0 -0.7500 1.0000
>> A1(2,:)=A1(2,:)/A1(2,2)
A1 =
1.0000 0 0 0 -2.0000 1.0000 0 0
0 1.0000 0 0 1.5000 -0.5000 0 0
0 0 4.0000 0 0 0 1.6000 -0.8000
0 0 0 1.2500 0 0 -0.7500 1.0000
>> A1(3,:)=A1(3,:)/A1(3,3)
A1 =
1.0000 0 0 0 -2.0000 1.0000 0 0
0 1.0000 0 0 1.5000 -0.5000 0 0
0 0 1.0000 0 0 0 0.4000 -0.2000
0 0 0 1.2500 0 0 -0.7500 1.0000
>> A1(4,:)=A1(4,:)/A1(4,4)
A1 =
1.0000 0 0 0 -2.0000 1.0000 0 0
0 1.0000 0 0 1.5000 -0.5000 0 0
0 0 1.0000 0 0 0 0.4000 -0.2000
0 0 0 1.0000 0 0 -0.6000 0.8000
%求出逆矩阵%
>> [A1(:,5) A1(:,6) A1(:,7) A1(:,8)]
ans =
-2.0000 1.0000 0 0
1.5000 -0.5000 0 0
0 0 0.4000 -0.2000
0 0 -0.6000 0.8000
%验证计算结果%
>> A=[1 2 0 0;3 4 0 0;0 0 4 1;0 0 3 2]
A =
1 2 0 0
3 4 0 0
0 0 4 1
0 0 3 2
>> A^-1
ans =
-2.0000 1.0000 0 0
1.5000 -0.5000 0 0
0 0 0.4000 -0.2000
0 0 -0.6000 0.8000
四’、计算公式
第二题
(2)
·······
五’、结构程序设计
%矩阵赋值%
>> A=[1 0 -1 2 1;3 2 -3 5 -3;2 2 1 4 -2;0 4 3 3 1;1 0 8 -11 4]
A =
1 0 -1 2 1
3 2 -3 5 -3
2 2 1 4 -2
0 4 3 3 1
1 0 8 -11 4
>> E=eye(5)
E =
1 0 0 0 0
0 1 0 0 0
0 0 1 0 0
0 0 0 1 0
0 0 0 0 1
>> A1=[A E]
A1 =
1 0 -1 2 1 1 0 0 0 0
3 2 -3 5 -3 0 1 0 0 0
2 2 1 4 -2 0 0 1 0 0
0 4 3 3 1 0 0 0 1 0
1 0 8 -11 4 0 0 0 0 1
%交换第一行和最大值的行%
>> t=A1(1,:)
t =
1 0 -1 2 1 1 0 0 0 0
>> A1(1,:)=A1(2,:)
A1 =
3 2 -3 5 -3 0 1 0 0 0
3 2 -3 5 -3 0 1 0 0 0
2 2 1 4 -2 0 0 1 0 0
0 4 3 3 1 0 0 0 1 0
1 0 8 -11 4 0 0 0 0 1
>> A1(2,:)=t
A1 =
3 2 -3 5 -3 0 1 0 0 0
1 0 -1 2 1 1 0 0 0 0
2 2 1 4 -2 0 0 1 0 0
0 4 3 3 1 0 0 0 1 0
1 0 8 -11 4 0 0 0 0 1
%通过行变换消去非主元%
>> v=A1(2,1)/A1(1,1)
v =
0.3333
>> A1(2,:)=A1(2,:)-A1(1,:)*v
A1 =
Columns 1 through 9
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0
2.0000 2.0000 1.0000 4.0000 -2.0000 0 0 1.0000 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000
1.0000 0 8.0000 -11.0000 4.0000 0 0 0 0
Column 10
0
0
0
0
1.0000
>> v=A1(3,1)/A1(1,1)
v =
0.6667
>> A1(3,:)=A1(3,:)-A1(1,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0 0
0 0.6667 3.0000 0.6667 0 0 -0.6667 1.0000 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
1.0000 0 8.0000 -11.0000 4.0000 0 0 0 0 1.0000
>> v=A1(4,1)/A1(1,1)
v =
0
>> v=A1(5,1)/A1(1,1)
v =
0.3333
>> A1(5,:)=A1(5,:)-A1(1,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0 0
0 0.6667 3.0000 0.6667 0 0 -0.6667 1.0000 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 -0.6667 9.0000 -12.6667 5.0000 0 -0.3333 0 0 1.0000
>> t=A1(4,:)
t =
0 4 3 3 1 0 0 0 1 0
>> A1(4,:)=A1(2,:)
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0 0
0 0.6667 3.0000 0.6667 0 0 -0.6667 1.0000 0 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0 0
0 -0.6667 9.0000 -12.6667 5.0000 0 -0.3333 0 0 1.0000
>> A1(2,:)=t
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0.6667 3.0000 0.6667 0 0 -0.6667 1.0000 0 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0 0
0 -0.6667 9.0000 -12.6667 5.0000 0 -0.3333 0 0 1.0000
>> v=A1(3,2)/A1(2,2)
v =
0.1667
>> A1(3,:)=A1(3,:)-A1(2,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
0 -0.6667 0 0.3333 2.0000 1.0000 -0.3333 0 0 0
0 -0.6667 9.0000 -12.6667 5.0000 0 -0.3333 0 0 1.0000
>> v=A1(4,2)/A1(2,2)
v =
-0.1667
>> A1(4,:)=A1(4,:)-A1(2,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
0 0 0.5000 0.8333 2.1667 1.0000 -0.3333 0 0.1667 0
0 -0.6667 9.0000 -12.6667 5.0000 0 -0.3333 0 0 1.0000
>> v=A1(5,2)/A1(2,2)
v =
-0.1667
>> A1(5,:)=A1(5,:)-A1(2,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
0 0 0.5000 0.8333 2.1667 1.0000 -0.3333 0 0.1667 0
0 0 9.5000 -12.1667 5.1667 0 -0.3333 0 0.1667 1.0000
>> A1(5,:)=A1(3,:)
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
0 0 0.5000 0.8333 2.1667 1.0000 -0.3333 0 0.1667 0
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
>> A1(3,:)=t
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 9.5000 -12.1667 5.1667 0 -0.3333 0 0.1667 1.0000
0 0 0.5000 0.8333 2.1667 1.0000 -0.3333 0 0.1667 0
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
>> v=A1(4,3)/A1(3,3)
v =
0.0526
>> A1(4,:)=A1(4,:)-A1(3,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 9.5000 -12.1667 5.1667 0 -0.3333 0 0.1667 1.0000
0 0 0 1.4737 1.8947 1.0000 -0.3158 0 0.1579 -0.0526
0 0 2.5000 0.1667 -0.1667 0 -0.6667 1.0000 -0.1667 0
>> v=A1(5,3)/A1(3,3)
v =
0.2632
>> A1(5,:)=A1(5,:)-A1(3,:)*v
A1 =
3.0000 2.0000 -3.0000 5.0000 -3.0000 0 1.0000 0 0 0
0 4.0000 3.0000 3.0000 1.0000 0 0 0 1.0000 0
0 0 9.5000 -12.1667 5.1667 0 -0.3333 0 0.1667 1.0000
0 0 0 1.4737 1.8947 1.0000 -0.3158 0 0.1579 -0.0526
0 0 0 3.3684 -1.5263 0 -0.5789 1.0000 -0.2105 -0.2632
>> t=A1(5,:)
t =