Accelerated convergence treatment of R−n lattice sums

Williams, D. E.

doi:10.1107/97809553602060000562

International
Tables for
Crystallography
Volume B
Reciprocal space
Edited by U. Shmueli

pdf | chapter contents | chapter index | related articles

International Tables for Crystallography (2006). Vol. B. ch. 3.4, pp. 385-397 | 1 | 2 |
https://doi.org/10.1107/97809553602060000562

Chapter 3.4. Accelerated convergence treatment of R⁻ⁿ lattice sums¹

D. E. Williams^a ^‡

^aDepartment of Chemistry, University of Louisville, Louisville, Kentucky 40292, USA

Accelerated-convergence treatment of the Coulombic lattice sum of ionic crystals and the dispersion-energy sum of ionic or molecular crystals is described. The Coulombic energy is very long-range and convergence of the Coulombic lattice-energy sum is extremely slow. Accelerated-convergence treatment is indispensable to achieve accuracy with a reasonable amount of computational effort. The dispersion-energy sum has somewhat better convergence properties than the Coulombic sum, but accelerated-convergence treatment of the dispersion sum is also strongly recommended since its use can yield at least an order of magnitude improvement in accuracy for a given calculation effort.

Keywords: accelerated convergence; lattice sums; direct-space sums; Coulombic lattice energy; dispersion energy; composite lattices; incomplete gamma function; intramolecular energy terms.

3.4.1. Introduction

| top | pdf |

The electrostatic energy of an ionic crystal is often represented by taking a pairwise sum between charge sites interacting via Coulomb's law (the [n = 1] sum). The individual terms may be positive or negative, depending on whether the pair of sites have charges of the same or different signs. The Coulombic energy is very long-range, and it is well known that convergence of the Coulombic lattice-energy sum is extremely slow. For simple structure types Madelung constants have been calculated which represent the Coulombic energy in terms of the cubic lattice constant or a nearest-neighbour distance. Glasser & Zucker (1980 ) give tables of Madelung constants and review the subject giving references dating back to 1884. If the ionic crystal structure is not of a simple type usually no Madelung constant will be available and the Coulombic energy must be obtained for the specific crystal structure being considered. In carrying out this calculation, accelerated-convergence treatment of the Coulombic lattice sum is indispensable to achieve accuracy with a reasonable amount of computational effort. A model of a molecular crystal may include partial net atomic charges or other charge sites such as lone-pair electrons. The [(n = 1)] sum also applies between these site charges.

The dispersion energy of ionic or molecular crystals may be represented by an [(n = 6)] sum over atomic sites, with possible inclusion of $[(n = 8, 10, \ldots)]$ terms for higher accuracy. The dispersion-energy sum has somewhat better convergence properties than the Coulombic sum. Nevertheless, accelerated-convergence treatment of the dispersion sum is strongly recommended since its use can yield at least an order of magnitude improvement in accuracy for a given calculation effort. The repulsion energy between non-bonded atoms in a crystal may be represented by an exponential function of short range, or possibly by an [(n = 12)] function of short range. The convergence of the repulsion energy is fast and no accelerated-convergence treatment is normally required.

3.4.2. Definition and behaviour of the direct-space sum

| top | pdf |

