用Matlab学习线性代数线性方程组与矩阵代数.docx

用Matlab学习线性代数线性方程组与矩阵代数.docx

《用Matlab学习线性代数线性方程组与矩阵代数.docx》由会员分享，可在线阅读，更多相关《用Matlab学习线性代数线性方程组与矩阵代数.docx（50页珍藏版）》请在冰豆网上搜索。

用Matlab学习线性代数线性方程组与矩阵代数

用Matlab学习线性代数

线性方程组与矩阵代数

实验目的：

熟悉线性方程组的解法和矩阵的基本运算及性质验证。

Matlab命令：

本练习中用到的Matlab命令有：

inv，floor，rand，tic，toc，rref，abs，max，round，sum，eye，triu，ones，zeros。

本练习引入的运算有：

+，-，*，’，，\。

其中+和-表示通常标量及矩阵的加法和减法运算；*表示标量或矩阵的乘法；对所有元素为实数的矩阵，’运算对应于转置运算。

若

为一个

非奇异矩阵（det!

=0）且

为一个

矩阵，则运算

等价于

。

实验内容：

1．用Matlab随机生成

的矩阵

和

。

求下列指定的

，并确定那些矩阵是相等的。

你可以利用Matlab计算两个矩阵的差来测试两个矩阵是否相等。

（1）C=A*B,D=B*A,G=（A’*B’）’,H=（B’*A’）’C=H;D=G；

（2）C=A’*B’,D=（A*B）’,G=B’*A’,H=（B*A）’C=H;D=G；

（3）C=inv（A*B）,D=inv（A）*inv（B）,G=inv（B*A）,H=inv（B）*inv（A）

（4）C=inv（（A*B）’）,D=inv（A’*B’）,G=inv（A’）*inv（B’）,H=（inv（A）*inv（B））’

（3）（4）中无相等的

2．令n=200，并使用命令

A=floor（10*rand（n））;

b=sum（A’）’

z=ones（n,1）;注释：

（n行一列全为1的矩阵）

生成一个

矩阵和两个

中的向量，它们的元素均为整数。

（因为矩阵和向量都很大，我们添加分号来控制输出。

（1）方程组

的真解应为

。

为什么？

【A中的每一行的元素之和正好等于对应b的每一列，故z为其一解，又det不等于0，RA=RAb=n，故z为其解】试说明，可在Matlab中利用”\”运算或计算

，然后用计算

来求解。

比较这两种计算方法的速度和精度。

我们将使用Matlab命令tic和toc来测量每一个计算过程消耗的时间。

只需要用下面的命令：

tic,x=A\b；toc

tic,y=inv（A）*b;toc

哪一种方法更快？

tic,x=A\b；更快！

为了比较这两种方法的精度，可以测量求得的解x和y与真解z接近的程度。

利用下面的命令：

max（abs（x-z））

max（abs（y-z））

哪种方法的到的解更精确？

>>max（abs（x-z））=4.0168e-013更精确！

>>max（abs（y-z））=6.1107e-013

（2）用n=500和n=1000替换

（1）中的

。

如

（1）结果一样！

3．令A=floor（10*rand（6））。

根据构造，矩阵A将有整数元。

将矩阵A的第六列更改，使得矩阵A为奇异的。

令

B=A’，A（:

6）=-sum（B（1:

5,:

））’

（1）设x=ones（6,1）,并利用Matlab计算Ax。

为什么我们知道A必为奇异的？

【因化简列，————>列成比例】试说明。

通过化为行最简形来判断A是奇异的。

（2）令B=x*[1:

6]，乘积AB应为零矩阵。

为什么？

【因A的每一行的前五个元素之和等于第六个元素的相反数，且在A上的每一行的元素同乘以相同的数，则仍等于0】试说明。

用Matlab的*运算计算AB进行验证。

（3）令C=floor（10*rand（6））和D=B+C，尽管

但乘积AC和AD是相等的。

