常微分方程（十）分岔理论
Chen Kai BOSS
  2019-05-24 09:15:00 2019-05-24 09:15      6.4k 字  29 分钟  

为什么湖泊会突然从清澈变成浑浊？为什么股市会毫无征兆地崩盘？为什么气候系统似乎在某些"临界点"会发生剧变？这些看似无关的现象背后，隐藏着一个共同的数学原理——分岔。当系统参数缓慢变化时，系统的定性行为可能突然改变：稳定的平衡点消失、周期解出现、混沌爆发。分岔理论正是研究这种"质变"的数学工具，它让我们能够预测和理解自然界中各种突变现象。

什么是分岔？

从生活中的例子说起

想象你在煮水。当温度低于 100 ° C 时，水保持液态；当温度达到 100 ° C，水开始沸腾变成气态。这就是一种"分岔"——系统的定性状态在临界参数值处发生根本改变。

再想象一座桥。当承重在设计范围内时，桥是稳定的；但当重量超过临界值，桥会突然坍塌。这种从稳定到不稳定的转变，正是分岔的体现。

数学定义

分岔（ bifurcation）是指：当系统的参数越过某个临界值时，系统的拓扑结构发生定性改变。

考虑参数化的动力系统：

其中 是参数。当 变化时： - 平衡点的数量可能改变 - 平衡点的稳定性可能改变 - 周期解可能出现或消失 - 混沌可能产生

发生这种改变的参数值 叫做分岔值或临界点。

分岔图：全景可视化

分岔图（ bifurcation diagram）是理解分岔最直观的工具。它以参数 为横轴，以系统的稳态（平衡点、周期解的振幅等）为纵轴，展示系统行为如何随参数变化。

一维系统的分岔

一维系统 是研究分岔的最佳起点。虽然简单，它已经展示了分岔的核心特征。

鞍结点分岔（ Saddle-Node Bifurcation）

最基本的分岔类型——两个平衡点"碰撞"并消失。

标准形式：

分析： - 当 ：无实平衡点（ 无解） - 当 ：一个平衡点 （半稳定） - 当 ：两个平衡点 稳定性分析： - $ x = ： f'() = -2 < 0（稳定） x = -： f'(-) = 2 > 0$（不稳定）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
import numpy as np
import matplotlib.pyplot as plt

def saddle_node_bifurcation():
    """鞍结点分岔的可视化"""
    
    fig, axes = plt.subplots(1, 2, figsize=(14, 5))
    
    # 左图：不同 mu 值下的相线
    ax1 = axes[0]
    x = np.linspace(-2, 2, 500)
    
    mu_values = [-0.5, 0, 0.5, 1.0]
    colors = ['blue', 'green', 'orange', 'red']
    
    for mu, color in zip(mu_values, colors):
        f = mu - x**2
        ax1.plot(x, f, color=color, linewidth=2, label=f'μ = {mu}')
        
        # 标记平衡点
        if mu >= 0:
            x_eq = np.sqrt(mu)
            ax1.plot(x_eq, 0, 'o', color=color, markersize=10)
            ax1.plot(-x_eq, 0, 's', color=color, markersize=10)
    
    ax1.axhline(y=0, color='k', linestyle='-', linewidth=0.5)
    ax1.set_xlabel('x', fontsize=12)
    ax1.set_ylabel('dx/dt = μ - x ²', fontsize=12)
    ax1.set_title('Phase Lines for Different μ', fontsize=13, fontweight='bold')
    ax1.legend(fontsize=10)
    ax1.grid(True, alpha=0.3)
    ax1.set_xlim(-2, 2)
    ax1.set_ylim(-2, 2)
    
    # 右图：分岔图
    ax2 = axes[1]
    mu_range = np.linspace(-1, 2, 300)
    
    # 稳定分支 (x = sqrt(mu))
    mu_pos = mu_range[mu_range >= 0]
    x_stable = np.sqrt(mu_pos)
    ax2.plot(mu_pos, x_stable, 'b-', linewidth=2.5, label='Stable')
    ax2.plot(mu_pos, -x_stable, 'r--', linewidth=2.5, label='Unstable')
    
    # 分岔点
    ax2.plot(0, 0, 'ko', markersize=12, zorder=5)
    ax2.annotate('Bifurcation point', xy=(0, 0), xytext=(0.5, 0.8),
                fontsize=11, arrowprops=dict(arrowstyle='->', color='black'))
    
    ax2.axvline(x=0, color='gray', linestyle=':', alpha=0.7)
    ax2.set_xlabel('μ', fontsize=12)
    ax2.set_ylabel('x* (equilibria)', fontsize=12)
    ax2.set_title('Saddle-Node Bifurcation Diagram', fontsize=13, fontweight='bold')
    ax2.legend(fontsize=11)
    ax2.grid(True, alpha=0.3)
    ax2.set_xlim(-1, 2)
    ax2.set_ylim(-2, 2)
    
    plt.tight_layout()
    plt.savefig('saddle_node_bifurcation.png', dpi=150, bbox_inches='tight')
    plt.close()

saddle_node_bifurcation()