This pairwise sum is taken between atoms (or sites) in the reference unit cell and all other atoms (or sites) in the crystal, excluding the self terms. Thus, the second atom (or site) is taken to range over the entire crystal, with elimination of self-energy terms. If $[V_{n}]$ represents an energy, each atom is assigned one half of the pair energy. Therefore, the energy per unit cell is $[V_{n} = (1/2) {\textstyle\sum\limits_{j}^{\rm one\; cell}} \ {\textstyle\sum\limits_{k}^{\rm all\; cells}\!\!\!'}\ Q_{jk} R_{jk}^{-n},]$ where $[Q_{jk}]$ is a given coefficient, $[R_{jk}]$ is an interatomic distance, and the prime on the second sum indicates that self terms are omitted. In the case of the Coulombic sum, [n = 1] and $[Q_{jk} = q_{j}q_{k}]$ is the product of the site charges.

Table 3.4.2.1 gives an example of the convergence behaviour of the untreated [(n = 1)] Coulombic sum for sodium chloride. Even at the rather large summation limit of 20 Å the Coulombic lattice sum has not converged and is incorrect by about 8%. The 20 Å sum included 832 molecules and 2494 individual distances. At various smaller summation limits the truncation error fluctuates wildly and can be either positive or negative. Note that the results shown in the table always refer to summation over whole molecules, that is, over neutral charge units.

Table 3.4.2.1 | top | pdf |
Untreated lattice-sum results for the Coulombic energy ( [n = 1] ) of sodium chloride (kJ mol⁻¹, Å); the lattice constant is taken as 5.628 Å

Truncation limit	Number of molecules	Number of terms	Calculated energy
6.0	23	67	−696.933
8.0	59	175	−597.371
10.0	108	322	−915.152
12.0	201	601	−773.475
14.0	277	829	−796.248
16.0	426	1276	−826.502
18.0	587	1759	−658.995
20.0	832	2494	−794.619
Converged value			−862.825

If the Coulombic summation is not carried out over neutral charge units the truncation error is even larger. These considerations support the conclusion that accelerated-convergence treatment of the Coulombic lattice sum should be regarded as mandatory. Table 3.4.2.2 gives an example of the convergence behaviour of the untreated [(n = 6)] dispersion sum for benzene. In obtaining this sum it is not necessary to consider whole molecules as in the Coulombic case. The exclusion of atoms (or sites) in the portions of molecules outside the summation limit greatly reduces the number of terms to be considered. At the summation limit of 20 Å, 439 benzene molecules and 22 049 individual distances are considered; the dispersion-sum truncation error is 0.4%. Thus, if sufficient computer time is available it may be possible to obtain a moderately accurate dispersion sum without the use of accelerated convergence. However, as shown below, the use of accelerated convergence will greatly speed up the calculation, and is in practice necessary if higher accuracy is required.

Table 3.4.2.2 | top | pdf |
Untreated lattice-sum results for the dispersion energy ( [n = 6] ) of crystalline benzene (kJ mol⁻¹, Å)

Truncation limit	Number of molecules	Number of terms	Calculated energy
6.0	26	524	−69.227
8.0	51	1313	−76.007
10.0	77	2631	−78.179
12.0	126	4718	−79.241
14.0	177	7531	−79.726
16.0	265	11274	−80.013
18.0	344	15904	−80.178
20.0	439	22049	−80.295
Converged value			−80.589

3.4.3. Preliminary description of the method

| top | pdf |

Ewald (1921 ) developed a method which modified the mathematical representation of the Coulombic lattice sum to improve the rate of convergence. This method was based on partially transforming the lattice sum into reciprocal space. Bertaut (1952 ) presented another method for derivation of the Ewald result which used the concept of the crystallographic structure factor. His formula extended the Ewald treatment to a composite lattice with more than one atom per lattice point. Nijboer & DeWette (1957 ) developed a general Fourier transform method for the evaluation of $[R^{-n}]$ sums in simple lattices. Williams (1971 ) extended this treatment to a composite lattice and gave general formulae for the $[R^{-n}]$ sums for any crystal. A review article, on which this chapter is based, appeared later (Williams, 1989 ).

Consider a function, W(R), which is unity at [R = 0] and smoothly declines to zero as R approaches infinity. If each term of the lattice sum is multiplied by W(R), the rate of convergence is increased. However, the rate of convergence of the remainder of the original sum, which contains the difference terms, is not increased. $[\eqalign{ V_{n} &= (1/2) {\textstyle\sum\limits_{j}^{\rm one\; cell}} \ {\textstyle\sum\limits_{k}^{\rm all\; cells}\!\!\!'}\ Q_{jk} R_{jk}^{-n} W(R)\cr &\quad + (1/2) {\textstyle\sum\limits_{j}^{\rm one\; cell}} \ {\textstyle\sum\limits_{k}^{\rm all\; cells}\!\!\!'}\ Q_{jk} R_{jk}^{-n} [1 - W(R)].}]$

In the accelerated-convergence method the difference terms are expressed as an integral of the product of two functions. According to Parseval's theorem (described below) this integral is equal to an integral of the product of the two Fourier transforms of the functions. Finally, the integral over the Fourier transforms of the functions is converted to a sum in reciprocal (or Fourier-transform) space. The choice of the convergence function W(R) is not unique; an obvious requirement is that the relevant Fourier transforms must exist and have correct limiting behaviour. Nijboer and DeWette suggested using the incomplete gamma function for W(R). More recently, Fortuin (1977 ) showed that this choice of convergence function leads to optimal convergence of the sums in both direct and reciprocal space: $[W(R) = \Gamma (n/2, \pi w^{2} R^{2})/\Gamma (n/2),]$ where $[\Gamma (n/2)]$ and $[\Gamma (n/2, \pi w^{2} R^{2})]$ are the gamma function and the incomplete gamma function , respectively: $[\Gamma (n/2, \pi w^{2} R^{2}) = \textstyle\int\limits_{\pi w^{2} R^{2}}^{\infty} t^{(n/2)-1} \exp (-t)\;\hbox{d}t]$ and $[\Gamma (n/2) = \Gamma (n/2, 0).]$

The complement of the incomplete gamma function is $[\gamma (n/2, \pi w^{2} R^{2}) = \Gamma (n/2) - \Gamma (n/2, \pi w^{2} R^{2}).]$

3.4.4. Preliminary derivation to obtain a formula which accelerates the convergence of an $[{\bi R}^{\bf -{\bi n}}]$ sum over lattice points X(d)

| top | pdf |

The three-dimensional direct-space crystal lattice is specified by the origin vectors $[{\bf a}_{1}]$ , $[{\bf a}_{2}]$ and $[{\bf a}_{3}]$ . A general vector in direct space is defined as $[{\bf X(x)} = x_{1} {\bf a}_{1} + x_{2} {\bf a}_{2} + x_{3} {\bf a}_{3},]$ where $[x_{1}, x_{2}, x_{3}]$ are the fractional cell coordinates of X. A lattice vector in direct space is defined as $[{\bf X(d)} = d_{1} {\bf a}_{1} + d_{2} {\bf a}_{2} + d_{3} {\bf a}_{3},]$ where $[d_{1}, d_{2}, d_{3}]$ are integers (specifying particular values of $[x_{1}, x_{2}, x_{3}]$ ) designating a lattice point. $[V_{d}]$ is the direct-cell volume which is equal to $[{\bf a}_{1} \cdot {\bf a}_{2} \times {\bf a}_{3}]$ . A general point in the direct lattice is X(x); the contents of the lattice are by definition identical as the components of x are increased or decreased by integer amounts.

The reciprocal-lattice vectors are defined by the relations $[\eqalign{ {\bf a}_{j} \cdot {\bf b}_{k} &= 1 \qquad j = k\cr &=0 \qquad j \neq k.}]$ A general vector in reciprocal space H(r) is defined as $[{\bf H(r)} = r_{1} {\bf b}_{1} + r_{2} {\bf b}_{2} + r_{3} {\bf b}_{3}.]$ A reciprocal-lattice vector H(h) is defined by the integer triplet $[h_{1}, h_{2}, h_{3}]$ (specifying particular values of $[r_{1}, r_{2}, r_{3}]$ ) so that $[{\bf H(h)} = h_{1} {\bf b}_{1} + h_{2} {\bf b}_{2} + h_{3} {\bf b}_{3}.]$ In other sections of this volume a shortened notation h is used for the reciprocal-lattice vector. In this section the symbol H(h) is used to indicate that it is a particular value of H(r).

The three-dimensional Fourier transform $[g({\bf t})]$ of a function $[f({\bf x})]$ is defined by $[g({\bf t}) = FT_{3} [\;f({\bf x})] = \textstyle\int f({\bf x}) \exp (2 \pi i {\bf x} \cdot {\bf t})\; \hbox{d}{\bf x}.]$ The Fourier transform of the set of points defining the direct lattice is the set of points defining the reciprocal lattice, scaled by the direct-cell volume. It is useful for our purpose to express the lattice transform in terms of the Dirac delta function $[\delta (x - x_{o})]$ which is defined so that for any function $[f({\bf x})]$ $[f({\bf x}_{o}) = \textstyle\int \delta ({\bf x} - {\bf x}_{o}) f({\bf x})\; \hbox{d}{\bf x}.]$ We then write $[FT_{3} \{{\textstyle\sum\limits_{\bf d}} \delta [{\bf X(x)} - {\bf X(d)}]\} = V_{d}^{-1} {\textstyle\sum\limits_{\bf h}} \delta [{\bf H(r)} - {\bf H(h)}].]$ First consider the lattice sum over the direct-lattice points X(d), relative to a particular point $[{\bf X(x)} = {\bf R}]$ , with omission of the origin lattice point. $[S'(n, {\bf R}) = {\textstyle\sum\limits_{{\bf d} \neq 0}} |{\bf X(d)} - {\bf R}|^{-n}.]$ The special case with $[{\bf R} = 0]$ will also be needed: $[S'(n, 0) = {\textstyle\sum\limits_{{\bf d} \neq 0}} |{\bf X(d)}|^{-n}.]$ Now define a sum of Dirac delta functions $[f'[{\bf X(d)}] = {\textstyle\sum\limits_{d \neq 0}} \delta [{\bf X(x)} - {\bf X(d)}].]$ Then S′ can be represented as an integral $[S'(n, {\bf R}) = \textstyle\int f'[{\bf X(d)}] |{\bf X - R}|^{-n}\; \hbox{d}{\bf X},]$ in which a term is contributed to S′ whenever the direct-space vector X coincides with the lattice vector X(d), except for $[{\bf d} = 0]$ . Now apply the convergence function to S′: $[\eqalign{ S'(n, {\bf R}) &= [\Gamma (n/2)]^{-1} \textstyle\int f'[{\bf X(d)}]|{\bf X - R}|^{-n}\cr &\quad \times \Gamma (n/2, \pi w^{2} |{\bf X - R}|^{2})\; \hbox{d}{\bf X}\cr &\quad + [\Gamma (n/2)]^{-1} \textstyle\int f'[{\bf X (d)}]|{\bf X - R}|^{-n}\cr &\quad \times \gamma (n/2, \pi w^{2} |{\bf X - R}|^{2})\; \hbox{d}{\bf X}.}]$

The first integral is shown here only for the purpose of giving a consistent representation of S′; in fact, the first integral will be reconverted back into a sum and evaluated in direct space. The second integral will be transformed to reciprocal space using Parseval's theorem [see, for example, Arfken (1970 )], which states that $[\textstyle\int f ({\bf X}) g^{*} ({\bf X})\; \hbox{d}{\bf X} = \textstyle\int FT_{3} [\; f({\bf X})] FT_{3} [g^{*} ({\bf X})]\; \hbox{d}{\bf H}.]$ Considering only the second integral in the formula for S′ and explicitly introducing the $[{\bf d} = 0]$ term we have $[\eqalign{ &[\Gamma (n/2)]^{-1} \textstyle\int f[{\bf X(d)}]| {\bf X (d) - R}|^{-n} \gamma (n/2, \pi w^{2}|{\bf X - R}|^{2})\; \hbox{d}{\bf X}\cr &\quad - [\Gamma (n/2)]^{-1} \textstyle\int \delta ({\bf X})|{\bf R}|^{-n} \gamma (n/2, \pi w^{2} |{\bf R}|^{2})\; \hbox{d}{\bf X},}]$ where the unprimed f includes the $[{\bf h} = 0]$ term which was earlier omitted from f′: $[f({\bf X}) = {\textstyle\sum\limits_{\bf d}} \delta [{\bf X (x) - X(d)}].]$ Using Parseval's theorem, and evaluating the origin term, we have $[\eqalign{ &[\Gamma (n/2)]^{-1} \textstyle\int FT_{3} \{f [{\bf X(d)}]\} FT_{3}\cr &\quad \times [|{\bf X (d) - R}|^{-n} \gamma (n/2, \pi w^{2}|{\bf X - R}|^{2})]\; \hbox{d}{\bf H}\cr &\quad - [\Gamma (n/2)]^{-1} |{\bf R}|^{-n} \gamma (n/2, \pi w^{2} |{\bf R}|^{2}).}]$ The Fourier transform of the complement of the incomplete gamma function divided by $[|{\bf X}|^{n}]$ is (Nijboer & DeWette, 1957 ) $[\eqalign{ &FT_{3} [\gamma (n/2, \pi w^{2} |{\bf X}|^{2}) |{\bf X}|^{-n}]\cr &\quad = \pi^{n - (3/2)} |{\bf H}|^{n-3} \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}].}]$ If there is a change of origin and the point $[({\bf X} - {\bf R})]$ is used instead of X the transform is $[\eqalign{ &FT_{3} [\gamma (n/2, \pi w^{2} |{\bf X - R}|^{2}) |{\bf X - R}|^{-n}]\cr &\quad = \pi^{n-(3/2)} |{\bf H}|^{n-3} \Gamma [(-n/2)\cr &\qquad + (3/2), \pi w^{-2} |{\bf H}|^{2}] \exp (2 \pi i{\bf H \cdot R}).}]$ Evaluation of the two Fourier transforms in the first term gives $[\eqalign{ &[\Gamma (n/2)]^{-1} \textstyle\int V_{d}^{-1} {\textstyle\sum\limits_{\bf h}} \delta [{\bf H (h) - H}] \pi^{n-(3/2)} |{\bf H}|^{n-3}\cr &\quad\times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}] \exp (2 \pi i{\bf H \cdot R})\; \hbox{d}{\bf H}.}]$ Because of the presence of the Dirac delta function in each integral, we can convert the integrals with h unequal to zero into a sum $[\eqalign{ &[\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{h \neq 0}} |{\bf H (h)}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H (h)}|^{2}] \exp [2 \pi i{\bf H(h) \cdot R}].}]$ The $[{\bf h} = 0]$ term needs to be evaluated in the limit. Clearly, the complex exponential goes to unity. If n is greater than 3 the limit of the indeterminate form infinity/infinity is needed: $[\eqalign{ &\lim\limits_{|{\bf H}| \rightarrow 0} {\Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}]\over |{\bf H}|^{3-n}}\cr &\quad = \lim\limits_{|{\bf H}| \rightarrow 0} {\int_{\pi w^{-2} |{\bf H}|^{2}}^{\infty} t^{(-n/2)+(1/2)} \exp (-t)\; \hbox{d}t\over |{\bf H}|^{3-n}}.}]$ The limit can be found by L'Hospital's rule [see, for example, Widder (1961 )] which states that if [f(x)] and [g(x)] both approach infinity as x approaches a constant, c, and the limit of the ratio of the first derivatives [f'(x)] and [g'(x)] exists, that limit is also true for the limit of the ratio of the functions: $[\lim\limits_{x \rightarrow c} \ {f(x)\over g(x)} = \lim\limits_{x \rightarrow c} \ {f'(x)\over g'(x)}.]$ To differentiate the definite integral function, Leibnitz's formula may be used [see, for example, Arfken (1970 )]. This formula states that $[\eqalign{ {\hbox{d}\over\; \hbox{d}x} \int\limits_{g(x)}^{h(x)} f(t, x)\; \hbox{d}t &\ = \int\limits_{g(x)}^{h(x)} {\hbox{d}f (t, x)\over\; \hbox{d}x}\; \hbox{d}t\cr &\quad + f[h (x)] {\hbox{d}h(x)\over \hbox{d}x} - f[g(x)] {\hbox{d}g(x)\over\; \hbox{d}x}.}]$ In our case, x becomes $[|{\bf H}|]$ ; f becomes $[t^{(-n/2)+(1/2)} \exp (-t)]$ which is independent of $[|{\bf H}|]$ ; g becomes $[\pi w^{-2} |{\bf H}|^{2}]$ ; and h is infinite. Thus only the last term of Leibnitz's formula is nonzero and we obtain for the ratio of the first derivatives $[\eqalign{ &\lim\limits_{|{\bf H}| \rightarrow 0} {-(\pi w^{-2} |{\bf H}|^{2})^{(-n/2)+(1/2)} \exp (-\pi w^{2} |{\bf H}|^{2}) 2 \pi w^{-2} |{\bf H}|\over (3 - n) |{\bf H}|^{2-n}}\cr &\quad = \pi^{(-n/2)+(3/2)} w^{n-3} [2/(n-3)],}]$ so that the limiting value for the $[{\bf h} = 0]$ term for n greater than 3 is $[+ [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n/2} w^{n-3} [2/(n-3)].]$ The final result for S′ is $[\eqalign{ S'(n, {\bf R}) &= [\Gamma (n/2)]^{-1} {\textstyle\sum\limits_{{\bf d} \neq 0}} |{\bf X(d) - R}|^{-n}\cr &\quad \times \Gamma (n/2, \pi w^{2} |{\bf X(d) - R}|^{2})\cr &\quad - [\Gamma (n/2)]^{-1} |{\bf R}|^{-n} \gamma (n/2, \pi w^{2} |{\bf R}|^{2})\cr &\quad + [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n - (3/2)} {\textstyle\sum\limits_{{\bf h} \neq 0}} |{\bf H (h)}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}]\cr &\quad \times \exp [2 \pi i{\bf H(h) \cdot R}]\cr &\quad + [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n/2} w^{n-3} 2(n - 3)^{-1}.}]$

The significance of the terms is as follows. The first term represents the convergence-accelerated direct sum, which does not include the origin term; the next term, also in direct space, corrects for the remainder resulting from the subtraction of the origin term; the third term comes from Parseval's theorem and is a sum over the nonzero h reciprocal-lattice points; and the last term is the reciprocal-lattice $[{\bf h} = 0]$ term.

If $[{\bf R} = 0]$ the second term becomes an indeterminate form 0/0. The limit can be found with use of L'Hospital's rule again, this time for the 0/0 form. We need the limit of [f(x)/g(x)] , where $[f(R) = \gamma (n/2, \pi w^{2} R^{2})]$ and $[g(R) = R^{n}]$ . To differentiate the incomplete gamma function, we can again use Leibnitz's formula. In this case only the second term of the formula is nonzero and we obtain for the ratio of the first derivatives $[{2 \pi^{n/2} w^{n} |{\bf R}|^{n-1} \exp (- \pi w^{2}|{\bf R}|^{2})\over n |{\bf R}|^{n-1}},]$ so that the limiting value for this term as $[|{\bf R}|]$ approaches zero is $[-[\Gamma (n/2)]^{-1} 2 \pi^{n/2} w^{n} n^{-1}.]$ Therefore, the value of the sum when $[{\bf R} = 0]$ is $[\eqalign{ S' (n,0) &= [\Gamma (n/2)]^{-1} {\textstyle\sum\limits_{{\bf d} \neq 0}} |{\bf X(d)}|^{-n} \Gamma (n/2, \pi w^{2} |{\bf X(d)}|^{2})\cr &\quad - [\Gamma (n/2)]^{-1} 2 \pi^{n/2} w^{n} n^{-1}\cr &\quad + [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{{\bf h} \neq 0}} |{\bf H(h)}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}]\cr &\quad + [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n/2} w^{n-3} 2(n - 3)^{-1}.}]$

