常微分方程（二）：一阶微分方程的求解方法

银行账户复利、药物代谢、盐水稀释、电容充电——这些看似无关的现象，背后都遵循同一种数学规律：一阶常微分方程（ODE）。关键在于快速识别方程属于哪一种典型形式，因为每种形式都有对应的解析解法。读完本章，你将能在几秒内对陌生的一阶方程完成模式匹配，并精准调用相应的求解工具。

总结#

如何识别可分离方程并直接积分求解
利用积分因子将线性方程转化为一个“完美导数”
恰当方程的判定条件及其几何意义：解曲线是某个势函数的等高线
通过伯努利代换 $v = y^{1-n}$ 将一大类非线性方程线性化
五个真实应用场景：药物半衰期、RC 电路充电、盐水混合、Logistic 增长、指数复利——每个都是可迁移的小型建模范例

前置知识#

第一章内容：ODE 的定义、初值问题的概念、方向场（slope field）的直观形态
基础积分技巧：换元法、分部积分、部分分式分解

可分离方程：最自然的形式#

形式特征#

\frac{dy}{dx} \;=\; g(x)\,h(y).

\frac{dy}{h(y)} \;=\; g(x)\,dx \qquad\Longrightarrow\qquad \int \frac{dy}{h(y)} \;=\; \int g(x)\,dx + C.

常数 $$C$$ 是解族的参数——每个取值对应一条解曲线。例如，对方程 $\frac{dy}{dx} = -x/y$ 进行分离，得到 $y\,dy + x\,dx = 0$ ，积分后即为 $$x^2 + y^2 = C$$ ，表示一族同心圆。方向场中的小线段处处与这些圆相切——这正是“解”的几何含义。

分离变量将二维方向场坍缩为一族一维曲线 $$x^2 + y^2 = C$$ 。

例题：指数增长与衰减#

求解 $\dfrac{dy}{dx} = ky$ ，其中 $$k$$ 为常数。

分离变量： $\dfrac{dy}{y} = k\,dx$
两边积分： $\ln|y| = kx + C_1$
指数化简： $y = A e^{kx}$ ，其中 $A = \pm e^{C_1}$ 同时吸收了符号和常数

这里发生了什么？增长率 $$k$$ 控制的是时间尺度，而非曲线形状。两种典型情形：

$$k > 0$$ ：每 $\ln 2 / k$ 时间翻倍——适用于无资源限制的种群增长、连续复利计息、失控的核链式反应
$$k < 0$$ ：半衰期为 $\ln 2 / |k|$ ——如放射性衰变、电容放电、肝脏清除药物

应用：药物代谢与半衰期法则#

\frac{dC}{dt} = -kC, \qquad C(t) = C_0 e^{-kt}.

t_{1/2} = \frac{\ln 2}{k} \approx \frac{0.693}{k}.

以布洛芬为例，其 $t_{1/2} \approx 2$ 小时。若初始剂量为 400 mg：

时间	残留量	含义
0 h	400 mg	峰浓度
2 h	200 mg	经历一个半衰期
4 h	100 mg	两个半衰期
6 h	50 mg	已低于治疗窗口

这解释了为何布洛芬需每 4–6 小时服用一次：间隔过长则药效不足，过短则药物蓄积。“五个半衰期基本清除”（约 99%）的经验法则也源于同一公式。

应用：Logistic 方程#

\frac{dP}{dt} = r P\left(1 - \frac{P}{K}\right).

\int \!\left(\frac{1}{P} + \frac{1/K}{1 - P/K}\right)\,dP \;=\; \int r\,dt,

P(t) \;=\; \frac{K}{1 + \left(\dfrac{K}{P_0} - 1\right)e^{-rt}}.

有趣的是，即使不解方程，也能从结构中直接看出三个关键性质：

当 $P \ll K$ 时，括号内近似为 1，增长接近指数形式
当 $$P = K/2$$ 时，右侧达到最大值 $$rK/4$$ ，即最大可能增长率
当 $P \to K$ 时，括号趋于零，种群趋于饱和

下图右侧展示了增长率 $\dot P$ 随 $$P$$ 的变化：它是一条开口向下的抛物线，零点位于 $$P=0$$ 和 $$P=K$$ ，峰值恰好在 $$P=K/2$$ 处。

不同初值下的 Logistic 曲线，以及增长率 dP/dt 关于 P 的抛物线。