物理意义：鞍结点分岔描述了"阈值现象"——当参数低于临界值时没有平衡态，越过临界值后突然出现。

生活例子： - 激光器：低于阈值电流无激光输出，超过阈值突然发光 - 神经元：低于阈值电压不放电，超过阈值产生动作电位

跨临界分岔（ Transcritical Bifurcation）

两个平衡点交换稳定性，但都不消失。

标准形式：

分析： - 平衡点始终是 和 - 稳定性： - 当 ： 稳定， 不稳定 - 当 ： 不稳定， 稳定

在 处，两个平衡点"交换"稳定性。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
import numpy as np
import matplotlib.pyplot as plt

def transcritical_bifurcation():
    """跨临界分岔的可视化"""
    
    fig, ax = plt.subplots(figsize=(10, 7))
    
    mu_range = np.linspace(-2, 2, 400)
    
    # x* = 0 分支
    ax.plot(mu_range[mu_range <= 0], np.zeros(sum(mu_range <= 0)), 
            'b-', linewidth=2.5, label='x*=0 (stable when μ<0)')
    ax.plot(mu_range[mu_range >= 0], np.zeros(sum(mu_range >= 0)), 
            'b--', linewidth=2.5, label='x*=0 (unstable when μ>0)')
    
    # x* = μ 分支
    ax.plot(mu_range[mu_range <= 0], mu_range[mu_range <= 0], 
            'r--', linewidth=2.5, label='x*=μ (unstable when μ<0)')
    ax.plot(mu_range[mu_range >= 0], mu_range[mu_range >= 0], 
            'r-', linewidth=2.5, label='x*=μ (stable when μ>0)')
    
    # 分岔点
    ax.plot(0, 0, 'ko', markersize=15, zorder=5)
    ax.annotate('Transcritical\nbifurcation', xy=(0, 0), xytext=(0.8, -1.2),
                fontsize=11, arrowprops=dict(arrowstyle='->', color='black'))
    
    ax.set_xlabel('μ', fontsize=12)
    ax.set_ylabel('x*', fontsize=12)
    ax.set_title('Transcritical Bifurcation: Exchange of Stability', 
                fontsize=13, fontweight='bold')
    ax.legend(fontsize=10, loc='upper left')
    ax.grid(True, alpha=0.3)
    ax.set_xlim(-2, 2)
    ax.set_ylim(-2, 2)
    ax.axhline(y=0, color='gray', linestyle='-', linewidth=0.5)
    ax.axvline(x=0, color='gray', linestyle='-', linewidth=0.5)
    
    plt.tight_layout()
    plt.savefig('transcritical_bifurcation.png', dpi=150, bbox_inches='tight')
    plt.close()

transcritical_bifurcation()

物理意义：跨临界分岔常见于带有"零状态"的系统，零状态始终是平衡点，但其稳定性随参数改变。

生活例子： - 传染病模型： 是临界点 - ：无病平衡态稳定（疾病消失） - ：无病平衡态不稳定（疾病爆发） - 种群增长：低于某环境质量物种灭绝，高于则生存

叉形分岔（ Pitchfork Bifurcation）

一个平衡点分裂成三个——对称系统的典型分岔。

超临界叉形分岔（ supercritical）：

分析： - 当 ：只有 （稳定） - 当 ：三个平衡点 $ x = 0（不稳定）、 x = $（稳定）

亚临界叉形分岔（ subcritical）： - 当 ：只有 （不稳定） - 当 ：三个平衡点 $ x = 0（稳定）、 x = $（不稳定）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
import numpy as np
import matplotlib.pyplot as plt

def pitchfork_bifurcations():
    """超临界和亚临界叉形分岔"""
    
    fig, axes = plt.subplots(1, 2, figsize=(14, 6))
    
    mu_range = np.linspace(-2, 2, 400)
    mu_pos = mu_range[mu_range >= 0]
    mu_neg = mu_range[mu_range <= 0]
    
    # 左图：超临界
    ax1 = axes[0]
    
    # x = 0 分支
    ax1.plot(mu_range[mu_range <= 0], np.zeros(sum(mu_range <= 0)), 
             'b-', linewidth=2.5)
    ax1.plot(mu_range[mu_range >= 0], np.zeros(sum(mu_range >= 0)), 
             'r--', linewidth=2.5)
    
    # x = ± sqrt(mu) 分支
    ax1.plot(mu_pos, np.sqrt(mu_pos), 'b-', linewidth=2.5)
    ax1.plot(mu_pos, -np.sqrt(mu_pos), 'b-', linewidth=2.5)
    
    ax1.plot(0, 0, 'ko', markersize=12)
    ax1.set_xlabel('μ', fontsize=12)
    ax1.set_ylabel('x*', fontsize=12)
    ax1.set_title('Supercritical Pitchfork Bifurcation\n(dx/dt = μ x - x ³)', 
                  fontsize=12, fontweight='bold')
    ax1.grid(True, alpha=0.3)
    ax1.set_xlim(-2, 2)
    ax1.set_ylim(-2, 2)
    ax1.legend(['Stable', 'Unstable'], fontsize=10, loc='upper left')
    
    # 右图：亚临界
    ax2 = axes[1]
    
    # x = 0 分支
    ax2.plot(mu_range[mu_range >= 0], np.zeros(sum(mu_range >= 0)), 
             'r--', linewidth=2.5)
    ax2.plot(mu_range[mu_range <= 0], np.zeros(sum(mu_range <= 0)), 
             'b-', linewidth=2.5)
    
    # x = ± sqrt(-mu) 分支
    ax2.plot(mu_neg, np.sqrt(-mu_neg), 'r--', linewidth=2.5)
    ax2.plot(mu_neg, -np.sqrt(-mu_neg), 'r--', linewidth=2.5)
    
    ax2.plot(0, 0, 'ko', markersize=12)
    ax2.set_xlabel('μ', fontsize=12)
    ax2.set_ylabel('x*', fontsize=12)
    ax2.set_title('Subcritical Pitchfork Bifurcation\n(dx/dt = μ x + x ³)', 
                  fontsize=12, fontweight='bold')
    ax2.grid(True, alpha=0.3)
    ax2.set_xlim(-2, 2)
    ax2.set_ylim(-2, 2)
    ax2.legend(['Unstable', 'Stable'], fontsize=10, loc='upper left')
    
    plt.tight_layout()
    plt.savefig('pitchfork_bifurcations.png', dpi=150, bbox_inches='tight')
    plt.close()

