Solving already solved equations takes a long time

[drhagen]

When solving a system of nonlinear equations, if the equations are already solved for all unknowns, then it should be very fast to return the solution. However, this is not the case. For example, this script:

import time
from sympy import *

var('kn v ah ldkd_1 ldkd_2 eh_1 lc_1 lc_2 rec_1 rec_2')
var('lf_1 lf_2 rh_1 rh_2 sc_1 sc_2 sh_1 sh_2 lrkd_1 lrkd_2')
var('dv_1 dv_2 d m')

dependent_vars = (var('ka kl_D1 kn2_T1 kf_D1T1 kf_L1R1 kl_R1 kl_L1R1 kl_L1 kl_S1 t_R1')
+ var('L1_0 t_S1 L1_R1_0 R1_0 L1_S1_0 S1_0 kS_1 ky_L1 ky_R1 a_1_0')
+ var('kn2_T2 kf_D1T2 kf_L2R2 kl_R2 kl_L2R2 kl_L2 kl_S2 t_R2')
+ var('L2_0 t_S2 L2_R2_0 R2_0 L2_S2_0 S2_0 kS_2 ky_L2 ky_R2 a_2_0'))


eqs = [
    Eq(ka, log(2) / (ah * d)),
    Eq(kl_D1, log(2) / (eh_1 * d)),
    Eq(kn2_T1, floor(dv_1 / 2) * kn),
    Eq(kf_D1T1, ldkd_1 * kn),
    Eq(kf_L1R1, lrkd_1 * kn),
    Eq(kl_R1, log(2) / (rh_1 * m)),
    Eq(kl_L1R1, kl_R1),
    Eq(kl_L1, log(2) / (lf_1 * m)),
    Eq(kl_S1, log(2) / (sh_1 * m)),
    Eq(t_R1, rec_1 * v),
    Eq(L1_0, lc_1 * v),
    Eq(t_S1, sc_1 * v),
    Eq(L1_R1_0, (kn * L1_0 * t_R1) / (v * (kf_L1R1 + kl_L1R1 + kn * L1_0 / v))),
    Eq(R1_0, t_R1 - L1_R1_0),
    Eq(L1_S1_0, (kn * L1_0 * t_S1) / (v * (kf_L1R1 + kl_L1 + kn * L1_0 / v))),
    Eq(S1_0, t_S1 - L1_S1_0),
    Eq(kS_1, (kl_S1 * S1_0 + kn * L1_0 * S1_0 / v - kf_L1R1 * L1_S1_0) / R1_0),
    Eq(ky_L1, kn * L1_0 * R1_0 / v - kf_L1R1 * L1_R1_0 + kn * L1_0 * S1_0 / v - kf_L1R1 * L1_S1_0 + kl_L1 * L1_0),
    Eq(ky_R1, kn * L1_0 * R1_0 / v - kf_L1R1 * L1_R1_0 + kl_R1 * R1_0 + kS_1 * R1_0),
    Eq(a_1_0, L1_R1_0),
    Eq(kn2_T2, floor(dv_2 / 2) * kn),
    Eq(kf_D1T2, ldkd_2 * kn),
    Eq(kf_L2R2, lrkd_2 * kn),
    Eq(kl_R2, log(2) / (rh_2 * m)),
    Eq(kl_L2R2, kl_R2),
    Eq(kl_L2, log(2) / (lf_2 * m)),
    Eq(kl_S2, log(2) / (sh_2 * m)),
    Eq(t_R2, rec_2 * v),
    Eq(L2_0, lc_2 * v),
    Eq(t_S2, sc_2 * v),
    Eq(L2_R2_0, (kn * L2_0 * t_R2) / (v * (kf_L2R2 + kl_L2R2 + kn * L2_0 / v))),
    Eq(R2_0, t_R2 - L2_R2_0),
    Eq(L2_S2_0, (kn * L2_0 * t_S2) / (v * (kf_L2R2 + kl_L2 + kn * L2_0 / v))),
    Eq(S2_0, t_S2 - L2_S2_0),
    Eq(kS_2, (kl_S2 * S2_0 + kn * L2_0 * S2_0 / v - kf_L2R2 * L2_S2_0) / R2_0),
    Eq(ky_L2, kn * L2_0 * R2_0 / v - kf_L2R2 * L2_R2_0 + kn * L2_0 * S2_0 / v - kf_L2R2 * L2_S2_0 + kl_L2 * L2_0),
    Eq(ky_R2, kn * L2_0 * R2_0 / v - kf_L2R2 * L2_R2_0 + kl_R2 * R2_0 + kS_2 * R2_0),
    Eq(a_2_0, L2_R2_0),
]

start_time = time.time()
print(solve(eqs, dependent_vars, rational=False, manual=True, dict=True))
print("--- {} seconds ---".format((time.time() - start_time)))

This script has many equations, but they are already solved for all unknowns. All that is needed should be a little substitution. However, it take 17 seconds to run on my machine, which seems excessive.

[asmeurer]

If you add simplify=False, check=False, and remove manual=True, the time goes down to 6 seconds.

It looks like the simplify time is spent in cancel, which is slow because it has so many variables. The check code is slow in subs. For some reason if you do simplify=False, check=False, manual=True, it hangs.

With check=False, manual=False, simplify=False, it is spending its time in the algorithm itself, specifically _solve_reduced_system. We can probably be doing better checking for trivial cases like this so that the solution is returned immediately.

[asmeurer]

I tried throwing it at nonlinsolve. There are some duplicated variables. After removing them, it hangs. I think it is also stuck doing checking and simplification.

[drhagen]

I removed the duplicate variables from the code in the original post.

[asmeurer]

I think this is also a place where splitting the connected components of the system before solving would help a lot.

[asmeurer]

It looks like _solve_reduced_system doesn't do any preprocessing at all. If it gets a polynomial that is already in one variable, it shouldn't pass it through to the groebner basis algorithm.

[oscarbenjamin]

I tried with manual=False, check=False, simplify=False and compare master with #16550 and it does get faster. On master it's 10 seconds and with #16550 it's 4.8 seconds.

Note that #16550 only splits fully disconnected components so any equations that have any unknowns in common will still be solved together.

[smichr]

One step of preprocessing that is especially applicable to systems passed as Eq objects with a Symbol on the left is to apply repsort from #6257. A successful return indicates that the system is already implicitly defined for the Symbols. Backsubstitution will then resolve the Symbols fully to give the rhs independent of the lhs:

    # using eqs and dependent_vars defined in OP
    # and repsort from #6257
    r = repsort(*[i.args for i in eqs])
    sol = {}
    for i in range(len(r)-1,-1,-1):
        k, v = r[i]
        sol[k] = v.xreplace(sol)
        #assert not sol[k].has(*dependent_vars)

So this is not "solving" as much as it is fully resolving the system fully to give variables in terms of independent variables. The given system of equations takes less than a second to resolve with the above code snippet.