我在Python中有一个ODE系统要解决，它分为两部分：过去时间解和未来时间解：
以下是ODE系统：

# defining system
def sys_to_solve(t, funcs):
    return [
        funcs[0] * funcs[1],
        (
            3
            / (3 + 2 * omega_BD)
            * psi_x
            / funcs[0]
            * H0 ** 2
            * Omega_m
            * funcs[2] ** -3
            - (funcs[1]) ** 2
            - 3 * H(funcs[2], funcs[1], funcs[0]) * funcs[1]
        ),
        funcs[2] * H(funcs[2], funcs[1], funcs[0]),
    ]

我知道3个变量的现值（在t_now=0时）。
我正在寻找初始值（t=t_开始），一旦解决，将匹配当前值。
当前时间的已知值为x_0, x_1, x_2
以下是过去的解决方案实现：

sol_past = solve_ivp(
        sys_to_solve,
        [t_begin, t_now],
        [z1, z2 , z3],
        method="RK45",
        t_eval=t_past,
        rtol=1e-12,
        )

为了找到这些初始值，有人告诉我，要找到合适的z1,z2,z3以上，方法如下：

def objective(x):
    sol_past = solve_ivp(
        sys_to_solve,
        [t_begin, t_now],
        [x[0], x[1], x[2]],
        method="RK45",
        t_eval=t_past,
        rtol=1e-12,
    )
    y1 = sol_past.y[0,:]
    y2 = sol_past.y[1,:]
    y3 = sol_past.y[2,:]

    return np.abs(y1[-1] - x_0), np.abs(y2[-1] - x_1), np.abs(y3[-1] - x_2)

z1,z2,z3 = fsolve(objective, [6.45, 1.317, 0.3] , xtol=1e-12)

但是6.45, 1.317 and 0.3的值是通过随机或逐步的方法找到的。
现在，我想证明并找到这些值与数值方法建立良好，而不仅仅是一个伪随机微调。
我听说这种方法被称为“射击法”。
我考虑最小化方法，如Chi 2的最小化：谁能帮我找到与当前值x_0,x_1,x_3最匹配的类似值6.45, 1.317 and 0.3？

编辑：@Lutz Lehmann给出的解决方案在前两个参数（Psi = G * phi和F0 = d_Phi_0 / dt /Phi 0）上运行良好，但我对比例因子H（t）= da/dt / a没有太多的数学知识，因为我只能访问a（t）而不能访问da/dt /a。

这里是H（t）上的坏匹配：
Bad match witht expected H(z) in blue point
然而，对于比例因子与Psi(t)和F0 = d_Phi_0 / Phi0的曲线，我得到了一致的曲线
Well match with a(t
作为信息，我们得到funcs[0] = Psi（t），funcs 1 =F（t）和funcs 2 =a（t）
因此，可以肯定的是，射击方法可以解决H（t）当前值过高的问题。