3.4.5. Extension of the method to a composite lattice

| top | pdf |

Define a general lattice sum over direct-space points $[{\bf R}_{j}]$ which interact with pairwise coefficients $[Q_{jk}]$ , where $[Q_{jk} = Q_{kj}]$ : $[V (n, {\bf R}_{j}) = (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} | {\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{-n},]$ where the prime indicates that when $[{\bf d} = 0]$ the self-terms with [j = k] are omitted. For convenience the terms may be divided into three groups: the first group of terms has $[{\bf d} = 0]$ , where j is unequal to k; the second group has d not zero and j not equal to k; and the third group had d not zero and [j = k] . (A possible fourth group with $[{\bf d} = 0]$ and [j = k] is omitted, as defined.) $[\eqalign{ V (n, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} |{\bf R}_{k} - {\bf R}_{j}|^{-n}\cr &\quad + (1/2) {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} S' (n, |{\bf R}_{j} - {\bf R}_{k}|) + (1/2) {\textstyle\sum\limits_{j}}\ Q_{jj} S' (n, 0).}]$ By expanding this expression we obtain $[\eqalignno{ V (n, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} |{\bf R}_{k} - {\bf R}_{j}|^{-n} &(1)\cr &\quad + \left\{[1/2 \Gamma (n/2)] {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} {\textstyle\sum\limits_{{\bf d} \neq 0}} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{-n}\right.\cr &\quad \left.\times\ \Gamma (n/2, \pi w^{2} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{2})\vphantom{\sum_{d}}\right\} &(2)\cr &\quad - \left\{[1/2 \Gamma (n/2)] {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} |{\bf R}_{k} - {\bf R}_{j}|^{-n}\right.\cr &\quad \left. \times\ \gamma (n/2, \pi w^{2} |{\bf R}_{k} - {\bf R}_{j}|^{2})\vphantom{\sum_{d}}\right\} &(3)\cr &\quad + \left\{[1/2 \Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} {\textstyle\sum\limits_{{\bf h} \neq 0}} |{\bf H(h)}|^{n-3}\right.\cr &\quad \times\ \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}]\cr &\quad \left.\times \exp [2 \pi i {\bf H(h)} \cdot ({\bf R}_{k} - {\bf R}_{j})]\vphantom{\sum_{d}}\right\} &(4)\cr&\quad + [1/2 \Gamma (n/2)] V_{d}^{-1} \pi^{n/2} w^{n-3} 2(n - 3)^{-1} {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} &(5)\cr &\quad + \left\{[1/2 \Gamma (n/2)] {\textstyle\sum\limits_{j}}\ Q_{jj} {\textstyle\sum\limits_{{\bf d} \neq 0}} |{\bf X(d)}|^{-n}\right.\cr &\quad \left.\times\ \Gamma (n/2, \pi w^{2} |{\bf X(d)}|^{2})\vphantom{\sum_{d}}\right\} &(6)\cr &\quad - [1 / \Gamma (n/2)] \pi^{n/2} w^{n} n^{-1} {\textstyle\sum\limits_{j}}\ Q_{jj} &(7)\cr &\quad + \left\{[1/2 \Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{j}}\ Q_{jj} {\textstyle\sum\limits_{{\bf h} \neq 0}} |{\bf H(h)}|^{n-3}\right.\cr &\quad \left.\times\ \Gamma [(- n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}]\vphantom{\sum_{d}}\right\} &(8)\cr &\quad + [1/2 \Gamma (n/2)] V_{d}^{-1} \pi^{n/2} w^{n-3} 2(n - 3)^{-1} {\textstyle\sum\limits_{j}}\ Q_{jj}. &(9)}]$ This expression for V has nine terms, which are numbered on the right-hand side. Term (3) can be expressed in terms of Γ rather than γ: $[\eqalign{ (3) &= -(1/2) {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} |{\bf R}_{k} - {\bf R}_{j}|^{-n}\cr &\quad + [1 / \Gamma (n/2)] {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} |{\bf R}_{k} - {\bf R}_{j}|^{-n} \Gamma (n/2, \pi w^{2}| {\bf R}_{k} - {\bf R}_{j}|^{2}).}]$ It is seen that cancellation occurs with term (1) so that $[\eqalign{ (1) + (3) &= [1 / \Gamma (n/2)] {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} |{\bf R}_{k} - {\bf R}_{j}|^{-n}\cr &\quad \times \Gamma (n/2, \pi w^{2} |{\bf R}_{k} - {\bf R}_{j}|^{2}),}]$ which is the $[{\bf d} = 0]$ , j unequal to k portion of the treated direct-lattice sum. The d unequal to 0, j unequal to k portion corresponds to term (2) and the d unequal to 0, [j = k] portion corresponds to term (6). The direct-lattice terms may be consolidated as $[\eqalign{ (1)\! +\! (2)\! +\! (3)\! +\! (6) &= [1/2 \Gamma (n/2)] {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} |{\bf R}_{k} + {\bf X(d)}- {\bf R}_{j}|^{-n}\cr &\quad \times \Gamma [n/2, \pi w^{2} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{2}].}]$ Now let us combine terms (4) and (8), carrying out the h summation first: $[\eqalign{ (4) + (8) &= [1/2 \Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{\bf h}} |{\bf H(h)}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}]\cr &\quad \times {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}\ Q_{jk} \exp [2 \pi i {\bf H(h)} \cdot ({\bf R}_{k} - {\bf R}_{j})].}]$ Terms (5) and (9) may be combined: $[(5) + (9) = [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n/2} w^{n-3} (n-3)^{-1} \left({\textstyle\sum\limits_{j}}\ Q_{ij} + {\textstyle\sum\limits_{j \neq k}}\ Q_{jk}\right).]$ The final formula is shown below. The significance of the four terms is: (1) the treated direct-lattice sum; (2) a correction for the difference resulting from the removal of the origin term in direct space; (3) the reciprocal-lattice sum, except $[{\bf h} = 0]$ ; and (4) the $[{{\bf h} = 0}]$ term of the reciprocal-lattice sum. $[\eqalignno{ &V (n, {\bf R}_{j})\cr &\quad = [1/2 \Gamma (n/2)] {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{-n}\cr &\qquad \times \Gamma (n/2, \pi w^{2} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{2}) &(1)\cr &\qquad - [1 / \Gamma (n/2)] \pi^{n/2} w^{n} n^{-1} {\textstyle\sum\limits_{j}}\ Q_{jj} &(2)\cr &\qquad + [1/2 \Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{\bf h}} |{\bf H(h)}|^{n-3}\cr &\qquad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}]\cr &\qquad \times {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}\ Q_{jk} \exp [2 \pi i {\bf H(h)} \cdot ({\bf R}_{k} - {\bf R}_{j})] &(3)\cr &\qquad + [\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n/2} w^{n-3} (n - 3)^{-1} \left({\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}\ Q_{jk}\right). &(4)}]$

