python用四阶Runge-Kutta公式编写求解常微分方程数值解的通用程序，并输出误差和误差阶

好的，下面提供一个常微分方程数值解的通用程序，使用四阶Runge-Kutta公式进行求解：

import numpy as np

def rk4(f, y0, t0, tf, h):
    """
    使用四阶Runge-Kutta方法求解常微分方程初值问题的数值解

    Args:
        f: 函数，表示常微分方程 dy/dt = f(t, y) 中的右侧函数
        y0: 数组，表示初始时刻 t0 时的 y 值
        t0: 实数，表示初始时刻
        tf: 实数，表示终止时刻
        h: 实数，表示步长

    Returns:
        ys: 数组，表示求解出的数值解 y 的序列
        ts: 数组，表示求解出的时间 t 的序列
        error: 实数，表示求解误差
        order: 实数，表示误差阶
    """
    n = int((tf - t0) / h)
    ts = np.linspace(t0, tf, n+1)
    ys = np.zeros((n+1, len(y0)))
    ys[0] = y0
    for i in range(n):
        k1 = f(ts[i], ys[i])
        k2 = f(ts[i] + 0.5 * h, ys[i] + 0.5 * h * k1)
        k3 = f(ts[i] + 0.5 * h, ys[i] + 0.5 * h * k2)
        k4 = f(ts[i] + h, ys[i] + h * k3)
        ys[i+1] = ys[i] + h * (k1 + 2*k2 + 2*k3 + k4) / 6
    # 计算误差和误差阶
    y_exact = np.exp(-ts)
    error = np.max(np.abs(ys - y_exact))
    order = np.log2(error / (h**4))
    return ys, ts, error, order

使用时，可以调用上述函数并传入相应的参数，例如：

def f(t, y):
    return -y

y0 = np.array([1.0])
t0 = 0.0
tf = 5.0
h = 0.1

ys, ts, error, order = rk4(f, y0, t0, tf, h)
print("数值解：", ys)
print("误差：", error)
print("误差阶：", order)

其中，`f` 表示常微分方程 dy/dt = f(t, y) 中的右侧函数，`y0` 表示初始时刻 t0 时的 y 值，`t0` 和 `tf` 分别表示初始时刻和终止时刻，`h` 表示步长。函数将返回求解出的数值解 `ys` 和时间序列 `ts`，以及误差 `error` 和误差阶 `order`。

以上是一个简单的例子，可以根据具体的常微分方程和初始条件进行修改。

需要注意的是，误差和误差阶的计算是基于一个已知的精确解（例如上述例子中的 y=exp(-t)），如果没有提供精确解，就无法计算误差和误差阶。