The series $ \sum\limits_{k=1}^{\infty} \frac1{\sqrt{{k}{(k^2+1)}}}$

Given the series

$$\sum_{k=1}^{\infty} \frac1{\sqrt{{k}{(k^2+1)}}}$$

How can I calculate its exact limit (if that is possible)?

Solution 1:

Rewrite the sum in terms of the zeta function. The resulting sum is alternating, $$\begin{eqnarray*} \sum_{k=1}^{\infty} \frac{1}{\sqrt{{k}{(k^2+1)}}} &=& \sum_{k=1}^{\infty} \frac{1}{k^{3/2}} \frac{1}{\left(1+\frac{1}{k^2}\right)^{1/2}} \\ &=& \sum_{k=1}^{\infty} \frac{1}{k^{3/2}} \sum_{j=0}^\infty (-1)^j \frac{\Gamma(j+\frac{1}{2})}{j!\Gamma(\frac{1}{2})} \left(\frac{1}{k^2}\right)^j \\ &=& \sum_{j=0}^\infty (-1)^j \frac{\Gamma(j+\frac{1}{2})}{j!\Gamma(\frac{1}{2})} \sum_{k=1}^\infty \frac{1}{k^{2j+3/2}} \\ &=& \sum_{j=0}^\infty (-1)^j \underbrace{\frac{(2j)!}{(j!)^2 2^{2j}} \textstyle\zeta(2j+\frac{3}{2})}_{a_j}. \end{eqnarray*}$$ In the second step we expanded $(1+1/k^2)^{-1/2}$ in a binomial series. The last line suggests that there is no nice closed form for the sum. Unfortunately, the convergence properties of this sum are not much better than the original. For large $j$, $a_j \sim 1/\sqrt{j}$, so the sum goes like $\sum_j(-1)^j/\sqrt{j}$.

Let's accelerate the sum using the Euler transform. Then, $$\begin{eqnarray*} \sum_{k=1}^{\infty} \frac{1}{\sqrt{{k}{(k^2+1)}}} &=& \sum_{n=0}^\infty \underbrace{(-1)^n \frac{1}{2^{n+1}} \sum_{k=0}^n (-1)^k {n\choose k} a_{n-k}}_{b_n}. \end{eqnarray*}$$ The sum $\sum_{n=0}^\infty b_n$ is not alternating, but has much better convergence properties than the original. Numerically we find $b_n \sim 1/2^n$.

Below we give the partial sums to $25$ digits. $$\begin{array}{ll} N & \sum_{n=0}^N b_n \\ \hline 1 & 1.818439778099484176707949 \\ 2 & 2.052201427955914030757063 \\ 4 & 2.214850676717845386847939 \\ 8 & 2.261274490824097920144642 \\ 16 & 2.264052732086094056012626 \\ 32 & 2.264062399012245290957591 \\ 64 & 2.264062399141221028566328 \\ 128 & 2.264062399141221028592305 \end{array}$$

Another good approach suggested by @qoqosz in the comments is to estimate the sum with an integral. That estimate can be improved somewhat by applying the Euler-Maclaurin formula. The resulting asymptotic series gets us four digits before it starts to diverge.

Solution 2:

Using the binomial theorem to expand $$ \begin{align} {\small\frac1{\sqrt{k(k^2+1)}}} &=k^{-3/2}\left(1+\frac1{k^2}\right)^{-1/2}\\ &\sim k^{-3/2}-\frac12k^{-7/2}+\frac38k^{-11/2}-\frac5{16}k^{-15/2}+\frac{35}{128}k^{-19/2}-\frac{63}{256}k^{-23/2}\tag{1} \end{align} $$ with an error of approximately $\frac14k^{-27/2}$.

Applying the Euler Maclaurin Sum Formula to $(1)$ yields $$ \begin{align} \sum_{k=1}^n{\small\frac1{\sqrt{k(k^2+1)}}} &\sim C-2n^{-1/2}+\frac12n^{-3/2}+\frac3{40}n^{-5/2}-\frac14n^{-7/2}+\frac{31}{384}n^{-9/2}+\frac3{16}n^{-11/2}\\ &-\frac{12959}{66560}n^{-13/2}-\frac5{32}n^{-15/2}+\frac{214861}{557056}n^{-17/2}+\frac{35}{256}n^{-19/2}\\[6pt] &-\frac{24553421}{27525120}n^{-21/2}-\frac{63}{512}n^{-23/2}\tag{2} \end{align} $$ with an error of approximately $3n^{-25/2}$. Since the series is convergent, $C$ represents the sum.

Computing the partial sum to $n=100$ and using $(2)$ to compute $C$, we get $$ C\stackrel{.}{=}2.264062399141221028592305 $$ Just as a check, with $n=1000$, $(2)$ gives $$ C\stackrel{.}{=}2.264062399141221028592304976971183194 $$

Solution 3:

Here is another possible approach for numerically summing this series. Through Laplace transform methods, we find that

$$\frac1{\sqrt{k(k^2+1)}}=\sqrt{\frac{\pi}{2}}\int_0^\infty \exp(-ku)\sqrt{u}\,J_{-\frac14}\left(\frac{u}{2}\right)J_{\frac14}\left(\frac{u}{2}\right)\mathrm du$$

where $J_\nu(u)$ is the Bessel function of the first kind. Substituting this into the summation expression, interchanging summation and integration, and some finagling yields the integral expression

$$\sum_{k=1}^\infty \frac1{\sqrt{k(k^2+1)}}=2\sqrt\pi\int_0^\infty \frac{\sqrt{t}}{\exp(2t)-1} J_{-\frac14}(t) J_{\frac14}(t)\mathrm dt$$

Oscillatory integrals such as the one given above can be treated as alternating series; if you look at a plot of the integrand, you can see that the integral across the zeroes of the function alternate between positive and negative signs. Thus, letting $t_k$ be the $k$-th zero of the function $J_{-\frac14}(t) J_{\frac14}(t)$, we consider the series

$$2\sqrt\pi\left(\int_0^{t_1} \frac{\sqrt{t}}{\exp(2t)-1} J_{-\frac14}(t) J_{\frac14}(t)\mathrm dt+\sum_{k=1}^\infty \int_{t_k}^{t_{k+1}} \frac{\sqrt{t}}{\exp(2t)-1} J_{-\frac14}(t) J_{\frac14}(t)\mathrm dt\right)$$

(This strategy for numerically evaluating oscillatory integrals is due to Longman.)

We can use a computing environment like Mathematica to generate the $t_k$:

\begin{array}{c|c}k&t_k\\\hline 1&2.006299672\\2&2.780887724\\3&5.123062743\\ 4&5.906142699\\5&8.257951176\\6&9.042383664\end{array}

One could use the Euler transform as in oen's answer to sum the alternating series, but the terms decrease rather quickly in magnitude that even direct summation works fine. Letting $S_k$ be the $k$-th partial sum, we have

\begin{array}{c|c}k&S_k\\\hline 1&2.2642181990016\\2&2.2632428613646\\3&2.2640625078551\\ 4&2.2640612519867\\5&2.2640623992549\\6&2.2640623973568\end{array}

The last approximation is already good to eight digits. Taking more terms of the series, we have, to twenty digits

$$\sum_{k=1}^\infty \frac1{\sqrt{k(k^2+1)}}\approx 2.26406239914122102859$$

Thirteen terms were needed to obtain this approximation.