为什么？

试说明。

计算A*C和A*D，并验证它们确实相等。

【此处B为令B=x*[1:

6]；A为A（:

6）=-sum（B（1:

5,:

））’】

由于A*B=0；故AC=AD；A（B+C）=AB+AC;

4．采用如下方式构造一个矩阵。

令

B=eye（10）-triu（ones（10）,1），参见最后附表二：

为什么我们知道B必为非奇异的？

【上三角矩阵的行列式的值等于对角线上的元素相乘】

令C=inv（B）且x=C（:

10），

现在用B（10,1）=-1/256将B进行微小改变。

利用Matlab计算乘积Bx。

由这个计算结果，你可以得出关于新矩阵B的什么结论？

【化简此时B，得行最简式，RB=9<10，可以得出B的第10列（从1—9行）与x互为相反数，且都是2的指数幂数，且第十行为0，】它是否为奇异的？

【是】试说明。

用Matlab计算它的行最简形。

5．生成一个矩阵A：

A=floor（20*rand（6））

并生成一个向量b：

B=floor（20*rand（6,1））-10

（1）因为A是随机生成的，我们可以认为它是非奇异的。

那么方程组

应有唯一解。

用运算“\”求解。

用Matlab计算[Ab]的行最简形U。

比较U的最后一列和解x，结果是什么？

【相等】在精确算术运算时，它们应当是相等的。

为什么？

【行最简式中可写出对应元素的实际含义，对应处的未知元就等于最后的数】试说明。

为比较他们两个，计算差U（:

7）-x或用formatlong考虑它们。

（2）现在改变A，试它成为奇异的。

令A（:

3）=A（:

1:

2）*[43]’【第一列乘以4加上第二列乘以3替换到第三列上】，利用Matlab计算rref（[Ab]）。

方程组

有多少组解？

【无解】试说明。

【RA

（3）令y=floor（20*rand（6,1））-10且c=A*y，为什么我们知道方程组Ax=c必为相容？

的？

【x此时必有一解y，故为相容的】试说明。

计算[Ac]的行最简形U。

方程组

有多少组解？

【无穷多解】试说明。

【RA=RAc<6】

（4）由行最简形确定的自由变量应为

。

通过考察矩阵U对应的方程组，可以求得

时所对应的解。

将这个解作为列向量

输入Matlab中。

为检验

，计算剩余向量

。

（5）令

。

矩阵U应对应于

的行最简形。

用

求自变量

时齐次线性方程组的解（手工计算），并将你的结果输入为向量Z。

用A*Z检验你的结论。

（6）令

。

向量v应为方程组

的解。

为什么？

试说明。

用Matlab计算剩余向量来验证v为方程组的解。

在这个解中，自由变量

的取值是什么？

【

=3】如何使用向量w和z来求所有可能的方程组的解？

【v=w+n*z,其中n为任意实数】试说明。

6．考虑下图：

（1）确定图的邻接矩阵A，将其输入Matlab；

（2）计算A2并确定长度为2的路的条数【72】，其起止点分别为：

【A^2+A中的数值之和，数字表示有几种路径，具体看程序】

（3）计算A4、A6、A8并回答

（2）中各种情况长度为4、【368】6、【2362】8、【15800】的路的条数。

试推测什么时候从顶点Vi到Vj没有长度为偶数【即为0】的路。

【i=1，j=6;i=2,j=5;i=3,j=6或8；i=4,j=7;i=5,j=8;i=6,j=1或3;i=7,j=4;i=8,j=3或6;】

（4）计算A3、A5、A7并回答

（2）中各情况长度为3、【154】5、【922】7【6098】的路的条数。

你由（3）得到的推测对长度为奇数的路是否成立？

【不成立】，试说明【见程序】。

推测根据i+j+k的奇偶性，是否存在长度为k的路。

【若i+j+k为偶数，不存在；相反，则存在】【路径见程序】

