Ewald Summation Rewrites the Interaction Potential As the Sum of Two Terms

Ewald summation

Ewald summation is a method for computing the interaction energies of periodic systems (e.g. crystals), particularly electrostatic energies. Ewald summation is a special case of the Poisson summation formula, replacing the summation of interaction energies in real space with an equivalent summation in Fourier space. The advantage of this approach is the rapid convergence of the Fourier-space summation compared to its real-space equivalent when the real-space interactions are long-range. Because electrostatic energies consist of both short- and long-range interactions, it is maximally efficient to decompose the interaction potential into a short-range component summed in real space and a long-range component summed in Fourier space.

Derivation

Ewald summation rewrites the interaction potential as the sum of two terms

where represents the short-range term that sums quickly in real space and represents the long-range term that sums quickly in Fourier space. The long-ranged part should be finite for all arguments (most notably r = 0) but may have any convenient mathematical form, most typically a Gaussian distribution. The method assumes that the short-range part can be summed easily; hence, the problem becomes the summation of the long-range term. Due to the use of the Fourier sum, the method implicitly assumes that the system under study is infinitely periodic (a sensible assumption for the interiors of crystals). One repeating unit of this hypothetical periodic system is called a unit cell. One such cell is chosen as the "central cell" for reference and the remaining cells are called images.

The long-range interaction energy is the sum of interaction energies between the charges of a central unit cell and all the charges of the lattice. Hence, it can be represented as a double integral over two charge density fields representing the fields of the unit cell and the crystal lattice

where the unit-cell charge density field is a sum over the positions of the charges qk in the central unit cell

and the total charge density field is the same sum over the unit-cell charges qk and their periodic images

Here, is the Dirac delta function, , and are the lattice vectors and n1, n2 and n3 range over all integers. The total field can be represented as a convolution of with a lattice function

Since this is a convolution, the Fourier transformation of is a product

where the Fourier transform of the lattice function is another sum over delta functions

where the reciprocal space vectors are defined (and cyclic permutations) where is the volume of the central unit cell (if it is geometrically a parallelepiped, which is often but not necessarily the case). Note that both and are real, even functions.

For brevity, define an effective single-particle potential

Since this is also a convolution, the Fourier transformation of the same equation is a product

where the Fourier transform is defined

The energy can now be written as a single field integral

Using Parseval's theorem, the energy can also be summed in Fourier space

where in the final summation.

This is the essential result. Once is calculated, the summation/integration over is straightforward and should converge quickly. The most common reason for lack of convergence is a poorly defined unit cell, which must be charge neutral to avoid infinite sums.

[edit]Particle mesh Ewald (PME) method

Ewald summation was developed as a method of theoretical physics, long before the advent of computers. However, the Ewald method has enjoyed widespread use since the 1970's in computer simulations of particle systems, especially those interacting via an inverse squareforce law such as gravity or electrostatics. Applications include simulations of plasmas, galaxies and biomolecules.

As in normal Ewald summation, a generic interaction potential is separated into two terms - a short-ranged part that sums quickly in real space and a long-ranged part that sums quickly in Fourier space. The basic idea of particle mesh Ewald summation is to replace the direct summation of interaction energies between point particles

with two summations, a direct sum Esr of the short-ranged potential in real space

