线性代数二三事

第6章:线性代数与微分方程

6.2 常数系数下的齐次线性系统

第6章 线性代数与微分方程

在本节里面, 我们将对常数系数下的齐次线性系统进行求解. 假设y(t)y(t)是关于变量tt的函数, 满足

a0y(t)+a1y′(t)+⋯+any(n)(t)=0,(6.1)a_0 y(t)+a_1y'(t)+\cdots+ a_n y^{(n)}(t) = 0,\tag{6.1}

且a0,⋯ ,ana_0,\cdots,a_n均为与tt无关的非零常数, y(n)(t)y^{(n)}(t)为yy关于tt的nn阶导数.形如上式的方程便是我们本节所讨论的重点. 按照上一节提出的方法, 我们可以将等式(6.1)写成是一个线性方程组的形式: 我们不妨设y1(t)=y(t),y2(t)=y′(t),⋯ ,yn(t)=y(n−1)(t)y_1(t) = y(t), y_2(t) = y'(t),\cdots, y_n(t) = y^{(n-1)}(t), 那么这样一来我们便有y1′(t)=y2(t),y2′(t)=y3(t),⋯y_1'(t) = y_2(t), y_2'(t) = y_3(t),\cdots, 以及

yn′(t)=y(n)(t)=−a0any(t)−a1any′(t)−⋯−an−1any(n−1)(t)=−a0any1(t)−a1any2(t)−⋯−an−1anyn(t).\begin{aligned} y_n'(t) = y^{(n)}(t) &= -\frac{a_0}{a_n} y(t)-\frac{a_1}{a_n}y'(t)-\cdots-\frac{a_{n-1}}{a_{n}}y^{(n-1)}(t)\\ &=-\frac{a_0}{a_n} y_1(t) - \frac{a_1}{a_n}y_2(t)-\cdots-\frac{a_{n-1}}{a_{n}} y_{n}(t). \end{aligned}

因此我们便有

(y1′(t)y2′(t)⋮yn′(t))=(y2(t)y3(t)⋮−a0any1(t)−a1any2(t)−⋯−an−1anyn(t))=(010⋯0001⋯0⋮⋮⋮⋱⋮−a0an−a1an⋯⋯−an−1an)(y1(t)y2(t)⋮yn(t)).\begin{aligned} \begin{pmatrix} y_1'(t) \\ y_2'(t) \\ \vdots \\ y_n'(t) \end{pmatrix} &= \begin{pmatrix} y_2(t)\\ y_3(t)\\ \vdots\\ -\frac{a_0}{a_n} y_1(t) - \frac{a_1}{a_n}y_2(t)-\cdots-\frac{a_{n-1}}{a_n} y_{n}(t) \end{pmatrix}\\ &=\begin{pmatrix} 0 & 1 & 0 & \cdots & 0 \\ 0 & 0 & 1 & \cdots & 0\\ \vdots & \vdots & \vdots & \ddots & \vdots \\ -\frac{a_0}{a_n} & -\frac{a_1}{a_n} & \cdots & \cdots & -\frac{a_{n-1}}{a_n} \end{pmatrix} \begin{pmatrix} y_1(t) \\ y_2(t) \\ \vdots \\ y_n(t) \end{pmatrix}. \end{aligned}

因此我们便有形如

y′(t)=Ay(t)\vec y'(t) = \vec A \vec y(t)

的常微分方程组. 形如这样的方程组也被称作自治系统 (Autonomous System). 即矩阵A\vec A中的元素与tt无关. 我们在今后的篇幅中只研究自治系统.

例题 6.4