（5）如果我们在图中增加边{V3,V6}，{V5,V8}，新图的邻接矩阵B可首先令B=A，然后令B（3,6）=144,B（6,3）=1,B（5,8）=1,B（8,5）=1，对k=2,3,4,5计算Bk。

（4）中的推测在新的图形中是否还是成立的？

【不成立】见程序】

（6）在图中增加{V6,V8}，并构造得到的图的邻接矩阵C，计算C的幂次，并验证你在（4）中的推测对这个新图是否仍然成立。

【不成立】【见程序】

7．令A=magic（8），然后计算其行最简形。

使得首1对应于前三个变量

，且其余的五个变量均为自由的。

（1）令c=[1:

8]’，通过计算矩阵[Ac]的行最简形确定方程组Ax=c是否相容。

方程组是相容的吗？

【不相容】试说明。

【RA

（2）令b=[8-8-888-8-88]’;

并考虑方程组Ax=b。

该方程组应为相容的。

通过U=rref（[Ab]）验证。

对五个自由变量的任一组取值，我们都应可以得到一组解。

事实上，令

x2=floor（10*rand（5,1）），若x2表示方程组解的最后5个坐标，则我们由x2求得x1=（x1,x2,x3）’。

要这样做，只需要令U=rref（[Ab]）。

U的非零行对应于分块形式的线性方程组

为解此方程组，令V=U（1:

3,4:

8）,c=U（1:

3,9）

并利用Matlab，根据x2，c和V计算x1。

令x=[x1;x2]，验证x是方程组的解。

8．令B=[-1,-1;1,1]和A=[zeros

（2）,eye

（2）;eye

（2）,B]验证B2=0。

（1）用Matlab计算A2，A4，A6，A8。

猜想用子矩阵E，O和B如何表示分块形式的A2k。

用数学归纳法证明你的猜想对任何正整数k都是成立的。

（2）用Matlab计算A3，A5，A7和A9。

猜想用子矩阵E，O和B如何表示分块形式的A2k-1。

用数学归纳法证明你的猜想对任何正整数k都是成立的。

9．

（1）Matlab命令

A=floor（10*rand（6））,B=A’*A

将得到元素为整数的对称矩阵。

为什么？

试说明。

【第i行第j列的数等于第i列的数分别乘以第j列的数之和；第j行第i列的数等于第j列的数分别乘以第i列的数之和，故为对称矩阵】

用这种方法计算B来验证结论，然后将B划分成四个3x3的子矩阵。

在Matlab中求子矩阵，令

B11=B（1:

3,1:

3）,B12=B（1:

3,4:

6）

并用B的第四行到第6行类似定义B21和B22。

（2）令C=inv（B11）。

应有CT=C和B21T=B12。

为什么？

【对称阵的逆矩阵与该逆矩阵的转置是相等的，B12的第i行的数等于B21的第i列的数】

试说明。

用Matlab运算符’计算转置，并验证结论。

然后，令

G=B21*C和H=B22-B21*C*B21’

利用Matlab函数eye和zeros构造

计算W=L*D*L’，并通过计算W-B与B进行比较。

证明：

若用算术运算精确计算LDLT，它应准确等于B。

附表：

第一题：

（1）

>>A=rand（4）;

>>B=rand（4）;

>>C=A*B;

>>D=B*A;

>>G=（A'*B'）';

>>H=（B'*A'）';

>>C-D

ans=

2.2376e-0014.7289e-0011.3979e+0001.3204e+000

-6.3633e-001-3.0354e-0012.2485e-002-1.5056e-001

-1.7227e-001-1.1938e-0012.9484e-0012.3624e-001

-8.7955e-001-6.5016e-0018.0370e-002-2.1506e-001

>>C-G

ans=

2.2376e-0014.7289e-0011.3979e+0001.3204e+000

-6.3633e-001-3.0354e-0012.2485e-002-1.5056e-001