pitchfork_bifurcations()

物理意义： - 超临界：分岔后系统平稳过渡到新状态（"软"分岔） - 亚临界：分岔前系统可能突然跳跃到远处状态（"硬"分岔，危险！）

生活例子： - 欧拉压杆失稳：垂直细杆受压，超过临界载荷后弯曲（左或右，对称性破缺） - 激光：超过阈值后，光场的相位打破对称性

分岔的分类小结

分岔类型 标准形式 特征 危险程度
鞍结点 平衡点创生/消失 中等
跨临界 稳定性交换 低
超临界叉形 对称分裂，软过渡 低
亚临界叉形 对称分裂，硬跳跃 高

二维系统的分岔

二维系统可以出现更丰富的分岔，特别是与周期解相关的分岔。

Hopf 分岔：周期解的诞生

Hopf 分岔是最重要的二维分岔——稳定焦点变成不稳定焦点，同时产生极限环。

标准形式（极坐标）：

分析 方程（与超临界叉形分岔相同）： - 当 ： 稳定（焦点） - 当 ： 不稳定， 稳定（极限环）

极限环的半径： 直角坐标形式： 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

def hopf_system(state, t, mu, omega):
    """Hopf 分岔标准形式"""
    x, y = state
    r_sq = x**2 + y**2
    dxdt = mu * x - omega * y - x * r_sq
    dydt = omega * x + mu * y - y * r_sq
    return [dxdt, dydt]

fig, axes = plt.subplots(1, 3, figsize=(15, 5))
omega = 1.0
t = np.linspace(0, 50, 5000)

mu_values = [-0.5, 0, 0.5]
titles = ['μ = -0.5 (Stable Focus)', 'μ = 0 (Bifurcation Point)', 
          'μ = 0.5 (Limit Cycle)']

for ax, mu, title in zip(axes, mu_values, titles):
    # 多个初始条件
    initial_conditions = [(0.1, 0), (0.5, 0.5), (1.5, 0), (0, 1.5), (-1, -1)]
    
    for x0, y0 in initial_conditions:
        sol = odeint(hopf_system, [x0, y0], t, args=(mu, omega))
        ax.plot(sol[:, 0], sol[:, 1], linewidth=0.7, alpha=0.7)
        ax.plot(x0, y0, 'go', markersize=5)
    
    # 如果有极限环，画出来
    if mu > 0:
        theta = np.linspace(0, 2*np.pi, 100)
        r = np.sqrt(mu)
        ax.plot(r * np.cos(theta), r * np.sin(theta), 'r-', 
                linewidth=2.5, label=f'Limit cycle (r = {r:.2f})')
        ax.legend(fontsize=9)
    
    ax.plot(0, 0, 'ko', markersize=8)
    ax.set_xlabel('x', fontsize=11)
    ax.set_ylabel('y', fontsize=11)
    ax.set_title(title, fontsize=12, fontweight='bold')
    ax.grid(True, alpha=0.3)
    ax.set_xlim(-2, 2)
    ax.set_ylim(-2, 2)
    ax.set_aspect('equal')

plt.tight_layout()
plt.savefig('hopf_bifurcation.png', dpi=150, bbox_inches='tight')
plt.close()

超临界 vs 亚临界 Hopf 分岔

超临界 Hopf（如上例）： - 分岔后极限环半径从零连续增长 - 系统平稳过渡，相对安全

亚临界 Hopf： - 分岔前就存在不稳定极限环 - 分岔后系统可能跳跃到远处 - 存在滞后（ hysteresis）现象

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
import numpy as np
import matplotlib.pyplot as plt

def hopf_bifurcation_diagram():
    """Hopf 分岔图对比"""
    
    fig, axes = plt.subplots(1, 2, figsize=(14, 5))
    
    mu_range = np.linspace(-1, 1, 200)
    mu_pos = mu_range[mu_range >= 0]
    mu_neg = mu_range[mu_range <= 0]
    
    # 左图：超临界
    ax1 = axes[0]
    
    # 原点稳定性
    ax1.plot(mu_neg, np.zeros(len(mu_neg)), 'b-', linewidth=2.5, label='Stable equilibrium')
    ax1.plot(mu_pos, np.zeros(len(mu_pos)), 'r--', linewidth=2.5, label='Unstable equilibrium')
    
    # 极限环
    ax1.plot(mu_pos, np.sqrt(mu_pos), 'b-', linewidth=2.5, label='Stable limit cycle')
    ax1.plot(mu_pos, -np.sqrt(mu_pos), 'b-', linewidth=2.5)
    
    ax1.fill_between(mu_pos, -np.sqrt(mu_pos), np.sqrt(mu_pos), alpha=0.2, color='blue')
    
    ax1.plot(0, 0, 'ko', markersize=12)
    ax1.set_xlabel('μ', fontsize=12)
    ax1.set_ylabel('Amplitude', fontsize=12)
    ax1.set_title('Supercritical Hopf Bifurcation', fontsize=13, fontweight='bold')
    ax1.legend(fontsize=9, loc='upper left')
    ax1.grid(True, alpha=0.3)
    ax1.set_xlim(-1, 1)
    ax1.set_ylim(-1.5, 1.5)
    
    # 右图：亚临界（带滞后）
    ax2 = axes[1]
    
    # 原点稳定性
    ax2.plot(mu_neg, np.zeros(len(mu_neg)), 'b-', linewidth=2.5)
    ax2.plot(mu_pos, np.zeros(len(mu_pos)), 'r--', linewidth=2.5)
    
    # 不稳定极限环（存在于 mu<0 区域）
    mu_unstable = mu_range[(mu_range >= -0.5) & (mu_range <= 0)]
    r_unstable = np.sqrt(-mu_unstable + 0.5)
    ax2.plot(mu_unstable, r_unstable, 'r--', linewidth=2.5, label='Unstable limit cycle')
    ax2.plot(mu_unstable, -r_unstable, 'r--', linewidth=2.5)
    
    # 稳定极限环（大振幅）
    mu_stable = mu_range[mu_range >= -0.5]
    r_stable = np.sqrt(mu_stable + 0.5)
    ax2.plot(mu_stable, r_stable, 'b-', linewidth=2.5, label='Stable limit cycle')
    ax2.plot(mu_stable, -r_stable, 'b-', linewidth=2.5)
    
    # 滞后区域
    ax2.fill_between(mu_neg[mu_neg >= -0.5], 
                     -np.sqrt(-mu_neg[mu_neg >= -0.5] + 0.5),
                     np.sqrt(-mu_neg[mu_neg >= -0.5] + 0.5),
                     alpha=0.2, color='red', label='Bistable region')
    
    ax2.plot(0, 0, 'ko', markersize=12)
    ax2.plot(-0.5, 0, 'ko', markersize=12)
    ax2.set_xlabel('μ', fontsize=12)
    ax2.set_ylabel('Amplitude', fontsize=12)
    ax2.set_title('Subcritical Hopf Bifurcation\n(with hysteresis)', 
                  fontsize=13, fontweight='bold')
    ax2.legend(fontsize=9, loc='upper left')
    ax2.grid(True, alpha=0.3)
    ax2.set_xlim(-1, 1)
    ax2.set_ylim(-1.5, 1.5)
    
    plt.tight_layout()
    plt.savefig('hopf_comparison.png', dpi=150, bbox_inches='tight')
    plt.close()

hopf_bifurcation_diagram()

Hopf 分岔定理

定理（ Hopf, 1942）：设 是系统 在 处的平衡点，雅可比矩阵有一对复共轭特征值 。如果：

  1. （实部过零）
  2. （虚部非零）
  3. （横截条件）

则在 处发生Hopf 分岔，存在周期解分支，周期约为 。

全局分岔

前面讨论的都是局部分岔——发生在平衡点或极限环附近。全局分岔涉及相空间中更大范围的结构改变。

同宿轨道分岔

同宿轨道（ homoclinic orbit）是从鞍点出发又回到同一鞍点的轨道。

当参数变化使同宿轨道出现或消失时，发生同宿分岔。

这种分岔的特殊之处： - 轨道的周期趋向无穷大（因为在鞍点附近停留无穷长时间） - 可能导致混沌（ Shilnikov 定理）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

def homoclinic_example(state, t, mu):
    """展示同宿轨道的示例系统"""
    x, y = state
    dxdt = y
    dydt = x - x**3 + mu * y  # 阻尼可调
    return [dxdt, dydt]

fig, axes = plt.subplots(1, 3, figsize=(15, 5))
t = np.linspace(0, 50, 5000)

# 不同 mu 值
mu_values = [-0.3, -0.05, 0.3]
titles = ['μ = -0.3: Stable limit cycles', 
          'μ ≈ 0: Near homoclinic bifurcation',
          'μ = 0.3: Trajectories escape to infinity']

