Solving a quadratic symbolically (Sympy vs. Mathematica)

[Justin Carden]

Hi all,

I'm trying to make the switch from Mathematica to Python in the lab, but I'm running into a small problem binding the result space to a boundary. 

Specifically, I'm trying to assign equality to each polynomial equation in a system to a total t1,t2 and t3 to solve the system symbolically. Since double equals (==)

tests equality exactly and not symbolically in Sympy, setting the equations to t1,t2 or t3 does not work. This works in Mathematica with no problem. 

How would I go about doing this with Sympy ? 

(Complete code for the system in both mathematica and Sympy below)

Mathematica:

f1 = x + s1 + s2 + s5 == t1
f2 = y + s2 + s3 + s6 == t2
f3 = z + s4 + s5 + s6 == t3

Current Sympy:

f_1 = x + s1 + s2 + s5

f_2 = y + s2 + s3 + s6

f_3 = z + s4 + s5 + s6

=====================

Sympy

=====================

k1 = kdAA = 0.11

k2 = kdBB = 0.0

k3 = kdCC = 0.0

k4 = kdAB = 0.25

k5 = kdAC = 0.5

k6 = kdBC = 0.33

kineticRates = [k1,k2,k3,k4,k5,k6]

t1 = totA = 1000

t2 = totB = 500

t3 = totC = 100

totals = [t1,t2,t3]

freeA = x

freeB = y

freeC = z

free = [freeA,freeB,freeC]

s1 = kdAA*x**2

s2 = kdAB*x*y

s3 = kdBB*y**2

s4 = kdCC*z**2

s5 = kdAC*x*z

s6 = kdBC*y*z

states = [s1,s2,s3,s4,s5,s6]

f_1 = x + s1 + s2 + s5

f_2 = y + s2 + s3 + s6

f_3 = z + s4 + s5 + s6

system = [f_1,f_2,f_3]

results = solve(system, x,y,z, dict=True, set=True)

--------------------------------------

[{x: -9.09090909090909, z: 0, y: 0}, 

{x: -4.00000000000000, z: 0, y: -2.24000000000000}, 

{x: -2.00000000000000, z: -1.56000000000000, y: 0}, 

{x: -1.96537201684605, z: -1.54138483572269, y: -0.0524666411423548}, 

{x: 0, z: -3.03030303030303, y: -3.03030303030303}, 

{x: 0, z: 0, y: 0}]

=====================

Mathematica

=====================

eq1 = kdAA = 0.11
eq2 = kdBB = 0.0
eq3 = kdCC = 0.0
eq4 = kdAB = 0.25
eq5 = kdAC = 0.5
eq6 = kdBC = 0.33

t1 = totA = 1000
t2 = totB = 500
t3 = totC = 100

freeA = x
freeB = y
freeC = z

s1 = kdAA*x^2
s2 = kdAB*x*y
s3 = kdBB*y^2
s4 = kdCC*z^2
s5 = kdAC*x*z
s6 = kdBC*y*z

f1 = x + s1 + s2 + s5 == t1
f2 = y + s2 + s3 + s6 == t2
f3 = z + s4 + s5 + s6 == t3

NSolve[{f1, f2, f3}, {x, y, z}]

--------------------------------------

Out[22] = 

{{x -> -1.4282, y -> -1.08259, z -> -1401.51}, 

{x -> -6.67923, y -> -598.932,  z -> -0.500032}, 

{x -> -2.17279 - 53.7467 I, y -> -1.88296 + 63.8985 I, z -> -2.08239 + 17.0259 I},

 {x -> -2.17279 + 53.7467 I,  y -> -1.88296 - 63.8985 I, z -> -2.08239 - 17.0259 I}, 

{x -> -66.1313, y -> -30.6463, z -> -2.37085},

 {x -> 61.528, y -> 29.1041, z -> 2.41731}}

[Stefan Krastanov]

f_1 = x + s1 + s2 + s5 - t1

this is all that you need to do if I understand your question correctly

Or if you wish, you can create 'Eq(right_hand, left_hand)' instances.

> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to sympy+un...@googlegroups.com.
> To post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> For more options, visit https://groups.google.com/groups/opt_out.
>
>

[Justin Carden]

Hi Stefan,

Both of these solutions return an empty array. Where I should get an array of possible solutions

similar to the output of the Mathematica version. Such as:

{{x -> -1.4282, y -> -1.08259, z -> -1401.51}, 

{x -> -6.67923, y -> -598.932,  z -> -0.500032}, 

{x -> -2.17279 - 53.7467 I, y -> -1.88296 + 63.8985 I, z -> -2.08239 + 17.0259 I},

 {x -> -2.17279 + 53.7467 I,  y -> -1.88296 - 63.8985 I, z -> -2.08239 - 17.0259 I}, 

{x -> -66.1313, y -> -30.6463, z -> -2.37085},

 {x -> 61.528, y -> 29.1041, z -> 2.41731}}

[Stefan Krastanov]

The equations are nonlinear and it seems that sympy is unable to solve
them. Maybe there is a way to instruct sympy how to deal with them but
I do not know how.

However, there is a very big difference between your mathematica code
and your sympy code. In mathematica you are calling `NSolve`, which is
a numeric solver, not a symbolic one. In sympy you are calling
`solve`, which is searching for a general symbolic solution. You
should either use `sympy.nsolve` or something from `scipy`.
`sympy.nsolve` returns only one solution and needs a starting point.

[Justin Carden]

So, I switched to the 'sympy.nsolve' method and got the approximate solutions similar to Mathematica.

Yes, the system is nonlinear. While the method in the Mathematica code was 'NSolve', the solution set

returned was the same for the 'Solve' method.

Regardless, this helped a bunch. Thanks!

[Stefan Krastanov]

For the moment mathematica is more advanced than sympy both in finding
symbolic solutions and in finding _multiple_ numeric solutions. This
explains most of the differences in results.

[Aaron Meurer]

It's also possible that Mathematica automatically returns numeric
solutions when explicit ones can't be found.

Aaron Meurer