带容量约束的 SIR 疫情干预¶

背景¶

疫情干预策略需要权衡两方面影响：减少传染性接触可以降低患病人数，但长期限制本身也有代价。易感者-感染者-移出者（SIR）模型可以简洁地描述这一权衡。本例进一步为感染比例设置硬上限，用医疗系统容量表示优化轨迹在任何时刻都不能超过的水平。

这是一个用于讲解最优控制的简化模型，不应作为疫情预测或公共卫生决策工具。

问题定义¶

令 \(S(t)\)、\(I(t)\) 和 \(R(t)\) 分别为易感、感染和移出人群比例，\(u(t)\) 为潜在传染性接触的减少比例。受控 SIR 动力学为

\[ \begin{aligned} \dot S &= -\beta(1-u)SI, \\ \dot I &= \beta(1-u)SI-\gamma I, \\ \dot R &= \gamma I, \end{aligned} \]

其中 \(\beta=0.32\ \mathrm{day}^{-1}\)、\(\gamma=0.10\ \mathrm{day}^{-1}\)，初始条件为

\[ (S(0),I(0),R(0))=(0.99,0.01,0). \]

在固定的 100 天时域内，目标函数为

\[ \min_u J=\int_0^{100}\left[I(t)^2+0.05u(t)^2\right]\,\mathrm dt. \]

第一项惩罚感染负担，第二项惩罚干预强度。路径约束为

\[ 0\le u(t)\le0.8, \qquad 0\le I(t)\le0.08. \]

由于三个状态导数之和为零，且初始人群比例之和为一，因此整条轨迹始终满足 \(S+I+R=1\)。

变量与单位¶

符号	含义	单位
\(t\)	时间	天
\(S,I,R\)	易感、感染和移出人群比例	-
\(u\)	接触减少比例	-
\(\beta\)	基准传播率	天\(^{-1}\)
\(\gamma\)	移出率	天\(^{-1}\)

建模选择¶

该问题使用单个固定时间的 Lobatto 阶段，包含 80 个网格区间，每个区间有两个插值点。干预上限和医疗容量都作为路径约束施加。由于重构后的状态和控制是分段线性的，节点边界在每个完整区间内都成立。终端人群状态不固定，而是由动力学和积分代价共同决定。优化完成后，程序还会通过独立的 4,001 点检查验证重构路径和人口守恒。

本例特意选择了使容量约束真正参与最优决策的干预权重。若权重过小，优化器会把感染比例压到远低于容量上限，无法展示活跃路径约束。

运行示例¶

在仓库根目录运行：

python -m examples.sir_epidemic_control

保存图片且不打开窗口：

python -m examples.sir_epidemic_control --save sir-epidemic-control.png --no-show

关键实现¶

system = System(0)
phase = system.new_phase(
    ["susceptible", "infected", "removed"], ["intervention"]
)
susceptible, infected, _ = phase.x
(intervention,) = phase.u

incidence = BETA * (1 - intervention) * susceptible * infected
phase.set_dynamics([-incidence, incidence - GAMMA * infected,
                    GAMMA * infected])
phase.set_integral([infected**2 + EFFORT_WEIGHT * intervention**2])
phase.set_phase_constraint(
    [intervention, infected], [0.0, 0.0], [INTERVENTION_MAX, CAPACITY]
)
phase.set_boundary_condition([0.99, 0.01, 0.0], [None] * 3,
                             0.0, HORIZON)
system.set_phase([phase])
system.set_objective(phase.I[0])

已验证结果¶

Ipopt 成功终止，目标函数值为 \(J=1.29856860\)。密集网格上的感染比例峰值为 \(0.080000\)；相对容量上限的超限量为 \(7.375\times10^{-11}\)，处于求解器容差内。最优策略先允许感染人数上升，再改变接触减少强度，使感染比例在一段时间内贴近容量，最后逐步解除干预。终端移出比例为 \(0.710622\)，最大人口守恒误差小于 \(3\times10^{-16}\)。

优化后的 SIR 人群比例、活跃容量约束与接触减少强度

源码¶

完整可运行示例见 examples/sir_epidemic_control.py。