电磁场中微积分方程的仿真实验(欧拉算法)
一、Euler公式
以当前点的值为初值,当前点的导数为导数,计算下一个步长的点的值。
代码实现如下:
from matplotlib import pyplot as plt
import numpy as np
h=0.01
x = 0
xx=np.zeros(100)
yy=np.zeros(100)
y = 1
for i in range (100):
y=h*(y-2.0*x/y)+y
x=x+h
xx[i]=x
yy[i]=y
print(xx[i],yy[i])
plt.plot(xx,yy, ls='-', lw=2, label='OULA', color='blue')
plt.legend()
plt.show()
实验运行截图如下
由于是用前一点的导数来计算,因此明显大于精确解,且之间的差越来越大
二、改进欧拉法
以当前点的值为初值,当前点的导数和下一个步长的点的导数平均数为导数,计算下一个步长的点的值,不需要进行迭代计算,只用计算一步即可
from matplotlib import pyplot as plt
import numpy as np
t=0.1
x0 = 0
y0 = 1
yb=1
x=np.zeros(100)
y=np.zeros(100)
xx=np.zeros(100)
yy=np.zeros(100)
for i in range (100):
y0 = 1.1 * y0 - 0.2 * x0 / y0
x0 = x0 + t
x1 = x0 + t
k0 = y0 - 2 * x0 / y0
yb = y0 + t * k0
k1 = yb - 2 * x1 / yb
y1 = y0 + t * (k0 + k1) / 2
print(x0, y1)
x[i] = x0
y[i] = y0
xx[i] = x1
yy[i] = y1
plt.plot(xx, yy, color='blue') #改进后用蓝色表示
plt.scatter(xx,yy)
plt.plot(x,y,color='orange')
plt.legend()
plt.show()
由于综合了两点的导数,因此求出来的解跟精确解相差不大,结果如截图所示:



