6.5 非线性系统与线性化

在本节里面, 我们将学习非线性系统的稳定性以及相图. 我们前面提到的Lorenz方程便是一个很好的非线性系统的例子:

\begin{pmatrix} x'(t) \\ y'(t) \\ z'(t) \end{pmatrix} = \begin{pmatrix} \sigma(y-x) \\ x(\rho-z) - y \\ xy-\beta z(t)\end{pmatrix},\quad \sigma,\rho,\beta\in\mathbb{R}.

由于 $xz, xy$ 项的出现, 使得这个系统成为了非线性系统. 这样一来最直接的结果便是矩阵 $\vec A$ 的结构就变复杂了:

\begin{pmatrix} x'(t) \\ y'(t) \\ z'(t) \end{pmatrix} = \begin{pmatrix} -1 & \sigma & 0 \\ \rho & -1 & -x(t) \\ y(t) & 0 & -\beta \end{pmatrix} \begin{pmatrix} x(t) \\ y(t) \\ z(t) \end{pmatrix}.

因此我们为了简化数学符号, 我们随后会将 $t$ 略去. 类似地, 为了方便理解, 我们不妨先考虑一元函数的情况. 假设我们有微分方程 $x'(t)=f(x)$ , 且 $x^*$ 为平衡点. 因此 $f(x^*)=0$ . 我们令 $u(t) = x(t)-x^*$ , 那么 $u'(t)=x'(t)$ , 且在系统 $u'(t)=f(u)$ 中 $u=0$ 便是一个平衡点. 因此我们考虑 $f(x)$ 在 $x^*$ 附近的Taylor展开式:

f(x)=f(u+x^*) = f(x^*)+uf'(x^*)+\frac{1}{2}u^2f''(x^*)+\frac{1}{6}u^3f'''(x^*)+\cdots

结合 $u'(t)=x'(t)=f(x), f(x^*)=0$ , 我们有

u'(t) = uf'(x^*)+O(u^2).

当 $x$ 很接近 $x^*$ 时, 我们可以类比函数在某一点处的切线. 此时我们可以忽略二次项的影响, 因而我们便有 $u'(x) = uf'(x^*)$ . 这样的过程也被称作是线性化 (Linearization). 这样的过程对于一个多元函数同样适用. 如果我们设 $\vec x'(t) = f(\vec x)$ , $\vec x \in \mathbb{R}^n, f=(f_1,\cdots,f_n)^\top$ , 那么其线性化结果为

\vec u'(t) = \mathcal{J} \vec u(t),

其中 $\mathcal{J}$ 即为函数 $f$ 在 $x^*$ 处的Jacobi矩阵:

\mathcal{J} = \mathcal{J}_f(\vec x^*) = \begin{pmatrix} \frac{\pl f_1}{\pl x_1} & \frac{\pl f_1}{\pl x_2} & \cdots & \frac{\pl f_1}{\pl x_n} \\ \vdots & \vdots & \ddots & \vdots\\ \frac{\pl f_n}{\pl x_1} & \frac{\pl f_n}{\pl x_2} & \cdots & \frac{\pl f_n}{\pl x_n} \end{pmatrix}_{\vec x=\vec x^*}.

这样一来, 对于一个非线性系统而言, 在每一个平衡点附近我们便有了形如 $\vec y'(t)=\vec A\vec y(t)$ 的常数微分方程组. 因此, 我们此时便想知道 $\mathcal{J}$ 在每一个平衡点的特征向量与特征值. 通过我们以前的学习, 我们自然会猜测: 若特征值的实部大于零, 则平衡点不稳定; 若特征值的实部小于零, 则平衡点渐近稳定.

定理 6.10

在非线性系统 $\vec x'(t) = f(\vec x(t))$ 中, 设 $\vec x^*$ 为一平衡点, $\mathcal{J}$ 为系统在该点的线性化矩阵. 若 $\mathcal{J}$ 存在实部大于零的特征值, 则 $\vec x^*$ 不稳定; 若 $\mathcal{J}$ 的所有特征向量的实部均小于零, 则 $\vec x^*$ 渐近稳定.

