Merge pull request #19432 from sachin-4099/gsoc#4

[GSoC] Series: Using is_meromorphic() for limit evaluations

sympy/core/add.py

@@ -901,7 +901,7 @@ def _eval_as_leading_term(self, x):
                 if not min or order not in min:
                     min = order
                     new_expr = term
-                elif order == min:
+                elif min in order:
                     new_expr += term
 
         except TypeError:

sympy/core/expr.py

@@ -2717,8 +2717,8 @@ def is_meromorphic(self, x, a):
         False
 
         """
-        if not x.is_Symbol:
-            raise TypeError("{} should be of type Symbol".format(x))
+        if not x.is_symbol:
+            raise TypeError("{} should be of symbol type".format(x))
         a = sympify(a)
 
         return self._eval_is_meromorphic(x, a)
@@ -3869,7 +3869,8 @@ def _eval_is_rational_function(self, syms):
         return True
 
     def _eval_is_meromorphic(self, x, a):
-        return True
+        from sympy.calculus.util import AccumBounds
+        return (not self.is_Number or self.is_finite) and not isinstance(self, AccumBounds)
 
     def _eval_is_algebraic_expr(self, syms):
         return True

sympy/core/power.py

@@ -1340,7 +1340,7 @@ def _eval_is_meromorphic(self, x, a):
         # E**(log(f)*g) is meromorphic if log(f)*g is meromorphic
         # and finite.
         base_merom = self.base._eval_is_meromorphic(x, a)
-        exp_integer = self.exp.is_integer
+        exp_integer = self.exp.is_Integer
         if exp_integer:
             return base_merom
 

sympy/core/tests/test_expr.py

@@ -626,15 +626,15 @@ def test_is_meromorphic():
     assert f.is_meromorphic(x, zoo) is True
 
     g = 3 + 2*x**(log(3)/log(2) - 1)
-    assert g.is_meromorphic(x, 0) is None
+    assert g.is_meromorphic(x, 0) is False
     assert g.is_meromorphic(x, 1) is True
-    assert g.is_meromorphic(x, zoo) is None
+    assert g.is_meromorphic(x, zoo) is False
 
     n = Symbol('n', integer=True)
     h = sin(1/x)**n*x
     assert h.is_meromorphic(x, 0) is False
     assert h.is_meromorphic(x, 1) is True
-    assert h.is_meromorphic(x, zoo) is True
+    assert h.is_meromorphic(x, zoo) is False
 
     e = log(x)**pi
     assert e.is_meromorphic(x, 0) is False
@@ -643,7 +643,7 @@ def test_is_meromorphic():
     assert e.is_meromorphic(x, zoo) is False
 
     assert (log(x)**a).is_meromorphic(x, 0) is False
-    assert (log(x)**a).is_meromorphic(x, 1) is None
+    assert (log(x)**a).is_meromorphic(x, 1) is False
     assert (a**log(x)).is_meromorphic(x, 0) is None
     assert (3**log(x)).is_meromorphic(x, 0) is False
     assert (3**log(x)).is_meromorphic(x, 1) is True

sympy/functions/elementary/trigonometric.py

@@ -468,7 +468,10 @@ def _eval_as_leading_term(self, x):
         if x in arg.free_symbols and Order(1, x).contains(arg):
             return arg
         else:
-            return self.func(arg)
+            if not arg.subs(x, 0).is_finite:
+                return self
+            else:
+                return self.func(arg)
 
     def _eval_is_extended_real(self):
         if self.args[0].is_extended_real:
@@ -908,7 +911,10 @@ def _eval_as_leading_term(self, x):
         if x in arg.free_symbols and Order(1, x).contains(arg):
             return S.One
         else:
-            return self.func(arg)
+            if not arg.subs(x, 0).is_finite:
+                return self
+            else:
+                return self.func(arg)
 
     def _eval_is_extended_real(self):
         if self.args[0].is_extended_real:

sympy/functions/special/gamma_functions.py

@@ -202,13 +202,15 @@ def _sage_(self):
         return sage.gamma(self.args[0]._sage_())
 
     def _eval_as_leading_term(self, x):
-        from sympy import Order
         arg = self.args[0]
-        arg_1 = arg.as_leading_term(x)
-        if Order(x, x).contains(arg_1):
-            return S(1) / arg_1
-        if Order(1, x).contains(arg_1):
-            return self.func(arg_1)
+        x0 = arg.subs(x, 0)
+
+        if x0.is_integer and x0.is_nonpositive:
+            n = -x0
+            res = (-1)**n/self.func(n + 1)
+            return res/(arg + n).as_leading_term(x)
+        elif x0.is_finite:
+            return self.func(x0)
         ####################################################
         # The correct result here should be 'None'.        #
         # Indeed arg in not bounded as x tends to 0.       #

sympy/series/limits.py

@@ -2,10 +2,8 @@
 
 from sympy.core import S, Symbol, Add, sympify, Expr, PoleError, Mul
 from sympy.core.exprtools import factor_terms
-from sympy.core.numbers import GoldenRatio
 from sympy.core.symbol import Dummy
 from sympy.functions.combinatorial.factorials import factorial
-from sympy.functions.combinatorial.numbers import fibonacci
 from sympy.functions.special.gamma_functions import gamma
 from sympy.polys import PolynomialError, factor
 from sympy.series.order import Order
@@ -47,10 +45,7 @@ def limit(e, z, z0, dir="+"):
     >>> limit(1/x, x, 0, dir="-")
     -oo
     >>> limit(1/x, x, 0, dir='+-')
-    Traceback (most recent call last):
-        ...
-    ValueError: The limit does not exist since left hand limit = -oo and right hand limit = oo
-
+    zoo
     >>> limit(1/x, x, oo)
     0
 
@@ -193,6 +188,7 @@ def doit(self, **hints):
         hints : optional keyword arguments
             To be passed to ``doit`` methods; only used if deep is True.
         """
