def rk3(f, y0, t0, tf, h):

    t_values = np.arange(t0, tf + h, h)

    y_values = np.zeros_like(t_values)

    y_values[0] = y0

    for i in range(1, len(t_values)):

        k1 = f(t_values[i-1], y_values[i-1])

        k2 = f(t_values[i-1] + h/2, y_values[i-1] + h*k1/2)

        k3 = f(t_values[i-1] + h, y_values[i-1] - 2*h*k1 + h*k2)

        y_values[i] = y_values[i-1] + h*(k1 + 2*k2 + k3)/6

    return t_values, y_values

def rk3(f, y0, t0, tf, h):

    t_values = np.arange(t0, tf + h, h)

    y_values = np.zeros_like(t_values)

    y_values[0] = y0

    for i in range(1, len(t_values)):

        k1 = f(t_values[i-1], y_values[i-1])

        k2 = f(t_values[i-1] + h/2, y_values[i-1] + h*k1/2)

        k3 = f(t_values[i-1] + h, y_values[i-1] - h*k1 + 2*h*k2)

        y_values[i] = y_values[i-1] + h*(k1 + 4*k2 + k3)/6

    return t_values, y_values

def rk3(f, y0, t0, tf, h):

    t_values = np.arange(t0, tf + h, h)

    y_values = np.zeros_like(t_values)

    y_values[0] = y0

    for i in range(1, len(t_values)):

        k1 = f(t_values[i-1], y_values[i-1])

        k2 = f(t_values[i-1] + h, y_values[i-1] + h*k1)

        k3 = f(t_values[i-1] + h, y_values[i-1] + h*k2)

        y_values[i] = y_values[i-1] + h*(k1 + k2 + k3)/3

    return t_values, y_values

def rk3(f, y0, t0, tf, h):

    t_values = np.arange(t0, tf + h, h)

    y_values = np.zeros_like(t_values)

    y_values[0] = y0

    for i in range(1, len(t_values)):

        k1 = f(t_values[i-1], y_values[i-1])

        k2 = f(t_values[i-1] + h/3, y_values[i-1] + h*k1/3)

        k3 = f(t_values[i-1] + 2*h/3, y_values[i-1] + 2*h*k2/3)

        y_values[i] = y_values[i-1] + h*(k1 + 3*k2 + 3*k3)/7

    return t_values, y_values

def rk3(f, y0, t0, tf, h):

    t_values = np.arange(t0, tf + h, h)

    y_values = np.zeros_like(t_values)

    y_values[0] = y0

    for i in range(1, len(t_values)):

        k1 = f(t_values[i-1], y_values[i-1])

        k2 = f(t_values[i-1] + h/2, y_values[i-1] + h*k1/2)

        k3 = f(t_values[i-1] + h, y_values[i-1] + h*k2)

        y_values[i] = y_values[i-1] + h*(k1 + 4*k2 + k3)/6

    return t_values, y_values