使用 RK4 的二阶 ODE答案

【问题标题】：Second order ODE with RK4使用 RK4 的二阶 ODE
【发布时间】：2017-06-21 21:49:49
【问题描述】：

我需要绘制一个显示 2 个变量的图形，以及带有 RK4 的二阶 ODE，到目前为止我已经完成了这个

from numpy import arange
from pylab import plot,xlabel,ylabel,show
Qger = 400
K = 20
T1 = 150
T2 = 60
N = 1000
h = (T2-T1)/N
rpoints = arange(6.0,8.0,h)
xpoints = []
x = 423
def df(s,t):
    dTdt = -Qger*t/(2*K) + 172.8/t
    return dTdt

for r in rpoints:
    xpoints.append(x)
    k1 = h*df(x,r)
    k2 = h*df(x+0.5*k1,r+0.5*h)
    k3 = h*df(x+0.5*k2,r+0.5*h)
    k4 = h*df(x+k3,r+h)
    x += (k1+2*k2+2*k3+k4)/6
pylab.plot(rpoints,xpoints)
pylab.xlabel("Raio")
pylab.ylabel("Temperatura")
pylab.show

但对于一阶 ODE，这是 RK4，因为我不知道并由手，但我不能这样做，也不使用 scipy，所以谁能向我解释如何集成这个功能或使用 RK4 和二阶 ODE。函数如下。

这是函数，只有T和r是变量，其余都是0

【问题讨论】：

scipy documentation for odeint 有一个二阶 ODE 解的详细示例。研究该示例并返回（也许在另一个问题中......）以获取有关特定问题的帮助。正如你的问题是一个糟糕的问题......

标签： python computer-science scientific-computing runge-kutta

【解决方案1】：

您应该能够将上述内容以“半离散”形式表示，也就是说 dT/dt 仅表示相对于 r 的偏导数。如果您可以找到等效于 dT/dt 的术语的数值或其他近似值，即 dT/dt= df(r,...) 的 RHS，则可以应用显式 RK4。

在这种方法中，您的时间步进方法 (RK4) 仅适用于温度相对于时间的一阶导数。

【讨论】：