常微分方程模型

常微分方程

数学建模中的微分方程模型是用来描述系统动态行为的重要工具，通过建立微分方程模型，可以刻画系统随时间或其他变量变化的规律。微分方程模型在自然科学、工程、经济学、生物学等领域有广泛应用，能够帮助人们理解和预测复杂系统的行为。

微分方程模型的作用

1. 描述动态系统：微分方程可以描述物理、化学、生物等系统的动态变化，例如人口增长、化学反应速率、传染病传播等。
2. 预测系统行为：通过求解微分方程，可以预测系统在未来的行为，为决策提供依据。
3. 优化与控制：在工程中，通过微分方程模型可以设计和优化控制策略，例如自动控制系统中的PID控制器设计。
4. 解释现象：微分方程模型可以帮助解释观察到的现象，通过拟合实验数据建立模型。

应用问题示例

1. 人口增长模型：

   • 描述人口随时间变化的规律。
   • 常见的模型如Malthus模型和Logistic模型。

2. 传染病模型：

   • 描述传染病在群体中的传播规律。
   • 常见的模型如SIR模型和SEIR模型。

3. 物理运动模型：

   • 描述物体在力作用下的运动。
   • 常见的模型如自由落体运动、阻尼振动等。

4. 生态系统模型：

   • 描述不同物种之间的相互作用和演化。
   • 常见的模型如捕食-被捕食模型（Lotka-Volterra模型）。

5. 经济学模型：

   • 描述经济指标如GDP、通货膨胀、失业率等随时间的变化。
   • 常见的模型如Solow经济增长模型。

示例：SIR传染病模型

SIR模型是一种经典的传染病模型，用于描述易感者(Susceptible)、感染者(Infected)和康复者(Recovered)在时间上的变化。

模型的微分方程为：

[ \frac{dS}{dt} = -\beta SI ] [ \frac{dI}{dt} = \beta SI - \gamma I ] [ \frac{dR}{dt} = \gamma I ]

其中， ( S ) 是易感者数量， ( I ) 是感染者数量， ( R ) 是康复者数量， ( \beta ) 是传播率， ( \gamma ) 是康复率。

Python代码实现

以下代码展示了如何使用Python中的scipy.integrate.solve_ivp求解SIR模型，并绘制结果：

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

# 定义SIR模型
def sir_model(t, y, beta, gamma):
    S, I, R = y
    dSdt = -beta * S * I
    dIdt = beta * S * I - gamma * I
    dRdt = gamma * I
    return [dSdt, dIdt, dRdt]

# 初始条件
S0, I0, R0 = 0.99, 0.01, 0.0  # 初始易感者、感染者、康复者比例
y0 = [S0, I0, R0]

# 参数
beta = 0.3  # 传播率
gamma = 0.1  # 康复率

# 定义时间范围
t_span = (0, 160)
t_eval = np.linspace(0, 160, 400)

# 求解SIR模型
solution = solve_ivp(sir_model, t_span, y0, args=(beta, gamma), t_eval=t_eval)

# 绘制结果
plt.plot(solution.t, solution.y[0], label='Susceptible')
plt.plot(solution.t, solution.y[1], label='Infected')
plt.plot(solution.t, solution.y[2], label='Recovered')
plt.xlabel('Time')
plt.ylabel('Proportion')
plt.title('SIR Model')
plt.legend()
plt.grid()
plt.show()

代码解释

1. 定义SIR模型：

   • sir_model(t, y, beta, gamma)：定义SIR模型的微分方程。
   • dSdt、dIdt、dRdt：分别表示易感者、感染者和康复者的变化率。

2. 初始条件：

   • S0, I0, R0：初始易感者、感染者和康复者的比例。

3. 参数：

   • beta：传播率。
   • gamma：康复率。

4. 定义时间范围：

   • t_span：定义时间范围。
   • t_eval：定义评估的时间点。