+        from sympy import sign
         from sympy.functions import RisingFactorial
 
         e, z, z0, dir = self.args
@@ -212,36 +208,51 @@ def doit(self, **hints):
         if not e.has(z):
             return e
 
+        if e.is_meromorphic(z, z0):
+            if abs(z0) is S.Infinity:
+                newe = e.subs(z, -1/z)
+            else:
+                newe = e.subs(z, z + z0)
+            try:
+                coeff, exp = newe.leadterm(z)
+                if coeff is not S.ComplexInfinity:
+                    if exp > 0:
+                        return S.Zero
+                    elif exp == 0:
+                        return coeff
+                    if str(dir) == "+" or not(int(exp) & 1):
+                        return S.Infinity*sign(coeff)
+                    elif str(dir) == "-":
+                        return S.NegativeInfinity*sign(coeff)
+                    else:
+                        return S.ComplexInfinity
+            except ValueError:
+                pass
+
         # gruntz fails on factorials but works with the gamma function
         # If no factorial term is present, e should remain unchanged.
         # factorial is defined to be zero for negative inputs (which
         # differs from gamma) so only rewrite for positive z0.
         if z0.is_extended_positive:
             e = e.rewrite([factorial, RisingFactorial], gamma)
 
-        if e.is_Mul:
-            if abs(z0) is S.Infinity:
-                e = factor_terms(e)
-                e = e.rewrite(fibonacci, GoldenRatio)
-                ok = lambda w: (z in w.free_symbols and
-                                any(a.is_polynomial(z) or
-                                    any(z in m.free_symbols and m.is_polynomial(z)
-                                        for m in Mul.make_args(a))
-                                    for a in Add.make_args(w)))
-                if all(ok(w) for w in e.as_numer_denom()):
-                    u = Dummy(positive=True)
-                    if z0 is S.NegativeInfinity:
-                        inve = e.subs(z, -1/u)
+        if e.is_Mul and abs(z0) is S.Infinity:
+            e = factor_terms(e)
+            u = Dummy('u', positive=True)
+            if z0 is S.NegativeInfinity:
+                inve = e.subs(z, -1/u)
+            else:
+                inve = e.subs(z, 1/u)
+            try:
+                f = inve.as_leading_term(u)
+                if f.is_meromorphic(u, S.Zero):
+                    r = limit(f, u, S.Zero, "+")
+                    if isinstance(r, Limit):
+                        return self
                     else:
-                        inve = e.subs(z, 1/u)
-                    try:
-                        r = limit(inve.as_leading_term(u), u, S.Zero, "+")
-                        if isinstance(r, Limit):
-                            return self
-                        else:
-                            return r
-                    except ValueError:
-                        pass
+                        return r
+            except (ValueError, NotImplementedError, PoleError):
+                pass
 
         if e.is_Order:
             return Order(limit(e.expr, z, z0), *e.args[1:])

sympy/series/tests/test_limits.py

@@ -60,7 +60,7 @@ def test_basic1():
     assert limit(1/x**2, x, 0, dir="+-") is oo
 
     # test failing bi-directional limits
-    raises(ValueError, lambda: limit(1/x, x, 0, dir="+-"))
+    assert limit(1/x, x, 0, dir="+-") is zoo
     # approaching 0
     # from dir="+"
     assert limit(1 + 1/x, x, 0) is oo
@@ -111,7 +111,7 @@ def eval(cls, arg):
 
 
 def test_issue_3885():
-    assert limit(x*y + x*z, z, 2) == x*y + 2*x
+    assert limit(x*y + x*z, z, 2) == x*(y + 2)
 
 
 def test_Limit():
@@ -424,7 +424,7 @@ def test_issue_6560():
                              35*x**4/8 - 15*x**2/4 + Rational(3, 8))/(2*(y + 1)))
     assert limit(e, y, oo) == (5*x**3 + 3*x**2 - 3*x - 1)/4
 
-
+@XFAIL
 def test_issue_5172():
     n = Symbol('n')
     r = Symbol('r', positive=True)
@@ -475,6 +475,10 @@ def test_issue_8730():
     assert limit(subfactorial(x), x, oo) is oo
 
 
+def test_issue_9558():
+    assert limit(sin(x)**15, x, 0, '-') == 0
+
+
 def test_issue_10801():
     # make sure limits work with binomial
     assert limit(16**k / (k * binomial(2*k, k)**2), k, oo) == pi
@@ -596,7 +600,7 @@ def test_issue_17325():
     assert Limit(sin(x)/x, x, 0, dir="+-").doit() == 1
     assert Limit(x**2, x, 0, dir="+-").doit() == 0
     assert Limit(1/x**2, x, 0, dir="+-").doit() is oo
-    raises(ValueError, lambda: Limit(1/x, x, 0, dir="+-").doit())
+    assert Limit(1/x, x, 0, dir="+-").doit() is zoo
 
 
 def test_issue_10978():
@@ -698,6 +702,11 @@ def test_issue_19026():
     assert limit(Abs(log(x) + 1)/log(x), x, oo) == 1
 
 
+def test_issue_19067():
+    x = Symbol('x')
+    assert limit(gamma(x)/(gamma(x - 1)*gamma(x + 2)), x, 0) == -1
+
+
 def test_issue_13715():
     n = Symbol('n')
     p = Symbol('p', zero=True)