Sum of exponential terms and binomial

Computing the binomial coefficients efficiently

If by "gets stuck" you mean that the computation is slow, I would guess that you are computing the binomial term inefficiently. Indeed, you shouldn't recompute the binomial term for every summand, but instead use the fact that $$ \binom{m}{q}=\frac{m!}{q!\left(m-q\right)!}=\frac{m-\left(q-1\right)}{q}\frac{m!}{\left(q-1\right)!\left(m-\left(q-1\right)\right)!}=\frac{m-q+1}{q}\binom{m}{q-1}. $$ Defining $$ C_{q}=\frac{m-q+1}{q}C_{q-1}\text{ if }q\geq1\qquad\text{and}\qquad C_{0}=1, $$ it follows from the previous claim that $C_{q}=\binom{m}{q}$. Therefore, you can rewrite the sum you are interested as $$ S\equiv \sum_{q=1}^{m}\frac{\left(-1\right)^{q+1}}{q+1}C_{q}\exp\left(-\frac{q}{q+1}\Gamma\right). $$

Removing some terms by symmetry

We can use the fact that $C_{q}=C_{m-q}$ to reduce the number of terms. Note that $$ S-1=\sum_{q=0}^{m}\frac{\left(-1\right)^{q+1}}{q+1}\exp\left(-\frac{q}{q+1}\Gamma\right)C_{q}. $$ Assuming $m=2j+1$ is odd, we get $$ S-1=\sum_{q=0}^{j}\left(-1\right)^{q+1}\left(\frac{1}{q+1}\exp\left(-\frac{q}{q+1}\Gamma\right)-\frac{1}{m-q+1}\exp\left(-\frac{m-q}{m-q+1}\Gamma\right)\right)C_{q}. $$ Assuming $m=2j$ is even, we get \begin{multline*} S-1=\frac{\left(-1\right)^{j+1}}{j+1}\exp\left(-\frac{j}{j+1}\Gamma\right)C_{j}\\ +\sum_{q=0}^{j}\left(-1\right)^{q+1}\left(\frac{1}{q+1}\exp\left(-\frac{q}{q+1}\Gamma\right)+\frac{1}{m-q+1}\exp\left(-\frac{m-q}{m-q+1}\Gamma\right)\right)C_{q}. \end{multline*}


This computation can be done numerically without a problem up to $m=2^{13}$. It took about 15 seconds of cpu time using Mathematica on my 2 year old Mac laptop for $m=2^{13}$.

Here is an example of the Mathematica code:

s[m_, g_] := Sum[ (-1)^(q+1)/(q + 1) Binomial[ m , q] Exp[ - q g/(q + 1)], {q, 1, m}];
Print[ Timing[ s[10000, N[1/4, 10000]]//N]; 

The output for the program above is {27.7445,0.999574} indicating that it took 27 seconds to compute the answer. Note that ${1000\choose 500}$ has about 3000 digits, so the program used 10000 digits of precision. The running time is order $m^3$.

The answer is usually close to 1 when $0<q<1$ and $m> 2^{10}$.


I wrote the code in Python and got the same result for $m=2^{13}$ and $q=1/4$.

from mpmath import mp

mp.dps =5000;

m = 2**13;

mp.pretty = True

rS = mp.mpf(0);

g = mp.mpf(1)/mp.mpf(4);

for q in range(1, m+1):
    rS = rS + mp.mpf((-1)**(q+1))* mp.binomial(m, q)/(q+1)*mp.exp(-q*g/(q+1));


mp.nprint(rS, 10)