Fix digamma integral

.mailmap

@@ -1225,6 +1225,7 @@ Robin Richard <[email protected]> Nornort <[email protected]>
 Robin Richard <[email protected]> robinechuca <[email protected]>
 Robin Richard <[email protected]> robinechuca <[email protected]>
 Rodrigo Luger <[email protected]>
+Roger Balsach <[email protected]> Roger Balsach <[email protected]>
 Rohit Jain <[email protected]>
 Rohit Rango <[email protected]>
 Rohit Sharma <[email protected]>

sympy/functions/special/tests/test_zeta_functions.py

@@ -84,6 +84,12 @@ def test_zeta_series():
  zeta(x, z) - x*(a-z)*zeta(x+1, z) + O((a-z)**2, (a, z))
 
 
+def test_zeta_leadterm():
+ assert zeta(2+x, 2-x)._eval_as_leading_term(x) == zeta(2, 2)
+ assert zeta(1+x, 4)._eval_as_leading_term(x) == 1 / x
+ assert zeta(1-I*x**3, 4)._eval_as_leading_term(x) == I / x**3
+
+
 def test_dirichlet_eta_eval():
  assert dirichlet_eta(0) == S.Half
  assert dirichlet_eta(-1) == Rational(1, 4)
@@ -229,6 +235,47 @@ def test_lerchphi_expansion():
  assert myexpand(lerchphi(exp(I*pi*Rational(2, 5)), s, a), None)
 
 
+def test_lerchphi_finite():
+ assert lerchphi(1, 1, a).is_finite is False
+ assert lerchphi(0, -3, 0).is_finite is True
+ assert lerchphi(2, s, 0).is_finite is None
+ assert lerchphi(2, 1, 0).is_finite is False
+ assert lerchphi(2, -1, 0).is_finite is True
+ assert lerchphi(2, 1j, 0).is_finite is None
+ assert lerchphi(z, s, -1).is_finite is None
+ assert lerchphi(0, s, -1).is_finite is True
+ assert lerchphi(z, -3, 0).is_finite is True
+ assert lerchphi(0, s, 0).is_finite is None
+ assert lerchphi(0, 1, 0).is_finite is False
+ assert lerchphi(0, -1, 0).is_finite is True
+ assert lerchphi(0.5, -1, 0).is_finite is True
+ assert lerchphi(0.5+0.5j, 0, 0).is_finite is True
+ assert lerchphi(1, 3, 0).is_finite is False
+ assert lerchphi(0, 1, a).is_finite is None
+ assert lerchphi(1, (1+I)/2, 0).is_finite is True
+
+
+def test_lerchphi_leadterm():
+ from sympy.functions.special.gamma_functions import digamma
+ h = S.Half
+ assert (lerchphi(1-x, 1, a)._eval_as_leading_term(x)
+ == -log(x) - S.EulerGamma - digamma(a))
+ assert lerchphi(1+x, 3, 4)._eval_as_leading_term(x) == zeta(3, 4)
+ assert lerchphi(1, 1+x, 3)._eval_as_leading_term(x) == 1 / x
+ assert lerchphi(1, 1+x**2, 3)._eval_as_leading_term(x) == 1 / x**2
+ assert lerchphi(S.One/10, h, x)._eval_as_leading_term(x) == 1 / sqrt(x)
+ assert (lerchphi(S.One/10, h, x-2)._eval_as_leading_term(x)
+ == 1 / (100*sqrt(x)))
+ assert lerchphi(z+x, h/2, x-1)._eval_as_leading_term(x) == z / x**(h/2)
+ assert lerchphi(1-x, 2+x**2, a)._eval_as_leading_term(x) == zeta(2, a)
+
+
+def test_lerchphi_nseries():
+ from sympy.core.function import PoleError
+ raises(PoleError,
+ lambda: lerchphi(1-x, 1, z)._eval_nseries(x, n=1, logx=None))
+
+
 def test_stieltjes():
  assert isinstance(stieltjes(x), stieltjes)
  assert isinstance(stieltjes(x, a), stieltjes)

sympy/functions/special/zeta_functions.py