3.4.6. The case of $[{\bf{\bi n} = 1}]$ (Coulombic lattice energy )

| top | pdf |

As taken above, the limit of the reciprocal-lattice $[{\bf h} = 0]$ term of $[S'(n, {\bf R})]$ or [S'(n, 0)] existed only if n was greater than 3. The corresponding contributions to $[V(n, {\bf R}_{j})]$ were terms (5) and (9) of Section 3.4.5 . To extend the method to [n = 1] we will show in this section that these $[{\bf h} = 0]$ terms vanish if conditions of unit-cell neutrality and zero dipole moment are satisfied.

The integral representation of the term (5) is $[\eqalign{ &[1/2\Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{j \neq k}}\ Q_{jk} \textstyle\int \delta (0) |{\bf H}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}]\cr &\quad \times \exp [2 \pi i {\bf H} \cdot ({\bf R}_{k} - {\bf R}_{j})]\; \hbox{d}{\bf H}}]$ and for term (9) is $[\eqalign{ &[1/2\Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} {\textstyle\sum\limits_{j}}\ Q_{jj} \textstyle\int \delta (0) |{\bf H}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}]\; \hbox{d}{\bf H}.}]$ Combining these two sums of integrals into one integral sum gives $[\eqalign{ &[1/2\Gamma (n/2)] V_{d}^{-1} \pi^{n-(3/2)} \textstyle\int \delta (0) |{\bf H}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}] {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}\ Q_{jk}\cr &\quad \times \exp [2\pi i {\bf H}\cdot ({\bf R}_{k} - {\bf R}_{j})] \hbox{ d}{\bf H}.}]$

