In the chatroom we discussed about the asymptotic of $\displaystyle \sum_{i = 1}^n \sum_{j = 1}^n \frac1{i^2+j^2}$, and if we think of the inverse tangent integral, it's easy to see that $\displaystyle \sum_{i = 1}^n \sum_{j = 1}^n \frac1{i^2+j^2}\approx \operatorname{Ti}_2(n)\approx \frac{\pi}{2} \log(n)$ where I used
the well-known relation $\displaystyle \operatorname{Ti}_2(x)-\operatorname{Ti}_2(1/x)=\frac{\pi}{2}\operatorname{sgn}(x) \log|x|$. At this point, @robjohn posed the following limit

$$\lim_{n\to\infty} \left(\displaystyle \sum_{i = 1}^n \sum_{j = 1}^n \frac1{i^2+j^2}-\frac{\pi}{2} \log(n)\right)$$ that looks like a pretty tough limit. Using $\coth(z)$ one can see the limit is approximately $-\dfrac{\pi^2}{12}$,
but how about finding a way to precisely compute the limit?

  • $\begingroup$ The limit appears to be $-\pi^2/12$. Won't let me correct a two character typo. $\endgroup$ – amcalde Sep 22 '14 at 17:21
  • $\begingroup$ I can see why it's close to that limit, but I cannot prove that it IS that limit. $\endgroup$ – amcalde Sep 22 '14 at 17:23
  • 1
    $\begingroup$ Extrapolating from $n=100$, $n=1000$, and $n=10000$, this looks to be asymptotic to $-0.825861+1.785/n$. $\endgroup$ – robjohn Sep 22 '14 at 18:11

I happened to struggle on the same problem 3 years ago. Here's another approach.

Start from $$\sum_{n=1}^\infty\frac1{x^2+n^2}=-\frac1{2x^2}+\frac{\pi}{2x}+\frac{\pi}{x(e^{2\pi x}-1)}$$ From this we have \begin{array} 1\sum_{m\le N}\sum_{n\le N}\frac{1}{m^2+n^2}&=\sum_{m\le N}\left(-\frac{1}{2m^2}+\frac{\pi}{2m}+\frac{\pi}{m(e^{2\pi m}-1)}-\sum_{n>N}\frac{1}{m^2+n^2}\right)\\ &=-\frac12\zeta(2)+\frac\pi2H_N+\pi\sum_{m\le N}\frac1{m(e^{2\pi m}-1)}-\sum_{m\le N}\sum_{n>N}\frac1{m^2+n^2}+O\left(\frac1N\right)\\ \end{array} Now $$\lim_{N\rightarrow\infty}\sum_{m\le N}\sum_{n>N}\frac1{m^2+n^2}=\lim_{N\rightarrow\infty}\int_N^\infty\int_0^N\frac{\mathrm{d}x\,\mathrm{d}y}{x^2+y^2}=\int_1^\infty\frac1t\tan^{-1}\frac1t\,\mathrm{d}t=G$$ and $$\sum_{m=1}^\infty\frac1{m(e^{2\pi m}-1)}=\sum_{m=1}^\infty\frac1m\sum_{k=1}^\infty e^{-2\pi mk}=-\sum_{k=1}^\infty \log(1-e^{-2\pi k})=-\log(\eta(i)e^{\pi/12})$$ The value of $\eta(i)$ is known to be $\Gamma(1/4)/(2\pi^{3/4})$.

