3
$\begingroup$

I am looking a an efficient algorithm to compute in C derivatives of the Hurwitz zeta function $$ \zeta^{(n)}(s,a) = (-1)^{n} \sum_{k=0}^\infty \frac{\log^{n}(k+a)}{(k+a)^{s}} . $$ for $1<s$ and $0<q$.

(This question is related to the previous question Numerical evaluation of Hurwitz zeta function .)

$\endgroup$

3 Answers 3

2
$\begingroup$

$$\sum_{k\geq 0}\frac{1}{(k+a)^s} \stackrel{\mathcal{L}^{-1}}{=}\frac{1}{\Gamma(s)}\int_{0}^{+\infty}\frac{e^{(1-a)u}u^{s-1}}{e^u-1}\,du $$ so $$ (-1)^n \sum_{k\geq 0}\frac{\log(k+a)^n}{(k+a)^s} = \frac{1}{\Gamma(s)}\int_{0}^{+\infty}\frac{e^{(1-a)u}u^{s-1}\log(u)^n}{e^u-1}\,du $$ and any accurate algorithm for numerical integration can do the job.

$\endgroup$
2
  • $\begingroup$ did I get your answer wrong or is it really missing something? See my answer. $\endgroup$ Feb 6, 2019 at 15:38
  • $\begingroup$ It seems to me that this answer is wrong, because we have to take multiple derivatives of a product $u^s/ \Gamma(s)$ and not just $u^s$. Could you please check? $\endgroup$
    – Yuriy S
    Jul 15, 2019 at 21:21
2
$\begingroup$

I tried commenting on @JackD'Aurizio's answer, but I'm three reputation points shy of being able to. While the idea of using the integral representation is great, I believe the answer given by @JackD'Aurizio is in error, at least for the first derivative.

In what follows I will provide what I believe to be the correct expression for the first derivative. I believe the correct expression is $$ \zeta^{(1)}(s, a) = \frac{1}{\Gamma(s)}\int_0^\infty \frac{e^{(1-a)x}x^{s-1}\left(\log(x) - \psi_0(s) \right)}{e^x -1}dx. $$ where $\psi_0$ is the polygamma function. I obtained this by differentiating $ \zeta(s) = \frac{1}{\Gamma(s)}\int_{0}^\infty \frac{e^{(1-a)x}x^{s-1}}{e^x -1}dx$ with respect to $s$ and then applying the product and reciprocal rules to get $$ \zeta^{(1)}(s, a) = -\frac{\psi_0(s)}{\Gamma(s)}\int_0^\infty \frac{e^{(1-a)x}x^{s-1}}{e^x -1}dx + \frac{1}{\Gamma(s)}\int_0^\infty \frac{e^{(1-a)x}x^{s-1}\log(x)}{e^x -1}dx, $$ from which the result follows. Accompanying R code (for $a = 1$):

zeta_int <- function(s){
  func <- function(x){
    x^(s-1)/(exp(x) - 1)
  }
  func <- Vectorize(func)
  return(integrate(func, 0, Inf)$value/gamma(s))
}
zeta_deriv <- function(s){
  k <- digamma(s)
  func <- function(x){
    (x^(s-1) *(log(x) - k))/(exp(x) -1)
  }
  func <- Vectorize(func)
  ans <- integrate(func, 0, Inf)$value/gamma(s)
  return(ans)
}
####
# Testing
ss <- 3/2
zeta_int(ss)
VGAM::zeta(ss)
zeta_deriv_alt(ss)
VGAM::zeta(ss, deriv = 1)

I know this does not answer the question as posed, but I think it's a valid contribution nonetheless.

$\endgroup$
1
$\begingroup$

If you can handle working in C++ rather than C, use the Chebyshev transform in Boost.Math. You will need a way of evaluate the Hurwitz zeta function to compute its value at the Chebyshev nodes, but once you do, you'll get a stable method of numerical differentiation via differentiating the associated Chebyshev series.

$\endgroup$

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .