python - 为什么我执行的四阶龙格库塔方法不起作用？

Question

在这个网站和我的引用书中搜索之后，我发现我不知道为什么我的代码不能工作。

正如我的教授在类里面向我们展示的那样，我为质量 Spring 系统(带摊销)制作了四阶龙格库塔实现。但是，如您所见，生成的图形非常奇怪。

我最后写的代码是这样的:

#! /usr/bin/env python3
#-*- coding: utf-8 -*-

def f(data, t, x1, v1):
    from math import cos

    F = data["F"]
    c = data["c"]
    k = data["k"]
    m = data["m"]
    omega = data["omega"]

    return( [v1, (F*cos(omega*t) - c*v1 - k*x1)/m] )

def run(data = {}):
    xi, vi, ti = [data["u1"]], [data["v1"]], [data["t_ini"]]
    dt = data["dt"]

    while ti[-1] <= data["t_end"]:
        xn = xi[-1]
        vn = vi[-1]
        tn = ti[-1]

        K1 = f(data, t = tn, x1 = xn, v1 = vn)
        K1 = [dt*K1[i] for i in range(len(K1))]

        K2 = f(data, t = tn + 0.5*dt, x1 = xn + 0.5*K1[0], v1 = vn + 0.5*K1[1])
        K2 = [dt*K2[i] for i in range(len(K2))]

        K3 = f(data, t = tn + 0.5*dt, x1 = xn + 0.5*K2[0], v1 = vn + 0.5*K2[1])
        K3 = [dt*K3[i] for i in range(len(K3))]

        K4 = f(data, t = tn + dt, x1 = xn + K3[0], v1 = vn + K3[1])
        K4 = [dt*K4[i] for i in range(len(K4))]

        xn = xn + (K1[0] + 2*K2[0] + 2*K3[0] + K4[0])/6
        vn = xn + (K1[1] + 2*K2[1] + 2*K3[1] + K4[1])/6

        ti.append(tn+dt)
        xi.append(xn)
        vi.append(vn)

    return(ti, xi, vi)

这是由main.py文件引入的，它只包含程序的GUI和plot部分，功能是在类中推导的，所以我认为错误在Runge-Kutta本身。 (可能是我搞砸了一些愚蠢的事情。)

我尝试切换“xn”和“vn”中的 K，在 f() 中强制使用“F”和“c”值，重写所有内容并手动编写每个 K 的每个元素(如 K11， K12、K21 等)，但它仅给出指数结果。此外，将 f() 的返回切换为 numpy 数组也无济于事。

我发现了一些关于RK4方法的问题，但是我无法解决这个问题，也不明白哪里出了问题。我对这个方法有一些了解，但这实际上是我第一次实现它，所以欢迎任何帮助。

如果重要的话，我正在使用 python3 的 Anaconda 发行版。

Answer 1

难道是这一行？

vn = xn + (K1[1] + 2*K2[1] + 2*K3[1] + K4[1])/6

应该是:

vn = vn + (K1[1] + 2*K2[1] + 2*K3[1] + K4[1])/6