Evaluating the Madelung Constant


The Madelung constant is used in determining the electrostatic potential of a single ion in a crystal by approximating the ions by point charges:

$$ V_i = \frac{e}{4\pi\epsilon_0r_0}\sum_j\frac{z_jr_0}{r_{ij}} = \frac{e}{4\pi\epsilon_0r_0}M_i $$

where $r_0$ is the nearest neighbour distance, and the Madelung constant of the ith ion is given by the sum over all the other ions in the crystal:

$$ M_i = \sum_j\frac{z_j}{r_{ij}/r_0}. $$

For the cubic crystal NaCl, the summation may be performed over three orthogonal co-ordinates:

$$ M_\mathrm{Na} = \sum_{j,k,l=-\infty}^{\infty}' \frac{(-1)^{j+k+l}}{\sqrt{j^2+k^2+l^2}}, $$

where the prime indicates that the term $(0,0,0)$ is excluded. Benson's formula provides a practical and efficient way of evaluating this sum as:

$$ M_\mathrm{Na^+} = -12\pi\sum_{m,n=1,3,...}^{\infty}\mathrm{sech}^2\left(\frac{1}{2}\pi\sqrt{m^2+n^2}\right), $$

where the summation is performed over all positive odd integers $m$ and $n$.

Write a program to calculate the Madelung constant for $\mathrm{Na^+}$ ions in NaCl, -1.74756...