左：所有解最终都收敛到 $$K$$ 。右：增长率 $\dot P = rP(1-P/K)$ 在 $$P=K/2$$ 处取得最大值。

一阶线性方程：积分因子法#

标准形式#

\frac{dy}{dx} + P(x)\,y \;=\; Q(x).

这是应用数学中最常用的方程形式。当某物理量的变化率由“自身比例项”和“外部驱动项”共同决定时，通常会归结为此类方程。

核心技巧#

\mu(x) \;=\; e^{\int P(x)\,dx}.

\frac{d}{dx}\bigl[\mu(x)\,y\bigr] \;=\; \mu(x)\,Q(x).

\boxed{\,y(x) \;=\; \frac{1}{\mu(x)} \left[\int \mu(x)\,Q(x)\,dx + C\right].\,}

整个方法至此结束。为何有效？因为 $\mu$ 被精心构造为满足 $\mu' = \mu P$ ，这正是乘积法则所要求的。

例题： $$y' + 2xy = x$$ #

步骤	计算过程	结果
1. 识别系数	$$P = 2x$$ ， $$Q = x$$	—
2. 积分因子	$\mu = \exp\!\int 2x\,dx$	$\mu = e^{x^2}$
3. 两边同乘	左边变为 $\frac{d}{dx}[e^{x^2}y]$	右边为 $x e^{x^2}$
4. 积分	$\int x e^{x^2} dx = \tfrac12 e^{x^2}$	$e^{x^2}y = \tfrac12 e^{x^2} + C$
5. 解出 $$y$$	两边除以 $e^{x^2}$	$y = \tfrac12 + C e^{-x^2}$

当 $x \to \pm\infty$ 时， $C e^{-x^2} \to 0$ ，因此所有解最终都趋近于平衡值 $$y = 1/2$$ 。积分常数 $$C$$ 仅在原点附近影响瞬态行为。

积分因子的几何作用#

三幅图能清晰说明其效果。以 $$y' + 2y = 4$$ 为例（此时 $$P = 2$$ ， $\mu = e^{2x}$ ，平衡解为 $$y=2$$ ）：

y’ + 2y = 4 的方向场与解族、积分因子 e^{2x}、以及左端坍缩为完整导数。

左：所有解 $y = 2 + Ce^{-2x}$ 都迅速被吸引至平衡点 $$y=2$$ 。中：积分因子 $e^{2x}$ 是一条光滑的指数曲线。右：乘上 $\mu$ 后， $\mu(x)y(x)$ 与 $\int 4\mu\,dx$ 完全重合——整个左边已化为一个原函数。

应用：RC 电路充电#

RC\,\frac{dV_c}{dt} + V_c \;=\; V_s.

V_c(t) \;=\; V_s\bigl(1 - e^{-t/\tau}\bigr), \qquad \tau = RC.

时间常数 $\tau$ 是每位电气工程师必须牢记的参数：

时间	$$V_c/V_s$$	工程解读
$\tau$	63.2%	“一个时间常数”
$3\tau$	95.0%	大多数数字逻辑电路中可视为充满
$5\tau$	99.3%	工程上认为“完全充满”

若将符号反转并设 $$V_s = 0$$ ，同一公式描述放电过程——这也解释了为何相机闪光灯在两次拍摄之间需要等待约一秒才能重新充电。

恰当方程#

几何思想#

M(x,y)\,dx + N(x,y)\,dy = 0.

\frac{\partial F}{\partial x} = M, \qquad \frac{\partial F}{\partial y} = N.

F(x,y) \;=\; C.

解曲线正是 $$F$$ 的等高线。这一观点极为深刻：求解 ODE 被转化为寻找一个标量函数 $$F$$ ，其等高线图本身就构成了所有解。

恰当性判据#

\boxed{\,\dfrac{\partial M}{\partial y} \;=\; \dfrac{\partial N}{\partial x}.\,}

若不满足，则方程不是恰当的（但仍有补救办法，见 §3.4）。

求解步骤示例#

考虑方程 $(2x + y)\,dx + (x + 2y)\,dy = 0$ 。先验证恰当性： $$M_y = 1 = N_x$$ ，满足条件。接下来重构 $$F$$ ：