-1.7227e-001-1.1938e-0012.9484e-0012.3624e-001

-8.7955e-001-6.5016e-0018.0370e-002-2.1506e-001

>>C-H

ans=

0000

0000

0000

0000

>>D-G

ans=

0000

0000

0000

0000

>>D-H

ans=

-2.2376e-001-4.7289e-001-1.3979e+000-1.3204e+000

6.3633e-0013.0354e-001-2.2485e-0021.5056e-001

1.7227e-0011.1938e-001-2.9484e-001-2.3624e-001

8.7955e-0016.5016e-001-8.0370e-0022.1506e-001

>>G-H

ans=

-2.2376e-001-4.7289e-001-1.3979e+000-1.3204e+000

6.3633e-0013.0354e-001-2.2485e-0021.5056e-001

1.7227e-0011.1938e-001-2.9484e-001-2.3624e-001

8.7955e-0016.5016e-001-8.0370e-0022.1506e-001

>>

（2）

>>C=A'*B';

>>D=（A*B）';

>>G=B'*A';

>>H=（B*A）';

>>C-D

ans=

-2.2376e-0016.3633e-0011.7227e-0018.7955e-001

-4.7289e-0013.0354e-0011.1938e-0016.5016e-001

-1.3979e+000-2.2485e-002-2.9484e-001-8.0370e-002

-1.3204e+0001.5056e-001-2.3624e-0012.1506e-001

>>C-G

ans=

-2.2376e-0016.3633e-0011.7227e-0018.7955e-001

-4.7289e-0013.0354e-0011.1938e-0016.5016e-001

-1.3979e+000-2.2485e-002-2.9484e-001-8.0370e-002

-1.3204e+0001.5056e-001-2.3624e-0012.1506e-001

>>C-H

ans=

0000

0000

0000

0000

>>D-G

ans=

0000

0000

0000

0000

>>D-H

ans=

2.2376e-001-6.3633e-001-1.7227e-001-8.7955e-001

4.7289e-001-3.0354e-001-1.1938e-001-6.5016e-001

1.3979e+0002.2485e-0022.9484e-0018.0370e-002

1.3204e+000-1.5056e-0012.3624e-001-2.1506e-001

>>G-H

ans=

2.2376e-001-6.3633e-001-1.7227e-001-8.7955e-001

4.7289e-001-3.0354e-001-1.1938e-001-6.5016e-001

1.3979e+0002.2485e-0022.9484e-0018.0370e-002

1.3204e+000-1.5056e-0012.3624e-001-2.1506e-001

>>

（3）

>>C=inv（A*B）;

>>D=inv（A）*inv（B）;

>>G=inv（B*A）;

>>H=inv（B）*inv（A）;

>>C-D

ans=

-3.9602e+001-1.4016e+0011.4537e+0012.2261e+001

1.5266e+0011.5778e+001-1.9398e+001-3.9304e+001

1.0821e+0011.4313e+000-2.7296e+001-4.8956e+001

1.3845e+001-5.5182e-0012.6289e+0015.1120e+001

>>C-G

ans=

-3.9602e+001-1.4016e+0011.4537e+0012.2261e+001

1.5266e+0011.5778e+001-1.9398e+001-3.9304e+001

1.0821e+0011.4313e+000-2.7296e+001-4.8956e+001

1.3845e+001-5.5182e-0012.6289e+0015.1120e+001

>>C-H

ans=

-5.6843e-014-1.2879e-0143.0198e-0147.1054e-014

-6.5370e-013-1.4744e-0133.3396e-0138.2423e-013

-1.5774e-012-3.5527e-0137.8870e-0131.9895e-012

1.8758e-0124.2988e-013-9.4502e-013-2.4016e-012

>>D-G

ans=

4.9738e-0131.1013e-013-8.3489e-014-3.1264e-013

1.7053e-0133.7303e-014-2.4869e-014-1.0747e-013

