CODEWITHSUNDEEP
pythonpython-2.7sympy
A_Pollo
A_Pollo

Reputation: 41

sympy solve() unable to find solution

I have a fairly basic code right now that's supposed to solve the equation shown below but it's not converging to a solution. It just hangs on a pulsing "_"

from sympy import *

gamma = 1.4
M_a = 1.0

y = Symbol('y', real=True) 
eqn = Eq((1.0/((gamma/2.0)*(y**2.0))) * ((((1.0 + ((gamma - 1.0)/2.0)*(y**2.0))/(1.0 + ((gamma - 1.0)/2.0)*(M_a**2.0)))**(gamma/(gamma - 1.0))) - 1.0) ,(-0.5704/((1.0 - (y**2.0))**(1.0/2.0))))
print solve(eqn, y)

simply printing the equation yields

1.42857142857143*y**(-2.0)*((0.166666666666667*y**2.0 + 0.833333333333333)**3.5- 1.0) == -0.5704*(-y**2.0 + 1.0)**(-0.5)

which when plugged into wolfram or maple yields the correct solution. ~= 0.696256

So I'm trying to figure out why sympy is unable to solve the equation.

The equation should look like this Picture

If sympy can not be used to solve this equation what can I use instead?

Thanks

Phil

Upvotes: 3

Views: 3832

Answers (1)

asmeurer
asmeurer

Reputation: 91680

If you only care about numeric solutions, don't use solve. It is trying really hard to find a symbolic solution. But symbolic solutions to polynomials can often be quite complicated, and in general they do not exist.

If you just want a numeric solution, use nsolve:

In [60]: nsolve(eqn, 0)
Out[60]: mpc(real='0.69625557901731519', imag='-3.4211388289180104e-49')

(the 0 is a guess for the solution). The imaginary part there is small enough to ignore. You can remove it with N(solution, chop=True).

If you do care about symbolic solutions, another piece of advice is to avoid floating point exponents. Use rational or integer exponents. In order to find the roots of a polynomial symbolically, it must first be converted to an actual polynomial, i.e., with integer coefficients. This means that floating point coefficients must be converted to rational numbers. But this can often result in huge numbers, and hence very large degrees (this I suspect is why this ends up hanging for you, consider:

In [61]: gamma
Out[61]: 1.4

In [62]: Rational(gamma)
Out[62]:
3152519739159347
────────────────
2251799813685248

Upvotes: 4

Related Questions