15

I solved a quadratic equation using sympy:

import sympy as sp
q,qm,k,c0,c,vt,vm = sp.symbols('q qm k c0 c vt vm')
c = ( c0 * vt - q * vm) / vt
eq1 = sp.Eq(qm * k * c / (1 + k * c) ,q)
q_solve = sp.solve(eq1,q)

Based on some testing I figured out that only q_solve[0] makes physical sense. Will sympy always put (b - sqrt(b**2 - 4*a*c))/2a in the first place ? I guess, it might change with an upgrade ?

Moritz
  • 5,130
  • 10
  • 40
  • 81
  • Which of the above defined parameters (q, qm, k,...) are supposed to be positive, how do you define 'physical sense' and which kind of testing did you do, which kind of values did you use,...? – Cleb Aug 17 '15 at 10:56
  • I do know that q_solve should be positive and q_solve should be smaller qm. – Moritz Aug 17 '15 at 11:04
  • And there are no fixed values (c0 is probably the speed of light)? Or other variables that are greater than 0? Which values do you take for testing? – Cleb Aug 17 '15 at 11:12
  • c is a concentration, q is a capacity (concentration on solid), qm describes the maximum adsorption capacity, k the ratio of speed of dissoziation to adsorption, vt and vm are volumes. All values are in the range of 1e-5 and 1 – Moritz Aug 17 '15 at 12:57
  • I know it is possible to include a feasibility calculation at each step but I was wondering how sympy orders the results. Furhtermore, it would be faster without these checks. – Moritz Aug 17 '15 at 12:59
  • 2
    There is no order, at least to my knowledge. But you can use e.g. q=sp.symbols('q', positive=True) to make sure that q is larger than 0. Unfortunately, this might only help if you have some more assumptions about your parameters. – Cleb Aug 17 '15 at 13:24

1 Answers1

5

A simple test to answer your question is to symbolically solve the quadratic equation using sympy per below:

import sympy as sp
a, b, c, x = sp.symbols('a b c x')
solve( a*x**2 + b*x + c, x)

this gives you the result:

[(-b + sqrt(-4*a*c + b**2))/(2*a), -(b + sqrt(-4*a*c + b**2))/(2*a)]

which leads me to believe that in general the order is first the + sqrt() solution and then the - sqrt() solution.

For your program q_solve[0] gives you:

(c0*k*vt + k*qm*vm + vt - sqrt(c0**2*k**2*vt**2 - 2*c0*k**2*qm*vm*vt + 2*c0*k*vt**2 + k**2*qm**2*vm**2 + 2*k*qm*vm*vt + vt**2))/(2*k*vm)

this is still the x= (-b + sqrt(b**2-4*a*c))/(2*a) answer, the negative sign from the b term goes away as a result of the distribution of the signs of the variables within the solution

bern
  • 366
  • 1
  • 10