show/hide this revision's text 3 Tweak to fix MathJax pathology

Reid Barton's approximation

[Revised and expanded to give the answer for all $k>1$ and incorporatefurther terms of an asymptotic expansion as $f(n,2) n \simeq R(nrightarrow \infty$]

Fix $k>1$, and write $a_1=f(1,k)=1$ and$$a_n = f(n,k) :=\frac12 \sum_{k\in{\bf Z}} frac1{1-q^{-n}} \phantom. 2^k sum_{r=1}^{n-1} {n \exp(-2^k n) choose r} (1/k)^{n-r} (1/q)^r a_r$$is corroborated by further numerical computation where $q := k/(k-1)$, so $(1/k) + (see plot below for 1/q) = 1$. Set$2^6 $a_\infty := \leq n frac1{k \leq 2^{13}$)log q}.$$For example, if $k=2$ then $a_\infty = 1 / \log 4 = 0.72134752\ldots$,which $a_n$ seems to approach for large $n$, and accounts likewise forboth $k=6$ (the averaged limit of dice-throwing case) with $.72134752\ldots$, which arises a_\infty = 1/(6 \log 1.2) =0.9141358\ldots$. Indeed as $1/(2 n \log 2)$rightarrow \infty$ we have"$a_n \rightarrow a_\infty$ on average", and in the residual oscillationsense that (for instance)$\sum_{n=1}^N (a_n/n) \sim a_\infty \phantom. \sum_{n=1}^N (1/n)$as $N \rightarrow \infty$. But, whose limiting form as suggested by earlierposted answers to Tim Chow's question, $a_n$ does not converge,though it stays quite close to $a_\infty$: we have$$a_n = a_\infty + \epsilon^{\phantom.}_0(\log_q n) + O(1/n)$$as $n \rightarrow \infty$, where $\epsilon^{\phantom.}_0$ is a smoothfunction of period $t := \log_2 n$ 1$ whose average over ${\bf R} / {\bf Z}$ vanishesbut is not identically zero; for large $k$ (already $k=2$ is largeenough), $\epsilon^{\phantom.}_0$ is a nearly perfect sine waveof period $1$ and with a tiny amplitude $\exp(-\pi^2 k + O(\log k))$, namely\frac1{\log 2}\left|\phantom.\Gamma\bigl(1 frac2{k\log q}\left|\phantom.\Gamma\bigl(1 + \frac{2\pi i}{\log 2}\bigr)\right| q}\bigr)\right|\frac1{\log 2}\left[\frac{(2\pi^2/ frac2{k \log2)}{\sinh(2\pi^2/ log q} \log2)}\right]^{1/2} = left[\frac{(2\pi^2/ \phantom. log q)}{\sinh(2\pi^2/ \log q)}\right]^{1/2}.$$For example, for $k=2$ the amplitude is $7.130117\ldots \cdot 10^{-6}10^{-6}$,in accordance with numerical observation (see previously posted answersand the plot below). For $k=6$ the amplitude is only$Period 8.3206735\ldots \cdot 10^{-23}$ so one must compute well beyondthe usual "double precision" to see the oscillations.

More precisely, there is an asymptotic expansion$1$ $a_n \sim a_\infty + \epsilon^{\phantom.}_0(\log_q n) + n^{-1} \epsilon^{\phantom.}_1(\log_q n) + n^{-2} \epsilon^{\phantom.}_2(\log_q n) + n^{-3} \epsilon^{\phantom.}_3(\log_q n) + \cdots,$$where each $\epsilon^{\phantom.}_j$ is clear smooth function of period $1$whose average over ${\bf R} / {\bf Z}$ vanishes, and— while the series need not converge —truncating it before the term $n^{-j} \epsilon^{\phantom.}_j(\log_q n)$yields an approximation good to within $O(n^{-j})$. The first few$\epsilon^{\phantom.}_j$ still have exponentially small amplitudes,but larger that of $\epsilon^{\phantom.}_0$ by a factor$\sim C_j k^{2j}$ for some $C_j > 0$; for instance, the amplitude of$\epsilon^{\phantom.}_1$ exceeds that of $\epsilon^{\phantom.}_0$by about $2(\pi / \log q)^2 \sim 2 \pi^2 k^2$. So $a_n$ must becomputed up to a somewhat large multiple of $k^2$ before it becomesexperimentally plausible that the residual oscillation $a_n - a_\infty$won't tend to zero in the limit as $n \rightarrow \infty$.