将微分方程y′′′(t)−2y′′(t)+4y′(t)−5y(t)=0y'''(t) -2y''(t)+4y'(t)-5y(t)=0写成y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)的形式.

解答 6.4

设y1(t)=y(t),y2(t)=y′(t),y3(t)=y′′(t)y_1(t) = y(t), y_2(t) = y'(t), y_3(t) = y''(t). 则

y1′(t)=y2(t),y2′(t)=y3(t),y3′(t)=2y′′(t)−4y′(t)+5y(t)=2y3(t)−4y2(t)+5y1(t).y_1'(t) = y_2(t),y_2'(t) = y_3(t), y_3'(t) = 2y''(t)-4y'(t)+5y(t) = 2y_3(t)-4y_2(t)+5y_1(t).

由此我们可知

(y1′(t)y2′(t)y3′(t))=(0100015−42)(y1(t)y2(t)y3(t)).\begin{pmatrix} y_1'(t)\\ y_2'(t)\\ y_3'(t) \end{pmatrix} = \begin{pmatrix} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 5 & -4 & 2 \end{pmatrix}\begin{pmatrix} y_1(t) \\ y_2(t) \\ y_3(t) \end{pmatrix}.

随后我们来学习如何求解这样的常微分方程组. 在求解之前, 我们需要先给出有关微分方程组解的存在唯一性定理:

定理 6.2

设y′(t)=A(t)y(t)+b(t)\vec y'(t) = \vec A(t) \vec y(t) + \vec b(t), A(t)∈Mn(R),b(t)∈Rn\vec A(t) \in M_n(\mathbb{R}), \vec b(t) \in \mathbb{R}^n. 若A(t),b(t)\vec A(t), \vec b(t)中的每一个元素均为定义在II上关于tt的连续函数, 那么对任意的t0∈It_0\in I, 初值问题y′(t)=A(t)y(t)+b(t),y(t0)=y0\vec y'(t) = \vec A(t) \vec y(t) + \vec b(t), y(t_0) = y_0在II上存在唯一解y(t)\vec y(t).

该定理的证明和上一节中的定理5.15.1的证明十分相似.

在求解这样的方程时, 我们不妨先尝试将问题简化: 假设矩阵为1×11\times 1, 那么这样一来, 我们可以得到形如

y′(t)=λy(t)(6.2)y'(t) = \lambda y(t)\tag{6.2}

的方程, 其中λ\lambda为常数. 我们将y′(t)y'(t)改写作dydt\frac{dy}{dt}, 于是我们有

dydt=λy.\frac{dy}{dt} = \lambda y.

形式上分离变量可得

1ydy=λdt  ⟹  ∫1ydy=∫λdt.\frac{1}{y} dy = \lambda dt \implies \int \frac{1}{y} dy = \int \lambda dt.

由此不难看出(6.2)的解为y(t)=Ceλty(t) = Ce^{\lambda t}, 其中CC为常数. 那么对于形如y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)的方程组, 我们是不是也可以类比这种方法呢? 我们可以尝试将该方程的解写成若干个形如eλte^{\lambda t}的线性组合. 设A∈Mn(R)\vec A \in M_n(\mathbb{R}), 在方程组y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)里面, 我们假设y1(t)=eλtv\vec y_1(t) = e^{\lambda t}\vec v为满足条件的一个解. 因此我们有 y1′(t)=λeλtv\vec y'_1(t) = \lambda e^{\lambda t}\vec v, 那么

y1′(t)=Ay1(t)  ⟹  λeλtv=Aeλtv.y'_1(t) = \vec A \vec y_1(t) \implies \lambda e^{\lambda t} \vec v = \vec A e^{\lambda t}\vec v.

如果我们将eλtve^{\lambda t}\vec v看作一个整体u\vec u, 我们便会得到Au=λu\vec A\vec u = \lambda\vec u. 也就是说, 在y1(t)=eλtv\vec y_1(t) = e^{\lambda t}\vec v中, 若y1′(t)=Ay1(t)\vec y_1'(t) = \vec A\vec y_1(t), 则λ\lambda为A\vec A的一个特征值, v\vec v为与之相关的特征向量. 如果此时A\vec A有多个特征值, 这些特征值与特征向量所构成的解yi(t)=eλitvi\vec y_i(t) = e^{\lambda_i t}\vec v_i的线性组合是否也为满足条件的解呢?

定理 6.3

在y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)中, 若y1(t),y2(t)\vec y_1(t),\vec y_2(t)均为原方程的解, 那么∀c1,c2∈R\forall c_1,c_2 \in \mathbb{R}, c1y1(t)+c2y2(t)c_1\vec y_1(t) + c_2 \vec y_2(t)也为原方程的解.

证明

显然.

∎

因此, 我们知道微分方程组的解即可看作是若干个线性无关的“基础解”的线性组合.

定义 6.4

我们若将y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)的解y1(t),⋯ ,yn(t)\vec y_1(t),\cdots,\vec y_n(t)当作列向量写到一个矩阵Y(t)=(y1(t)⋯yn(t))\vec Y(t) = \begin{pmatrix} \vec y_1(t) & \cdots & \vec y_n(t) \end{pmatrix}, 我们则称Y(t)\vec Y(t)为一个矩阵解 (Matrix Solution).

根据线性代数的知识, y1(t),⋯ ,yn(t)\vec y_1(t),\cdots,\vec y_n(t)彼此线性无关的充要条件为对任意的tt, det⁡(Y(t))≠0\det(\vec Y(t)) \neq 0. 此时我们称y1(t),⋯ ,yn(t)\vec y_1(t),\cdots,\vec y_n(t) 为y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)的一个基础解系 (Fundamental Solution). 如果A∈Mn(R)\vec A \in M_n(\mathbb{R})有nn个不同的实数特征值λ1,⋯ ,λn\lambda_1,\cdots,\lambda_n, 那么我们便可以找到nn个与特征值所对应的彼此线性无关的特征向量v1,⋯ ,vn\vec v_1,\cdots,\vec v_n. 由此, 方程组y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)的解便可以看作是由yi(t)=eλitvi\vec y_i(t) = e^{\lambda_i t}\vec v_i所构成的线性组合:

y(t)=c1eλ1tv1+⋯+cneλntvn.\vec y(t) = c_1e^{\lambda_1 t}\vec v_1+\cdots+c_n e^{\lambda_n t} \vec v_n.

定理 6.4

设常数系数下的线性齐次常微分方程组y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t), A∈Mn(R)\vec A \in M_n(\mathbb{R}). 若A\vec A有nn个不同的实数特征值λ1,⋯ ,λn\lambda_1,\cdots,\lambda_n, 其对应的特征向量分别为v1,⋯ ,vn\vec v_1,\cdots,\vec v_n. 则y′(t)=Ay(t)\vec y'(t) = \vec A \vec y(t)的解集为

y(t)=c1eλ1tv1+⋯+cneλntvn,c1,⋯ ,cn∈R.\vec y(t) = c_1e^{\lambda_1 t}\vec v_1+\cdots+c_n e^{\lambda_n t} \vec v_n, \quad c_1,\cdots,c_n \in \mathbb{R}.

对于存在nn个不同的实数特征值的矩阵而言, 线性齐次常微分方程组y′(t)=Ay(t)\vec y'(t) = \vec A \vec y(t)便很好求得. 我们不禁要问: 如果A\vec A含有复数特征值呢? 如果A\vec A无法对角化呢? 我们将在本节讨论特征值为复数的情况.

推论 6.4

设A∈Mn(R)\vec A \in M_n(\mathbb{R}), 且λ=a+bi\lambda = a+bi为A\vec A的一个特征值, 那么λˉ=a−bi\bar{\lambda} = a-bi也为A\vec A的一个特征值. 特别地, 如果u+iv\vec u+i\vec v为λ=a+bi\lambda=a+bi所对应的特征向量, 那么u−iv\vec u-i\vec v为λˉ=a−bi\bar{\lambda}=a-bi所对应的特征向量.

设A∈Mn(R)\vec A \in M_n(\mathbb{R})存在复数特征值. 我们取λ=a+bi,(a,b∈R)\lambda=a+bi, (a,b \in \mathbb{R})和特征向量u+iv,(u,v∈Rn)\vec u+i\vec v, (\vec u,\vec v \in \mathbb{R}^n), 那么根据定义,

y(t)=e(a+bi)t(u+iv)(6.3)\vec y(t) = e^{(a+bi)t}(\vec u+i\vec v)\tag{6.3}

即为满足条件的一个解. 但由于A∈Mn(R)\vec A \in M_n(\mathbb{R}), 所以我们不能就这样把这一个复数解直接放上去. 我们可以利用Euler公式对式子(6.3)进行变形. 根据Euler公式, eix=cos⁡(x)+isin⁡(x)e^{ix} = \cos(x) + i\sin(x). 那么

y(t)=eλt⋅(u+iv)=eat⋅(ebitu+ebit(iv))=eat⋅((cos⁡(bt)+isin⁡(bt))u+i⋅(cos⁡(bt)+isin⁡(bt))v)=eat⋅(cos⁡(bt)u−sin⁡(bt)v)+i⋅eat(sin⁡(bt)u+cos⁡(bt)v).\begin{aligned} \vec y(t) &= e^{\lambda t} \cdot (\vec u+i\vec v)\\ &=e^{at}\cdot\Big( e^{bit}\vec u + e^{bit} (i\vec v)\Big)\\ &=e^{at}\cdot\Big( (\cos(bt)+i\sin(bt)) \vec u + i\cdot(\cos(bt) + i\sin(bt))\vec v \Big)\\ &=e^{at}\cdot\Big(\cos(bt) \vec u - \sin(bt)\vec v \Big) + i\cdot e^{at}\Big(\sin(bt)\vec u+\cos(bt)\vec v\Big). \end{aligned}

推论 6.5

设y1(t)=eat⋅(cos⁡(bt)u−sin⁡(bt)v),y2(t)=eat(sin⁡(bt)u+cos⁡(bt)v)\vec y_1(t) = e^{at}\cdot\Big(\cos(bt) \vec u - \sin(bt)\vec v \Big), \vec y_2(t) = e^{at}\Big(\sin(bt)\vec u+\cos(bt)\vec v\Big). 则y1(t),y2(t)\vec y_1(t), \vec y_2(t)均为原方程组y′(t)=Ay(t)\vec y'(t) = \vec A\vec y(t)的解.

证明

由于y(t)=y1(t)+iy2(t)\vec y(t) = \vec y_1(t) + i \vec y_2(t) and y′(t)=Ay(t)y'(t) = \vec A\vec y(t), 那么

y′(t)=y1′(t)+iy2′(t)=Ay1(t)+iAy2(t).\vec y'(t) = \vec y_1'(t) + i \vec y'_2(t) = \vec A\vec y_1(t) + i \vec A\vec y_2(t).

由于b≠0b \neq 0, 则y1(t),y2(t)\vec y_1(t),\vec y_2(t)线性无关. 因此y1′(t)=Ay1(t);y2′(t)=Ay2(t)\vec y_1'(t) = \vec A\vec y_1(t); \vec y_2'(t) = \vec A\vec y_2(t).

∎

例题 6.5

设α,β,γ∈R\alpha,\beta,\gamma \in \mathbb{R}, 矩阵A=(α−β0βα000γ)A = \begin{pmatrix} \alpha &-\beta &0 \\ \beta & \alpha & 0 \\ 0 & 0 & \gamma \end{pmatrix}. 求出 y′(t)=Ay(t)\vec y'(t)=A\vec y(t)的解集.

解答 6.5

A\vec A的特征多项式为

p(λ)=det⁡(A−λI)=det⁡(α−λ−β0βα−λ000γ−λ)=(γ−λ)[(α−λ)2+β2],p(\lambda) = \det(A-\lambda I) = \det\begin{pmatrix} \alpha-\lambda & - \beta & 0 \\ \beta & \alpha-\lambda & 0 \\ 0 & 0 & \gamma-\lambda \end{pmatrix} = (\gamma-\lambda)[(\alpha-\lambda)^2+\beta^2],

因此A\vec A的特征值为

λ1=α+iβ,λ2=α−iβ,λ3=γ,\lambda_1=\alpha+i\beta, \hspace{0.2cm} \lambda_2=\alpha-i\beta, \hspace{0.2cm} \lambda_3=\gamma,

与特征值所对应的特征向量为

u1∈ker⁡(A−λ1I)=ker⁡(−iβ−β0β−iβ000γ−α−iβ)∈Span(1−i0)  ⟹  令u1=(1−i0).\vec u_1 \in \ker(A-\lambda_1 I) = \ker\begin{pmatrix} -i\beta & -\beta & 0 \\ \beta & -i\beta & 0 \\ 0 & 0 & \gamma-\alpha-i\beta \end{pmatrix} \in \Span\begin{pmatrix} 1 \\ -i \\ 0 \end{pmatrix} \implies \text{令} \vec u_1 = \begin{pmatrix} 1 \\ -i \\ 0 \end{pmatrix}.

根据推论6.46.4, 我们很容易得到 u2=(1i0)\vec u_2 = \begin{pmatrix} 1 \\ i \\ 0 \end{pmatrix}.

最后,

u3∈ker⁡(A−λ3I)=Span(001)  ⟹  令u3=(001)\vec u_3 \in \ker(A-\lambda_3 I) = \Span\begin{pmatrix} 0 \\ 0\\ 1 \end{pmatrix} \implies \text{令} \vec u_3 = \begin{pmatrix} 0 \\ 0\\ 1 \end{pmatrix}

对于特征值α±βi\alpha\pm\beta i而言, 其对应的解满足

y(t)=eαt(cos⁡(βt)+isin⁡(βt))⋅((100)+i(0−10))=eαt[cos⁡(βt)(100)−sin⁡(βt)(0−10)]+ieαt[sin⁡(βt)(100)+cos⁡(βt)(0−10)].\begin{aligned} \vec y(t) &= e^{\alpha t}\left(\cos(\beta t) + i\sin(\beta t) \right)\cdot \left( \begin{pmatrix} 1 \\ 0 \\ 0 \end{pmatrix} +i\begin{pmatrix} 0\\-1 \\ 0 \end{pmatrix} \right)\\ &=e^{\alpha t} \left[ \cos(\beta t)\begin{pmatrix} 1 \\ 0 \\ 0\end{pmatrix} -\sin(\beta t) \begin{pmatrix} 0 \\ -1 \\0 \end{pmatrix} \right]+ie^{\alpha t}\left[ \sin(\beta t)\begin{pmatrix} 1\\0\\0\end{pmatrix} + \cos(\beta t) \begin{pmatrix} 0 \\ -1 \\0 \end{pmatrix} \right]. \end{aligned}

则

y1(t)=eαt(cos⁡(βt)sin⁡(βt)0),y2(t)=eαt(sin⁡(βt)−cos⁡(βt)0).\vec y_1(t) = e^{\alpha t} \begin{pmatrix} \cos(\beta t) \\ \sin(\beta t) \\ 0 \end{pmatrix}, \quad\vec y_2(t) = e^{\alpha t} \begin{pmatrix} \sin(\beta t) \\ -\cos(\beta t) \\ 0 \end{pmatrix}.

对于λ3=γ\lambda_3=\gamma而言, 显然y3=eγtu3\vec y_3 = e^{\gamma t} \vec u_3. 则原方程的解集为

y(t)=c1eαt(cos⁡(βt)sin⁡(βt)0)+c2eαt(sin⁡(βt)−cos⁡(βt)0)+c3eγt(001),c1,c2,c3∈R.\vec y(t) = c_1e^{\alpha t} \begin{pmatrix} \cos(\beta t) \\ \sin(\beta t) \\ 0 \end{pmatrix}+c_2 e^{\alpha t} \begin{pmatrix} \sin(\beta t) \\ -\cos(\beta t) \\ 0 \end{pmatrix}+c_3e^{\gamma t}\begin{pmatrix} 0 \\0\\1\end{pmatrix},\quad c_1,c_2,c_3\in\mathbb{R}.

{6.2 练习}

1. [带阻尼的LRCLRC电路] 考虑串联LRCLRC电路: 其中LL为电感, RR为电阻, CC为电容. I(t)=I(t)= 为电路中的电流. 已知微分方程

LI′′(t)+RI′(t)+1CI(t)=0.LI''(t)+RI'(t)+\frac{1}{C}I(t)=0.

(i) [过阻尼] 设L=1,R=3,C=0.5L=1, R=3,C=0.5. 求出I(t)I(t)随时间的变化关系. 并大致画出其图像;

(ii) [欠阻尼] 设L=1,R=1,C=1L=1, R=1, C=1. 求出I(t)I(t)随时间的变化关系. 并大致画出其图像.

2. [溶液溶质问题] 如图, 有A,B,CA,B,C三个密封良好且装满水的水缸, 且A,B;B,C;A,CA,B; B,C; A,C水缸之间通过单向阀连接. A,CA,C水缸容积均为60L60L, BB水缸容积为30L30L. 已知A,B,CA,B,C缸中均匀且充分地溶解了0kg,1kg,4kg0kg, 1kg, 4kg的盐, 且三个单向阀均处于闭合状态. 在某一时刻 (t=0t=0), 三个阀门突然同时开启. 随后每单位时间内有10L10L的盐水分别由AA流向BB, 由BB流向CC, 以及由CC流向AA. 设A,B,CA,B,C缸中盐水的浓度在tt时刻分别为y1(t),y2(t),y3(t)y_1(t), y_2(t), y_3(t) (单位: kg/Lkg/L). 建立合适的微分方程组, 求出y1(t),y2(t),y3(t)y_1(t), y_2(t), y_3(t)随时间的变化关系.

3. [弹簧振子系统] 如图,物块m1,m2m_1,m_2被三根劲度系数分别为k1,k2,k3k_1,k_2,k_3的弹簧连接. 设在tt时刻对m1,m2m_1,m_2施加的力分别为F1(t),F2(t)F_1(t),F_2(t),物块m1,m2m_1,m_2相较初始位置的位移分别为x1,x2x_1,x_2. 取右边为正方向,不计一切摩擦,根据牛顿运动定律和胡克定律我们有

m1x1′′(t)=k2(x2−x1)−k1x1+F1(t);m_1x_1''(t) = k_2(x_2-x_1)-k_1x_1+F_1(t); m2x2′′(t)=−k3x2−k2(x2−x1)+F2(t).m_2x''_2(t) = -k_3x_2 - k_2(x_2-x_1)+F_2(t).

设m1=m2=1,k1=k3=1,k2=2,F1=F2=0\displaystyle{m_1=m_2=1,k_1=k_3=1,k_2=2, F_1=F_2=0}. 求出x1(t),x2(t)x_1(t),x_2(t)的通解, 并大致描述不同的初始条件x1(0),x2(0),x1′(0),x2′(0)x_1(0),x_2(0),x_1'(0),x_2'(0)对系统可能产生的影响.