for ax, mu, title in zip(axes, mu_values, titles):
    # 平衡点
    ax.plot(0, 0, 'rs', markersize=10, label='Saddle at origin')
    ax.plot(1, 0, 'bo', markersize=8, label='Stable/unstable node')
    ax.plot(-1, 0, 'bo', markersize=8)
    
    # 多条轨迹
    initial_conditions = [
        (0.5, 0.2), (-0.5, 0.2), (0.5, -0.2), (-0.5, -0.2),
        (0.1, 0.5), (-0.1, 0.5), (0.1, -0.5), (-0.1, -0.5),
        (1.5, 0), (-1.5, 0)
    ]
    
    for x0, y0 in initial_conditions:
        try:
            sol = odeint(homoclinic_example, [x0, y0], t, args=(mu,))
            # 只画有限部分
            mask = (np.abs(sol[:, 0]) < 3) & (np.abs(sol[:, 1]) < 3)
            if np.any(mask):
                ax.plot(sol[mask, 0], sol[mask, 1], 'g-', linewidth=0.7, alpha=0.7)
        except:
            pass
        ax.plot(x0, y0, 'ko', markersize=4)
    
    ax.set_xlabel('x', fontsize=11)
    ax.set_ylabel('y', fontsize=11)
    ax.set_title(title, fontsize=11, fontweight='bold')
    ax.grid(True, alpha=0.3)
    ax.set_xlim(-2, 2)
    ax.set_ylim(-1.5, 1.5)
    ax.legend(fontsize=8, loc='upper right')

plt.tight_layout()
plt.savefig('homoclinic_bifurcation.png', dpi=150, bbox_inches='tight')
plt.close()

异宿轨道分岔

异宿轨道（ heteroclinic orbit）连接两个不同的鞍点。

当存在多个鞍点时，它们之间的连接关系可能随参数改变，导致异宿分岔。

SNIC 分岔

SNIC（ Saddle-Node on Invariant Circle）——鞍结点分岔发生在极限环上。

结果：极限环通过一个鞍结点"打开"成为不变曲线。

周期趋向无穷大的特点：周期 发散。

余维数与普适性

余维数的概念

余维数（ codimension）是指使分岔发生所需调节的独立参数的最小数目。

  • 余维数 1 分岔：调节一个参数即可发生（鞍结点、跨临界、 Hopf 等）
  • 余维数 2 分岔：需要同时调节两个参数（ Bogdanov-Takens 、尖点等）

余维数越高，分岔越"不稳定"——参数空间中越难遇到。

普适展开与标准形式

标准形式（ normal form）是分岔的"最简代表"。

通过坐标变换，任何满足某种分岔条件的系统都可以化为标准形式（在分岔点附近）。

这就是分岔理论的强大之处：虽然具体系统千差万别，但分岔的类型是有限的！

分岔在自然科学中的应用

生态学：物种灭绝与恢复

考虑带 Allee 效应的种群模型：

其中 是 Allee 阈值——低于此值种群会走向灭绝。

当环境恶化（ 减小）时，可能发生鞍结点分岔：正平衡点消失，种群突然崩溃！

这解释了为什么渔业资源可能"突然"枯竭——不是线性衰减，而是越过临界点后崩溃。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
import numpy as np
import matplotlib.pyplot as plt

def allee_effect_bifurcation():
    """Allee 效应导致的分岔"""
    
    fig, axes = plt.subplots(1, 2, figsize=(14, 5))
    
    # 参数
    r = 1.0
    A = 10.0  # Allee 阈值
    
    K_values = [50, 30, 20, 15]  # 不同承载力
    colors = ['green', 'blue', 'orange', 'red']
    
    N = np.linspace(0, 60, 500)
    
    # 左图：相图
    ax1 = axes[0]
    
    for K, color in zip(K_values, colors):
        dNdt = r * N * (1 - N/K) * (N/A - 1)
        ax1.plot(N, dNdt, color=color, linewidth=2, label=f'K = {K}')
        
        # 标记平衡点
        from scipy.optimize import brentq
        try:
            # 找正平衡点
            if K > A:
                N_stable = brentq(lambda n: r*n*(1-n/K)*(n/A-1), A+1, K-1)
                N_unstable = brentq(lambda n: r*n*(1-n/K)*(n/A-1), 1, A-1)
                ax1.plot(N_stable, 0, 'o', color=color, markersize=8)
                ax1.plot(N_unstable, 0, 's', color=color, markersize=8)
        except:
            pass
    
    ax1.axhline(y=0, color='k', linestyle='-', linewidth=0.5)
    ax1.set_xlabel('Population N', fontsize=12)
    ax1.set_ylabel('dN/dt', fontsize=12)
    ax1.set_title('Allee Effect: Phase Lines for Different K', 
                  fontsize=13, fontweight='bold')
    ax1.legend(fontsize=10)
    ax1.grid(True, alpha=0.3)
    ax1.set_xlim(0, 60)
    ax1.set_ylim(-5, 10)
    
    # 右图：分岔图
    ax2 = axes[1]
    
    K_range = np.linspace(10, 60, 200)
    
    # 稳定高密度分支
    N_stable = []
    N_unstable = []
    K_valid_s = []
    K_valid_u = []
    
    for K in K_range:
        # 稳定：介于 A 和 K 之间的根
        if K > A:
            try:
                ns = brentq(lambda n: r*n*(1-n/K)*(n/A-1), A+0.1, K-0.1)
                N_stable.append(ns)
                K_valid_s.append(K)
            except:
                pass
            try:
                nu = brentq(lambda n: r*n*(1-n/K)*(n/A-1), 0.1, A-0.1)
                N_unstable.append(nu)
                K_valid_u.append(K)
            except:
                pass
    
    ax2.plot(K_valid_s, N_stable, 'b-', linewidth=2.5, label='Stable (high density)')
    ax2.plot(K_valid_u, N_unstable, 'r--', linewidth=2.5, label='Unstable (Allee threshold)')
    ax2.plot(K_range, np.zeros(len(K_range)), 'b-', linewidth=2.5, label='Stable (extinction)')
    
    # 分岔点
    ax2.axvline(x=A, color='gray', linestyle=':', linewidth=2)
    ax2.text(A+1, 50, 'Saddle-node\nbifurcation', fontsize=10)
    
    ax2.set_xlabel('Carrying capacity K', fontsize=12)
    ax2.set_ylabel('Equilibrium population N*', fontsize=12)
    ax2.set_title('Bifurcation Diagram: Population Collapse', 
                  fontsize=13, fontweight='bold')
    ax2.legend(fontsize=10)
    ax2.grid(True, alpha=0.3)
    ax2.set_xlim(10, 60)
    ax2.set_ylim(-5, 60)
    
    plt.tight_layout()
    plt.savefig('allee_bifurcation.png', dpi=150, bbox_inches='tight')
    plt.close()