对 $$M$$ 关于 $$x$$ 积分： $$F = x^2 + xy + g(y)$$
对 $$y$$ 求偏导并与 $$N$$ 对比： $x + g'(y) = x + 2y \Rightarrow g'(y) = 2y \Rightarrow g(y) = y^2$
因此 $$F(x,y) = x^2 + xy + y^2$$ ，解族为 $$x^2 + xy + y^2 = C$$

几何上，这些是倾斜的椭圆。下图将梯度场 $\nabla F = (M, N)$ 叠加在等高线上；梯度处处与等高线垂直——这正是“沿解曲线 $$dF = 0$$ ”的几何体现：梯度指向“上坡”方向，而解曲线始终保持在同一高度。

F(x,y) = x^2 + xy + y^2 的等高线与梯度场，二者处处垂直。

每条等高线对应一条解。灰色箭头表示 $\nabla F$ ，始终垂直于等高线——这正是“沿解 $$dF = 0$$ ”的几何内涵。

补救策略：使非恰当方程变为恰当#

若 $M_y \neq N_x$ ，可尝试乘以一个积分因子 $\mu(x,y)$ 来修正。在两种常见情形下， $\mu$ 有显式表达式：

条件	积分因子选择
$$(M_y - N_x)/N$$ 仅依赖 $$x$$	$\mu(x) = \exp\!\int \dfrac{M_y - N_x}{N}\,dx$
$$(N_x - M_y)/M$$ 仅依赖 $$y$$	$\mu(y) = \exp\!\int \dfrac{N_x - M_y}{M}\,dy$

这实际上与线性方程的积分因子构造一致——事实上，线性方程的解法正是此方法的一个特例。

伯努利方程#

方程形式#

\frac{dy}{dx} + P(x)\,y \;=\; Q(x)\,y^{n}, \qquad n \neq 0, 1.

当 $$n = 0$$ 时，方程已是线性的；当 $$n = 1$$ 时，它是可分离的。真正有趣的情形是其他所有 $$n$$ 值，此时 $$y^n$$ 项破坏了线性结构。

代换技巧#

\frac{dv}{dx} + (1-n)P(x)\,v \;=\; (1-n)Q(x).

这是一个关于 $$v$$ 的线性方程。可用 §2 中的积分因子法求解 $$v$$ ，再通过 $y = v^{1/(1-n)}$ 回代还原。

代换的几何意义#

以 $$y' + y = y^2$$ 为例（此时 $$n = 2$$ ，故 $$v = 1/y$$ ）。代换后变为线性方程 $$v' - v = -1$$ ，其解是以 $$v = 1$$ 为中心的指数族。在 $$v$$ -空间中，轨迹整齐有序；而在原始 $$y$$ -空间中，轨迹弯曲，且可能出现有限时间内的爆破（blow-up）。

y’ + y = y^2 在原 (x, y) 平面与替换后 (x, v) 平面（v = 1/y）的对比。

相同的五条轨迹。左图在 $$(x, y)$$ 平面：非线性， $$y=1$$ 是不稳定平衡点。右图在 $$(x, v)$$ 平面：线性， $$v=1$$ 是稳定平衡点，对应 $$y$$ 从下方趋近于 1。

同样的代换还能还原 Logistic 方程：它本质上是 $$n = 2$$ 、 $$P(t) = -r$$ 、 $$Q(t) = -r/K$$ 的伯努利方程。因此，§1.4 中看似复杂的部分分式积分，其实只是伯努利代换加线性求解的另一种表现形式。

综合应用：水箱中的盐#

这是经典的“房室模型”（compartment model），其骨架广泛应用于流行病学、药代动力学和化学工程。

问题设定：一个 1000 L 的水箱初始装有纯水。浓度为 2 g/L 的盐水以 5 L/min 的速率流入，水箱充分搅拌后以相同速率流出，因此体积恒定。求盐的质量 $$Q(t)$$ 。

\underbrace{\frac{dQ}{dt}}_{\text{变化率}} = \underbrace{5 \cdot 2}_{\text{流入：g/min}} \,-\, \underbrace{5 \cdot \tfrac{Q}{1000}}_{\text{流出：g/min}} = 10 - \tfrac{Q}{200}, \qquad Q(0) = 0.

Q(t) \;=\; 2000\bigl(1 - e^{-t/200}\bigr).

结果解读：平衡盐量 $Q^\ast = 2000$ g，对应浓度 2 g/L——恰好等于流入浓度。时间常数 $\tau = 200$ 分钟控制趋近速度。若初始盐量 $Q(0) > Q^\ast$ ，则描述的是用清水冲洗偏咸水箱的过程——下图展示了两种情形，最终都收敛到同一平衡点。