For [n = 1] , suppose $[q_{j}]$ are net atomic charges so that the geometric combining law holds for $[Q_{jk} = q_{j}q_{k}]$ . Then the double sum over j and k can be factored so that the limit that needs to be considered is $[\lim\limits_{|{\bf H}| \rightarrow 0} {{\left[{\textstyle\sum_{k{\phantom j}}} q_{k} \exp (2\pi i {\bf H} \cdot {\bf R}_{k})\right] \left[{\textstyle\sum_{j}}\ q_{j} \exp (-2\pi i {\bf H} \cdot {\bf R}_{j})\right]} \over |{\bf H}|^{2}}.]$ If the unit cell does not have a net charge, the sum over the q's goes to zero in the limit and this is a 0/0 indeterminate form. Let $[|{\bf H}|]$ approach zero along the polar axis so that $[{\bf H}\cdot {\bf R}_{k} = H_{3} R_{3k}]$ , where subscript 3 indicates components along the polar axis. To find the limit with L'Hospital's rule the numerator and denominator are differentiated twice with respect to $[H_{3}]$ . Represent the numerator of the limit by the product [(uv)] and note that $[{\hbox{d}^{2} (uv)\over\; \hbox{d}x^{2}} = u {\hbox{d}^{2} v\over\; \hbox{d}x^{2}} + v {\hbox{d}^{2} u\over\; \hbox{d}x^{2}} + 2 {\hbox{d}u\over\; \hbox{d}x} {\hbox{d}v\over\; \hbox{d}x}.]$ It is seen that in addition to cell neutrality the product of the first derivatives of the sums must exist. These sums are $[\left[2\pi i {\textstyle\sum\limits_{k}}\ q_{k} R_{3k} \exp (2\pi i H_{3} R_{3k})\right]]$ and $[\left[-2\pi i {\textstyle\sum\limits_{j}}\ q_{j} R_{3j} \exp (-2\pi i H_{3} R_{3j})\right],]$ which vanish if the unit cell has no dipole moment in the polar direction, that is, if $[\sum q_{j} R_{3j} = 0]$ . Since the second derivative of the denominator is a constant, the desired limit is zero under the specified conditions. Now the polar direction can be chosen arbitrarily, so the unit cell must not have a dipole moment in any direction for the limit of the numerator to be zero. Thus we have the formula for the Coulombic lattice sum $[\eqalign{ V(1, {\bf R}_{j}) &= [ 1/2\Gamma (1/2)] {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{-1}\cr &\quad \times \Gamma (1/2, \pi w^{2} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{2})\cr &\quad + [1/2 \Gamma (1/2) V_{d}^{-1} \pi^{-1/2} {\textstyle\sum\limits_{\bf h}} |{\bf H(h)}|^{-2}\cr &\quad \times \Gamma (1/2, \pi w^{-2} |{\bf H(h)}|^{2}) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}\ Q_{jk}\cr &\quad \times \exp [2\pi i {\bf H(h)} \cdot ({\bf R}_{k} - {\bf R}_{j})]\cr &\quad - [1/\Gamma (1/2)] \pi^{1/2} w {\textstyle\sum\limits_{j}}\ q_{j}^{2},}]$ which holds on conditions that the unit cell be electrically neutral and have no dipole moment. If the unit cell has a dipole moment, the limiting value discussed above depends on the direction of H. For methods of obtaining the Coulombic lattice sum where the unit cell does have a dipole moment, the reader is referred to the literature (DeWette & Schacher, 1964 ; Cummins et al., 1976 ; Bertaut, 1978 ; Massidda, 1978 ).

3.4.7. The cases of $[{\bf{\bi n} = 2}]$ and $[{\bf{\bi n} = 3}]$

| top | pdf |

If [n = 2] the denominator considered for the limit in the preceding section is linear in |H| so that only one differentiation is needed to obtain the limit by L'Hospital's method. Since a term of the type $[\sum q_{j} \exp (2\pi i{\bf H}\cdot {\bf R}_{j})]$ is always a factor, the requirement that the unit cell have no dipole moment can be relaxed. For [n = 2] the zero-charge condition is still required: $[\sum q_{j} = 0]$ . When [n = 3] the expression becomes determinate and no differentiation is required to obtain a limit. In addition, factoring the $[Q_{jk}]$ sums into $[q_{j}]$ sums is not necessary so that the only remaining requirement for this term to be zero is $[\sum \sum Q_{jk} = 0]$ , which is a further relaxation beyond the requirement of cell neutrality.

3.4.8. Derivation of the accelerated convergence formula via the Patterson function

| top | pdf |

The structure factor with generalized coefficients $[q_{j}]$ is defined by $[F[{\bf H(h)}] = {\textstyle\sum\limits_{j}}\ q_{j} \exp [2\pi i {\bf H(h)} \cdot {\bf R}_{j}].]$ The corresponding Patterson function is defined by $[P({\bf X}) = V_{d}^{-1} {\textstyle\sum\limits_{\bf h}} |F[{\bf H(h)}]|^{2} \exp [2\pi i {\bf H(h)}\cdot {\bf X}].]$ The physical interpretation of the Patterson function is that it is nonzero only at the intersite vector points $[{\bf R}_{k} + {\bf X(h)} - {\bf R}_{j}]$ . If the origin point is removed, the lattice sum may be expressed as an integral over the Patterson function. This origin point in the Patterson function corresponds to intersite vectors with [j = k] and $[{\bf H(h)} = 0]$ : $[S_{n} = (1/2 V_{d}) \textstyle\int |{\bf X}|^{-n} [P({\bf X}) - P({\bf X}) \delta ({\bf X})] \hbox{ d}{\bf X}.]$ Using the incomplete gamma function as a convergence function, this formula expands into two integrals $[\eqalign{ S_{n} &= [1/2 V_{d} \Gamma (n/2)] \textstyle\int |{\bf X}|^{-n}\cr &\quad \times [P ({\bf X}) - P({\bf X}) \delta ({\bf X})] \Gamma (n/2, \pi w^{2} |{\bf X}|^{2})\; \hbox{d}{\bf X}\cr &\quad + [1/2 V_{d} \Gamma (n/2)] \textstyle\int |{\bf X}|^{-n}\cr &\quad \times [P({\bf X}) - P({\bf X}) \delta ({\bf X})] \gamma (n/2, \pi w^{2} |{\bf X}|^{2})\; \hbox{d}{\bf X}.}]$ The first integral is shown only for a consistent representation; actually it will be reconverted to a sum and evaluated in direct space. The first part of the second integral will be evaluated with Parseval's theorem and the second part in the limit as $[|{\bf X}|]$ approaches zero: $[\eqalign{ &[1/2 V_{d} \Gamma (n/2)] \textstyle\int FT_{3} [P({\bf X})]\cr &\quad \times FT_{3} [|{\bf X}|^{-n} \gamma (n/2, \pi w^{2} |{\bf X}|^{2})] \hbox{ d}{\bf H}\cr &\quad -\lim\limits_{{\bf X} \rightarrow 0} [1/2 V_{d} \Gamma (n/2)] [P(0) |{\bf X}|^{-n} \gamma (n/2, \pi w^{2} |{\bf X}|^{2})].}]$ The first Fourier transform (of the Patterson function) is the set of amplitudes of the structure factors and the second Fourier transform has already been discussed above; the method for obtaining the limit (for n equal to or greater than 1) was also discussed above. The result obtained is $[\eqalign{ &[1/2 V_{d} \Gamma (n/2)] \pi^{n-(3/2)} \textstyle\int |F [{\bf H (h)}]|^{2} |{\bf H}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H}|^{2}] \hbox{ d}{\bf H}\cr &\quad -[1/2 V_{d} \Gamma (n/2)] |F(0)|^{2} 2 \pi^{n/2}w^{n} n^{-1}.}]$ The integral can be converted into a sum, since $[|F[{\bf H(h)}]|]$ is nonzero only at the reciprocal-lattice points: $[\eqalign{ &[1/2 V_{d} \Gamma (n/2)] \pi^{n-(3/2)} {\textstyle\sum\limits_{\bf h}} |F[{\bf H(h)}]|^{2} |{\bf H(h)}|^{n-3}\cr &\quad \times \Gamma [(-n/2) + (3/2), \pi w^{-2} |{\bf H(h)}|^{2}].}]$ The term with $[{\bf H(h)} = 0]$ is evaluated in the limit, for n greater than 3, as $[[\Gamma (n/2)]^{-1} V_{d}^{-1} \pi^{n/2} w^{n-3} (n-3)^{-1} |F(0)|^{2}.]$ Since $[|F(0)|^{2} = \sum\sum q_{j} q_{k}]$ , this term is identical with the third term of $[V(n, {\bf R}_{j})]$ as derived earlier. The case of [n = 1] is handled in the same way as previously discussed, where the limit of this term is zero provided the unit cell has no net charge or dipole moment.