Here's a plot that shows $a_n$ for $k=2$ (so also $q=2$) and$2^6 \leq n \leq 2^{13}$, and compares with the periodic approximation$a_\infty + \epsilon^{\phantom.}_0(\log_q n)$ and the refined approximation$a_\infty + \sum_{j=0}^2 n^{-j} \epsilon^{\phantom.}_j(\log_q n)$.(See http://math.harvard.edu/~elkies/mo11255+.pdf forthe original PDF plot, which can be "zoomed in" to view details.) Thehorizontal coordinate is why $\log_2 n$; the above vertical coordinate is centered at$a_\infty = 1/(2 \log 2)$, showing also the lines $a_\infty \pm 2|a_1|$;black cross-hairs, eventually merging visually into a continuous curve,show the numerical values of $a_n$; and the red and green contours show the smooth approximations.

To obtain this asymptotic expansion, we start by generalizingR.Barton's formula for from $R(n)$ k=2$ to arbitrary $k>1$:$$a_n = \frac1k \sum_{r=0}^\infty \phantom. n q^{-r} (1-q^{-r})^{n-1}.$$[The proof is equivalent the same, but note the exponent $n$ has been correctedto Reid's even though he wrote $n-1$ since we want $n-1$ players eliminated at the $r$-th step,not all $n$; this does not affect the limiting behavior$a_\infty+\epsilon^{\phantom.}_0(\log_q n)$, but is needed to get$\epsilon^{\phantom.}_m$ right for $m>1$.] We would like to approximatethe sum by an integral, which can be evaluated by the change of variable$q^{-r} = z$:$$\frac1k \int_{r=0}^\infty \phantom. n q^{-r} (1-q^{-r})^{n-1}= \frac1{k \log q} \int_0^1 \phantom. n (1-z)^{n-1} dz= \left[-a_\infty(1-z)^n\right]_{z=0}^1 = a_\infty.$$But it takes some effort to get at the error in terms this approximation.

We start by comparing $(1-q^{-r})^{n-1}$ with $\exp(-nq^{-r})$:$$\exp(-nq^{-r}) \cdot \exp \phantom. [nq^{-r} + (n-1) \log(1-q^{-r})]\exp(-nq^{-r}) \cdot \exp \left[q^{-r} - (n-1) \left( \frac{q^{-2r}}2 + \frac{q^{-3r}}3 + \frac{q^{-4r}}4 + \cdots \right) \right].$$The next two steps require justification (as R.Barton noted forthe corresponding steps at the end of his analysis), but thefractional part justification should be straightforward. Expand the second factorin powers of $t$). The rest u := nq^{-r}$, and collect like powers of $n$, obtaining$$\exp(-nq^{-r}) \cdot \left(1 - \frac{u^2-2u}{2n} + \frac{3u^4-20u^3+24u^2}{24n^2}- \frac{u^6-14u^5+52u^4-48u^3}{48n^3} + - \cdots \right).$$Each term $n^{-j} \epsilon_j(\log_q(n))$ ($j=0,1,2,3,\ldots$)will arise from the $n^{-j}$ term in this expansion.