allee_effect_bifurcation()

气候科学：临界点与不可逆变化

气候系统可能存在多个稳态（如"温室地球"和"冰室地球"）。

当温室气体浓度越过临界点时，可能发生不可逆的状态转变——这就是气候科学家警告的"临界点"（ tipping points）。

例子： - 北极海冰：反射率-温度正反馈可能导致突然消失 - 亚马逊雨林：干旱-火灾-森林退化可能导致稀树草原化 - 大西洋环流：淡水注入可能导致环流崩溃

这些都是鞍结点分岔或亚临界 Hopf 分岔的表现。

神经科学：神经元放电模式

神经元可以在"静息"和"放电"两种状态之间切换。

Hodgkin-Huxley 模型（简化版 FitzHugh-Nagumo）展示了多种分岔： - 鞍结点分岔： Class I 神经元（频率可从零连续调节） - Hopf 分岔： Class II 神经元（频率有下限）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

def fitzhugh_nagumo(state, t, I, a, b, tau):
    """FitzHugh-Nagumo 模型"""
    v, w = state
    dvdt = v - v**3/3 - w + I
    dwdt = (v + a - b*w) / tau
    return [dvdt, dwdt]

# 不同外部电流下的行为
I_values = [0, 0.3, 0.5, 0.8]
a, b, tau = 0.7, 0.8, 12.5

fig, axes = plt.subplots(2, 2, figsize=(12, 10))
t = np.linspace(0, 200, 4000)

for ax, I in zip(axes.flatten(), I_values):
    # 相图
    v_range = np.linspace(-2.5, 2.5, 20)
    w_range = np.linspace(-1.5, 2, 20)
    V, W = np.meshgrid(v_range, w_range)
    dV = V - V**3/3 - W + I
    dW = (V + a - b*W) / tau
    
    # 向量场
    M = np.sqrt(dV**2 + dW**2)
    M[M == 0] = 1
    ax.quiver(V, W, dV/M, dW/M, M, cmap='Greys', alpha=0.3)
    
    # 零等倾线
    v_null = np.linspace(-2.5, 2.5, 100)
    w_vnull = v_null - v_null**3/3 + I
    w_wnull = (v_null + a) / b
    
    ax.plot(v_null, w_vnull, 'b-', linewidth=2, label='v-nullcline')
    ax.plot(v_null, w_wnull, 'r-', linewidth=2, label='w-nullcline')
    
    # 轨迹
    sol = odeint(fitzhugh_nagumo, [-1, 0], t, args=(I, a, b, tau))
    ax.plot(sol[:, 0], sol[:, 1], 'g-', linewidth=1.5, alpha=0.7)
    ax.plot(sol[0, 0], sol[0, 1], 'go', markersize=8)
    
    ax.set_xlabel('v (membrane potential)', fontsize=11)
    ax.set_ylabel('w (recovery variable)', fontsize=11)
    ax.set_title(f'I = {I}', fontsize=12, fontweight='bold')
    ax.legend(fontsize=9, loc='upper right')
    ax.grid(True, alpha=0.3)
    ax.set_xlim(-2.5, 2.5)
    ax.set_ylim(-1.5, 2)

plt.suptitle('FitzHugh-Nagumo Model: Hopf Bifurcation with Increasing Current', 
             fontsize=14, fontweight='bold', y=1.02)
plt.tight_layout()
plt.savefig('fhn_bifurcation.png', dpi=150, bbox_inches='tight')
plt.close()

