1

I encountered a problem when I using Sympy to solve problems Here is my code:

from math import pi, hypot
from sympy import solve, solveset, sqrt, Symbol
one_x=-0.08
one_y=1.28
second_x=0
second_y=0
second_r=7
one_r=7.3
slopes=-16.0000000000  (maybe more trailing 0s)
intercepts=0.0
x=Symbol('x')
solveset(sqrt((x-second_x)**2+(slope*x+intercept-second_y)**2)+second_r-one_r-sqrt((x-one_x)**2+(slope*x+intercept-one_y)**2),x)

That's only part of my code but it raises a lot of errors but instead, I replaced all of the variables with its value like

x=Symbol('x')

solveset(sqrt((x)**2+((-16)*x)**2)+7-7.3-sqrt((x+0.08)**2+((-16)*x-1.28)**2),x)

It works nicely and i can get an output {-0.0493567429232771}

I think It's because of the type of slopes(-16 compared with -16.000000), I really wanna know why an equation with float number cannot be calculated, and how I can fix it (cause I need it to be more precise so I cannot just ignore the number after the dot ) Thanks so much!

Jasper Zhou
  • 449
  • 1
  • 5
  • 13

1 Answers1

3

SymPy + algebraic equation + floating point numbers => trouble. Floating point math does not work like normal math, and SymPy is designed for the latter. Small things like 16 (integer) versus 16.0 (float) make a lot of difference in solving equations with SymPy: ideally, you would have no floating point numbers there, creating exact rational numbers instead, like this.

from sympy import S 
one_x = S('-0.08')

However, you have floating point data and are looking for a floating point solution. This makes SymPy the wrong tool for the job. SymPy is for doing math with symbols, not for crunching floating point numbers. The correct solution is to use an appropriate solve from SciPy, such as brentq. It takes a bracketing interval as an input (where the function has different signs at both ends). For example: 

from scipy.optimize import brentq
eq = lambda x: np.sqrt((x-second_x)**2 + (slope*x+intercept-second_y)**2) + second_r - one_r - np.sqrt((x-one_x)**2 + (slope*x + intercept - one_y)**2)
brentq(eq, -10, 10)   # returns -0.049356742923277075

If you stick with SymPy, that means your equation is going outsourced to mpmath library, which is much more limited in the numerical root finding and optimization. To get a solution to converge with its methods, you'll need a really good starting point: apparently, one_x/2 is such a point. 

from sympy import sqrt, Symbol, nsolve
# ... as in your code
nsolve(sqrt((x-second_x)**2+(slope*x+intercept-second_y)**2)+second_r-one_r-sqrt((x-one_x)**2+(slope*x+intercept-one_y)**2), one_x/2)

returns -0.0493567429232771.

By using sympy.solveset, which is intended for symbolic solution, you deprive yourself not only of SciPy's powerful numeric solvers, but also of an opportunity to set a good starting value for the numeric search which sympy.nsolve provides. Hence the lack of convergence in this numerically tricky problem. By the way, this is what makes it numerically tricky: the function is nearly constant most of the time, with one rapid change.

tricky

  • You 're right, now I know that Sympy is not a correct choice to solve my problem, and I used Scipy you mentioned above and solved my problems, and it works perfectly, your description is so clear and detailed, thanks so much!!! – Jasper Zhou Sep 14 '18 at 04:39