We start with the main term, for $j=0$,which is obtained by applying Poisson summation tothe only one that does not decay with $n$. Define\varphi(x) varphi_0(x) := 2^x q^x \exp(-2^x),and as $x \rightarrow -\infty$. We haveOur zeroth-order approximation to $a_n$ isR(n) = \frac12 frac1k \sum_{k sum_{r=0}^\infty \phantom. \varphi_0(\log_q(n)-r),$$which as $n \rightarrow \infty$ rapidly nears$$\frac1k \sum_{r=-\infty}^\infty \varphi_0(\log_q(n)-r).$$For $k=q=2$, this is equivalent with Reid's formula for $R(n)$,even though he wrote it in terms of the fractional part of $\log_2(n)$,because the sum is clearly invariant under translation of $\log_q(n)$by integers.

We next apply Poisson summation. Since$\sum_{r \in {\bf Z}} \phantom. \varphi(t+k),is a smooth ${\bf Z}$-periodic function of period $1$ in $t = \log_2 n$, which thus t$, it has a Fourier expansionR(n) = \sum_{m\in{\bf Z}} \phantom. a_m c_m e^{2\pi i m t}a_m c_m = \frac12 \int_0^1 R(2^t) \left[ \sum_{r \in {\bf Z}} \phantom. \varphi_0(t+r) \right] \phantom. e^{-2\pi i m t} dt= \frac12 \int_{-\infty}^\infty \varphi(t) varphi_0(t) \phantom. e^{-2\pi i m t} dtChanging the variable of integration from $t$ to $2^t$ q^t$ lets us recognize theFourier transform $\hat\varphi$ as $1/\log(2)$ 1/\log(q)$ times a Gamma integral:\hat\varphi(y) hat\varphi_0(y)= \frac1{\log 2q} \Gamma\Bigl(1 + \frac{2 \pi i y} {\log 2}\Bigr)q}\Bigr).$a_0 = 1 / (2\log 2)$ \log q)$ can again be interpreted as the approximation of theRiemann sum $R(n)$ \sum_{r \in {\bf Z}} \phantom. \varphi_0(t+r)$ by an integral;$|\Gamma(1+i\tau)| = (\pi\tau / \sinh(\pi\tau))^{1/2}$. So we have$$\frac1k \sum_{r \in \bf Z} \phantom. \varphi_0(\log_q(n)-r)= \frac1k \sum_{m \in \bf Z} \phantom. \hat\varphi_0(-m) e^{2\pi i \log_q(n)}= a_\infty + \epsilon_0(\log_q(n))$$where $a_\infty = a_0 / k = 1 / (k \log q)$ as above, and$\epsilon^{\phantom.}_0$, defined by$$\epsilon^{\phantom.}_0(t) =\left[ \sum_{r\in\bf Z} \phantom. \varphi_0(t+r) \right] - a_\infty,$$has the Fourier series$$= \frac1k \sum_{m \neq 0} \hat\varphi_0(-m) e^{2\pi i m t}.$$Taking $m = \pm 1$ recovers the amplitude $2|a_1|$ 2|a_1|/k$ exhibited above;the $m = \pm 2$ and further terms yield an oscillation faster but tinier oscillations,e.g. for $k=2$ the $m=\pm 2$ terms oscillate twice as fast but withmagnitude amplitude only $6.6033857\ldots \cdot 10^{-12}$, and .

The functions $\epsilon^{\phantom.}_j$ appearing in the further terms$n^{-j} \epsilon^{\phantom.}_j(\log_q(n))$of the asymptotic expansion of $a_n$ are smaller yet.