工程学：结构稳定性

欧拉失稳是经典的叉形分岔：

细长杆受轴向压力 ： - ：直立状态稳定 - ：分岔点（临界载荷） - ：直立不稳定，弯曲状态稳定

临界载荷 （: 弹性模量，: 截面惯性矩，: 长度）

这是超临界叉形分岔——杆会平稳地弯向一侧。

经济学：市场崩溃

金融市场可能存在多个"均衡"： - 高信心、高投资的繁荣状态 - 低信心、低投资的萧条状态

当某些指标越过临界值时，可能发生鞍结点分岔——市场从繁荣突然跳到萧条（崩盘）。

这解释了为什么市场崩溃往往是突然的，而不是渐进的。

周期倍增与通往混沌之路

周期倍增级联

当参数变化时，周期解可能经历周期倍增分岔： - 周期-1 → 周期-2 → 周期-4 → 周期-8 → ... → 混沌

Feigenbaum 常数：

这个常数是普适的——无论具体系统如何，只要经历周期倍增，比值都趋向这个常数！

Logistic 映射：混沌的简单模型

Logistic 映射： $$

x_{n+1} = rx_n(1 - x_n) $$

这是研究混沌最简单的模型。随着 增大：

范围 行为
所有轨道趋向 0
趋向非零不动点
周期-2
周期-4
... 周期倍增级联
混沌（有周期窗口） 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
import numpy as np
import matplotlib.pyplot as plt

def logistic_bifurcation_diagram():
    """Logistic 映射的分岔图"""
    
    r_values = np.linspace(2.5, 4.0, 2000)
    iterations = 1000
    last_points = 100
    
    fig, ax = plt.subplots(figsize=(14, 8))
    
    for r in r_values:
        x = 0.5  # 初始值
        
        # 迭代到稳态
        for _ in range(iterations - last_points):
            x = r * x * (1 - x)
        
        # 收集最后若干点
        x_values = []
        for _ in range(last_points):
            x = r * x * (1 - x)
            x_values.append(x)
        
        ax.scatter([r] * len(x_values), x_values, 
                   c='blue', s=0.1, alpha=0.5)
    
    ax.set_xlabel('r', fontsize=12)
    ax.set_ylabel('x*', fontsize=12)
    ax.set_title('Bifurcation Diagram of Logistic Map: $ x_{n+1} = rx_n(1-x_n)$', 
                 fontsize=14, fontweight='bold')
    ax.grid(True, alpha=0.3)
    
    # 标注关键点
    ax.axvline(x=3, color='red', linestyle='--', alpha=0.5)
    ax.text(3.02, 0.1, 'Period-2\nbifurcation', fontsize=10, color='red')
    
    ax.axvline(x=3.5699, color='green', linestyle='--', alpha=0.5)
    ax.text(3.59, 0.1, 'Onset of\nchaos', fontsize=10, color='green')
    
    plt.tight_layout()
    plt.savefig('logistic_bifurcation.png', dpi=150, bbox_inches='tight')
    plt.close()

logistic_bifurcation_diagram()

周期窗口与危机

即使在混沌区域，也有周期窗口——参数的某些区间内系统回到周期运动。

最著名的是 附近的周期-3 窗口。

Li-Yorke 定理："周期 3 蕴含混沌"——如果一维映射有周期-3 轨道，则有所有周期的轨道，且存在不可数多条混沌轨道。

数值方法：延拓与检测

延拓方法

如何数值地"跟踪"平衡点或周期解随参数的变化？

伪弧长延拓（ pseudo-arclength continuation）： 1. 沿着解曲线定义弧长参数 2. 在弧长方向上前进，同时调整参数和状态 3. 可以追踪"折回"的分支（鞍结点处）

常用软件： - AUTO（经典工具） - MATCONT（ MATLAB 工具箱） - PyDSTool（ Python） - XPPAut

分岔检测

检测鞍结点分岔：监测雅可比矩阵的特征值，当有特征值过零时发出警报。

检测 Hopf 分岔：监测特征值的实部，当一对共轭复特征值实部过零时发生 Hopf 分岔。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
import numpy as np
from scipy.linalg import eigvals

def detect_bifurcations(jacobian_func, param_range, x_eq_func):
    """
    沿参数路径检测分岔
    
    jacobian_func: 计算雅可比矩阵的函数 J(x, mu)
    param_range: 参数值序列
    x_eq_func: 计算平衡点的函数 x*(mu)
    """
    
    bifurcations = []
    prev_eigenvalues = None
    
    for mu in param_range:
        x_eq = x_eq_func(mu)
        J = jacobian_func(x_eq, mu)
        eigenvalues = eigvals(J)
        
        if prev_eigenvalues is not None:
            # 检测实特征值过零（鞍结点、跨临界、叉形）
            prev_real = np.sort(np.real(prev_eigenvalues))
            curr_real = np.sort(np.real(eigenvalues))
            
            for pr, cr in zip(prev_real, curr_real):
                if pr * cr < 0:  # 符号改变
                    bifurcations.append({
                        'type': 'Real eigenvalue crossing',
                        'mu': mu,
                        'eigenvalues': eigenvalues
                    })
                    break
            
            # 检测 Hopf 分岔（复特征值实部过零）
            complex_eigs = eigenvalues[np.abs(np.imag(eigenvalues)) > 1e-10]
            if len(complex_eigs) > 0:
                prev_complex = prev_eigenvalues[np.abs(np.imag(prev_eigenvalues)) > 1e-10]
                if len(prev_complex) > 0:
                    pr = np.real(prev_complex[0])
                    cr = np.real(complex_eigs[0])
                    if pr * cr < 0:
                        bifurcations.append({
                            'type': 'Hopf bifurcation',
                            'mu': mu,
                            'eigenvalues': eigenvalues
                        })
        
        prev_eigenvalues = eigenvalues
    
    return bifurcations