两种初始条件，一个共同平衡。从 $$Q(0)=0$$ 开始灌注，或从 $$Q(0)=3000$$ 开始冲洗，最终浓度都趋于 2 g/L；区别仅在于偏离平衡的方向。

数值方法预览#

一旦 $$f(x,y)$$ 足够复杂，解析技巧就会失效——而大多数实际问题恰恰如此。此时的解决方案是：沿着方向场一步步推进。

欧拉方法#

y_{n+1} \;=\; y_n + h\,f(t_n, y_n).

单步局部截断误差： $$O(h^2)$$
全局误差（经过 $$1/h$$ 步）： $$O(h)$$ ，因此被称为“一阶方法”

实现简单，但精度通常不足。

四阶 Runge–Kutta 方法（RK4）#

k_1 = f(t_n, y_n), \qquad k_2 = f\!\left(t_n + \tfrac{h}{2},\, y_n + \tfrac{h}{2}k_1\right),

k_3 = f\!\left(t_n + \tfrac{h}{2},\, y_n + \tfrac{h}{2}k_2\right), \qquad k_4 = f(t_n + h,\, y_n + h k_3),

y_{n+1} \;=\; y_n + \frac{h}{6}\bigl(k_1 + 2k_2 + 2k_3 + k_4\bigr).

全局误差为 $$O(h^4)$$ ——在相同步长下，精度比欧拉方法高出四个数量级，代价是多进行三次函数求值。

方向场揭示长期行为#

在调用数值求解器前，先观察右端函数。两个表面相似的方程，长期行为可能截然不同：

左（线性）：唯一直线解 $$y=x+1$$ ，其余解均指数发散。右（伯努利型）：三个平衡点 $$y=-1, 0, +1$$ ；无论初值正负，轨迹都被“挤压”到 $\pm 1$ 。同属一阶方程，几何结构却大相径庭。

SciPy 三行代码搞定#

以牛顿冷却定律为例： $T' = -k(T - T_{\text{env}})$ ，初始温度 $T(0) = 90^\circ\text{C}$ ，环境温度 $T_{\text{env}} = 20^\circ\text{C}$ ， $$k = 0.1$$ /min。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
from scipy.integrate import solve_ivp
import numpy as np

def cooling(t, T, k=0.1, T_env=20):
    return -k * (T - T_env)

sol = solve_ivp(cooling, [0, 60], [90.0], method="RK45", dense_output=True)

t = np.linspace(0, 60, 300)
T = sol.sol(t)[0]
print(f"t=10 min 时: {T[50]:.1f} C")
print(f"t=30 min 时: {T[150]:.1f} C")

solve_ivp 会自动调整步长。对大多数问题，你不需要自己实现 RK4；真正需要关注的是：标准求解器在什么情况下会遇到困难（例如刚性系统——详见第 11 章）。

总结：四问筛选法#

遇到一阶方程时，按以下顺序逐一判断——第一个“是”就指明了解法：

提问	若成立	解法
能否写成 $y' = g(x)\,h(y)$ ？	可分离	移项后两边积分
是否为 $$y' + P(x)y = Q(x)$$ ？	线性	使用积分因子 $\mu = e^{\int P\,dx}$
能否写成 $M\,dx + N\,dy = 0$ 且 $$M_y = N_x$$ ？	恰当	重构势函数 $$F$$ ；解为 $$F = C$$
是否为 $$y' + Py = Qy^n$$ （ $n \notin \{0,1\}$ ）？	伯努利	代换 $v = y^{1-n}$ ，再用线性方法
以上皆否？	—	使用数值方法（第 11 章）或寻找巧妙代换

这四种方法并非彼此孤立。线性方程的积分因子，其实是修正非恰当方程时所用积分因子的特例；伯努利代换则将非线性问题转化为线性问题。因此，掌握这四种方法，本质上是在反复练习同一个核心思想——找到合适的变量替换或乘子，将方程转化为一眼就能积分的形式——只是在四种不同场景下分别应用而已。

练习题#

基础题