Hence we have \begin{array} 1&\lim_{N\rightarrow\infty}\left(\sum_{m\le N}\sum_{n\le N}\frac{1}{m^2+n^2}-\frac\pi2\log N\right)\\ &=\pi\left(-\log\Gamma\left(\frac14\right)-\frac\pi{12}+\log2+\frac34\log\pi\right)-\frac{\zeta(2)}2+\frac{\pi\gamma}2-G\\ &=\frac14\pi\log\left(\frac{16\pi^3e^{2\gamma}}{\Gamma(\frac14)^4}\right)-G-\zeta(2) \end{array} same as Kirill's answer.

  • $\begingroup$ I also love this approach because I can understand each step of your answer. BTW, how to prove the first line:$$\sum_{n=1}^\infty\frac1{x^2+n^2}=-\frac1{2x^2}+\frac{\pi}{2x}+\frac{\pi}{x(e^{2\pi x}-1)}$$It seems I have seen this equation here, but I can't find it. +1 anyway $\endgroup$ – Anastasiya-Romanova 秀 Sep 25 '14 at 10:09
  • $\begingroup$ @Anastasiya-Romanova: I use this formula (altered slightly) in this answer as $(2)$ and give a link to a proof there. $\endgroup$ – robjohn Sep 25 '14 at 11:37
  • $\begingroup$ @Anastasiya-Romanova: Differentiate the natural logarithm of Euler's infinite product formula for the sine function. :-) $\endgroup$ – Lucian Sep 25 '14 at 11:49
  • $\begingroup$ @robjohn Ah!? I remember now, I saw it when I tried to look up Laila Podlesny's question. I've upvoted your answer anyway. Thanks. Sorry for the bad word in chatroom, I was not in a good mood $\endgroup$ – Anastasiya-Romanova 秀 Sep 25 '14 at 12:42
  • $\begingroup$ @Lucian Did you mean this: $\displaystyle \sin( x) = x \prod_{k=1}^\infty \left( 1 - \frac{x^2}{\pi^2k^2} \right)$? I take the log & then differentiate it? Edit: I get it, like this one: math.stackexchange.com/a/143157/133248 Thanks! $\endgroup$ – Anastasiya-Romanova 秀 Sep 25 '14 at 12:47

The limit in question is equal to $$\def\tfrac#1#2{{\textstyle\frac{#1}{#2}}} \tfrac14\pi\log\left(\frac{16\pi^3e^{2\gamma}}{\Gamma(\frac14)^4}\right) -G-\zeta(2) \\ = -0.82586\ 11759\ 78831\ 08201\ 02008\ 35613\ 80953\ 63017\ 94512\ 34066\ 96955\ 08772 \ldots $$ in terms of Euler's $\gamma$ and Catalan $G$.

Demonstration. The sum $$ S(N) = \sum_{1\leq m,n\leq N} \frac1{m^2+n^2} \sim \tfrac12\pi\log N $$ is in an inconvenient form, because it is taken over the square $[1,N]^2$ instead of the circle of radius $N$: $$ S_2(N) = \sum_{m,n\geq 1}\frac{[m^2+n^2\leq N^2]}{m^2+n^2} \sim \tfrac12\pi\log N. $$ This is so, because the circle-sum can be rewritten as a one-dimensional sum involving the function $r_2(k)$ that counts the number of ways to write $k$ as a sum of two squares of (signed) integers. So we can write the sum $S_2$ as a quarter of the sum over a circle in $\mathbb{Z}^2$ excluding the axes, like so: $$\begin{eqnarray} S_2(N) &=& \frac14 \sum_{m,n\in\mathbb{Z}}' \frac{[m^2+n^2\leq N^2]}{m^2+n^2} - H_N^{(2)} \\&=& \frac14\sum_{1\leq n\leq N^2} \frac{r_2(n)}{n} - H_N^{(2)}. \end{eqnarray} $$ Here $H_N^{(2)}$ has the limit $\zeta(2)=\frac{\pi^2}{6}$, and the prime above the sum means the sum omits $(m,n)=(0,0)$.

The sum $$ g(n) = \sum_{1\leq k\leq n} \frac{r_2(k)}{k} = \pi\log N + \pi S + O(N^{-1/2}) $$ has a standard asymptotic form in terms of the Sierpinski constant $$ S = \log \left(\frac{4\pi^3e^{2\gamma}}{\Gamma(\frac14)^4}\right). $$

