zeta_functions.py 24 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786
  1. """ Riemann zeta and related function. """
  2. from sympy.core.add import Add
  3. from sympy.core.cache import cacheit
  4. from sympy.core.function import ArgumentIndexError, expand_mul, DefinedFunction
  5. from sympy.core.logic import fuzzy_not
  6. from sympy.core.numbers import pi, I, Integer
  7. from sympy.core.relational import Eq
  8. from sympy.core.singleton import S
  9. from sympy.core.symbol import Dummy
  10. from sympy.core.sympify import sympify
  11. from sympy.functions.combinatorial.numbers import bernoulli, factorial, genocchi, harmonic
  12. from sympy.functions.elementary.complexes import re, unpolarify, Abs, polar_lift
  13. from sympy.functions.elementary.exponential import log, exp_polar, exp
  14. from sympy.functions.elementary.integers import ceiling, floor
  15. from sympy.functions.elementary.miscellaneous import sqrt
  16. from sympy.functions.elementary.piecewise import Piecewise
  17. from sympy.polys.polytools import Poly
  18. ###############################################################################
  19. ###################### LERCH TRANSCENDENT #####################################
  20. ###############################################################################
  21. class lerchphi(DefinedFunction):
  22. r"""
  23. Lerch transcendent (Lerch phi function).
  24. Explanation
  25. ===========
  26. For $\operatorname{Re}(a) > 0$, $|z| < 1$ and $s \in \mathbb{C}$, the
  27. Lerch transcendent is defined as
  28. .. math :: \Phi(z, s, a) = \sum_{n=0}^\infty \frac{z^n}{(n + a)^s},
  29. where the standard branch of the argument is used for $n + a$,
  30. and by analytic continuation for other values of the parameters.
  31. A commonly used related function is the Lerch zeta function, defined by
  32. .. math:: L(q, s, a) = \Phi(e^{2\pi i q}, s, a).
  33. **Analytic Continuation and Branching Behavior**
  34. It can be shown that
  35. .. math:: \Phi(z, s, a) = z\Phi(z, s, a+1) + a^{-s}.
  36. This provides the analytic continuation to $\operatorname{Re}(a) \le 0$.
  37. Assume now $\operatorname{Re}(a) > 0$. The integral representation
  38. .. math:: \Phi_0(z, s, a) = \int_0^\infty \frac{t^{s-1} e^{-at}}{1 - ze^{-t}}
  39. \frac{\mathrm{d}t}{\Gamma(s)}
  40. provides an analytic continuation to $\mathbb{C} - [1, \infty)$.
  41. Finally, for $x \in (1, \infty)$ we find
  42. .. math:: \lim_{\epsilon \to 0^+} \Phi_0(x + i\epsilon, s, a)
  43. -\lim_{\epsilon \to 0^+} \Phi_0(x - i\epsilon, s, a)
  44. = \frac{2\pi i \log^{s-1}{x}}{x^a \Gamma(s)},
  45. using the standard branch for both $\log{x}$ and
  46. $\log{\log{x}}$ (a branch of $\log{\log{x}}$ is needed to
  47. evaluate $\log{x}^{s-1}$).
  48. This concludes the analytic continuation. The Lerch transcendent is thus
  49. branched at $z \in \{0, 1, \infty\}$ and
  50. $a \in \mathbb{Z}_{\le 0}$. For fixed $z, a$ outside these
  51. branch points, it is an entire function of $s$.
  52. Examples
  53. ========
  54. The Lerch transcendent is a fairly general function, for this reason it does
  55. not automatically evaluate to simpler functions. Use ``expand_func()`` to
  56. achieve this.
  57. If $z=1$, the Lerch transcendent reduces to the Hurwitz zeta function:
  58. >>> from sympy import lerchphi, expand_func
  59. >>> from sympy.abc import z, s, a
  60. >>> expand_func(lerchphi(1, s, a))
  61. zeta(s, a)
  62. More generally, if $z$ is a root of unity, the Lerch transcendent
  63. reduces to a sum of Hurwitz zeta functions:
  64. >>> expand_func(lerchphi(-1, s, a))
  65. zeta(s, a/2)/2**s - zeta(s, a/2 + 1/2)/2**s
  66. If $a=1$, the Lerch transcendent reduces to the polylogarithm:
  67. >>> expand_func(lerchphi(z, s, 1))
  68. polylog(s, z)/z
  69. More generally, if $a$ is rational, the Lerch transcendent reduces
  70. to a sum of polylogarithms:
  71. >>> from sympy import S
  72. >>> expand_func(lerchphi(z, s, S(1)/2))
  73. 2**(s - 1)*(polylog(s, sqrt(z))/sqrt(z) -
  74. polylog(s, sqrt(z)*exp_polar(I*pi))/sqrt(z))
  75. >>> expand_func(lerchphi(z, s, S(3)/2))
  76. -2**s/z + 2**(s - 1)*(polylog(s, sqrt(z))/sqrt(z) -
  77. polylog(s, sqrt(z)*exp_polar(I*pi))/sqrt(z))/z
  78. The derivatives with respect to $z$ and $a$ can be computed in
  79. closed form:
  80. >>> lerchphi(z, s, a).diff(z)
  81. (-a*lerchphi(z, s, a) + lerchphi(z, s - 1, a))/z
  82. >>> lerchphi(z, s, a).diff(a)
  83. -s*lerchphi(z, s + 1, a)
  84. See Also
  85. ========
  86. polylog, zeta
  87. References
  88. ==========
  89. .. [1] Bateman, H.; Erdelyi, A. (1953), Higher Transcendental Functions,
  90. Vol. I, New York: McGraw-Hill. Section 1.11.
  91. .. [2] https://dlmf.nist.gov/25.14
  92. .. [3] https://en.wikipedia.org/wiki/Lerch_transcendent
  93. """
  94. def _eval_expand_func(self, **hints):
  95. z, s, a = self.args
  96. if z == 1:
  97. return zeta(s, a)
  98. if s.is_Integer and s <= 0:
  99. t = Dummy('t')
  100. p = Poly((t + a)**(-s), t)
  101. start = 1/(1 - t)
  102. res = S.Zero
  103. for c in reversed(p.all_coeffs()):
  104. res += c*start
  105. start = t*start.diff(t)
  106. return res.subs(t, z)
  107. if a.is_Rational:
  108. # See section 18 of
  109. # Kelly B. Roach. Hypergeometric Function Representations.
  110. # In: Proceedings of the 1997 International Symposium on Symbolic and
  111. # Algebraic Computation, pages 205-211, New York, 1997. ACM.
  112. # TODO should something be polarified here?
  113. add = S.Zero
  114. mul = S.One
  115. # First reduce a to the interaval (0, 1]
  116. if a > 1:
  117. n = floor(a)
  118. if n == a:
  119. n -= 1
  120. a -= n
  121. mul = z**(-n)
  122. add = Add(*[-z**(k - n)/(a + k)**s for k in range(n)])
  123. elif a <= 0:
  124. n = floor(-a) + 1
  125. a += n
  126. mul = z**n
  127. add = Add(*[z**(n - 1 - k)/(a - k - 1)**s for k in range(n)])
  128. m, n = S([a.p, a.q])
  129. zet = exp_polar(2*pi*I/n)
  130. root = z**(1/n)
  131. up_zet = unpolarify(zet)
  132. addargs = []
  133. for k in range(n):
  134. p = polylog(s, zet**k*root)
  135. if isinstance(p, polylog):
  136. p = p._eval_expand_func(**hints)
  137. addargs.append(p/(up_zet**k*root)**m)
  138. return add + mul*n**(s - 1)*Add(*addargs)
  139. # TODO use minpoly instead of ad-hoc methods when issue 5888 is fixed
  140. if isinstance(z, exp) and (z.args[0]/(pi*I)).is_Rational or z in [-1, I, -I]:
  141. # TODO reference?
  142. if z == -1:
  143. p, q = S([1, 2])
  144. elif z == I:
  145. p, q = S([1, 4])
  146. elif z == -I:
  147. p, q = S([-1, 4])
  148. else:
  149. arg = z.args[0]/(2*pi*I)
  150. p, q = S([arg.p, arg.q])
  151. return Add(*[exp(2*pi*I*k*p/q)/q**s*zeta(s, (k + a)/q)
  152. for k in range(q)])
  153. return lerchphi(z, s, a)
  154. def fdiff(self, argindex=1):
  155. z, s, a = self.args
  156. if argindex == 3:
  157. return -s*lerchphi(z, s + 1, a)
  158. elif argindex == 1:
  159. return (lerchphi(z, s - 1, a) - a*lerchphi(z, s, a))/z
  160. else:
  161. raise ArgumentIndexError
  162. def _eval_rewrite_helper(self, target):
  163. res = self._eval_expand_func()
  164. if res.has(target):
  165. return res
  166. else:
  167. return self
  168. def _eval_rewrite_as_zeta(self, z, s, a, **kwargs):
  169. return self._eval_rewrite_helper(zeta)
  170. def _eval_rewrite_as_polylog(self, z, s, a, **kwargs):
  171. return self._eval_rewrite_helper(polylog)
  172. ###############################################################################
  173. ###################### POLYLOGARITHM ##########################################
  174. ###############################################################################
  175. class polylog(DefinedFunction):
  176. r"""
  177. Polylogarithm function.
  178. Explanation
  179. ===========
  180. For $|z| < 1$ and $s \in \mathbb{C}$, the polylogarithm is
  181. defined by
  182. .. math:: \operatorname{Li}_s(z) = \sum_{n=1}^\infty \frac{z^n}{n^s},
  183. where the standard branch of the argument is used for $n$. It admits
  184. an analytic continuation which is branched at $z=1$ (notably not on the
  185. sheet of initial definition), $z=0$ and $z=\infty$.
  186. The name polylogarithm comes from the fact that for $s=1$, the
  187. polylogarithm is related to the ordinary logarithm (see examples), and that
  188. .. math:: \operatorname{Li}_{s+1}(z) =
  189. \int_0^z \frac{\operatorname{Li}_s(t)}{t} \mathrm{d}t.
  190. The polylogarithm is a special case of the Lerch transcendent:
  191. .. math:: \operatorname{Li}_{s}(z) = z \Phi(z, s, 1).
  192. Examples
  193. ========
  194. For $z \in \{0, 1, -1\}$, the polylogarithm is automatically expressed
  195. using other functions:
  196. >>> from sympy import polylog
  197. >>> from sympy.abc import s
  198. >>> polylog(s, 0)
  199. 0
  200. >>> polylog(s, 1)
  201. zeta(s)
  202. >>> polylog(s, -1)
  203. -dirichlet_eta(s)
  204. If $s$ is a negative integer, $0$ or $1$, the polylogarithm can be
  205. expressed using elementary functions. This can be done using
  206. ``expand_func()``:
  207. >>> from sympy import expand_func
  208. >>> from sympy.abc import z
  209. >>> expand_func(polylog(1, z))
  210. -log(1 - z)
  211. >>> expand_func(polylog(0, z))
  212. z/(1 - z)
  213. The derivative with respect to $z$ can be computed in closed form:
  214. >>> polylog(s, z).diff(z)
  215. polylog(s - 1, z)/z
  216. The polylogarithm can be expressed in terms of the lerch transcendent:
  217. >>> from sympy import lerchphi
  218. >>> polylog(s, z).rewrite(lerchphi)
  219. z*lerchphi(z, s, 1)
  220. See Also
  221. ========
  222. zeta, lerchphi
  223. """
  224. @classmethod
  225. def eval(cls, s, z):
  226. if z.is_number:
  227. if z is S.One:
  228. return zeta(s)
  229. elif z is S.NegativeOne:
  230. return -dirichlet_eta(s)
  231. elif z is S.Zero:
  232. return S.Zero
  233. elif s == 2:
  234. dilogtable = _dilogtable()
  235. if z in dilogtable:
  236. return dilogtable[z]
  237. if z.is_zero:
  238. return S.Zero
  239. # Make an effort to determine if z is 1 to avoid replacing into
  240. # expression with singularity
  241. zone = z.equals(S.One)
  242. if zone:
  243. return zeta(s)
  244. elif zone is False:
  245. # For s = 0 or -1 use explicit formulas to evaluate, but
  246. # automatically expanding polylog(1, z) to -log(1-z) seems
  247. # undesirable for summation methods based on hypergeometric
  248. # functions
  249. if s is S.Zero:
  250. return z/(1 - z)
  251. elif s is S.NegativeOne:
  252. return z/(1 - z)**2
  253. if s.is_zero:
  254. return z/(1 - z)
  255. # polylog is branched, but not over the unit disk
  256. if z.has(exp_polar, polar_lift) and (zone or (Abs(z) <= S.One) == True):
  257. return cls(s, unpolarify(z))
  258. def fdiff(self, argindex=1):
  259. s, z = self.args
  260. if argindex == 2:
  261. return polylog(s - 1, z)/z
  262. raise ArgumentIndexError
  263. def _eval_rewrite_as_lerchphi(self, s, z, **kwargs):
  264. return z*lerchphi(z, s, 1)
  265. def _eval_expand_func(self, **hints):
  266. s, z = self.args
  267. if s == 1:
  268. return -log(1 - z)
  269. if s.is_Integer and s <= 0:
  270. u = Dummy('u')
  271. start = u/(1 - u)
  272. for _ in range(-s):
  273. start = u*start.diff(u)
  274. return expand_mul(start).subs(u, z)
  275. return polylog(s, z)
  276. def _eval_is_zero(self):
  277. z = self.args[1]
  278. if z.is_zero:
  279. return True
  280. def _eval_nseries(self, x, n, logx, cdir=0):
  281. from sympy.series.order import Order
  282. nu, z = self.args
  283. z0 = z.subs(x, 0)
  284. if z0 is S.NaN:
  285. z0 = z.limit(x, 0, dir='-' if re(cdir).is_negative else '+')
  286. if z0.is_zero:
  287. # In case of powers less than 1, number of terms need to be computed
  288. # separately to avoid repeated callings of _eval_nseries with wrong n
  289. try:
  290. _, exp = z.leadterm(x)
  291. except (ValueError, NotImplementedError):
  292. return self
  293. if exp.is_positive:
  294. newn = ceiling(n/exp)
  295. o = Order(x**n, x)
  296. r = z._eval_nseries(x, n, logx, cdir).removeO()
  297. if r is S.Zero:
  298. return o
  299. term = r
  300. s = [term]
  301. for k in range(2, newn):
  302. term *= r
  303. s.append(term/k**nu)
  304. return Add(*s) + o
  305. return super(polylog, self)._eval_nseries(x, n, logx, cdir)
  306. ###############################################################################
  307. ###################### HURWITZ GENERALIZED ZETA FUNCTION ######################
  308. ###############################################################################
  309. class zeta(DefinedFunction):
  310. r"""
  311. Hurwitz zeta function (or Riemann zeta function).
  312. Explanation
  313. ===========
  314. For $\operatorname{Re}(a) > 0$ and $\operatorname{Re}(s) > 1$, this
  315. function is defined as
  316. .. math:: \zeta(s, a) = \sum_{n=0}^\infty \frac{1}{(n + a)^s},
  317. where the standard choice of argument for $n + a$ is used. For fixed
  318. $a$ not a nonpositive integer the Hurwitz zeta function admits a
  319. meromorphic continuation to all of $\mathbb{C}$; it is an unbranched
  320. function with a simple pole at $s = 1$.
  321. The Hurwitz zeta function is a special case of the Lerch transcendent:
  322. .. math:: \zeta(s, a) = \Phi(1, s, a).
  323. This formula defines an analytic continuation for all possible values of
  324. $s$ and $a$ (also $\operatorname{Re}(a) < 0$), see the documentation of
  325. :class:`lerchphi` for a description of the branching behavior.
  326. If no value is passed for $a$ a default value of $a = 1$ is assumed,
  327. yielding the Riemann zeta function.
  328. Examples
  329. ========
  330. For $a = 1$ the Hurwitz zeta function reduces to the famous Riemann
  331. zeta function:
  332. .. math:: \zeta(s, 1) = \zeta(s) = \sum_{n=1}^\infty \frac{1}{n^s}.
  333. >>> from sympy import zeta
  334. >>> from sympy.abc import s
  335. >>> zeta(s, 1)
  336. zeta(s)
  337. >>> zeta(s)
  338. zeta(s)
  339. The Riemann zeta function can also be expressed using the Dirichlet eta
  340. function:
  341. >>> from sympy import dirichlet_eta
  342. >>> zeta(s).rewrite(dirichlet_eta)
  343. dirichlet_eta(s)/(1 - 2**(1 - s))
  344. The Riemann zeta function at nonnegative even and negative integer
  345. values is related to the Bernoulli numbers and polynomials:
  346. >>> zeta(2)
  347. pi**2/6
  348. >>> zeta(4)
  349. pi**4/90
  350. >>> zeta(0)
  351. -1/2
  352. >>> zeta(-1)
  353. -1/12
  354. >>> zeta(-4)
  355. 0
  356. The specific formulae are:
  357. .. math:: \zeta(2n) = -\frac{(2\pi i)^{2n} B_{2n}}{2(2n)!}
  358. .. math:: \zeta(-n,a) = -\frac{B_{n+1}(a)}{n+1}
  359. No closed-form expressions are known at positive odd integers, but
  360. numerical evaluation is possible:
  361. >>> zeta(3).n()
  362. 1.20205690315959
  363. The derivative of $\zeta(s, a)$ with respect to $a$ can be computed:
  364. >>> from sympy.abc import a
  365. >>> zeta(s, a).diff(a)
  366. -s*zeta(s + 1, a)
  367. However the derivative with respect to $s$ has no useful closed form
  368. expression:
  369. >>> zeta(s, a).diff(s)
  370. Derivative(zeta(s, a), s)
  371. The Hurwitz zeta function can be expressed in terms of the Lerch
  372. transcendent, :class:`~.lerchphi`:
  373. >>> from sympy import lerchphi
  374. >>> zeta(s, a).rewrite(lerchphi)
  375. lerchphi(1, s, a)
  376. See Also
  377. ========
  378. dirichlet_eta, lerchphi, polylog
  379. References
  380. ==========
  381. .. [1] https://dlmf.nist.gov/25.11
  382. .. [2] https://en.wikipedia.org/wiki/Hurwitz_zeta_function
  383. """
  384. @classmethod
  385. def eval(cls, s, a=None):
  386. if a is S.One:
  387. return cls(s)
  388. elif s is S.NaN or a is S.NaN:
  389. return S.NaN
  390. elif s is S.One:
  391. return S.ComplexInfinity
  392. elif s is S.Infinity:
  393. return S.One
  394. elif a is S.Infinity:
  395. return S.Zero
  396. sint = s.is_Integer
  397. if a is None:
  398. a = S.One
  399. if sint and s.is_nonpositive:
  400. return bernoulli(1-s, a) / (s-1)
  401. elif a is S.One:
  402. if sint and s.is_even:
  403. return -(2*pi*I)**s * bernoulli(s) / (2*factorial(s))
  404. elif sint and a.is_Integer and a.is_positive:
  405. return cls(s) - harmonic(a-1, s)
  406. elif a.is_Integer and a.is_nonpositive and \
  407. (s.is_integer is False or s.is_nonpositive is False):
  408. return S.NaN
  409. def _eval_rewrite_as_bernoulli(self, s, a=1, **kwargs):
  410. if a == 1 and s.is_integer and s.is_nonnegative and s.is_even:
  411. return -(2*pi*I)**s * bernoulli(s) / (2*factorial(s))
  412. return bernoulli(1-s, a) / (s-1)
  413. def _eval_rewrite_as_dirichlet_eta(self, s, a=1, **kwargs):
  414. if a != 1:
  415. return self
  416. s = self.args[0]
  417. return dirichlet_eta(s)/(1 - 2**(1 - s))
  418. def _eval_rewrite_as_lerchphi(self, s, a=1, **kwargs):
  419. return lerchphi(1, s, a)
  420. def _eval_is_finite(self):
  421. return fuzzy_not((self.args[0] - 1).is_zero)
  422. def _eval_expand_func(self, **hints):
  423. s = self.args[0]
  424. a = self.args[1] if len(self.args) > 1 else S.One
  425. if a.is_integer:
  426. if a.is_positive:
  427. return zeta(s) - harmonic(a-1, s)
  428. if a.is_nonpositive and (s.is_integer is False or
  429. s.is_nonpositive is False):
  430. return S.NaN
  431. return self
  432. def fdiff(self, argindex=1):
  433. if len(self.args) == 2:
  434. s, a = self.args
  435. else:
  436. s, a = self.args + (1,)
  437. if argindex == 2:
  438. return -s*zeta(s + 1, a)
  439. else:
  440. raise ArgumentIndexError
  441. def _eval_as_leading_term(self, x, logx, cdir):
  442. if len(self.args) == 2:
  443. s, a = self.args
  444. else:
  445. s, a = self.args + (S.One,)
  446. try:
  447. c, e = a.leadterm(x)
  448. except NotImplementedError:
  449. return self
  450. if e.is_negative and not s.is_positive:
  451. raise NotImplementedError
  452. return super(zeta, self)._eval_as_leading_term(x, logx=logx, cdir=cdir)
  453. class dirichlet_eta(DefinedFunction):
  454. r"""
  455. Dirichlet eta function.
  456. Explanation
  457. ===========
  458. For $\operatorname{Re}(s) > 0$ and $0 < x \le 1$, this function is defined as
  459. .. math:: \eta(s, a) = \sum_{n=0}^\infty \frac{(-1)^n}{(n+a)^s}.
  460. It admits a unique analytic continuation to all of $\mathbb{C}$ for any
  461. fixed $a$ not a nonpositive integer. It is an entire, unbranched function.
  462. It can be expressed using the Hurwitz zeta function as
  463. .. math:: \eta(s, a) = \zeta(s,a) - 2^{1-s} \zeta\left(s, \frac{a+1}{2}\right)
  464. and using the generalized Genocchi function as
  465. .. math:: \eta(s, a) = \frac{G(1-s, a)}{2(s-1)}.
  466. In both cases the limiting value of $\log2 - \psi(a) + \psi\left(\frac{a+1}{2}\right)$
  467. is used when $s = 1$.
  468. Examples
  469. ========
  470. >>> from sympy import dirichlet_eta, zeta
  471. >>> from sympy.abc import s
  472. >>> dirichlet_eta(s).rewrite(zeta)
  473. Piecewise((log(2), Eq(s, 1)), ((1 - 2**(1 - s))*zeta(s), True))
  474. See Also
  475. ========
  476. zeta
  477. References
  478. ==========
  479. .. [1] https://en.wikipedia.org/wiki/Dirichlet_eta_function
  480. .. [2] Peter Luschny, "An introduction to the Bernoulli function",
  481. https://arxiv.org/abs/2009.06743
  482. """
  483. @classmethod
  484. def eval(cls, s, a=None):
  485. if a is S.One:
  486. return cls(s)
  487. if a is None:
  488. if s == 1:
  489. return log(2)
  490. z = zeta(s)
  491. if not z.has(zeta):
  492. return (1 - 2**(1-s)) * z
  493. return
  494. elif s == 1:
  495. from sympy.functions.special.gamma_functions import digamma
  496. return log(2) - digamma(a) + digamma((a+1)/2)
  497. z1 = zeta(s, a)
  498. z2 = zeta(s, (a+1)/2)
  499. if not z1.has(zeta) and not z2.has(zeta):
  500. return z1 - 2**(1-s) * z2
  501. def _eval_rewrite_as_zeta(self, s, a=1, **kwargs):
  502. from sympy.functions.special.gamma_functions import digamma
  503. if a == 1:
  504. return Piecewise((log(2), Eq(s, 1)), ((1 - 2**(1-s)) * zeta(s), True))
  505. return Piecewise((log(2) - digamma(a) + digamma((a+1)/2), Eq(s, 1)),
  506. (zeta(s, a) - 2**(1-s) * zeta(s, (a+1)/2), True))
  507. def _eval_rewrite_as_genocchi(self, s, a=S.One, **kwargs):
  508. from sympy.functions.special.gamma_functions import digamma
  509. return Piecewise((log(2) - digamma(a) + digamma((a+1)/2), Eq(s, 1)),
  510. (genocchi(1-s, a) / (2 * (s-1)), True))
  511. def _eval_evalf(self, prec):
  512. if all(i.is_number for i in self.args):
  513. return self.rewrite(zeta)._eval_evalf(prec)
  514. class riemann_xi(DefinedFunction):
  515. r"""
  516. Riemann Xi function.
  517. Examples
  518. ========
  519. The Riemann Xi function is closely related to the Riemann zeta function.
  520. The zeros of Riemann Xi function are precisely the non-trivial zeros
  521. of the zeta function.
  522. >>> from sympy import riemann_xi, zeta
  523. >>> from sympy.abc import s
  524. >>> riemann_xi(s).rewrite(zeta)
  525. s*(s - 1)*gamma(s/2)*zeta(s)/(2*pi**(s/2))
  526. References
  527. ==========
  528. .. [1] https://en.wikipedia.org/wiki/Riemann_Xi_function
  529. """
  530. @classmethod
  531. def eval(cls, s):
  532. from sympy.functions.special.gamma_functions import gamma
  533. z = zeta(s)
  534. if s in (S.Zero, S.One):
  535. return S.Half
  536. if not isinstance(z, zeta):
  537. return s*(s - 1)*gamma(s/2)*z/(2*pi**(s/2))
  538. def _eval_rewrite_as_zeta(self, s, **kwargs):
  539. from sympy.functions.special.gamma_functions import gamma
  540. return s*(s - 1)*gamma(s/2)*zeta(s)/(2*pi**(s/2))
  541. class stieltjes(DefinedFunction):
  542. r"""
  543. Represents Stieltjes constants, $\gamma_{k}$ that occur in
  544. Laurent Series expansion of the Riemann zeta function.
  545. Examples
  546. ========
  547. >>> from sympy import stieltjes
  548. >>> from sympy.abc import n, m
  549. >>> stieltjes(n)
  550. stieltjes(n)
  551. The zero'th stieltjes constant:
  552. >>> stieltjes(0)
  553. EulerGamma
  554. >>> stieltjes(0, 1)
  555. EulerGamma
  556. For generalized stieltjes constants:
  557. >>> stieltjes(n, m)
  558. stieltjes(n, m)
  559. Constants are only defined for integers >= 0:
  560. >>> stieltjes(-1)
  561. zoo
  562. References
  563. ==========
  564. .. [1] https://en.wikipedia.org/wiki/Stieltjes_constants
  565. """
  566. @classmethod
  567. def eval(cls, n, a=None):
  568. if a is not None:
  569. a = sympify(a)
  570. if a is S.NaN:
  571. return S.NaN
  572. if a.is_Integer and a.is_nonpositive:
  573. return S.ComplexInfinity
  574. if n.is_Number:
  575. if n is S.NaN:
  576. return S.NaN
  577. elif n < 0:
  578. return S.ComplexInfinity
  579. elif not n.is_Integer:
  580. return S.ComplexInfinity
  581. elif n is S.Zero and a in [None, 1]:
  582. return S.EulerGamma
  583. if n.is_extended_negative:
  584. return S.ComplexInfinity
  585. if n.is_zero and a in [None, 1]:
  586. return S.EulerGamma
  587. if n.is_integer == False:
  588. return S.ComplexInfinity
  589. @cacheit
  590. def _dilogtable():
  591. return {
  592. S.Half: pi**2/12 - log(2)**2/2,
  593. Integer(2) : pi**2/4 - I*pi*log(2),
  594. -(sqrt(5) - 1)/2 : -pi**2/15 + log((sqrt(5)-1)/2)**2/2,
  595. -(sqrt(5) + 1)/2 : -pi**2/10 - log((sqrt(5)+1)/2)**2,
  596. (3 - sqrt(5))/2 : pi**2/15 - log((sqrt(5)-1)/2)**2,
  597. (sqrt(5) - 1)/2 : pi**2/10 - log((sqrt(5)-1)/2)**2,
  598. I : I*S.Catalan - pi**2/48,
  599. -I : -I*S.Catalan - pi**2/48,
  600. 1 - I : pi**2/16 - I*S.Catalan - pi*I/4*log(2),
  601. 1 + I : pi**2/16 + I*S.Catalan + pi*I/4*log(2),
  602. (1 - I)/2 : -log(2)**2/8 + pi*I*log(2)/8 + 5*pi**2/96 - I*S.Catalan
  603. }