求解 $y' = y^2 \sin x$ ，初值 $$y(0) = 1$$ 。解在何处发生爆破？
求解 $y' + y = e^{-x}$ 。指出其中的瞬态部分与稳态部分。
验证 $(2xy + 3)\,dx + (x^2 + 1)\,dy = 0$ 是恰当方程并求其通解。
一个 1000 L 水箱初始含有 50 kg 溶解盐。纯水以 10 L/min 流入，充分混合后的溶液以相同速率流出。求盐量 $$Q(t)$$ 及盐量减半所需时间。

进阶题

求解伯努利方程 $$y' + y = y^3$$ 。绘制相线并标出所有平衡点。
对方程 $y' = 3y^{2/3}$ ，初值 $$y(0) = 0$$ ，在区间 $[0, \infty)$ 上给出至少两个不同的解。为何这不违反存在唯一性定理？
某细菌种群满足 $P' = 0.5\,P(1 - P/1000)$ ，初值 $$P(0) = 50$$ 。显式求出 $$P(t)$$ ，并计算种群增长率最大的时刻。

编程题

从零实现欧拉法和 RK4。在方程 $$y' = -y$$ 、初值 $$y(0) = 1$$ 、区间 $$[0, 5]$$ 上，比较两者在步长 $$h = 0.5, 0.1, 0.01$$ 下的全局误差。绘制 $\log(\text{误差})$ 对 $\log(h)$ 的图像，验证斜率分别为 1 和 4。
复现 §6.3 中的方向场对比：使用 numpy.meshgrid 和 matplotlib.pyplot.quiver 绘制方向场，再用 scipy.integrate.solve_ivp 计算轨迹。尝试方程 $$y' = y - y^5$$ ，并标出所有平衡点。

下一步#

下一章把视野提升到二阶及以上。多自由度系统会引入’特征方程’这个核心工具——一个代数方程，它的根决定了原来微分方程解的结构（指数衰减/周期振荡/共振增长）。这是 ODE 第一次真正展示’用代数解几何’的力量。

参考文献#

Boyce & DiPrima，Elementary Differential Equations and Boundary Value Problems，Wiley (2012)。本科阶段的标准教材；第 2–3 章以经典严谨的方式覆盖了本章全部内容。
Zill，A First Course in Differential Equations with Modeling Applications，Cengage (2017)。更侧重应用；其中的混合槽与电路例题尤为清晰。
Strogatz，Nonlinear Dynamics and Chaos，Westview (2014)。第 2 章对 Logistic 方程的相线分析，是学完本章后的必读延伸。
SciPy 文档：scipy.integrate 。solve_ivp 的方法选择（RK45、Radau、BDF）及事件处理功能均在此处查阅。

常微分方程（二）：一阶微分方程的求解方法

总结#

前置知识#

可分离方程：最自然的形式#

形式特征#

例题：指数增长与衰减#

应用：药物代谢与半衰期法则#

应用：Logistic 方程#

一阶线性方程：积分因子法#

标准形式#

核心技巧#

例题： $$y' + 2xy = x$$ #

积分因子的几何作用#

应用：RC 电路充电#

恰当方程#

几何思想#

恰当性判据#

求解步骤示例#

补救策略：使非恰当方程变为恰当#

伯努利方程#

方程形式#

代换技巧#

代换的几何意义#

综合应用：水箱中的盐#

数值方法预览#

欧拉方法#

四阶 Runge–Kutta 方法（RK4）#

方向场揭示长期行为#

SciPy 三行代码搞定#

总结：四问筛选法#

练习题#

下一步#

参考文献#

ODE 入门精讲 18 篇

读有所得？

总结#

前置知识#

可分离方程：最自然的形式#

形式特征#

例题：指数增长与衰减#

应用：药物代谢与半衰期法则#

应用：Logistic 方程#

一阶线性方程：积分因子法#

标准形式#

核心技巧#

例题：$y' + 2xy = x$ #

积分因子的几何作用#

应用：RC 电路充电#

恰当方程#

几何思想#

恰当性判据#

求解步骤示例#

补救策略：使非恰当方程变为恰当#

伯努利方程#

方程形式#

代换技巧#

代换的几何意义#

综合应用：水箱中的盐#

数值方法预览#

欧拉方法#

四阶 Runge–Kutta 方法（RK4）#

方向场揭示长期行为#

SciPy 三行代码搞定#

总结：四问筛选法#

练习题#

下一步#

参考文献#

ODE 入门精讲 18 篇

读有所得？

例题： $$y' + 2xy = x$$ #