Sympy hangs on fairly simple expression

[sympy-issue-migrator]

Hi,

I'm trying to use sympy to create a symbolic Taylor series of a fairly simple expression:

from sympy import *
a, b, x, alpha = symbols('a b x alpha')
f = (a*x**2 +  b*x) ** (2 - alpha)
series(f, x)

returns 'nan'

[Wolfram alpha manages fine: ](http://www.wolframalpha.com/input/?i=taylor+series+%28a+x^2+%2B+bx%29^%282+-+2alpha%29 Help?)

In [49]: sympy.__version__
Out[49]: '0.7.3'

Original issue for #7129: http://code.google.com/p/sympy/issues/detail?id=4030
Original author: https://code.google.com/u/101140808581050767033/

[sympy-issue-migrator]

Sorry, I'm missing a factor of 2 in the exponent of the function. It still doesn't work though.

a, b, x, alpha = symbols('a b x alpha')
f = (a*x**2 +  b*x) ** (2 - 2*alpha)
series(f, x)

[skirpichev]

HEAD now returns O(log(x)**6)

**Status:** Valid  
**Labels:** Series  

Original comment: http://code.google.com/p/sympy/issues/detail?id=4030#c2
Original author: https://code.google.com/u/113518207426301808764/

[borntorock]

I am getting this NotImplementedError: unknown boundedness for [zoo*alpha]
We have to define it. Does sympy has the implementation for solving taylor series.If yes then could anyone specify the file and if not where it should go up. I would like to work on this, from where should i start off.

[avichaldayal]

I tried a simpler example such as series(x**y, x) and that also gives O(log(x)**6)

In exponential.py a condition checks if arg0 belongs in [-oo, oo]. However value of arg0 is oo*sign(y). Now since it does not know if sign(y) is +1 or -1, it evaluates the conditions as false.
However it is clear that oo*sign(y) is either -oo or oo.

Should a separate function be made to check if a given expression is in indeterminate form?

[asmeurer]

This is an interesting use-case for my work at #2508 . One could make ask(Q.positive(sign(y)) | Q.negative(sign(y)) work in that branch (one would additionally need to assume ~Q.zero(y), which I'm guessing the series code implicitly does anyway), and with a bit more work, make ask(Eq(sign(y)*oo, oo) | Eq(sign(y)*oo, -oo)) work as well. It's a lot of work, though, but it's what I'd like to see.

[anutosh491]

Just came across this issue now, this issue was fixed by a pr of mine earlier this year #22870 and the methodology falls around this line of thought as expalined by avichaldayal

I tried a simpler example such as series(x**y, x) and that also gives O(log(x)**6)
In exponential.py a condition checks if arg0 belongs in [-oo, oo]. However value of arg0 is oosign(y). Now since it does not know if sign(y) is +1 or -1, it evaluates the conditions as false.
However it is clear that oosign(y) is either -oo or oo.

I feel the implementation in the pr would suffice for solving this and related issues but not sure if the following approach has been explored yet.

One could make ask(Q.positive(sign(y)) | Q.negative(sign(y)) work in that branch (one would additionally need to assume ~Q.zero(y), which I'm guessing the series code implicitly does anyway), and with a bit more work, make ask(Eq(sign(y)*oo, oo) | Eq(sign(y)*oo, -oo)) work as well. It's a lot of work, though, but it's what I'd like to see.

`
Marking as could close as the issue is solved as of now .

[oscarbenjamin]

So now we have:

In [3]: f
Out[3]: 
            2 - α
⎛   2      ⎞     
⎝a⋅x  + b⋅x⎠     

In [4]: f.series(x)
Out[4]: 
            2 - α
⎛   2      ⎞     
⎝a⋅x  + b⋅x⎠ 

Is that the best that can be done here?

For particular values of alpha we have:

In [5]: f.subs(alpha, 1).series(x)
Out[5]: 
   2      
a⋅x  + b⋅x

In [6]: f.subs(alpha, 2).series(x)
Out[6]: 1

In [7]: f.subs(alpha, 3).series(x)
Out[7]: 
            2      3  2    4  3    5  4    6  5        
 1    a    a ⋅x   a ⋅x    a ⋅x    a ⋅x    a ⋅x     ⎛ 6⎞
─── - ── + ──── - ───── + ───── - ───── + ───── + O⎝x ⎠
b⋅x    2     3       4       5       6       7         
      b     b       b       b       b       b    

In [8]: f.subs(alpha, Rational(1, 2)).series(x)
Out[8]: 
                    5/2      2  7/2    3  9/2      4  11/2        
 3/2  3/2   3⋅a⋅√b⋅x      3⋅a ⋅x      a ⋅x      3⋅a ⋅x        ⎛ 6⎞
b   ⋅x    + ─────────── + ───────── - ─────── + ────────── + O⎝x ⎠
                 2           8⋅√b         3/2         5/2         
                                      16⋅b       128⋅b    

[anutosh491]

Is that the best that can be done here?

I have been coming across more of these issues recently . Though I was thinking it's best to generalize such expressions to self because at the end we can't calculate a limit for such a function at 0 , so it might be logical to just return self for a series at 0 or maybe a NotImplementedError just like for limits

>>> limit(f, x, 0)
Traceback (most recent call last):
.................................
NotImplementedError: Result depends on the sign of sign(alpha - 1)

But now I think , maybe we can try out series for these expressions and it helps overall . It won't involve a lot of changes to return something like the following

>>> f = (a*x**2 +  b*x) ** (2 - alpha)
>>> f.series(x, 0, 3)
(b*x)**(2 - alpha) + a*x*(b*x)**(2 - alpha)*(2 - alpha)/b + a**2*x**2*(b*x)**(2 - alpha)*(1 - alpha)*(2 - alpha)/(2*b**2) + O(x**3)

Though my only concern then would be the subs method. For eg ... the current way we deal with our _eval_nseries methods (in this case pow._eval_nseries) in most cases in we find terms seperately and then multiply then with x**k (0<= k <= n). This is what we also do in the series above

x**0 *   (b*x)**(2 - alpha)    + 
x**1 *   a*(b*x)**(2 - 2*alpha)*(2 - alpha)/b     + 
x**2 *   a**2*(b*x)**(2 - alpha)*(1 - alpha)*(2 - alpha)/(2*b**2)     +
O(x**3)

Everything looks fine .... .But in cases where the coefficient would have x or powers or x, the subs function there would mess up the number of terms being returned.

>>> _.subs(alpha, 3)
1/(b*x) - a/b**2 + a**2*x/b**3 + O(x**3)
>>> f.subs(alpha, 3).series(x, 0, 3)
-a**3*x**2/b**4 + a**2*x/b**3 - a/b**2 + 1/(b*x) + O(x**3)

We miss out on the last term -a**3*x**2/b**4 here because our last term has powers of x**2 *x**(2 - alpha) where alpha = 3 , hence it boils down to 1 . So currently I am not sure how we can address this , if we want to use the subs function after getting the series ..... if we stick to our methodology of calculating the term separately and multiplying it with power of x.

if coefficient doesn't have powers of x it look good though

>>> n, k = symbols('n k', positive=True)
>>> expr = (1 + 1/n)**k
>>> expr._eval_nseries(n, 3, None)
n**(-k) + k*n**(2 - k)*(k - 1)/2 + k*n**(1 - k) + O(n**3)