总结

本章我们系统学习了分岔理论的核心内容：

  1. 一维分岔：
    • 鞍结点分岔（平衡点创生/消失）
    • 跨临界分岔（稳定性交换）
    • 叉形分岔（对称性破缺）
  2. 二维分岔：
    • Hopf 分岔（周期解诞生）
    • 全局分岔（同宿、异宿）
  3. 余维数与普适性：
    • 分岔的分类取决于余维数
    • 标准形式提供统一描述
  4. 应用：
    • 生态学（物种灭绝）
    • 气候科学（临界点）
    • 神经科学（放电模式）
    • 工程学（结构失稳）
    • 经济学（市场崩溃）
  5. 通往混沌之路：
    • 周期倍增级联
    • Feigenbaum 常数

分岔理论为理解自然界中的"突变"提供了数学框架。它告诉我们：渐变的原因可以产生突变的结果——理解分岔，就是理解临界现象的钥匙。

练习题

概念题

  1. 什么是分岔？分岔值有什么特殊意义？

  2. 解释超临界和亚临界分岔的区别。为什么亚临界分岔更"危险"？

  3. Hopf 分岔与鞍结点分岔有什么本质区别？

  4. 什么是余维数？余维数 2 的分岔为什么比余维数 1 的分岔更难观察到？

计算题

  1. 分析系统 的分岔行为。找出所有分岔点并分类。

  2. 对于系统

  1. 找出平衡点
  2. 分析 处的分岔类型
  3. 计算 Hopf 分岔发生的 值（如果存在）
  1. 对于 Logistic 映射 ：
    1. 找出不动点并分析其稳定性
    2. 证明当 时发生周期倍增分岔
    3. 找出周期-2 轨道的显式表达式
  2. 考虑带时间延迟的 Logistic 方程：

定性分析延迟 如何影响系统的稳定性。

分析题

  1. 湖泊生态系统可能存在"清水态"和"浑水态"两种稳定状态。用分岔理论解释：
    1. 为什么湖泊可能"突然"从清变浑
    2. 为什么减少污染物不一定能使湖泊恢复清澈
    3. 什么是"滞后效应"
  2. 解释为什么周期倍增级联最终导致混沌。 Feigenbaum 常数的"普适性"意味着什么？

编程题

  1. 编写程序绘制以下系统的完整分岔图：

标出所有分岔点。

  1. 实现 AUTO 类型的延拓算法，追踪 的平衡点分支，包括折回部分。

  2. 对 FitzHugh-Nagumo 模型，数值计算 Hopf 分岔点的位置，并绘制分岔图（极限环振幅 vs 外部电流）。

  3. 验证 Feigenbaum 常数：对 Logistic 映射，数值计算周期倍增的分岔值 ，计算比值序列 ，观察其收敛到约 4.669 。

参考文献

  1. Strogatz, S. H. (2015). Nonlinear Dynamics and Chaos. CRC Press.

  2. Kuznetsov, Y. A. (2004). Elements of Applied Bifurcation Theory. Springer.

  3. Guckenheimer, J., & Holmes, P. (1983). Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields. Springer.

  4. Wiggins, S. (2003). Introduction to Applied Nonlinear Dynamical Systems and Chaos. Springer.

  5. Seydel, R. (2010). Practical Bifurcation and Stability Analysis. Springer.

  6. Feigenbaum, M. J. (1978). "Quantitative universality for a class of nonlinear transformations." Journal of Statistical Physics, 19(1), 25-52.

  7. May, R. M. (1976). "Simple mathematical models with very complicated dynamics." Nature, 261, 459-467.

  8. Scheffer, M., et al. (2009). "Early-warning signals for critical transitions." Nature, 461, 53-59.

本文是《常微分方程的世界》系列的第 10 章。

  • 本文标题：常微分方程（十）分岔理论
  • 本文作者：Chen Kai
  • 创建时间：2019-05-24 09:15:00
  • 本文链接：https://www.chenk.top/%E5%B8%B8%E5%BE%AE%E5%88%86%E6%96%B9%E7%A8%8B%EF%BC%88%E5%8D%81%EF%BC%89%E5%88%86%E5%B2%94%E7%90%86%E8%AE%BA/
  • 版权声明：本博客所有文章除特别声明外，均采用 BY-NC-SA 许可协议。转载请注明出处！
 评论