Something I do not understand with fractional decomposition (apart in symPy )

52 views
Skip to first unread message

Mikhael Myara

unread,
Aug 10, 2020, 4:28:56 AM8/10/20
to sympy
Dear all,

 I don't know if my problem is with the knowledge of fractionnal decomposition itself or with symPy's implementation.

 I start with the following code :

import sympy as sp
sp
.var('s')


Hstab = 1/(s**3+s**2+5*s+4)
Hstab = sp.apart(Hstab)
display
(Hstab)


The result is :

Capture d’écran 2020-08-10 à 10.18.18.png

Now I check another similar expression :

Hstab2 = 1/(s**3 + 2*s**2 + 5*s + 4)
Hstab2 = sp.apart(Hstab2,s)
display
(Hstab2)

and I get :

Capture d’écran 2020-08-10 à 10.21.14.png




If now I check for the roots of the denominator for each case :


print("\n\nFirst case :")
P1
=(1/Hstab).simplify()
display
(P1)
for sol in sp.solve(P1,s):
    display
(sol.evalf())


print("\n\nSecond case :")    
P2
=(1/Hstab2).simplify()
display
(P2)
for sol in sp.solve(P2,s):
    display
(sol.evalf())  

I get :

Capture d’écran 2020-08-10 à 10.26.18.png

Which is very similar except that the real root of the first case is not ass simple as the one of the second case. So I do not understand as I can't see any mathematical impossibility, the two cases are very close.
Is there a way to reach the fractional decomposition for the first case ?

Thanks a lot,
   Mike


Aaron Meurer

unread,
Aug 10, 2020, 2:48:04 PM8/10/20
to sympy
The roots are the key difference. The second polynomial has a rational root, -1, meaning it can be split out in a partial fraction decomposition. If you don't call evalf() on the roots you can also see that the roots of the first polynomial are more complicated, because they are coming from the general cubic formula. You can also see the difference in the two if you call factor(). 

>>> factor(s**3+s**2+5*s+4)
s**3 + s**2 + 5*s + 4
>>> factor(s**3+2*s**2+5*s+4)
(s + 1)*(s**2 + s + 4)

The first polynomial is irreducible over rationals, but the second factors. So if you do a partial fraction decomposition, it will not split because apart() only decomposes over rational numbers by default. If you want a full decomposition over all the roots, you can use something like apart(1/(s**3+s**2+5*s+4), full=True).doit().

Note that it's not uncommon for polynomials to work like this, where if you change a coefficient it changes the behavior of it. That's because it's easy for a polynomial to have a rational root with one coefficient but not with another close coefficient, like in this case.

Aaron Meurer

--
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 view this discussion on the web visit https://groups.google.com/d/msgid/sympy/4ccc519f-8a28-4478-bc2e-adbc302a9647o%40googlegroups.com.

Mikhael Myara

unread,
Aug 10, 2020, 3:53:29 PM8/10/20
to sympy
Thanks a lot Aaron. 
However : for my usage, an approximate root would be sufficient. Is there a way to allow this in Sympy ?
Thanks again,
   Mikhaël
To unsubscribe from this group and stop receiving emails from it, send an email to sy...@googlegroups.com.

Aaron Meurer

unread,
Aug 10, 2020, 3:58:23 PM8/10/20
to sympy
If you replace the coefficients with floating point numbers, you can use full=True and get an answer with floats (although they aren't fully simplified for some reason, so you may want to use apart_list() instead). I'm not sure if it is possible to get an answer with just the real roots, so that there are no complex numbers added.

Aaron Meurer

To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/b4c883c0-0840-4cd2-b964-cb32a27c8aaao%40googlegroups.com.

Mikhael Myara

unread,
Aug 10, 2020, 5:03:23 PM8/10/20
to sympy
Thanks again Aaron,
  here is what I did :

import sympy as sp

sp
.var('s t')


Hstab = 1/(s**3+s**2+5*s+4)
Hstab = sp.apart(Hstab,full=True).doit().evalf()
display
(Hstab)
And I get :

Capture d’écran 2020-08-10 à 23.00.45.png



Which is great. However I would have liked the conjugate denominator fractions combined. Is this possible ou should I write some code for this ?
Moreover : I am not sure I understood how is appart_list() interesting ?

Thanks again,
  Mikhaël

Aaron Meurer

unread,
Aug 10, 2020, 5:05:42 PM8/10/20
to sympy
Yes, that is what I mean by only getting the real roots. I don't think it's possible now, though it shouldn't be hard to do it manually with apart_list().

Aaron Meurer

To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/9a1d5c28-b8f7-4d66-b197-02d615c7af18o%40googlegroups.com.

Mikhael Myara

unread,
Aug 11, 2020, 11:23:38 AM8/11/20
to sympy
apart_list() does not seem to support the « full=True » mode.
So I guess I have to loop over the summed terms. Do you agree ?
thanks again,
  Mike

Aaron Meurer

unread,
Aug 11, 2020, 3:26:37 PM8/11/20
to sympy
Ah that's interesting. I would consider that to be a bug. I would expect apart() to just be assemble_partfrac_list(apart_list()), but apparently there is a lot of duplicated logic.  Can you open a bug report about this https://github.com/sympy/sympy/issues

Aaron Meurer

To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/64ac532d-bf89-41fa-b4e6-9a2d77d25940o%40googlegroups.com.

Mikhael Myara

unread,
Aug 12, 2020, 4:53:26 AM8/12/20
to sympy
Dear Aaron,

  I just did it : bug report #19955

Thanks for your help,
   Mike
Reply all
Reply to author
Forward
0 new messages