时滞微分方程

时滞微分方程（delay differential equation）把过去的状态放进当前的变化规律。它适合描述信号传输、运输、成熟期等需要时间的过程。与普通微分方程相比，确定未来通常需要一段已知历史。

温控器收到的消息晚了两分钟

设一个理想化温控系统的目标温度固定。令 $x(t)$ 表示实际温度减目标温度，单位为摄氏度，时间 $t$ 用分钟计。正值表示偏热，负值表示偏冷。为单独研究反馈延迟，假设调节产生的净温度变化率正比于检测到的温差，方向相反；这个线性模型只描述目标附近的小偏差。

如果测量立即到达，反馈规律是 $x^{\prime }=-ax(t)$，其中 $a>0$ 的单位为每分钟。从 $x(0)=x_0$ 出发，解为 $x_0{e}^{-at}$，温差单调衰减。

若检测与执行总共落后 $\tau$ 分钟，控制器使用的是 $x(t-\tau )$，于是 $$x^{\prime }(t)=-ax(t-\tau ).$$ 取教学参数 $a=0.4$ 每分钟、$\tau =2$ 分钟。在时刻0，如果控制器收到两分钟前的温差为0.8度，则当前变化率是 $-0.32$ 度每分钟。即使当前已经不热，它仍可能继续按旧信息降温。

图中向右的箭头表示时间经过，反馈箭头把过去的测量送到当前。$\tau$ 是模型中的延迟，不是数值求解时随意选择的步长。

反馈作用在现在，依据却来自过去。

同一个当前值，为何会有不同未来

在最初两分钟内，方程不断读取区间 $[-2,0]$ 上的温差。因此除了 $x(0)$，还要指定 $$x(s)=\varphi (s),-\tau \le s\le 0.$$ 函数 $\varphi$ 称为初始历史。比较两种历史： $${\varphi }_1(s)=0.8,{\varphi }_2(s)=0.8(1+s/2),-2\le s\le 0.$$ 两者在当前都等于0.8，但两分钟前分别是0.8与0。故 ${x_1}^{\prime }(0+)=-0.32$，而 ${x_2}^{\prime }(0+)=0$。当前值相同，当前变化率已经不同。

对第一种历史，$0\le t\le 2$ 时有 $x_1(t)=0.8-0.32t$。对第二种历史，${\varphi }_2(t-2)=0.4t$，所以 ${x_2}^{\prime }=-0.16t$，积分得 $x_2(t)=0.8-0.08t^2$。两分钟后分别为0.16与0.48度。这个差别由历史决定。

虚线表示事先给定的历史，实线表示方程决定的未来。

把过去算好的部分接着使用

对于常值历史，前两分钟已经得到 $x(t)=0.8-0.32t$。在 $2\le t\le 4$ 内，延迟项读取的正是这段解： $$x^{\prime }(t)=-0.4[0.8-0.32(t-2)]=-0.32+0.128(t-2).$$ 从 $x(2)=0.16$ 积分，得到 $$x(t)=0.16-0.32(t-2)+0.064(t-2{)}^2,2\le t\le 4.$$ 于是 $x(4)=-0.224$。温差已经越过零，这和没有延迟的单调回归不同。继续计算 $[4,6]$ 时，再将刚求出的二次多项式作为已知延迟输入。

这种分步法以一个延迟长度为区间，把已经知道的过去用于下一段。在更一般的方程 $$x^{\prime }(t)=F(t,x(t),x(t-\tau ))$$ 中，下一段仍可能需要求一个关于当前状态的普通微分方程；已知的是延迟项，并非整个右端。Marc R. Roussel 的时滞方程讲义讨论了历史、分步法和线性化。

如果 $F$ 连续并对状态变量局部 Lipschitz 连续，连续历史可以在这些条件成立的范围内逐段确定唯一解。全局存在还需排除有限时间发散。历史在0处的左导数不必符合方程给出的右导数，因此求解器还需要处理由接点传播的光滑性变化。

真正控制响应的是反应速度乘以延迟

令 $u=t/\tau$、$z(u)=x(\tau u)/0.8$，则链式法则给出 $$z^{\prime }(u)=-qz(u-1),q=a\tau .$$ $q$ 没有单位。它衡量控制器在“等待一个延迟”的时间内会作出多强的修正。增大反应强度或增加延迟，都可能把同一类反馈推向振荡。

常值历史 $z(u)=1$（$-1\le u\le 0$）下，可比较 $q=0.4,1.2,1.8$。前两组的偏差最终衰减，但衰减过程形状不同；第三组出现增长的振荡。下图中的水平时间单位是一个延迟长度。

振荡可以来自模型自身的延迟；判断是否由数值误差造成，还需缩小计算步长。

临界值为什么与圆周率有关

寻找指数形式的解 $z(u)=e^{\lambda u}$，代入得 $$\lambda +q{e}^{-\lambda }=0.$$ 这称为特征方程。指数项使它通常具有无穷多个复根。根的实部决定指数模式增长还是衰减，虚部决定振荡。

先证明 $0<q<\pi /2$ 时没有实部非负的根。设 $\lambda =\alpha +i\beta$ 且 $\alpha \ge 0$，分离实虚部： $$\alpha =-q{e}^{-\alpha }\cos \beta ,\beta =q{e}^{-\alpha }\sin \beta .$$ 第二式给出 $|\beta |\le q<\pi /2$，故 $\cos \beta >0$。第一式却要求 $\alpha <0$，矛盾。对这个常系数滞后型方程，特征根全在左半平面的谱判据因而给出零解的渐近稳定性。

临界时令 $\lambda =i\omega$、$\omega >0$，得到 $$q\cos \omega =0,\omega =q\sin \omega .$$ 第一对虚轴根出现在 $\omega =q=\pi /2$。此时有不衰减的正弦模式，因此零解不再渐近稳定。沿根分支对 $\lambda e^{\lambda }=-q$ 求导， $$\frac{d\lambda }{dq}=\frac{\lambda }{q(1+\lambda )}.$$ 在 $\lambda =i\pi /2$ 处，实部为正；随 $q$ 增大，这对根进入右半平面。其余虚轴穿越发生在 $q=\pi /2+2k\pi$（非负整数 $k$），穿越方向也相同。因此本方程的渐近稳定范围恰为 $0<a\tau <\pi /2$，超过它则不稳定。这是所写线性延迟方程的结论，不是所有控制系统的统一阈值。

如何判断计算出的振荡可信

求解 $x(t-\tau )$ 通常需要保存过去的计算点，并对网格间的值插值。简单地把它换成“上一个点”，会把延迟错误地改成求解步长。

本例的头两个区间有明确多项式，可以先用它们检查程序；更长时间则比较步长减半后的曲线。还可以把 $\tau =0$ 的无延迟模型单独求解，检查是否得到指数衰减。有限区间的数值对照支持计算精度，长期稳定性仍由前面的理论判断。

温控系统真实的散热、功率上限和执行器饱和会改变方程。当计算温差增长到不再“小”时，应回到建模假设。延迟还可出现在捕食者-猎物模型的成熟期、库存补货和神经信号中；它们需要各自的机制方程。

来源与继续阅读

• Marc R. Roussel，Delay-differential equations：历史函数、分步法、特征方程与根的位置。
• 先修：一阶线性微分方程、平衡点与稳定性、欧拉公式；后续：分岔、微分方程数值解、动态系统建模。