3.4.9. Evaluation of the incomplete gamma function

| top | pdf |

The incomplete gamma function may be expressed in terms of commonly available functions such as the exponential integral and the complement of the error function. The definition of the exponential integral is $[E_{1} (x^{2}) = \textstyle\int\limits_{x^{2}}^{\infty} t^{-1} \exp (-t)\; \hbox{d}t = \Gamma (0, x^{2}).]$ The definition of the complement of the error function is $[\hbox{erfc} (x) = \textstyle\int\limits_{x}^{\infty} \exp (-t^{2})\; \hbox{d}t = \pi^{-1/2} \Gamma (1/2, x^{2}).]$ Numerical approximations to these functions are given, for example, by Hastings (1955 ). The recursion formula for the incomplete gamma function (Davis, 1972 ) $[\Gamma (n+1, x^{2}) = n \Gamma (n, x^{2}) + x^{2n} \exp (-x^{2})]$ may be used to obtain working formulae starting from the special values of $[\Gamma (0, x^{2})]$ and $[\Gamma (1/2, x^{2})]$ which are defined above. Also we note that $[\Gamma (1, x^{2}) = \exp (-x^{2})]$ .

3.4.10. Summation over the asymmetric unit and elimination of intramolecular energy terms

| top | pdf |

Let us consider the case where the unit cell contains Z molecules which are related by Z symmetry operations, and it is desired to include only intermolecular distances in the summation. In the direct sum (1) the indices j and k will then run only over the asymmetric unit, and all terms with $[{\bf d} = 0]$ are eliminated. The calculated energy refers then to one molecule (or mole) rather than to one unit cell. The correction term (2) also refers to one molecule according to the range of j and k. Since the reciprocal-lattice sum refers to the entire unit cell, terms (3) and (4) need to be divided by Z to refer the energy to one molecule.

Both the direct and reciprocal sums must be corrected for the elimination of intramolecular terms. Using the convergence function [W(R)] , we have $[\eqalign{ V(n, {\bf R}_{j}) &= {\textstyle\sum\limits_{\rm inter}} |{\bf R}|^{-n} W(R) + {\textstyle\sum\limits_{\rm intra}} |{\bf R}|^{-n} W(R)\cr &\quad + {\textstyle\sum\limits_{\rm inter}} |{\bf R}|^{-n} [1 - W(R)] + {\textstyle\sum\limits_{\rm intra}} |{\bf R}|^{-n} [1 - W(R)].}]$ As mentioned above, the second summation term, which is the intramolecular term in direct space, is simply left out of the calculation. When using the accelerated-convergence method the third and fourth summation terms are always obtained, evaluated in reciprocal space. The undesired inclusion of the intramolecular term (fourth term above) in the reciprocal-space sum may be compensated for by explicit subtraction of this term from the sum.

3.4.11. Reference formulae for particular values of n

| top | pdf |