The plot below shows defined similarly by$R(n)$ for $2^6 \leq n epsilon^{\phantom.}_j(t) =\leq 2^{13}$ frac1k \sum_{r\in\bf Z} \phantom. \varphi_j(t+r),$$$$\varphi_j(x) = P_j(q^x) \varphi_0(x) = P_j(q^x) q^x \exp(-q^x)$$and compares it with actual values $f(n,2)$. (See http://math.harvard.edu/~elkies/mo11255.pdf for P_j$ is the original PDF plot, which can be "zoomed coefficient of $n^{-j}$ in " to view details.our power series$$(1-q^r)^{n-1} = \exp(-nq^{-r}) \phantom.\sum_{j=0}^\infty \frac{P_j(nq^{-r})}{n^j}.$$P_0(u)=1, \phantom+P_1(u) The horizontal coordinate is = -(u^2-2u)/2, \phantom+P_2(u) = (3u^4-20u^3+24u^2)/24$\log_2 n$; the vertical coordinate is centered at , etc.Again we apply Poisson to expand $1/(2 \epsilon^{\phantom.}_j(\log_q(n))$in a Fourier series:$$\log 2)$epsilon^{\phantom.}_j(t)= \frac1k \sum_{m \in \bf Z} \hat\varphi_j(-m) e^{2\pi i m t},showing also $$and evaluate the lines Fourier transform $1/(2 \hat\varphi_j$ by integratingwith respect to $q^t$. This yields a linear combination of Gammaintegrals evaluated at $1 + (2\pi i y / \log 2q) + j'$ forintegers $j' \pm 2|a_1|$; the red contour shows in [j,2j]$, giving $R(n)$; and the black cross-hairs, eventually merging visually into \hat\varphi_j$ as acontinuous curvedegree-$2j$ polynomial multiple of $\hat\varphi_0$. The first case is$$&=& \frac1{\log q} \left[ \Gamma\Bigl(2 + \frac{2 \pi i y} {\log q}\Bigr) - \frac12 \Gamma\Bigl(3 + \frac{2 \pi i y} {\log q}\Bigr) \right]&=& \frac1{\log q} \frac{\pi y}{\log q} \left(\frac{2 \pi y}{\log q} - i\right) \phantom. \Gamma\Bigl(1 + \frac{2 \pi i y} {\log q}\Bigr)&=& \frac{\pi y}{\log q} \left(\frac{2 \pi y}{\log q} - i\right) \phantom. \hat\varphi_0(y).$$Note that $\varphi_1(0) = 0$, show so the numerical values constant coefficient of theFourier series for $f(n,2)$.

Numerical computation \epsilon^{\phantom.}_1(t)$ vanishes;this is indeed true of $\epsilon^{\phantom.}_j(t)$ for each $j>0$,because $\hat\varphi_j(0) = \int_{-\infty}^\infty \phi_j(x) \phantom. dx$is the difference $f(n,2) - R(n)$ suggests n^{-j}$ coefficient of a power series in $n^{-1}$ that it oscillates we'vealready identified with the same period but with decreasing amplitude asymptotic to $C/n$ for some constant $C$; a_\infty$. Hence (as can alsobe seen in the blue curves show this with plot above) none of the numerical estimate decaying corrections$2.947 n^{-j} \cdot 10^{-4}$ for $C$.

Much epsilon^{\phantom.}_j(\log_q n)$ biases the same behavior should hold for other average of $k$, with a_n$away from $a_\infty$, even tighter oscillations; e.gwhen $n$ is small enough thatthose corrections are a substantial fraction of theresidual oscillation $\epsilon_0(\log_q n)$. for This leaves$k=6$ \hat\varphi_j(\mp1) e^{\pm 2 \pi i t} / k$ as the Gamma factor leading termsin the expansion of each $a_1$ would be \epsilon^{\phantom.}_j(t)$, so we see as promisedthat $\Gamma(1 + \frac{2\pi i}{\log 1.2})$ which has magnitude about \epsilon^{\phantom.}_j$ is still exponentially small but withan extra factor whose leading term is a multiple of $4.55 (2\pi / \cdot 10^{-23}$log q)^{2j}$.
show/hide this revision's text 2 Correct a few typos ($\varphi$ for $\hat\varphi$), add sentence on $s = \lfloor t \rfloor$

Reid Barton's approximation $$ f(n,2) \simeq R(n) := \frac12 \sum_{k\in{\bf Z}} \phantom. 2^k n \exp(-2^k n) \phantom{\infty\infty}(n\rightarrow\infty) $$ is corroborated by further numerical computation (see plot below for $2^6 \leq n \leq 2^{13}$), and accounts for both the averaged limit of $.72134752\ldots$, which arises as $1/(2 \log 2)$, and the residual oscillation, whose limiting form as a function of $t := \log_2 n$ is a nearly perfect sine wave of period $1$ and amplitude $$ \frac1{\log 2}\left|\phantom.\Gamma\bigl(1 + \frac{2\pi i}{\log 2}\bigr)\right| \phantom.=\phantom. \frac1{\log 2}\left[\frac{(2\pi^2/ \log2)}{\sinh(2\pi^2/ \log2)}\right]^{1/2} = \phantom. 7.130117\ldots \cdot 10^{-6}. $$ This Period $1$ is clear (which is why the above formula for $R(n)$ is equivalent to Reid's even though he wrote it in terms of the fractional part of $t$). The rest is obtained by applying Poisson summation to $$ \varphi(x) := 2^x \exp(-2^x), $$ which as Reid observed decays rapidly both as $x \rightarrow \infty$ and as $x \rightarrow -\infty$. We have $$ R(n) = \frac12 \sum_{k \in {\bf Z}} \phantom. \varphi(t+k), $$ a smooth function of period $1$ in $t = \log_2 n$, which thus has a Fourier expansion $$ R(n) = \sum_{m\in{\bf Z}} \phantom. a_m e^{2\pi i m t} $$ where $$ a_m = \frac12 \int_0^1 R(2^t) \phantom. e^{-2\pi i m t} dt = \frac12 \int_{-\infty}^\infty \varphi(t) \phantom. e^{-2\pi i m t} dt = \frac12\hat\varphi(-m). $$ Changing the variable of integration from $t$ to $2^t$ lets us recognize the Fourier transform $\varphi$ \hat\varphi$ as $1/\log(2)$ times a Gamma integral: $$ \varphi(y) hat\varphi(y) = \frac1{\log 2} \Gamma\Bigl(1 + \frac{2 \pi i y} {\log 2}\Bigr). $$ This gives us the coefficients $a_m$ in closed form. The constant coefficient $a_0 = 1 / (2\log 2)$ can be interpreted as the approximation of the Riemann sum $R(n)$ by an integral; the oscillating terms $a_m e^{2\pi i m t}$ for $m \neq 0$ are the corrections to this approximation, and are small due to the exponential decay of the Gamma function on vertical strips ; — indeed we can compute the magnitude $|a_m|$ in elementary closed form using the formula $|\Gamma(1+i\tau)| = (\pi\tau / \sinh(\pi\tau))^{1/2}$. Taking $m = \pm 1$ recovers the amplitude $2|a_1|$ exhibited above; the $m = \pm 2$ terms yield an oscillation twice as fast but with magnitude $6.6033857\ldots \cdot 10^{-12}$, and further terms are smaller yet.

