Why SymPy can't solve quadratic equation with complicated coefficients

Question

SymPy can easily solve quadratic equations with short simple coefficients. For example:

from pprint import pprint
from sympy import *
x,b,f,Lb,z = symbols('x b f Lb z')
eq31 = Eq((x*b + f)**2, 4*Lb**2*z**2*(1 - x**2))
pprint(eq31)
sol = solve(eq31, x)
pprint(sol)

But with a little bit larger coefficients - it can't:

from pprint import pprint
from sympy import *
c3,b,f,Lb,z = symbols('c3 b f Lb z')
phi,Lf,r = symbols('phi Lf r')
eq23 = Eq(
    (
        c3 * (2*Lb*b - 2*Lb*f + 2*Lb*r*cos(phi + pi/6))
        + (Lb**2 - Lf**2 + b**2 - 2*b*f + 2*b*r*cos(phi + pi/6) + f**2 - 2*f*r*cos(phi + pi/6) + r**2 + z**2)
    )**2,
    4*Lb**2*z**2*(1 - c3**2)
    )
pprint(eq23)
print("\n\nSolve (23) for c3:")
solutions_23 = solve(eq23, c3)
pprint(solutions_23)

Why?

Answer 1

This is not specific to Sympy - other programs like Maple or Mathematica suffer from same the problem: When solving an equation, solve needs to choose a proper solution strategy (see e.g. Sympy's Solvers) based on assumptions about the variables and the structure of the equation. These are choices are normally heuristic and often incorrect (hence no solution, or false strategies are tried first). Furthermore, the assumptions of variables is often to broad (e.g., complex instead of reals).

Thus, for complex equations the solution strategy often has to be given by the user. For your example, you could use:

sol23 = roots(eq23.lhs - eq23.rhs, c3)

Answer 2

Since symbolic solutions are supported, one thing you can do is solve the generic quadratic and substitute in your specific coefficients:

>>> eq = eq23.lhs-eq23.rhs
>>> a,b,c = Poly(eq,c3).all_coeffs()
>>> var('A:C')
(A, B, C)
>>> ans=[i.xreplace({A:a,B:b,C:c}) for i in solve(A*x**2 + B*x + C,x)]
>>> print filldedent(ans)
...

But you can get the same result if you just shut of simplification and checking:

>>> ans=solve(eq23,c3,simplify=False,check=False)

(Those are the really expensive parts of the call to solve.)