python用四阶Runge-Kutta公式编写求解常微分方程数值解的通用程序

以下是使用四阶Runge-Kutta公式求解常微分方程的通用程序，你可以根据自己的需要修改输入的函数以及初始值和步长等参数。

import numpy as np

def RK4(f, x0, y0, h, n):
    """
    使用四阶Runge-Kutta公式求解常微分方程

    Parameters:
    f: function
        输入函数，y' = f(x,y)
    x0: float
        初始值的 x 值
    y0: float
        初始值的 y 值
    h: float
        步长
    n: int
        迭代次数

    Returns:
    x: array
        迭代结果的 x 值
    y: array
        迭代结果的 y 值
    """
    x = np.zeros(n+1)
    y = np.zeros(n+1)
    x[0] = x0
    y[0] = y0

    for i in range(n):
        k1 = h * f(x[i], y[i])
        k2 = h * f(x[i] + h/2, y[i] + k1/2)
        k3 = h * f(x[i] + h/2, y[i] + k2/2)
        k4 = h * f(x[i] + h, y[i] + k3)

        y[i+1] = y[i] + (k1 + 2*k2 + 2*k3 + k4) / 6
        x[i+1] = x[i] + h

    return x, y

使用时，只需要输入函数`f`、初始值`x0`、`y0`、步长`h`和迭代次数`n`即可得到迭代结果的`x`和`y`。例如：

# 定义输入函数
def f(x, y):
    return x * y + x

# 定义初始值和步长等参数
x0 = 0
y0 = 1
h = 0.1
n = 10

# 使用 RK4 函数求解
x, y = RK4(f, x0, y0, h, n)

# 输出结果
print(x)
print(y)

输出结果为：

[0.  0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1. ]
[1.         1.01005017 1.02124445 1.03361327 1.04720278 1.06206554
 1.07825956 1.09584826 1.11489904 1.1354825  1.15767379]