The difference $S(N)-S_2(N)$ between the circle-sum and the square-sum does not diverge as $N\to\infty$, and we may sandwich it above and below with two integrals that converge to the same value of $$ \int_1^N dx\int_1^N dy \frac{[x^2+y^2>N^2]}{x^2+y^2}. $$ This integral can be computed in polar coordinates, by writing $$ \begin{eqnarray} 2\int_{N}^{N\sqrt{2}}\frac{dr}{r} \int_0^{\pi/4}d\phi\,[r\cos\phi < N] &=& 2\int_1^{\sqrt{2}}\frac{dr}{r}\left(\frac\pi4 - \arccos\frac1r\right) \\ &=& \tfrac12\pi\log2 - G. \end{eqnarray}$$

Putting everything together leads to the expression given above.

  • $\begingroup$ (+1) I like very much your approach, it looks like a gem. Do you know how we can get that asymptotic involving the Sierpinski constant? $\endgroup$ – user 1591719 Sep 25 '14 at 7:10

Here is the beginning of an asymptotic expansion. First note that $$ \begin{align} \arctan\left(\frac kn\right)-\arctan\left(\frac {k-1}n\right) &=\arctan\left(\frac{\frac1n}{1+\frac kn\frac {k-1}n}\right)\\ &=\arctan\left(\frac{n}{n^2+k^2-k}\right)\tag{1} \end{align} $$ Then $$ \begin{align}\hspace{-1cm} \sum_{k=1}^n\frac1{n^2+k^2}-\frac\pi{4n} &=\sum_{k=1}^n\left[\frac1{n^2+k^2}-\frac1n\arctan\left(\frac{n}{n^2+k^2-k}\right)\right]\\ &=\sum_{k=1}^n\left[\frac1{n^2+k^2}-\frac1{n^2+k^2-k}+\frac{n^2}{3(n^2+k^2-k)^3}-\dots\right]\\ &=\sum_{k=1}^n\left[-\frac{k}{(n^2+k^2)(n^2+k^2-k)}+\frac{n^2}{3(n^2+k^2-k)^3}-\dots\right]\\ &\sim\int_0^1\left[-\frac1{n^2}\frac{x}{(1+x^2)^2}+\frac1{n^3}\frac1{3(1+x^2)^3}-\dots\right]\mathrm{d}x\\ &=-\frac1{4n^2}+O\left(\frac1{n^3}\right)\tag{2} \end{align} $$ Therefore, $$ \sum_{k=1}^n\frac1{n^2+k^2}=\frac\pi{4n}-\frac1{4n^2}+O\left(\frac1{n^3}\right)\tag{3} $$ By computing the edges of the square in index space, we get $$ \begin{align} \sum_{j,k=1}^n\frac1{j^2+k^2}-\sum_{j,k=1}^{n-1}\frac1{j^2+k^2} &=2\sum_{k=1}^n\frac1{n^2+k^2}-\frac1{2n^2}\\ &=\frac\pi{2n}-\frac1{n^2}+O\left(\frac1{n^3}\right)\tag{4} \end{align} $$ The Euler-Maclaurin Sum Formula applied to $(4)$ gives $$ \sum_{j,k=1}^n\frac1{j^2+k^2}=\frac\pi2\log(n)+C+\frac{\pi+4}{4n}+O\left(\frac1{n^2}\right)\tag{5} $$ for some constant $C$.

Equation $(5)$ supports my comment that computationally, $$ \sum_{j,k=1}^n\frac1{j^2+k^2}-\frac\pi2\log(n)\sim C+\frac{\pi+4}{4n}\tag{6} $$ where $C\approx-0.82586118$ and $\frac{\pi+4}{4}\approx1.785398$.

  • $\begingroup$ Because of the limited precision, ISC says this might be $\frac15(5\cdot9^{2/3}-21^{1/2})^{1/2}$ $\endgroup$ – robjohn Sep 23 '14 at 2:00