In this section let $[a^{2} = \pi w^{2} |{\bf R}_{k} + {\bf X(d)} - {\bf R}_{j}|^{2}]$ and $[b^{2} = \pi w^{-2} |{\bf H(h)}|^{2}]$ . Let $[T_{0} = \sum Q_{jj} = \sum q_{j}^{2}]$ ; $[T_{1} = \sum \sum Q_{jk} =]$ $[T_{0} + 2 \sum \sum_{j\;\gt\; k} Q_{jk}]$ . If the geometric mean combining law holds, $[T_{1} = ({\textstyle\sum\limits_{j}} q_{j})^{2}]$ ; let $[\eqalign{ T_{2} (h) &= {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}\ Q_{jk} \exp [2 \pi i {\bf H(h)} \cdot ({\bf R}_{k} - {\bf R}_{j})]\cr &= T_{0} + 2 {\lower3pt\hbox{${\sum\sum\atop \scriptstyle j\;\gt\; k}$}}\ Q_{jk} \cos [2 \pi {\bf H(h)} \cdot ({\bf R}_{k} - {\bf R}_{j})].}]$ Then $[T_{2} ({\bf h}) = |F {\bf (h)}|^{2} = |{\textstyle\sum\limits_{j}}\ q_{j} \exp [2 \pi i{\bf H(h)} \cdot {\bf R}_{j}]|^{2} = A {\bf (h)}^{2} + B {\bf (h)}^{2},]$ where $[A {\bf (h)} = {\textstyle\sum\limits_{j}}\ q_{j} \cos [2 \pi {\bf H(h)} \cdot {\bf R}_{j}]]$ and $[B {\bf (h)} = {\textstyle\sum\limits_{j}}\ q_{j} \sin [2 \pi {\bf H(h)} \cdot {\bf R}_{j}].]$ The formulae below describe $[V (n, {\bf R}_{j})]$ in terms of $[T_{0},\ T_{1}]$ and $[T_{2}]$ ; the distance $[|{\bf R}_{k} + {\bf X (d)} - {\bf R}_{j}|]$ is simply represented by $[R_{jkd}]$ . $[\eqalign{ V (1, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ q_{j}\ q_{k} {\textstyle\sum\limits_{\bf d}} {R}_{jkd}^{-1} \ \hbox{erfc} (a)\cr &\quad + (1/2 \pi V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2} {\bf (h)} |{\bf H(h)}|^{-2} \exp (-b^{2}) - w T_{0}\cr V (2, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} {R}_{jkd}^{-2} \exp (-a^{2})\cr &\quad + (\pi / 2 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2} {\bf (h)}| {\bf H(h)}|^{-1} \hbox{erfc} (b) - (\pi / 2) w^{2} T_{0}\cr V (3, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} {R}_{jkd}^{-3} [\hbox{erfc} (a) + 2 a \pi^{-1/2} \exp (-a^{2})]\cr &\quad + (\pi / V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2} {\bf (h)} E_{1} (b^{2}) - (2 \pi/3) w^{3} T_{0}\cr\noalign{\vskip6pt} V (4, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} {R}_{jkd}^{-4} (1 + a^{2}) \exp (-a^{2})\cr &\quad + (\pi^{5/2} / V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2} {\bf (h)} |{\bf H(h)}|\cr &\quad \times [-\pi^{1/2} \hbox{erfc} (b) + b^{-1} \exp (-b^{2})]\cr &\quad + (\pi^{2} / V_{d}) w T_{1} - (\pi^{2}/4) w^{4} T_{0}\cr\noalign{\vskip6pt} {\bi V}(5, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} R_{jkd}^{-5}\cr &\quad \times [\hbox{erfc} (a) + 2\pi^{-1/2} a(1 + 2a^{2}/3)\exp (-a^{2})]\cr &\quad + (2\pi^{3}/3 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2}({\bf h})|{\bf H(h)}|^{2}\cr &\quad \times [b^{-2}\exp (-b^{2}) - E_{1}(b^{2})]\cr &\quad + (2\pi^{2}/3 V_{d}) w^{2}T_{1} - (4\pi^{2}/15) w^{5}T_{0}\cr\noalign{\vskip6pt} {\bi V}(6, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} R_{jkd}^{-6} [1 + a^{2} + (a^{4}/2)]\exp (-a^{2})\cr &\quad + (\pi^{9/2}/3 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2}({\bf h})|{\bf H(h)}|^{3}\cr &\quad \times [\pi^{1/2}\hbox{erfc}(b) + [(1/2b^{3}) - (1/b)]\exp (-b^{2})]\cr &\quad + (\pi^{3}/6 V_{d}) w^{3}T_{1} - (\pi^{3}/12) w^{6}T_{0}\cr}]$ $[\eqalign{ {\bi V}(7, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} R_{jkd}^{-7}\cr &\quad \times [\hbox{erfc}(a) + 2\pi^{-1/2} a[1 + (2/3)a^{2} + (4/15)a^{4}]\cr &\quad \times \exp (-a^{2})]\cr &\quad + (2\pi^{5}/15 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2}({\bf h})|{\bf H(h)}|^{4}\cr &\quad \times [(-b^{-2} + b^{-4})\exp (-b^{2}) + E_{1}(b^{2})]\cr &\quad + (2\pi^{3}/15 V_{d}) w^{4}T_{1} - (8\pi^{3}/105) w^{7}T_{0}\cr\noalign{\vskip6pt} {\bi V}(8, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} R_{jkd}^{-8}\cr &\quad \times [1 + a^{2} + (a^{4}/2) + (a^{6}/6)]\exp (-a^{2})\cr &\quad + (2\pi^{13/2}/45 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2}({\bf h})|{\bf H(h)}|^{5}\cr &\quad \times [-\pi^{1/2}\hbox{erfc}(b) + [(1/b) - (1/2b^{3}) + (3/4b^{5})]\cr &\quad \times \exp (-b^{2})]\cr &\quad + (\pi^{4}/30 V_{d}) w^{5}T_{1} - (\pi^{4}/48) w^{8}T_{0}\cr\noalign{\vskip6pt}{\bi V}(9, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} R_{jkd}^{-9}\cr &\quad \times [\hbox{erfc}(a) + 2\pi^{-1/2} a[1 + (2/3)a^{2} + (4/15)a^{4}\cr &\quad + (8/105)a^{6}]\exp (-a^{2})]\cr &\quad + (4\pi^{7}/315 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2}({\bf h})|{\bf H(h)}|^{6}\cr &\quad \times [(b^{-2} - b^{-4} + 2b^{-6})\exp (-b^{2}) - E_{1}(b^{2})]\cr &\quad + (8\pi^{4}/315 V_{d}) w^{6}T_{1} - (16\pi^{4}/945) w^{9}T_{0}\cr {\bi V}(10, {\bf R}_{j}) &= (1/2) {\textstyle\sum\limits_{j}} {\textstyle\sum\limits_{k}}'\ Q_{jk} {\textstyle\sum\limits_{\bf d}} R_{jkd}^{-10}\cr &\quad \times [1 + a^{2} + (a^{4}/2) + (a^{6}/6) + (a^{8}/24)]\exp (-a^{2})\cr &\quad + (\pi^{17/2}/315 V_{d}) {\textstyle\sum\limits_{{\bf h} \neq 0}} T_{2}({\bf h})|{\bf H(h)}|^{7}\cr &\quad \times [\pi^{1/2}\hbox{erfc}(b) + [(-1/b) + (1/2b^{3}) - (3/4b^{5})\cr &\quad + (15/8b^{7})]\exp (-b^{2})]\cr &\quad + (\pi^{5}/168 V_{d}) w^{7}T_{1} - (\pi^{5}/240) w^{10}T_{0}.}]$

3.4.12. Numerical illustrations

| top | pdf |

Consider the case of the sodium chloride crystal structure (a face-centred cubic structure) as a simple example for evaluation of the Coulombic sum. The sodium ion can be taken at the origin, and the chloride ion halfway along an edge of the unit cell. The results can easily be generalized for this structure type by using the unit-cell edge length, a, as a scaling constant.

First, consider the nearest neighbours. Each sodium and each chloride ion is surrounded by six ions of opposite sign at a distance of [a/2] . The Coulombic energy for the first coordination sphere is $[-(1/2)(12)(2/a)(1389.3654) = -16672.385/a\;\hbox{kJ mol}^{-1}]$ . Table 3.4.2.1 shows that the converged value of the lattice energy is [-4855.979/a] . Thus the nearest-neighbour energy is over three times more negative than the total lattice energy. In the second coordination sphere each ion is surrounded by 12 similar ions at a distance of $[a/2^{1/2}]$ . The energy contribution of the second sphere is $[+(1/2)(24)(2^{1/2}/a)(1389.3654) = +23578.313/a]$ . Thus, major cancellation occurs and the net energy for the first two coordination spheres is [+6905.928/a] which actually has the wrong sign for a stable crystal. The third coordination sphere again makes a negative contribution. Each ion is surrounded by eight ions of opposite sign at a distance of $[a/3^{1/2}]$ . The energy contribution is $[-(1/2)(16)(3^{1/2}/a)(1389.3654) = -19251.612/a]$ , now giving a total so far of [-12345.684/a] . In the fourth coordination sphere each ion is surrounded by six others of the same sign at a distance of a. The energy contribution is [+(1/2)(12)(1/a)(1389.3654) = +8336.19/a] to yield a total of [-4009.491/a] .

It is seen immediately by examining the numbers that the Coulombic sum is converging very slowly in direct space. Madelung (1918 ) devised a method for accurate evaluation of the sodium chloride lattice sum. However, his method is not generally applicable for more complex lattice structures. Evjen (1932 ) emphasized the importance of summing over a neutral domain, and replaced the sum with an integral outside of the first few shells of nearest neighbours. But the method of Ewald remained as the only completely general and accurate method of evaluating the Coulombic sum for a general lattice. Although it was derived in a somewhat different way, Ewald's method is equivalent to accelerated convergence for the special case of [n = 1] .

In the reciprocal lattice of sodium chloride only points with indices (hkl) all even or all odd are permitted by the face-centred symmetry. The reciprocal cell has edge length [1/a] and the reciprocal-axis directions coincide with the direct-lattice axis directions. The closest points to the origin are the eight (111) forms at a distance of $[(1/a)/3^{1/2}]$ . For sodium chloride this distance is 0.3078 Å⁻¹.

Table 3.4.12.1 shows the effect of convergence acceleration on the direct-space portion of the [(n = 1)] sum for the sodium chloride structure. The constant term $[-w \sum q_{j}^{2}]$ is included in the values given. This constant term is always large if w is not zero; for instance, when [w = 0.1] this term is [-277.872] (Table 3.4.12.2 ). For the reciprocal-lattice sum is zero to six figures. Thus, only the direct sum (plus the constant term) is needed, evaluated out to 20 Å in direct space, to obtain six-figure accuracy. As shown in Table 3.4.2.1 above, the same summation effort without the use of accelerated convergence gave 8% error, or only slightly better than one-figure accuracy. The accelerated-convergence technique therefore yielded nearly five orders of magnitude improvement in accuracy, even without evaluation of the reciprocal-lattice sum.

Table 3.4.12.1 | top | pdf |
Accelerated-convergence results for the Coulombic sum ( [n = 1] ) of sodium chloride (kJ mol⁻¹, Å): the direct sum plus the constant term

Limit
6.0	−779.087	−838.145	−860.393	−924.275	−1125.372
8.0	−818.549	−860.194	−863.764	−924.282	−1125.372
10.0	−865.323	−862.818	−863.811	−924.282
12.0	−861.183	−862.824	−863.811
14.0	−862.717	−862.828
16.0	−862.792	−862.828
18.0	−862.810
20.0	−862.825

Table 3.4.12.2 | top | pdf |
The reciprocal-lattice results (kJ mol⁻¹, Å) for the Coulombic sum ( [n = 1] ) of sodium chloride

Limit
0.0	−277.872	−416.806	−555.742	−833.613	−1111.483
0.4	0.000	0.003	0.986	61.451	261.042
0.5		0.003	0.986	61.451	261.042
0.6				61.457	262.542
0.7				61.457	262.542
0.8					262.547
0.9					262.547

The column showing [w = 0.15] shows an example of how the reciprocal-lattice sum can also be neglected if lower accuracy is required. Table 3.4.12.2 shows that the reciprocal-lattice sum is still only 0.003. But now the direct-lattice sum only needs to be evaluated out to 14 Å, with further savings in calculation effort. For w values larger than 0.15 the reciprocal sum is needed. For [w = 0.4] this sum must be evaluated out to 0.8 Å⁻¹ to obtain six-figure accuracy.

Table 3.4.12.3 illustrates an application for the [(n = 6)] dispersion sum. When $[w = 0.1\;\hbox{\AA}^{-1}]$ five figures of accuracy can be obtained without consideration of the reciprocal sum. The direct sum is required out to 18 Å. If [w = 0.15] , better than four-figure accuracy can still be obtained without evaluating the reciprocal-lattice sum. In this case, the direct lattice needs to be summed only to 12 Å, and there is a saving of an order of magnitude in the length of the calculation. As with the Coulombic sum, if w is greater than 0.15 the reciprocal-lattice summation is needed; Table 3.4.12.4 shows the values.

Table 3.4.12.3 | top | pdf |
Accelerated-convergence results for the dispersion sum ( [n = 6] ) of crystalline benzene (kJ mol⁻¹, Å); the figures shown are the direct-lattice sum plus the two constant terms

Limit
6.0	−73.452	−77.761	−79.651	−61.866	76.645
8.0	−79.029	−80.374	−80.256	−61.870	76.645
10.0	−80.217	−80.571	−80.265	−61.870
12.0	−80.527	−80.585	−80.265
14.0	−80.578	−80.585
16.0	−80.588
18.0	−80.589
20.0	−80.589

Table 3.4.12.4 | top | pdf |
The reciprocal-lattice results (kJ mol⁻¹, Å) for the dispersion sum ( [n = 6] ) of crystalline benzene

Limit
0.0	−5.547	−16.706	−32.326	−43.681	80.947
0.3	0.000	−0.004	−0.321	−16.792	−117.106
0.4	0.000	−0.004	−0.324	−18.656	−152.651
0.5			−0.324	−18.716	−155.940
0.6				−18.719	−157.102
0.7				−18.719	−157.227
0.8					−157.233
0.9					−157.234
1.0					−157.234

The time required to obtain a lattice sum of given accuracy will vary depending on the particular structure considered and of course on the computer and program which are used. An example of timing for the benzene dispersion sum is given in Table 3.4.12.5 for the PCK83 program (Williams, 1984 ) running on a VAX-11/750 computer. In this particular case direct terms were evaluated at a rate of about 200 terms s⁻¹ and reciprocal terms, being a sum themselves, were evaluated at a slower rate of about 5 terms s⁻¹.

Table 3.4.12.5 | top | pdf |
Approximate time (s) required to evaluate the dispersion sum ( [n = 6] ) for crystalline benzene within 0.001 kJ mol⁻¹ truncation error

w	Direct terms	Reciprocal terms	Time, direct	Time, reciprocal	Total time
0.0	(not yet converged at 20 Å summation limit)				>107
0.1	15904	0	77	0	77
0.15	4718	34	23	6	29
0.2	2631	78	13	14	27
0.3	1313	246	7	46	53
0.4	524	804	3	149	152

Table 3.4.12.5 shows the time required for evaluation of the dispersion sum using various values of the convergence constant, w. The timing figures show that there is an optimum choice for w; for the PCK83 program the optimum value indicated is 0.15–0.2 Å⁻¹. In the program of Pietila & Rasmussen (1984 ) values in the range 0.15–0.2 Å⁻¹ are also suggested. For the WMIN program (Busing, 1981 ) a slightly higher value of 0.25 Å⁻¹ is suggested. Trial calculations can be used to determine the optimum value of w for the situation of a particular crystal structure, program and computer.

References

Arfken, G. (1970). Mathematical methods for physicists, 2nd ed. New York: Academic Press.Google Scholar

Bertaut, E. F. (1952). L'énergie électrostatique de réseaux ioniques. J. Phys. (Paris), 13, 499–505.Google Scholar

Bertaut, E. F. (1978). The equivalent charge concept and its application to the electrostatic energy of charges and multipoles. J. Phys. (Paris), 39, 1331–1348.Google Scholar

Busing, W. R. (1981). WMIN, a computer program to model molecules and crystals in terms of potential energy functions. Oak Ridge National Laboratory Report ORNL-5747. Oak Ridge, Tennessee 37830, USA.Google Scholar

Cummins, P. G., Dunmur, D. A., Munn, R. W. & Newham, R. J. (1976). Applications of the Ewald method. I. Calculation of multipole lattice sums. Acta Cryst. A32, 847–853.Google Scholar

Davis, P. J. (1972). Gamma function and related functions. Handbook of mathematical functions with formulas, graphs, and mathematical tables, edited by M. Abramowitz & I. A. Stegun, pp. 260–262. London, New York: John Wiley. [Reprint, with corrections of 1964 Natl Bur. Stand. publication.]Google Scholar

DeWette, F. W. & Schacher, G. E. (1964). Internal field in general dipole lattices. Phys. Rev. 137, A78–A91.Google Scholar

Evjen, H. M. (1932). The stability of certain heteropolar crystals. Phys. Rev. 39, 675–694.Google Scholar

Ewald, P. P. (1921). Die Berechnung optischer und elektrostatischer Gitterpotentiale. Ann. Phys. (Leipzig), 64, 253–287.Google Scholar

Fortuin, C. M. (1977). Note on the calculation of electrostatic lattice potentials. Physica (Utrecht), 86A, 574–586.Google Scholar

Glasser, M. L. & Zucker, I. J. (1980). Lattice sums. Theor. Chem. Adv. Perspect. 5, 67–139.Google Scholar

Hastings, C. Jr (1955). Approximations for digital computers. New Jersey: Princeton University Press.Google Scholar

Madelung, E. (1918). Das elektrische Feld in Systemen von regelmässig angeordneten Punktladungen. Phys. Z. 19, 524–532.Google Scholar

Massidda, V. (1978). Electrostatic energy in ionic crystals by the planewise summation method. Physica (Utrecht), 95B, 317–334.Google Scholar

Nijboer, B. R. A. & DeWette, F. W. (1957). On the calculation of lattice sums. Physica (Utrecht), 23, 309–321.Google Scholar

Pietila, L.-O. & Rasmussen, K. (1984). A program for calculation of crystal conformations of flexible molecules using convergence acceleration. J. Comput. Chem. 5, 252–260.Google Scholar

Widder, D. V. (1961). Advanced calculus, 2nd ed. New York: Prentice-Hall.Google Scholar

Williams, D. E. (1971). Accelerated convergence of crystal lattice potential sums. Acta Cryst. A27, 452–455.Google Scholar

Williams, D. E. (1984). PCK83, a crystal molecular packing analysis program. Quantum Chemistry Program Exchange, Department of Chemistry, Indiana University, Bloomington, Indiana 47405, USA.Google Scholar

Williams, D. E. (1989). Accelerated convergence treatment of $[{\bi R}^{-n}]$ lattice sums. Crystallogr. Rev. 2, 3–23. Corrections: Crystallogr. Rev. 2, 163–166.Google Scholar

International Tables for Crystallography (2006). Vol. B. ch. 3.4, pp. 385-397
https://doi.org/10.1107/97809553602060000562

Chapter 3.4. Accelerated convergence treatment of R−n lattice sums1

3.4.1. Introduction

3.4.2. Definition and behaviour of the direct-space sum

3.4.3. Preliminary description of the method

3.4.4. Preliminary derivation to obtain a formula which accelerates the convergence of an sum over lattice points X(d)

3.4.5. Extension of the method to a composite lattice

3.4.6. The case of (Coulombic lattice energy)

3.4.7. The cases of and

3.4.8. Derivation of the accelerated convergence formula via the Patterson function

3.4.9. Evaluation of the incomplete gamma function

3.4.10. Summation over the asymmetric unit and elimination of intramolecular energy terms

3.4.11. Reference formulae for particular values of n

3.4.12. Numerical illustrations

References

Chapter 3.4. Accelerated convergence treatment of R⁻ⁿ lattice sums¹

3.4.4. Preliminary derivation to obtain a formula which accelerates the convergence of an $[{\bi R}^{\bf -{\bi n}}]$ sum over lattice points X(d)

3.4.6. The case of $[{\bf{\bi n} = 1}]$ (Coulombic lattice energy )

3.4.7. The cases of $[{\bf{\bi n} = 2}]$ and $[{\bf{\bi n} = 3}]$