电子科技大学数学实验3.线性代数实验

3 线性代数实验

1.1 基础训练

1. 用左除运算符求解线性方程组

$$\{\begin{array}{c}7x_1+x_2+2x_3-x_4=6,\\ 2x_1+6x_2+x_3+3x_4=5,\\ x_1-4x_2+10x_3+3x_4=4,\\ 2x_1-x_2-3x_3+11x_4=3.\end{array}$$

​ 返回解向量（1个返回参数）。注：要求本问题不使用符号计算函数solve解线性方程组。

​ 首先我们根据方程定义系数矩阵A和常数矩阵B 设未知数所组成的列向量为X

​ 由AX=B X=A^-1*B

A=[7 1 2 -1
    2 6 1 3
    1 -4 10 3
    2 -1 -3 11];
B=[6;5;4;3];
X=A\B
X =

    0.7299
    0.3828
    0.3953
    0.2826

2. 基础知识：假定现在要在平面上旋转一个物体，需要绘制物体的旋转过程。由于绘制过程需要知道物体上各点的坐标，而旋转过程中物体上点的坐标是不断变化的。因此，绘制过程需要跟踪物体坐标的变化. 问题可以归结为已知旋转角度$\theta$，旋转前P点的坐标为$(x,y)$，计算旋转后得点$P^{\prime }$，设其坐标为$(x^{\prime },y^{\prime })$. 这里P点代表任意待旋转的点. 旋转变换就是一种线性变换，看如何构造一个矩阵完成这种旋转变换，完成坐标的计算. 问题进一步化为构造矩阵$M=\left(\begin{array}{cc}{m}_{11}& m_{12}\\ m_{21}& m_{22}\end{array}\right)$，满足

$$\left(\begin{array}{c}{x}^{\prime }\\ y^{\prime }\end{array}\right)=\left(\begin{array}{cc}{m}_{11}& m_{12}\\ m_{21}& m_{22}\end{array}\right)\left(\begin{array}{c}x\\ y\end{array}\right).$$

​ 通过建模可得旋转变换矩阵$M=\left(\begin{array}{cc}cos\theta & -sin\theta \\ sin\theta & cos\theta \end{array}\right).$

​ 请设计一个正六边形绕原点旋转一圈的动画程序. 正六边形的位置、大小自行设置。

N=6;
t=0:2*pi/N:2*pi;
x=2*sin(t)
y=2*cos(t)
z=[x;y];
plot(x,y)
hold on
for th=0:pi/5:2*pi
    M=[cos(th) -sin(th)
        sin(th) cos(th)];
    hh=M*z;
    a=hh(1,:);
    b=hh(2,:);
    plot(a,b)
    pause(0.5)
end
![在这里插入图片描述](https://img-blog.csdnimg.cn/008c578a86aa49e58bd825589b30e52c.png#pic_center)

1.2实验任务

一. 实验任务

已知某种生物分3个年龄组，$02$岁为第一组，$3-4$岁为第二组，$56$岁为第三组. 每两年为一个时段. 设$x_j^{(k)}$为第j个年龄组在第k个时段的数量. 记$X^{(k)}=(x_1^{(k)},x_2^{(k)},x_3^{(k)}{)}^T$. 根据生物繁殖率以及存活率，得到该种生物数量变化规律的数学模型$X^{(k+1)}=\left(\begin{array}{ccc}0& 3& 2\\ 0.4& 0& 0\\ 0& 0.15& 0\end{array}\right)X^{(k)}(k=0,1,2,...)$，已知$X^{(0)}=(10,20,15{)}^T$.

请编程计算2、4、6、8、10年后该种生物各年龄组数量，当k足够大时，生物各年龄组数量变化规律是什么，说明导致该规律的原因。

（2）如果使用某种杀虫剂后，生物数量变化规律由式子$X^{(k+1)}=\left(\begin{array}{ccc}0& 3r& 2r\\ 0.4& 0& 0\\ 0& 0.15& 0\end{array}\right)X^{(k)}$

决定，$r$反应使用杀虫剂的剂量$（0.1\le r\le 0.9）$，请编写程序确定$r$使得三个年龄组昆虫数量最后基本趋于稳定不变，并找出初始稳定值。

二. 实验目的

​ 熟悉Matlab矩阵计算函数。

三. 实验过程

（1）

数学模型已经为我们构建完毕

设变化矩阵为A 根据递推关系得出关系式：

$X^{(k+1)}=A^k{X}^{(0)}$

计算2-10年的数量，实现程序如下

A=[0 3 2
   0.4 0 0
   0 0.15 0];
x=[10;20;15];
for i=2:2:10
    X=[];
    X=A^i*x;
    fprintf('第%d年的数量分别为\n',i)
    disp(vpa(X,4))
end

输出结果：

第2年的数量分别为

18.0

36.0

0.6

第4年的数量分别为

32.4

43.68

1.08

第6年的数量分别为

51.98

53.28

1.944

第8年的数量分别为

78.36

65.49

3.119

第10年的数量分别为

113.7

81.08

4.702

通过改变循环的次数 观察输出结果我们发现 三个年龄段的种群数量均随着时间的推移而增长 但是相较奇偶年数又有所不同 我们绘制三个年龄阶段的1-10年的数量

可以直观的看出变化规律：

年龄组1 偶数年份相对于前一年的奇数年有所下降 整体呈现上升的趋势

年龄组2偶数年份相对于前一年的奇数年有所上升 整体呈现上升的趋势

年龄组3 偶数年份相对于前一年的奇数年有所下降 整体呈现上升的趋势

下面我们来探讨原因

通过[x,y]=eig(A) 求得A的特征值和特征向量

x =

   -0.9429   -0.9323    0.1393
   -0.3302   0.3580   -0.5525
   -0.0434   -0.0516    0.8218

y =

    1.1424         0         0
         0    1.0415         0
         0         0   -0.1009

为了方便起见我们可以将特征值改为正值 同时改变特征向量的正负

设λ1=1.1424，λ2=-1.0415，λ3=-0.1009 α1，α2，α3为对应于λ1，λ2，λ3的特征向量，则三个特征向量线性无关。

因此x可由α1，α2，α3线性表示，即x=[x1;x2;x3]= c1α1+ c2α2+c3α3

syms c1 c2 c3;
C=[c1 c2 c3];
[c1,c2,c3]=solve(C*x(:,1)==10,C*x(:,2)==20,C*x(:,3)==15);
disp(vpa([c1,c2,c3],4))

其中[c1 c2 c3]=[ -17.36, 15.17, 31.4]

由于λ2=-1.0415，λ3=-0.1009对应的特征向量都有正有负 所以奇偶年份有所不同
$$
X(k)\ =A^kX(0)\ =\mathrm{\ }Ak(c1α1+ c2α2+c3α3) \

=\mathrm{\ }c1Akα1+ c2Akα2 + c3Akα3 \

=\mathrm{\ }c1λ1kα1+ c2λ2kα2 + c3λ3kα3 \
$$

当k趋于无穷时
$$X\to -17.36\ast 1.14^k\left(\begin{array}{c}-0.9429\\ -0.3302\\ -0.0434\end{array}\right)+15.71\ast 1.04k\left(\begin{array}{c}-0.9323\\ 0.3580\\ -0.0516\end{array}\right)$$

所以随着时间的推移 种群数量不断增加

（2）

通过（1）的推导 我们得知只有在主特征值等于1的时候才会在时间趋于无穷的时候出现均衡

即问题转化为当r为何值时特征值为1

根据特征值的定义|λI-A|=0 带入特征值1 便可以求出r的值

syms r
A=[0 3*r 2*r
   0.4 0 0
   0 0.15 0];
r=solve(det((eye(3)-A)))

解得r =25/33

下面计算稳定值 程序仿照（1）所示

syms r
A=[0 3*r 2*r
   0.4 0 0
   0 0.15 0];
r=solve(det((eye(3)-A)));
A=[0 3*r 2*r
   0.4 0 0
   0 0.15 0];
x=[10;20;15];
for i=2:1:100
    X=[];
    X=A^i*x;
    fprintf('第%d年的数量分别为\n',i)
    disp(vpa(X,4))
end

由于我们采用的是四位有效数字 当k>46的时候基本上趋于稳定：

三个种群的数量分别为：39.56 15.83 2.373

实验自评与改进方向

自评：关于第一问完成情况较好，但是第二问并没有通过编程找到r,而是借助第一问推导出的结果，通过特征值为1的条件去解得r值，理想的解决办法是通过数量相等去解得r值。

改进方向：

syms r;
A=[0 3*r 2*r
   0.4 0 0
   0 0.15 0];
x=[10;20;15];
X=[];
for i=100:100:200
    X=[X,A^i*x];
    X=vpa(X,2);
end
disp(X)
r=solve(X(1,1)==X(1,2),X(2,1)==X(2,2),X(3,1)==X(3,2))

但是由于方程次数太高还有精确度的问题，无法得出有效的结果。但是这样改进的方向是没有问题的。想用fix（100*X）/100来确定保留的小数位数，但是由于数据类型的原因无法实现。

实验体会，收获及建议

题目主要涉及的是矩阵的相关知识，但是也用到了许多之前学到的内容。这就需要及时的回顾一下前边的内容。同时也结合了上学期的矩阵相关知识，用matlab去解决实际问题。将知识真正的应用到了实际生活中。

更好的了解了上课所讲的增长模型，回顾了特征值和特征向量去分析解决迭代问题的方法。更为重要的是了解了用matlab解决问题的思路和方法。

在编程中遇到的一些意想不到的问题也让我意识到，要在平时多加练习，不能眼高手低。