Using the divergent series $$ \sum_{j=0}^\infty(-1)^j\binom{2j}{n} =\mathrm{Re}\left[\frac{i^n}{(1-i)^{n+1}}\right]\tag{1} $$ which is gotten by computing $$ \lim_{\alpha\to1^-}\sum_{j=0}^\infty(-1)^j\binom{2j}{n}\alpha^{2j}\tag{2} $$ we get $$ \begin{align}\hspace{-1cm} &\sum_{k=1}^n\frac1{n^2+k^2}\\ &=\frac1{n^2}\sum_{k=1}^n\frac1{1+(k/n)^2}\\ &=\frac1{n^2}\sum_{k=1}^n\sum_{j=0}^\infty(-1)^j\frac{k^{2j}}{n^{2j}}\\ &=\sum_{j=0}^\infty\frac{(-1)^j}{n^{2j+2}}\sum_{k=1}^nk^{2j}\\ &=\sum_{j=0}^\infty\frac{(-1)^j}{n^{2j+2}}\left[\frac{n^{2j+1}}{2j+1}+\binom{2j}{0}\frac{n^{2j}}{2}+\binom{2j}{1}\frac{n^{2j-1}}{12}-\binom{2j}{3}\frac{n^{2j-3}}{120}+\dots\right]\\ &=\frac\pi{4n}-\frac1{4n^2}-\frac1{24n^3}+\frac1{2016n^7}-\frac1{4224n^{11}}+\frac1{1536n^{15}}+O\left(\frac1{n^{19}}\right)\tag{3} \end{align} $$ which is using $(1)$ to evaluate $$ \lim_{\alpha\to1^-}\sum_{k=1}^n\frac1{n^2+\alpha^2k^2}\tag{4} $$ Therefore $$ \begin{align} &\sum_{j,k=1}^n\frac1{j^2+k^2}-\sum_{j,k=1}^{n-1}\frac1{j^2+k^2}\\ &=\frac\pi{2n}-\frac1{n^2}-\frac1{12n^3}+\frac1{1008n^7}-\frac1{2112n^{11}}+\frac1{768n^{15}}+O\left(\frac1{n^{19}}\right)\tag{5} \end{align} $$ Applying the Euler-Maclaurin Sum Formula to $(5)$ we get $$ \begin{align} &\sum_{j,k=1}^n\frac1{j^2+k^2}\\ &=\frac\pi2\log(n)+C+\frac{4+\pi}{4n}-\frac{11+\pi}{24n^2}+\frac1{8n^3}+\frac{5+\pi}{240n^4}-\frac1{30n^5}-\frac{43+12\pi}{6048n^6}\\ &+\frac7{288n^7}+\frac{55+18\pi}{8640n^8}-\frac1{30n^9}-\frac{149+48\pi}{12672n^{10}}+\frac{29}{384n^{11}}+\frac{104845+33168\pi}{3144960n^{12}}\\ &-\frac{691}{2730n^{13}}-\frac{1411+448\pi}{10752n^{14}}+\frac{1793}{1536n^{15}}+\frac{1634635+520848\pi}{2350080n^{16}}-\frac{3617}{510n^{17}}\\ &+O\left(\frac1{n^{18}}\right)\tag{6} \end{align} $$ Computing the double sum for $n=10000$ and using $(6)$ we get that $C$ is approximately $$ -0.82586 11759 78831 08201 02008 35613 80953 63017 94512 34066 96955 08772 47540 16332 $$ which matches Kirill's answer.


Here is a way a physicist gets a quick estimation:

First we note that

$$\sum_{i = 1}^n \sum_{j = 1}^n \frac1{i^2+j^2}\sim\int_0^{n}\int_0^{n}\frac{dxdy}{x^2+y^2}\;(n \to\infty)$$

Now, the double integral $I(n)$, because the integration is elementary I write only the final result:

$$I(n)=\arctan(n)\ln(n)-\int_1^{n}\frac{\ln x}{1+x^2}dx-\int_{\frac{1}{n}}^{1}\frac{\arctan x}{x}dx$$

Because the integrals in the last expression are bounded when $n \to\infty$, then

$$I(n)\rightarrow \frac{\pi}{2}\ln n$$


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.