青春有我
有两个小细节:在 RK-2 循环中,您使用的是 z,但要存储您使用的值 i在初始代码中,您在更新 K2 时使用了 y+K1/K2,这是错误的。我看到你在第二个代码中修复它。所以,固定代码是:import matplotlib.pyplot as pltimport numpy as npfrom math import exp # exponential functiondy = lambda x, y: x * yf = lambda x: exp(x ** 2 / 2) # analytical solution functionx_final = 2# analytical solutionx_a = np.arange(0, x_final, 0.01)y_a = np.zeros(len(x_a))for i in range(len(x_a)): y_a[i] = f(x_a[i])plt.plot(x_a, y_a, label="analytical")# Container for step sizes dt /dtdt = 0.5x = 0y = 1print("dt = " + str(dt))print("x \t\t y (Euler) \t y (analytical)")print("%f \t %f \t %f" % (x, y, f(x)))n = int((x_final - x) / dt)x_e = np.zeros(n + 1)y_e = np.zeros(n + 1)x_e[0] = xy_e[0] = y#Plot Euler's methodfor i in range(n): y += dy(x, y) * dt x += dt print("%f \t %f \t %f" % (x, y, f(x))) x_e[i + 1] = x y_e[i + 1] = yplt.plot(x_e, y_e, "x-", label="Euler dt=" + str(dt))###################33# Runge-Kutta's method 2'nd order (RK2)x = 0y = 1print("dt = " + str(dt))print("x \t\t y (rk2) \t y (analytical)")print("%f \t %f \t %f" % (x, y, f(x)))n = int((x_final - x) / dt)x_r = np.zeros(n + 1)y_r = np.zeros(n + 1)x_r[0] = xy_r[0] = y# Plot the RK2for i in range(n): K1 = dt*dy(x,y) # Step 1 K2 = dt*dy(x+dt/2,y+K1/2) # Step 2 y += K2 # Step 3 x += dt print("%f \t %f \t %f" % (x, y, f(x))) x_r[i + 1] = x y_r[i + 1] = yplt.plot(x_r, y_r, "s-", label="RK2 dt=" + str(dt))plt.title("numerical differential equation problem")plt.xlabel("x")plt.ylabel("y")plt.legend()plt.show()