@@ -214,6 +214,60 @@ def _eval_rewrite_as_zeta(self, z, s, a, **kwargs):
  def _eval_rewrite_as_polylog(self, z, s, a, **kwargs):
  return self._eval_rewrite_helper(polylog)
 
+ def _eval_is_finite(self):
+ z, s, a = self.args
+ if z.is_infinite or s.is_infinite or a.is_infinite:
+ # TODO: Write a valid implementation for infinite numbers.
+ return
+ if z == 0:
+ if a.is_zero and s.is_real:
+ return s.is_nonpositive
+ if a.is_nonzero:
+ return True
+ return
+ if a.is_integer and a.is_nonpositive:
+ if s.is_real:
+ return s.is_nonpositive
+ # This part assumes analytic continuation
+ if z == 1:
+ arg_is_one = (s - 1).is_zero
+ if arg_is_one is not None:
+ return not arg_is_one
+
+ def _eval_as_leading_term(self, x, logx=None, cdir=0):
+ '''
+ Use the relation lerchphi(1, s, a) = zeta(s, a) for z = 1 and the
+ expansions given in [1]
+
+ References
+ ==========
+
+ .. [1] https://en.wikipedia.org/wiki/Lerch_zeta_function#Series_representations
+
+ '''
+ from sympy.functions.special.gamma_functions import digamma, gamma
+ z, s, a = self.args
+ z0, s0, a0 = (arg.subs(x, 0).cancel() for arg in self.args)
+ lt = 0
+ if z == 1:
+ lt = zeta(s, a)
+ elif z0 == 1:
+ if s == 1:
+ if not a.has(x):
+ # This could be generalised for arbitrary a if we check
+ # that this term doesn't vanish.
+ lt = -log(-log(z)) - S.EulerGamma - digamma(a)
+ elif s.is_integer and s.is_positive or re(s0 - 1).is_positive:
+ return zeta(s0, a0)
+ elif re(s0 - 1).is_negative:
+ lt = gamma(1-s)*(-log(z))**(s-1)
+ elif a0.is_integer and a0.is_nonpositive:
+ if s.is_real and s.is_positive:
+ lt = z**(-a0)/(a-a0)**s
+ if lt == 0:
+ return super()._eval_as_leading_term(x, logx=logx, cdir=cdir)
+ return lt._eval_as_leading_term(x, logx=logx, cdir=cdir)
+
 ###############################################################################
 ###################### POLYLOGARITHM ##########################################
 ###############################################################################
@@ -579,7 +633,14 @@ def _eval_as_leading_term(self, x, logx=None, cdir=0):
  if e.is_negative and not s.is_positive:
  raise NotImplementedError
 
- return super(zeta, self)._eval_as_leading_term(x, logx, cdir)
+ s0, a0 = (arg.subs(x, 0).cancel() for arg in self.args)
+ if (s0 - 1).is_nonzero:
+ if (lt := zeta(s0, a0)).is_nonzero:
+ return lt
+ elif s.has(x) and (s0 - 1).is_zero:
+ return (1 / (s - 1))._eval_as_leading_term(x, logx=logx, cdir=cdir)
+
+ return super(zeta, self)._eval_as_leading_term(x, logx=logx, cdir=cdir)
 
 
 class dirichlet_eta(Function):

sympy/integrals/tests/test_integrals.py

@@ -2139,3 +2139,7 @@ def test_old_issues():
  # https://github.com/sympy/sympy/issues/6278
  I3 = integrate(1/(cos(x)+2),(x,0,2*pi))
  assert I3 == 2*sqrt(3)*pi/3
+
+def test_issue_26523():
+ Int = Integral((1 - t**z) / (1 - t), (t, 0, 1))
+ assert Int.doit() == polygamma(0, z + 1) + EulerGamma

sympy/printing/mathematica.py

@@ -105,7 +105,7 @@
  "airyaiprime": [(lambda x: True, "AiryAiPrime")],
  "airybiprime": [(lambda x: True, "AiryBiPrime")],
  "polylog": [(lambda *x: True, "PolyLog")],
- "lerchphi": [(lambda *x: True, "LerchPhi")],
+ "lerchphi": [(lambda *x: True, "HurwitzLerchPhi")],
  "gcd": [(lambda *x: True, "GCD")],
  "lcm": [(lambda *x: True, "LCM")],
  "jn": [(lambda *x: True, "SphericalBesselJ")],