在平衡点 $\vec x^*$ 处, 我们定义系统关于 $x^*$ 的稳定流形 (Stable Manifold)为所有实部小于零的特征值对应的特征向量的线性组合, 记作 $E^s$ ; 不稳定流形 (Unstable Manifold)为所有实部大于零的特征值对应的特征向量的线性组合, 记作 $E^u$ ;以及中心流形(Centre Manifold)为所有实部等于零的特征值对应的特征向量的线性组合, 记作 $E^c$ . 对于非线性系统, 真正的稳定流形、不稳定流形和中心流形一般是弯曲的流形, 它们在平衡点处分别与 $E^s, E^u, E^c$ 相切.

为了引出接下来的定理, 我们先定义一种特殊的平衡点:

定义 6.8

若在平衡点 $\vec x^*$ 处的线性化矩阵 $\mathcal{J}$ 没有纯虚数特征值, 我们则称 $\vec x^*$ 为双曲平衡点 (Hyperbolic Equilibrium Point).

双曲平衡点的重要性在于:它排除了平衡点处的中心流形. 也就是说，线性化矩阵 $\mathcal{J}$ 的每一个特征方向都对应明确的指数增长或指数衰减(对应实部为正的特征值和实部为负的特征值). 因此,平衡点附近的动力学主要由线性化系统决定，而高阶非线性项只会改变轨线的具体形状，不会改变其基本拓扑结构. 一个自然的问题便是:非线性系统在双曲平衡点附近是否可以用其线性化系统来描述？

定理 6.11

设多元函数 $f(\vec x(t))$ 连续且可导, 非线性系统 $\vec x'(t) = f(\vec x(t))$ 在 $\vec x^*$ 处为双曲平衡点. 那么该系统在 $\vec x^*$ 附近与对应的线性化系统 $\vec x'(t) = \mathcal{J}\vec x(t)$ 在 $\vec x^*$ 附近局部拓扑共轭. 即二者有着相同的局部相图类型.

它告诉我们，在双曲平衡点附近，非线性系统的相图与线性化系统的相图在拓扑意义下相同. 因此我们对于这样的平衡点而言只需要研究其线性化系统即可.

例题 6.10

设非线性系统满足

\begin{pmatrix} y_1'(t) \\ y_2'(t) \end{pmatrix} = f(\vec y(t)) = \begin{pmatrix} y_1(1-y_1-y_2) \\ y_2(3-2y_1-4y_2) \end{pmatrix}.

求出该系统在平衡点附近的稳定性及相图类型.

解答 6.10

我们首先找出该系统的所有平衡点, 因此我们将求解

\begin{cases} y_1 - y_1^2 - y_1y_2 &= 0 \\ 3y_2 - 2y_1y_2 - 4y_2^2 &= 0\end{cases}.

在此我们分类讨论: 若 $y_1=0,y_2=0$ , 那么显然 $(0,0)$ 为一个稳定点; 其次我们设 $y_1=0,y_2\neq 0$ , 因此我们有 $3y_2-4y_2^2=0, y_2=\frac{3}{4}$ . 因此 $(0,\frac{3}{4})$ 也为一个平衡点; 我们随后设 $y_1\neq 0, y_2=0$ . 那么 $y_1(1-y_1)=0, y_1=1$ . 则 $(1,0)$ 也为一个平衡点. 最后我们设 $y_1\neq 0, y_2 \neq 0$ . 因此 $1-y_1-y_2=0, 3-2y_1-4y_2=0$ . 即 $y_1=y_2=\frac{1}{2}$ . 因此这个系统中一共有四个平衡点: $(0,0); (0,\frac{3}{4}); (1,0); (\frac{1}{2},\frac{1}{2})$ .

随后我们先求出该系统的Jacobi矩阵. 令 $f_1(y_1,y_2) = y_1(1-y_1-y_2)$ , $f_2(y_1,y_2) = y_2(3-2y_1-4y_2)$ . 那么