5. 求解SIR模型：

   • 使用solve_ivp函数，指定SIR模型、时间范围、初始条件和参数。

6. 绘制结果：

   • 使用Matplotlib绘制易感者、感染者和康复者随时间变化的曲线。

一个实例-物体冷却问题

要解决这个问题，我们需要基于尸体冷却的模型来判断嫌疑人张某是否有可能在案发现场。尸体冷却常用的模型是Newton冷却定律，其形式如下：

$$\frac{dT}{dt} = -k (T - T_{\text{env}})$$

其中：

• $( T(t) )$ 是尸体在时间 $( t )$ 的温度
• $( T_{\text{env}} )$ 是环境温度
• $( k )$ 是冷却常数
• $( t )$ 是时间

根据Newton冷却定律，可以通过求解微分方程来模拟尸体冷却过程。以下是Python代码，用于基于给定条件求解微分方程，并判断张某是否有可能在案发现场。

问题中的条件

• 尸体发现时间：晚上 7:30
• 法医到达现场时间：晚上 8:20
• 法医测得尸体温度：32.6℃
• 尸体搬走时间：法医到达1小时后，即晚上9:20
• 搬走时测得尸体温度：31.4℃
• 环境温度：21.1℃
• 张某离开办公室时间：下午5:00
• 从张某到受害者现场步行需5分钟

Python代码

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

# 定义冷却常数
def cooling_law(t, T, k, T_env):
    return -k * (T - T_env)

# 初始条件
T_env = 21.1  # 环境温度
T1 = 32.6     # 法医到达时尸体温度
T2 = 31.4     # 尸体搬走时温度
t1 = 0        # 法医到达时间，假设为t=0
t2 = 60       # 搬走时间，假设为t=60分钟

# 根据温度变化计算冷却常数k
k = -np.log((T2 - T_env) / (T1 - T_env)) / (t2 - t1)

# 定义时间范围
t_span = (0, 120)  # 模拟时间范围，0到120分钟
t_eval = np.linspace(0, 120, 1000)

# 初始条件，法医到达时的温度
y0 = [T1]

# 求解微分方程
solution = solve_ivp(cooling_law, t_span, y0, args=(k, T_env), t_eval=t_eval)

# 获取时间和温度数据
t = solution.t
T = solution.y[0]

# 找出7:30之前的时间点
t_offset = 50  # 法医到达时间为7:30，搬走时间为8:20
t_actual = t - t_offset  # 时间往前偏移50分钟，即7:30之前

# 绘制结果
plt.plot(t_actual, T, label='Body Temperature')
plt.axvline(x=-140, color='r', linestyle='--', label='5:00 PM')
plt.axvline(x=-90, color='b', linestyle='--', label='7:30 PM')
plt.xlabel('Time (minutes from 5:00 PM)')
plt.ylabel('Temperature (°C)')
plt.title('Body Cooling Curve')
plt.legend()
plt.grid()
plt.show()

# 判断是否可信
if T[np.abs(t_actual + 90).argmin()] <= T1:
    print("张某的证言是可信的")
else:
    print("张某的证言不可信")

代码解释

1. 定义冷却定律：根据Newton冷却定律定义微分方程。
2. 初始条件：包括环境温度、法医测得的尸体温度、搬走时的尸体温度等。
3. 计算冷却常数：利用法医测得的两次温度和时间差计算冷却常数 ( k )。
4. 求解微分方程：使用solve_ivp求解微分方程，得到温度随时间的变化。
5. 绘制结果：绘制尸体温度随时间变化的曲线，并标出关键时间点（如张某离开办公室的时间）。
6. 判断可信性：通过比较法医到达时间前尸体温度，判断张某的证言是否可信。

结果

通过运行上述代码，可以得到尸体冷却的温度曲线，并且判断张某是否有可能在案发现场。如果法医到达之前的尸体温度符合计算结果，则张某的证言可信；否则，不可信。