sympy/printing/tests/test_mathematica.py

@@ -10,7 +10,7 @@
  FallingFactorial, harmonic, atan2, sec, acsc,
  hermite, laguerre, assoc_laguerre, jacobi,
  gegenbauer, chebyshevt, chebyshevu, legendre,
- assoc_legendre, Li, LambertW)
+ assoc_legendre, Li, LambertW, lerchphi)
 
 from sympy.printing.mathematica import mathematica_code as mcode
 
@@ -81,6 +81,7 @@ def test_Function():
  assert mcode(LambertW(x)) == "ProductLog[x]"
  assert mcode(LambertW(x, -1)) == "ProductLog[-1, x]"
  assert mcode(LambertW(x, y)) == "ProductLog[y, x]"
+ assert mcode(lerchphi(x, y, z)) == "HurwitzLerchPhi[x, y, z]"
 
 
 def test_special_polynomials():

sympy/series/limits.py

@@ -93,6 +93,9 @@ def heuristics(e, z, z0, dir):
  return heuristics(m, z, z0, dir)
  return
  return
+ elif (not (l.has(S.Infinity) or l.has(S.NegativeInfinity))
+ and l.is_finite is False):
+ r.append(S.ComplexInfinity)
  elif isinstance(l, Limit):
  return
  elif l is S.NaN:
@@ -347,6 +350,8 @@ def set_signs(expr):
  return coeff
  if coeff.has(S.Infinity, S.NegativeInfinity, S.ComplexInfinity, S.NaN):
  return self
+ from sympy.simplify.gammasimp import gammasimp
+ coeff = gammasimp(coeff)
  if not coeff.has(z):
  if ex.is_positive:
  return S.Zero

sympy/series/tests/test_limits.py

@@ -19,7 +19,8 @@
  atan, cos, cot, csc, sec, sin, tan)
 from sympy.functions.special.bessel import (besseli, bessely, besselj, besselk)
 from sympy.functions.special.error_functions import (Ei, erf, erfc, erfi, fresnelc, fresnels)
-from sympy.functions.special.gamma_functions import (digamma, gamma, uppergamma)
+from sympy.functions.special.gamma_functions import (digamma, gamma, polygamma, uppergamma)
+from sympy.functions.special.zeta_functions import lerchphi
 from sympy.functions.special.hyper import meijerg
 from sympy.integrals.integrals import (Integral, integrate)
 from sympy.series.limits import (Limit, limit)
@@ -311,8 +312,11 @@ def test_set_signs():
 
 def test_heuristic():
  x = Symbol("x", real=True)
+ a = Symbol("a")
  assert heuristics(sin(1/x) + atan(x), x, 0, '+') == AccumBounds(-1, 1)
  assert limit(log(2 + sqrt(atan(x))*sqrt(sin(1/x))), x, 0) == log(2)
+ # https://github.com/sympy/sympy/issues/26523
+ assert heuristics(lerchphi(x, 1, a) + log(x-1), x, 1, '-') is None
 
 
 def test_issue_3871():
@@ -1412,3 +1416,14 @@ def test_issue_26250():
  e1 = ((1-3*x**2)*e**2/2 - (x**2-2*x+1)*e*k/2)
  e2 = pi**2*(x**8 - 2*x**7 - x**6 + 4*x**5 - x**4 - 2*x**3 + x**2)
  assert limit(e1/e2, x, 0) == -S(1)/8
+
+
+def test_issue_26523():
+ expr = -x**(z + 1)*lerchphi(x, 1, z + 1) - log(x - 1)
+ L = Limit(expr, x, 1, '-')
+ assert L.doit() == polygamma(0, z + 1) + EulerGamma - I*pi
+ expr = (-x*x**z*z*lerchphi(x, 1, z + 1)*gamma(z + 1)/gamma(z + 2)
+ - x*x**z*lerchphi(x, 1, z + 1)*gamma(z + 1)/gamma(z + 2)
+ - log(x - 1))
+ L = Limit(expr, x, 1, '-')
+ assert L.doit() == polygamma(0, z + 1) + EulerGamma - I*pi