含阶跃函数的微分方程求解：代码中u的条件设置疑问

我来帮你搞定这个阶跃函数的实现！你现在的问题是把u设成了固定值5，但实际需要的是t<5时u=0，t≥5时u=2的阶跃输入，下面一步步给你修正代码：

问题分析 #

你的微分方程组是：

2dx/dt = -x(t) + u(t)
5dy/dt = -y(t) + x(t)
u(t) = 2S(t-5) （S是阶跃函数，t≥5时为1，否则0）
初始条件x(0)=0，y(0)=0

原来的代码有两个核心问题：

• u被硬编码为固定值，没有随时间t动态切换
• 求解y时参数传递错误，且没有正确关联随时间变化的x(t)

方法1：耦合求解（推荐） #

把x和y放到同一个状态向量里，让odeint同时求解两个方程，这样能保证x(t)和y(t)在每个时间点的取值完全匹配，避免插值误差：

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

def combined_model(z, t):
    x, y = z
    # 实现阶跃函数u=2S(t-5)
    u = 2 if t >= 5 else 0
    # 计算两个微分方程
    dxdt = (-x + u) / 2
    dydt = (-y + x) / 5
    return [dxdt, dydt]

# 初始条件 x(0)=0, y(0)=0
z0 = [0, 0]
# 时间范围，增加点数让曲线更平滑
t = np.linspace(0, 40, 1000)
# 求解微分方程
z = odeint(combined_model, z0, t)
# 提取x和y的结果
x = z[:, 0]
y = z[:, 1]

# 绘图
plt.figure(figsize=(10,6))
plt.plot(t, x, 'r-', label='x(t)')
plt.plot(t, y, 'b-', label='y(t)')
plt.axvline(x=5, color='gray', linestyle='--', label='t=5 (阶跃点)')
plt.xlabel('Time')
plt.ylabel('Value')
plt.legend()
plt.grid(True)
plt.show()

方法2：先求x再求y #

如果你想分开求解，需要先得到x(t)的数组，然后在求解y时，根据当前时间t找到对应的x值（用插值实现）：

import numpy as np
from scipy.integrate import odeint
from scipy.interpolate import interp1d
import matplotlib.pyplot as plt

# 定义阶跃函数
def step_u(t):
    return 2 if t >=5 else 0

# x的微分方程模型
def model_x(x, t):
    u = step_u(t)
    dxdt = (-x + u)/2
    return dxdt

# y的微分方程模型，依赖x(t)
def model_y(y, t, x_interp):
    x_t = x_interp(t)  # 根据当前t插值得到x(t)
    dydt = (-y + x_t)/5
    return dydt

# 初始条件
x0 = 0
y0 = 0
t = np.linspace(0, 40, 1000)

# 先求解x(t)
x = odeint(model_x, x0, t)
# 创建x(t)的插值函数
x_interp = interp1d(t, x.flatten(), kind='linear')

# 再求解y(t)，传入x的插值函数
y = odeint(model_y, y0, t, args=(x_interp,))

# 绘图
plt.figure(figsize=(10,6))
plt.plot(t, x, 'r-', label='x(t)')
plt.plot(t, y, 'b-', label='y(t)')
plt.axvline(x=5, color='gray', linestyle='--', label='t=5 (阶跃点)')
plt.xlabel('Time')
plt.ylabel('Value')
plt.legend()
plt.grid(True)
plt.show()

关键说明 #

• 阶跃函数的实现：通过if t >=5 else 0直接在模型函数里判断时间，简单直观
• 耦合求解的优势：避免了分开求解时的时间点对齐问题，结果更准确
• 注意odeint的模型函数参数顺序：必须是(状态变量, 时间, 其他参数)，这也是你原来代码里model2参数顺序错误的地方

内容的提问来源于stack exchange，提问作者Amit