显示欧拉法求解 y' = y 问题

简介

最近 无聊看了一下 数值解法

已知 $\frac{dy}{dx} = y$ 我们知道其有一个解析解为 $y = e^x$
同时我们知道 其初值 $y(0) = 1$, x 的范围为 $0<= x <= 1$
我们需要求解其数值解，也就是$y$ 值在 $0<= x <= 1$ 的分布

我们猜想显示欧拉是由 泰勒展开式进行构建的

$$f(x) = f(x_0) +f'(x_0)(x-x_0)+\frac{f''(x_0)}{2!}(x-x_0)^2+···+\frac{f^{(n)}(x_0)}{n!}(x-x_0)^n+R_n(x)$$

也就是可以简单等同于求解

$$y(x) = y(0) + \frac{dy}{dx(x == 0)}(x - 0) + 余项（即误差）$$

但是关于 $\frac{dy}{dx(x == 0)}$ 比较难以求解，采用 $\frac{\mathbf{y}\left( k+1 \right)-\mathbf{y}\left( k \right)}{h}$ 来近似求解
其中 h 即步长步长越小带来的精度越来，同时带来的计算量越大
这里比较特殊由于题目中已知了$\frac{dy}{dx} = y$所以我们可以使用

$$y_{n+1} = y_{n} + h(y')$$

$$y_{n+1} = y_{n} + h y_{n}$$

那么我们对于其数值解就可以比较简单求解出来了

code

import math
from matplotlib import pyplot as plt

t_0 = 0
y_0 = 1
tau = 0.01
i = 1
solve = []
Euler = []
t = []
while i < 1000:
    if i == 1:
        y_n = y_0
        t_n = t_0
    Euler.append(y_n)
    solve.append(math.exp(t_n)) # 精确解
    t.append(t_n)
    func = y_n
    y_n = y_n + tau * func
    t_n = t_n + tau
    i += 1

plt.plot(t, Euler, c='green', label=' Euler method')
plt.plot(t, solve, c='red', label=' accuracy')
plt.fill_between(t, solve, Euler, facecolor='blue', alpha=0.2)
plt.title('Euler method', fontsize=19)
plt.xlabel('t', fontsize=19)
plt.ylabel('y', fontsize=19)
plt.legend()
plt.show()

img

其中紫色部分是误差

参考链接

https://zhuanlan.zhihu.com/p/131303043
https://zhuanlan.zhihu.com/p/422074759