Add topological_sort from SymPy
Add and use partial_lt()
Add max_sols function to get the maximal solutions
Fix max_sols to be more usable from solve2
Use maximal_sols from solve2

There are still issues.

conda create -n test -c defaults sympy python=2 numpy scipy --dry gives a
ValueError: cycle detected.