我们应该编写一个程序来使用 4 阶 Runge-Kutta 数值求解以下初始值问题。该算法不是问题,我可以在完成后发布我的解决方案。
问题是,将它干净地分离成我可以放入 Runge-Kutta 的东西。
e^(-x') = x' −x + exp(−t^3)
x(t=0) = 1
知道这叫什么类型的 ODE?或解决此问题的方法?与数学相比,我对 CS 技能和编程数值方法更有信心……所以对这个问题的任何见解都会有所帮助。
更新:如果有人对解决方案感兴趣,代码如下。我认为这是一个有趣的问题。
import numpy as np
import matplotlib.pyplot as plt
def Newton(fn, dfn, xp_guess, x, t, tolerance):
iterations = 0
value = 100.
max_iter = 100
xp = xp_guess
value = fn(t, x, xp)
while (abs(value) > tolerance and iterations < max_iter):
xp = xp - (value / dfn(t,x,xp))
value = fn(t,x,xp)
iterations += 1
root = xp
return root
tolerance = 0.00001
x_init = 1.
tmin = 0.0
tmax = 4.0
t = tmin
n = 1
y = 0.0
xp_init = 0.5
def fn(t,x,xp):
'''
0 = x' - x + e^(-t^3) - e^(-x')
'''
return (xp - x + np.e**(-t**3.) - np.e**(-xp))
def dfn(t,x,xp):
return 1 + np.e**(-xp)
i = 0
h = 0.0001
tarr = np.arange(tmin, tmax, h)
y = np.zeros((len(tarr)))
x = x_init
xp = xp_init
for t in tarr:
# RK4 with values coming from Newton's method
y[i] = x
f1 = Newton(fn, dfn, xp, x, t, tolerance)
K1 = h * f1
f2 = Newton(fn, dfn, f1, x+0.5*K1, t+0.5*h, tolerance)
K2 = h * f2
f3 = Newton(fn, dfn, f2, x+0.5*K2, t+0.5*h, tolerance)
K3 = h * f3
f4 = Newton(fn, dfn, f3, x+K3, t+h, tolerance)
K4 = h * f4
x = x + (K1+2.*K2+2.*K3+K4)/6.
xp = f4
i += 1
fig = plt.figure()
ax = fig.add_subplot(1,1,1)
ax.plot(tarr, y)
plt.show()