5.8265e-0131.3145e-013-9.4147e-014-3.8369e-013

-1.0516e-012-2.3448e-0131.7053e-0136.6791e-013

>>D-H

ans=

3.9602e+0011.4016e+001-1.4537e+001-2.2261e+001

-1.5266e+001-1.5778e+0011.9398e+0013.9304e+001

-1.0821e+001-1.4313e+0002.7296e+0014.8956e+001

-1.3845e+0015.5182e-001-2.6289e+001-5.1120e+001

>>G-H

ans=

3.9602e+0011.4016e+001-1.4537e+001-2.2261e+001

-1.5266e+001-1.5778e+0011.9398e+0013.9304e+001

-1.0821e+001-1.4313e+0002.7296e+0014.8956e+001

-1.3845e+0015.5182e-001-2.6289e+001-5.1120e+001

>>

（4）

>>c=inv（（A*B）'）;

>>d=inv（A'*B'）;

>>g=inv（A'）*inv（B'）;

>>h=（inv（A）*inv（B））';

>>c-d

ans=

-3.9602e+0011.5266e+0011.0821e+0011.3845e+001

-1.4016e+0011.5778e+0011.4313e+000-5.5182e-001

1.4537e+001-1.9398e+001-2.7296e+0012.6289e+001

2.2261e+001-3.9304e+001-4.8956e+0015.1120e+001

>>c-g

ans=

-1.6875e-014-5.4712e-013-1.3216e-0121.5774e-012

-2.8866e-015-1.3145e-013-3.1264e-0133.7659e-013

8.8818e-0152.6290e-0136.3949e-013-7.6028e-013

2.5757e-0147.1765e-0131.7195e-012-2.0606e-012

>>c-h

ans=

-3.9602e+0011.5266e+0011.0821e+0011.3845e+001

-1.4016e+0011.5778e+0011.4313e+000-5.5182e-001

1.4537e+001-1.9398e+001-2.7296e+0012.6289e+001

2.2261e+001-3.9304e+001-4.8956e+0015.1120e+001

>>d-g

ans=

3.9602e+001-1.5266e+001-1.0821e+001-1.3845e+001

1.4016e+001-1.5778e+001-1.4313e+0005.5182e-001

-1.4537e+0011.9398e+0012.7296e+001-2.6289e+001

-2.2261e+0013.9304e+0014.8956e+001-5.1120e+001

>>d-h

ans=

-2.4158e-013-1.1724e-013-2.7711e-0135.2580e-013

-5.6843e-014-1.8652e-014-5.3291e-0141.0658e-013

4.2633e-0141.7764e-0144.7962e-014-8.8818e-014

1.5987e-0136.7502e-0141.8474e-013-3.3396e-013

>>g-h

ans=

-3.9602e+0011.5266e+0011.0821e+0011.3845e+001

-1.4016e+0011.5778e+0011.4313e+000-5.5182e-001

1.4537e+001-1.9398e+001-2.7296e+0012.6289e+001

2.2261e+001-3.9304e+001-4.8956e+0015.1120e+001

>>

第二题：

（1）

>>n=200;

>>A=floor（10*rand（n））;

>>b=sum（A'）';

>>z=ones（n,1）;

>>c=linsolve（A,b）;

>>d=c-z；精度为1e-14-----1e-13;

tic,x=A\b,toc=Elapsedtimeis0.016000seconds.

tic,inv（A）*b,toc=Elapsedtimeis0.031000seconds.

（2）

n=500;

>>tic,x=A\b;toc

Elapsedtimeis0.187000seconds.更快！

>>tic,y=inv（A）*b;toc

Elapsedtimeis0.343000seconds.

>>max（abs（x-z））=4.3987e-013更精确！

>>max（abs（y-z））=2.2524e-012

>>

>>n=1000;

>>tic,x=A\b;toc

Elapsedtimeis0.920000seconds.更快！

>>tic,y=inv（A）*b;toc

Elapsed