\mathcal{J} = \begin{pmatrix} \frac{\pl f_1}{\pl y_1} & \frac{\pl f_1}{\pl y_2} \\ \frac{\pl f_2}{\pl y_1} & \frac{\pl f_2}{\pl y_2} \end{pmatrix} = \begin{pmatrix} 1-2y_1-y_2 & -y_1 \\ -2y_2 & 3-2y_1-8y_2 \end{pmatrix}.

我们分别将 $4$ 个平衡点带入Jacobi矩阵:

\mathcal{J}(0,0) = \begin{pmatrix} 1 & 0 \\ 0 & 3 \end{pmatrix}.

此时的矩阵有两个特征值 $\lambda_1=1,\lambda_2=3$ . 因此根据Hartman-Grobman定理, 该系统在 $(0,0)$ 附近不稳定;

\mathcal{J}(0,\frac{3}{4}) = \begin{pmatrix} \frac{1}{4} & 0 \\ -\frac{3}{2} & -3 \end{pmatrix}.

此时的矩阵有两个特征值 $\lambda_1=\frac{1}{4}$ , $\lambda_2=-3$ . 因此根据Hartman-Grobman定理, 该系统在 $(0,\frac{3}{4})$ 处为不稳定的鞍点;

\mathcal{J}(1,0) = \begin{pmatrix} -1 & -1 \\ 0 & 1 \end{pmatrix}.

此时的矩阵有两个特征值 $\lambda_1=-1,\lambda_2=1$ . 因此根据Hartman-Grobman定理,该系统在 $(1,0)$ 处为不稳定的鞍点;

\mathcal{J}(\frac{1}{2},\frac{1}{2}) = \begin{pmatrix} -\frac{1}{2} & -\frac{1}{2} \\ -1 & -2 \end{pmatrix}.

此时的矩阵有两个负数特征值 $\lambda = \frac{1}{4}(-5\pm \sqrt{17})$ , 因此根据Hartman-Grobman定理, 该系统在 $(\frac{1}{2},\frac{1}{2})$ 处为渐近稳定的节点.

如果我们把整个系统在 $\mathbb{R}^2$ 内的大致相图画出来, 便会得到下图:

{6.5 练习}

1. 求出下面四个非线性系统中的平衡点; 当平衡点为双曲平衡点时判断平衡点的稳定性, 并画出大致相图.

\begin{aligned} &\text{(a)}.\begin{pmatrix} y_1'(t) \\ y_2'(t) \end{pmatrix} = \begin{pmatrix} y_2 \\ 8y_1-2y_2^3 \end{pmatrix} & \text{(b)}.\begin{pmatrix} y_1'(t) \\ y_2'(t) \end{pmatrix} = \begin{pmatrix} 2y_2(2-y_1) \\ y_2^2 - 2y_1+2 \end{pmatrix}\\ &\text{(c)}.\begin{pmatrix} y_1'(t) \\ y_2'(t) \end{pmatrix} = \begin{pmatrix} y_2 \\ -\sin(y_1) - y_2 \end{pmatrix} & \text{(d)}.\begin{pmatrix} y_1'(t) \\ y_2'(t) \end{pmatrix} = \begin{pmatrix}y_1(3y_2^2+y_1^2-1) \\ y_2(1-3y_1^2-y_2^2) \end{pmatrix} \end{aligned}

2. [Lotka-Volterra 方程组] 设 $y_1(t), y_2(t)$ 满足微分方程组

\begin{cases} \displaystyle{y_1'(t) = ay_1 - \frac{a}{K}y_1^2-by_1y_2} \\ \\\displaystyle{y_2'(t) = -cy_2+dy_1y_2} \end{cases},

其中 $a,b,c,d,K$ 均为大于零的常数.

(i) 求出该系统的所有平衡点, 并分析所有双曲平衡点的稳定性.

(ii) 如果将该方程视为生态系统中种群 $y_1, y_2$ 随时间的变化关系, 不同的 $a,b,c,d,K$ 的取值以及种群 $y_1,y_2$ 的初始数量会对系统产生什么样的影响?