(that's the particle part of particle mesh Ewald) and a summation in Fourier space of the long-ranged part

where and represent the Fourier transforms of the potential and the charge density (that's the Ewald part). Since both summations converge quickly in their respective spaces (real and Fourier), they may be truncated with little loss of accuracy and great improvement in required computational time. To evaluate the Fourier transform of the charge density field efficiently, one uses the Fast Fourier transform, which requires that the density field be evaluated on a discrete lattice in space (that's the mesh part).

Due to the periodicity assumption implicit in Ewald summation, applications of the PME method to physical systems require the imposition of periodic symmetry. Thus, the method is best suited to systems that can be simulated as infinite in spatial extent. In molecular dynamics simulations this is normally accomplished by deliberately constructing a charge-neutral unit cell that can be infinitely "tiled" to form images; however, to properly account for the effects of this approximation, these images are reincorporated back into the original simulation cell. The overall effect is one of a three-dimensional version of the Asteroids game, in which each dimension "wraps around" on itself. This property of the cell is called a periodic boundary condition. To visualize this most clearly, think of a unit cube; the upper face is effectively in contact with the lower face, the right with the left face, and the front with the back face. As a result the unit cell size must be carefully chosen to be large enough to avoid improper motion correlations between two faces "in contact", but still small enough to be computationally feasible. The definition of the cutoff between short- and long-range interactions can also introduce artifacts.

The restriction of the density field to a mesh makes the PME method more efficient for systems with "smooth" variations in density, or continuous potential functions. Localized systems or those with large fluctuations in density may be treated more efficiently with the fast multipole method of Greengard and Rokhlin.

Dipole term

The electrostatic energy of a polar crystal (i.e., a crystal with a net dipole in the unit cell) is conditionally convergent, i.e., depends on the order of the summation. For example, if the dipole-dipole interactions of a central unit cell with unit cells located on an ever-increasing cube, the energy converges to a different value than if the interaction energies had been summed spherically. Roughly speaking, this conditional convergence arises because (1) the number of interacting dipoles on a shell of radius R grows like R2; (2) the strength of a single dipole-dipole interaction falls like ; and (3) the mathematical summation diverges.

This somewhat surprising result can be reconciled with the finite energy of real crystals because such crystals are not infinite, i.e., have a particular boundary. More specifically, the boundary of a polar crystal has an effective surface charge density on its surface where is the surface normal vector and represents the net dipole moment per volume. The interaction energy U of the dipole in a central unit cell with that surface charge density can be written

where and Vuc are the net dipole moment and volume of the unit cell, dS is an infinitesimal area on the crystal surface and is the vector from the central unit cell to the infinitesimal area. This formula results from integrating the energy where represents the infinitesimal electric field generated by an infinitesimal surface charge (Coulomb's law)

The negative sign derives from the definition of , which points towards the charge, not away from it.

[edit]History

The Ewald summation was developed by Paul Peter Ewald in 1921 (see References below) to determine the electrostatic energy (and, hence, the Madelung constant) of ionic crystals.

La sommation d'Ewald (ou parfois somme d'Ewald) est une méthode de calcul des énergies d'interaction de systèmes périodiques (et particulier des cristaux), et tout particulièrement les énergies électrostatiques. La sommation d'Ewald est un cas particulier de la formule sommatoire de Poisson, avec le remplacement de la sommation des énergies d'interaction dans l'espace réel par une sommation équivalente dans un espace de Fourier. L'avantage de cette approche est la convergence rapide de la sommation dans l'espace de Fourier comparée à son équivalent dans l'espace réel lorsque les interactions se font à longue distance. Les énergies électrostatiques comprenant à la fois des termes d'interactions de courtes et de longues portées, il est très intéressant de décomposer le potentiel d'interaction en termes de courte portée - dont la sommation se fait dans l'espace réel - et de longue portée - dont la sommation se fait dans l'espace de Fourier.
La méthode fut développée par Paul Peter Ewald en 1921 (voir références) afin de déterminer l'énergie électrostatique (et, par là, la constante de Madelung) des cristaux ioniques.

La sommation d'Ewald réécrit le potentiel d'interaction comme la somme de deux termes:

où représente le terme de courte portée dont la sommation est rapide dans l'espace réel et le terme de longue portée dont la sommation est rapide est rapide dans l'espace de Fourier. La partie à longue portée doit être finie pour tous les arguments (en particulier pour r = 0) mais peuvent avoir une forme mathématique adéquate quelconque, le plus généralement une distribution gaussienne. La méthode postule que la partie à courte portée peut être sommée facilement, ainsi, le problème est de traiter la sommation du terme à longue portée. En raison de l'utilisation d'une intégrale de Fourier, la méthode postule implicitement que le système étudié est infiniment périodique (postulat pertinent pour un cristal parfait - i.e. pour l'intérieur d'un cristal réel). L'unité de base reproduite par périodicité est appelée maille unitaire. Une telle maille est choisie comme «maille centrale» de référence et les cellules reproduites par périodicité sont appelées images.
L'énergie d'interaction à longue portée est la somme des énergies d'interaction entre les charges d'une maille unitaire centrale et toutes les charges du réseau. Ainsi, elle peut être décrite comme une double intégrale deux champs de densités de charges correspondant aux champs de la maille unitaire et de celle du réseau cristallin.

où le champ de densité de charge de la maille unitaire est une intégration sur les positions des charges qk dans la maille unitaire centrale:

et le champ de densité de charge totaleest la même intégration sur les charges de la maille unitaire qk et leurs images périodiques:

Ici, est la fonction δ de Dirac, , et sont les vecteurs de maille et n1, n2 et n3 sont des entiers. Le champ total peut être décrit comme une convolution de avec une fonction de réseau:

Puisque l'on a une convolution, la transformation de Fourier de est un produit:

dans laquelle la transformée de Fourier de la fonction de réseau est une autre intégration sur des fonctions δ:

où les vecteurs de l'espace réciproque sont définis comme il suit: (aux permutations cycliques près) où est le volume de la maille unitaire centrale (si c'est géométriquement un parallélépipède, ce qui est souvent mais pas nécessairement le cas). Notons que et sont des fonctions réelles paires.
Pour la brièveté, on définit un potentiel effectif à une seule particule:

Cette équation étant aussi une convolution, la transformation de Fourier de cette équation est aussi un produit:

dans laquelle la transformée de Fourier est définié de la manière suivante:

L'énergie peut maintenant être écrite comme une intégrale d'un unique champ:

En utilisant le théorème de Parseval, l'énergie peut aussi être sommée dans l'espace de Fourier:

où dans la dernière intégrale.

C'est le résultat principal. Une fois que est calculé, la sommation/intégration sur est simple et devrait rapidement converger. La raison la plus courante d'absence de convergence est une mauvaise définition de la maille unitaire, qui doit être électriquement neutre afin d'éviter des sommes infinies.

Méthode du maillage particulier d'Ewald

La sommation d'Ewald fut developpée comme méthode de physique théorique, bien avant la venue des ordinateurs et de l'informatique. Cependant, elle s'est très largement répandue depuis les années 1970 dans les simulations numériques de systèmes de particules, et tout particulièrement celles interagissant par des lois de forces en carré inverse, comme la gravité ou l'électrostatique. Les applications de ces simulations comprennent des objets aussi divers que les plasmas, les galaxies ou les biomolécules.
Comme dans une sommation d'Ewald «normale», un potentiel d'interaction générique est séparé en deux termes: - un terme d'interaction à courte portée à intégration rapide dans l'espace réel est un terme d'interaction à longue portée qui s'intègre rapidement dans l'espace de Fourier. L'idée de base de la sommation sur un maillage particulier d'Ewald est de remplacer la sommation directe des énergies d'interaction entre les particules ponctuelles:

avec deux sommations, une intégration directe Esr du potentiel à courte portée dans l'espace réel:

