欧拉法解微分方程

欧拉法解微分方程

本文介绍如何使用简单的欧拉法求解微分方程，大部分内容出自吴一东老师在他的B站个人空间发布的课程

方法介绍

对于一个一般的微分方程：

\[\begin{cases} \begin{aligned} \frac{\mathrm{d} y}{\mathrm{d} t} &= f(y(t), t)\\ y(0) &= y_0 \\ \end{aligned} \end{cases} \]

​ 假如我们很难得到他的解析解或者不存在解析解，那我们可以尝试使用欧拉法将微分方程利用泰勒展开构造成一个递推式，通过程序将所有的点求出，进而得到他的数值解。以下说明求法：

1. 区间划分

   对于一个微分方程我们只关心函数的一段区域内的数值解。

   \(t \in [0, T]\) 将区间划成\(N\)​段，则有\(\Delta t = \frac{T}{N}\)

   \(t_1 = 0, t_2 = \Delta t, ..., t_{n} = (n-1)\Delta t ,...,t_{N+1} = T\)

2. 泰勒展开

   \[\begin{aligned} y(t+\Delta t) &= y(t) + \frac{\mathrm{d} y} {\mathrm{d} t}(t)\Delta t + O(\Delta t^2) \\ &=y(t) + f(y(t), t)\Delta t + O(\Delta t^2) \end{aligned} \]

3. 得到递推公式

   \[\begin{aligned} y(t_{n+1}) &= y(t_n + \Delta t)\\ &=y(t_n) + f(y(t_n), t_n)\Delta t + O(\Delta t^2) \end{aligned} \]

4. 近似计算

   \[y(t_n) \to t_n\\ y_{n + 1} = y_n +f(y_n, t)\Delta t \]

示例

例1：

\[\begin{cases} \begin{aligned} f(y, t) &= y\\ \frac{\mathrm{d} y}{\mathrm{d} t} &= y\\ y(0) &= 1 \end{aligned} \end{cases} \]

由题目其实我们很容易得知，\(y = e^t\) 所以我们可以选择解析解对数值解方法进行验证

按上文求解方法我们求其数值解

\(y_{n + 1} = y_n + y_n\Delta t\)​

求解代码如下：

import numpy as np
import matplotlib.pyplot as plt

# 设置初始条件
T = 5
N = 10000
dt = T / N
t = np.linspace(0, T, N + 1)

# 边界条件
y = np.zeros(N + 1)
y[0] = 1

for i in range(0, N):
    y[i + 1] = y[i] + y[i] * dt

ex = np.exp(t)
plt.plot(t, y, color='red')
plt.plot(t, ex, color='blue')
plt.show()