3. [三体问题中的Lagrange点] 我们考虑一个简化版的三体问题: 假设存在两个天体 $M_1,M_2$ , 其质量分别为 $M_1=1-\mu, M_2=\mu, M_1>M_2>0$ . 在此基础上添加一质量不计的小型天体 $M_3$ . 我们建立旋转坐标系 $xyz$ , 使得 $M_1,M_2$ 的位置固定在点 $(-\mu,0,0), (1-\mu,0,0)$ 处. 我们设 $M_3$ 的轨迹方程为 $(x(t), y(t), z(t))=(x,y, z)$ , 其速度向量为 $(x',y',z')$ , 其加速度向量为 $(x'',y'',z'')$ . 定义

d=\sqrt{(x+\mu)^2+y^2+z^2},\quad r = \sqrt{(x-1+\mu)^2+y^2+z^2}.

那么根据物理学规律, 我们有如下的微分方程组:

\begin{cases} x'' &= 2y' + \Omega_x \\ y'' &= -2x' + \Omega_y \\ z'' &= \Omega_z \end{cases},\tag{6.4}

其中 $\Omega = \frac{1-\mu}{d} + \frac{\mu}{r} + \frac{x^2+y^2}{2}$ 为势能函数, $\Omega_x,\Omega_y,\Omega_z$ 为势能函数关于 $x,y,z$ 的偏导数. 若 $M_3$ 在某一点 $(x_0,y_0,z_0)$ 处的速度与加速度均为零(即 $x'=y'=z'=x''=y''=z''|_{x=x_0}=0$ ), 我们则称 $(x_0,y_0,z_0)$ 为系统的一个Lagrange点.

(i) 证明系统(6.4)中的所有Lagrange点的坐标均位于 $xy$ 平面内;

(ii) 求出两个满足 $y$ 坐标不为零的Lagrange点. 将这两个点记作 $L_4, L_5$ . 这两个点与 $M_1,M_2$ 的位置之间有着什么样的几何关系?

(iii) 当 $y=z=0$ 时, 系统同样存在另外 $3$ 个Lagrange点: $L_1,L_2,L_3$ . 写出该情况下 $x$ 满足的关系式 (此时 $x$ 没有通解, 我们只能借助Newton Raphson或Gradient Descent算法来求出数值解);

(iv) 设 $(x_0,y_0,0)$ 为 $L_i$ , $i=1,\cdots,5$ 的坐标. 为研究其稳定性, 我们在该点施加一微小扰动: $(x_0+\xi, y_0+\eta, 0)$ . 此时系统满足微分方程 $\vec y' = \vec A\vec y$ :

\begin{pmatrix} \xi' \\ \eta' \\ \xi'' \\ \eta'' \end{pmatrix}= \begin{pmatrix} 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ \Omega_{xx} & \Omega_{xy} & 0 & 2 \\ \Omega_{yx} & \Omega_{yy} & -2 & 0 \end{pmatrix} \begin{pmatrix} \xi \\ \eta \\ \xi' \\ \eta' \end{pmatrix}.

• 设 $i=1,2,3$ :

(a) 通过判断 $\vec A$ 的特征值实部的符号, 解释为什么 $L_1,L_2,L_3$ 为不稳定的节点.

• 设 $i=4,5$ :

(a) 通过求出 $\vec A$ 的特征值, 解释为什么 $L_4,L_5$ 为稳定的节点.
(b) 当 $\vec A$ 的特征值为纯虚数时, 证明 $M_1,M_2$ 的质量关系大约满足 $M_2 \leq 0.04M_1$ .

太阳--木星系统提供了 $L_4, L_5$ 稳定性的一个重要实际例子.由于木星质量远小于太阳质量, 该系统满足 $L_4,L_5$ 的稳定条件. 因此在木星轨道前方和后方约 $60^\circ$ 的区域便有大量小行星在此相对静止, 随木星一起绕太阳运动. 这些小行星被称为 Trojan 小行星(特洛伊群).