(partie particulière du maillage particulier d'Ewald) et une intégration dans l'espace de Fourier de la partie à longue portée:

où et représentent les transformées de Fourier du potentiel et la densité de charge (partie Ewald). En raison de leurs convergences rapides dans leurs espaces respectifs (réel et de Fourier), elles peuvent être tronquées sans perte réelle de précision et un gain important du temps de calcul requis. Afin d'évaluer efficacement la transformée de Fourier du champ de densité de charge, on peut utiliser la transformée de Fourier rapide, qui nécessite que le champ de densité soit évaluée sur une maille discrète de l'espace (partie maillage).
En raison du postulat de périodicité implicite de la sommation d'Ewald, les applications de la méthode du maillage particulier d'Ewald (particular mesh Ewald en anglais, abrégé par l'acronyme PME) aux systèmes physiques nécessitent d'imposer une symétrie par périodicité. Pour cette raison, la méthode est bien plus efficace pour les systèmes pouvant être simulé comme étant infinis dans l'espace. Dans les simulations de dynamique moléculaire, cette configuration est normalement atteinte en construisant de manière délibérée une maille unitaire de charge neutre qui peut être infiniment reproduite par pavage afin de former des images; cependant, afin de prendre correctement en compte les effets de cette approximation, ces images sont réincorporées dans la boîte de simulation originale. L'effet global est similaire à une version tridimensionnelle du jeu Asteroids, dans laquelle chaque dimension se «replie» sur elle-même. Cette propriété de la maille est appelée condition périodique aux limites. Afin de mieux visualiser cela, on peut imaginer un cube unitaire; la face du dessus est en contact avec la face du dessois, la face de droite avec celle de gauche, et celle de devant avec celle de derrière. Il en résulte que la taille de la boîte doit être choisie avec précaution pour éviter de générer des corrélations de mouvement impropres entre deux faces «en contact», mais également assez petite pour pouvoir être traitée numériquement. La définition du rayon de coupure, distance à laquelle on passe de la partie d'interaction à courte distance de la partie interaction à longue distance, peut aussi introduire des artéfacts.
La restriction du champ de densité à un maillage rend la méthode PME encore plus efficace pour des systèmes avec des variations de densité «douces», ou des fonctions de potentiels continues. Les systèmes localisés (par exemple des systèmes radicalaires) ou avec de fortes variations de densité pourront être traités plus efficacement par la méthode des multipôles rapides de Greengard et Rokhlin.

Terme dipolaire

L'énergie électrostatique d'un cristal polaire (c'est-à-dire un cristal avec un dipôle identifié dans la maille unitaire) est semi-convergente, i.e. elle dépend de l'ordre de la sommation. Par exemple, si les interactions dipôle-dipôle d'une maille unitaire centrale avec des mailles unitaires sont localisées sur un cube toujours croissant, l'énergie converge vers une valeur différente de celle obtenue pour une intégration sphérique. En d'autres termes, cette semi-convergence se produit car: premièrement, le nombre de dipôles interagissants dans une couche de rayon R augmente en R2; deuxièmement, la force d'une interaction dipôle-dipôle simple décroît en ; et troisièmement, la sommation mathématique diverge.
Ce résultat relativement surprenant peut être concilié avec le fait que l'énergie des cristaux réels est finie car ceux-ci ne sont pas infinis, et donc des limites particulières. Plus spécifiquement, la frontière d'un cristal polaire montre une densité de charge effective surfacique sur sa surface où est le vecteur normal à la surface et représente le moment dipolaire net volumique. L'énergie d'interaction du dipôle dans une maille unitaire centrale avec cette densité de charge surfacique peut être écrite:

où et Vuc sont les moments dipolaires nets et le volume de la maille unitaire, dS une aire infinitésimale sur la surface du cristal et est le vecteur de la maille unitaire à cette aire infinitésimale. Cette formule provient de l'intégration de l'énergie où représente le champ électrique infinitésimale généré par une charge surfacique infinitésimale (loi de Coulomb):

Le signe négatif provient de la définition de , qui pointe vers la charge, et non de la charge.

Scalabilité de la méthode[modifier]

De manière générale, différentes méthodes de sommation d'Ewald donnent des scalabilités différentes. Ainsi, les calculs directs donnent une scalabilité en O(N2), où N est le nombre d'atomes du système. La PME procure elle une scalabilité en O(N log N).