如何在Python中数值求解指定微分方程组？（仅熟悉solve_ivp）

核心逻辑 #

单个微分方程用 solve_ivp 时，传入的是单个导数函数；求解方程组只需把所有未知函数的导数打包成向量，让传入的函数返回这个导数向量即可。

实操示例 #

假设你要解的微分方程组为：

dx/dt = -x + y
dy/dt = x - 3y
初始条件：x(0)=1，y(0)=0，求解区间 t ∈ [0, 10]

1. 导入依赖库

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

2. 定义方程组的导数函数
   把未知函数打包成向量 z = [x, y]，函数返回对应的导数向量 [dx/dt, dy/dt]：

def system(t, z):
    x, y = z
    dxdt = -x + y
    dydt = x - 3*y
    return [dxdt, dydt]

3. 设置求解参数

# 初始条件向量
z0 = [1, 0]
# 求解的时间区间
t_span = [0, 10]
# 可选：指定输出的时间采样点
t_eval = np.linspace(t_span[0], t_span[1], 100)

4. 调用 solve_ivp 求解

sol = solve_ivp(system, t_span, z0, t_eval=t_eval)

5. 提取结果并可视化

# 提取时间点和解
t = sol.t
x_sol = sol.y[0]
y_sol = sol.y[1]

# 绘图展示
plt.figure(figsize=(10, 6))
plt.plot(t, x_sol, label='x(t)')
plt.plot(t, y_sol, label='y(t)')
plt.xlabel('t')
plt.ylabel('x, y')
plt.title('微分方程组的解')
plt.legend()
plt.grid(True)
plt.show()

关键提示 #

• 无论方程组有多少个未知函数，都要统一打包成向量处理，输入输出保持维度一致。
• solve_ivp 返回的 sol.y 是二维数组，每一行对应一个未知函数的解，每一列对应一个时间点的取值。
• 如果是高阶微分方程，需要先通过变量替换转化为一阶微分方程组再求解。

内容的提问来源于stack exchange，提问作者Hug de Roda