The plot below shows $R(n)$ for $2^6 \leq n \leq 2^{13}$ and compares it with actual values $f(n,2)$. (See http://math.harvard.edu/~elkies/mo11255.pdf for the original PDF plot, which can be "zoomed in" to view details.) The horizontal coordinate is $\log_2 n$; the vertical coordinate is centered at $1/(2 \log 2)$, showing also the lines $1/(2 \log 2) \pm 2|a_1|$; the red contour shows $R(n)$; and the black cross-hairs, eventually merging visually into a continuous curve, show the numerical values of $f(n,2)$.

Numerical computation of the difference $f(n,2) - R(n)$ suggests that it oscillates with the same period but with decreasing amplitude asymptotic to $C/n$ for some constant $C$; the blue curves show this with the numerical estimate $2.947 \cdot 10^{-4}$ for $C$.

Much the same behavior should hold for other $k$, with even tighter oscillations; e.g. for $k=6$ the Gamma factor in $a_1$ would be $\Gamma(1 + \frac{2\pi i}{\log 1.2})$ which has magnitude about $4.55 \cdot 10^{-23}$.
show/hide this revision's text 1

Reid Barton's approximation $$ f(n,2) \simeq R(n) := \frac12 \sum_{k\in{\bf Z}} \phantom. 2^k n \exp(-2^k n) \phantom{\infty\infty}(n\rightarrow\infty) $$ is corroborated by further numerical computation (see plot below for $2^6 \leq n \leq 2^{13}$), and accounts for both the averaged limit of $.72134752\ldots$, which arises as $1/(2 \log 2)$, and the residual oscillation, whose limiting form as a function of $t := \log_2 n$ is a nearly perfect sine wave of period $1$ and amplitude $$ \frac1{\log 2}\left|\phantom.\Gamma\bigl(1 + \frac{2\pi i}{\log 2}\bigr)\right| \phantom.=\phantom. \frac1{\log 2}\left[\frac{(2\pi^2/ \log2)}{\sinh(2\pi^2/ \log2)}\right]^{1/2} = \phantom. 7.130117\ldots \cdot 10^{-6}. $$ This is obtained by applying Poisson summation to $$ \varphi(x) := 2^x \exp(-2^x), $$ which as Reid observed decays rapidly both as $x \rightarrow \infty$ and as $x \rightarrow -\infty$. We have $$ R(n) = \frac12 \sum_{k \in {\bf Z}} \phantom. \varphi(t+k), $$ a smooth function of period $1$ in $t = \log_2 n$, which thus has a Fourier expansion $$ R(n) = \sum_{m\in{\bf Z}} \phantom. a_m e^{2\pi i m t} $$ where $$ a_m = \frac12 \int_0^1 R(2^t) \phantom. e^{-2\pi i m t} dt = \frac12 \int_{-\infty}^\infty \varphi(t) \phantom. e^{-2\pi i m t} dt = \frac12\hat\varphi(-m). $$ Changing the variable of integration from $t$ to $2^t$ lets us recognize the Fourier transform $\varphi$ as $1/\log(2)$ times a Gamma integral: $$ \varphi(y) = \frac1{\log 2} \Gamma\Bigl(1 + \frac{2 \pi i y} {\log 2}\Bigr). $$ This gives us the coefficients $a_m$ in closed form. The constant coefficient $a_0 = 1 / (2\log 2)$ can be interpreted as the approximation of the Riemann sum $R(n)$ by an integral; the oscillating terms $a_m e^{2\pi i m t}$ for $m \neq 0$ are the corrections to this approximation, and are small due to the exponential decay of the Gamma function on vertical strips; indeed we can compute the magnitude $|a_m|$ in elementary closed form using the formula $|\Gamma(1+i\tau)| = (\pi\tau / \sinh(\pi\tau))^{1/2}$. Taking $m = \pm 1$ recovers the amplitude $2|a_1|$ exhibited above; the $m = \pm 2$ terms yield an oscillation twice as fast but with magnitude $6.6033857\ldots \cdot 10^{-12}$, and further terms are smaller yet.

The plot below shows $R(n)$ for $2^6 \leq n \leq 2^{13}$ and compares it with actual values $f(n,2)$. (See http://math.harvard.edu/~elkies/mo11255.pdf for the original PDF plot, which can be "zoomed in" to view details.) The horizontal coordinate is $\log_2 n$; the vertical coordinate is centered at $1/(2 \log 2)$, showing also the lines $1/(2 \log 2) \pm 2|a_1|$; the red contour shows $R(n)$; and the black cross-hairs, eventually merging visually into a continuous curve, show the numerical values of $f(n,2)$.

Numerical computation of the difference $f(n,2) - R(n)$ suggests that it oscillates with the same period but with decreasing amplitude asymptotic to $C/n$ for some constant $C$; the blue curves show this with the numerical estimate $2.947 \cdot 10^{-4}$ for $C$.

Much the same behavior should hold for other $k$, with even tighter oscillations; e.g. for $k=6$ the Gamma factor in $a_1$ would be $\Gamma(1 + \frac{2\pi i}{\log 1.2})$ which has magnitude about $4.55 \cdot 10^{-23}$.