Ewald summering

Ewald summation , uppkallad efter Paul Peter Ewald , är en metod för att beräkna långdistansinteraktioner (t.ex. elektrostatiska interaktioner) i periodiska system. Det utvecklades först som metoden för att beräkna elektrostatiska energier av joniska kristaller , och används nu allmänt för att beräkna långdistansinteraktioner i beräkningskemi . Ewald-summering är ett specialfall av Poisson-summeringsformeln , som ersätter summeringen av interaktionsenergier i det verkliga rummet med en ekvivalent summering i Fourierrymden . I denna metod är interaktionen på lång räckvidd uppdelad i två delar: ett kortdistansbidrag och ett långdistansbidrag som inte har en singularitet . Kortdistansbidraget beräknas i reellt utrymme, medan långdistansbidraget beräknas med hjälp av en Fouriertransform . Fördelen med denna metod är den snabba konvergensen av energin jämfört med en direkt summering. Detta innebär att metoden har hög noggrannhet och rimlig hastighet vid beräkning av långdistansinteraktioner, och den är därmed den de facto standardmetoden för att beräkna långdistansinteraktioner i periodiska system. Metoden kräver laddningsneutralitet hos molekylsystemet för att exakt kunna beräkna den totala Coulombic interaktionen. En studie av trunkeringsfelen som introduceras i energi- och kraftberäkningarna av oordnade punktladdningssystem tillhandahålls av Kolafa och Perram.

Härledning

Ewald summation skriver om interaktionspotentialen som summan av två termer,

\varphi (\mathbf {r} )\ {\stackrel {\mathrm {def} }{=}}\ \ varphi _{sr}(\mathbf {r} )+\varphi _{\ell r}(\mathbf {r} )

,

där $\varphi _{sr}(\mathbf {r} )$ representerar kortdistanstermen vars summa snabbt konvergerar i reella rymden och $\varphi _ {\ell r}(\mathbf {r} )$ representerar långdistanstermen vars summa snabbt konvergerar i Fourierrymd (reciprok). Den långdistansdelen bör vara finit för alla argument (främst r = 0) men kan ha vilken lämplig matematisk form som helst, oftast en Gauss-fördelning . Metoden förutsätter att kortdistansdelen lätt kan summeras; därför blir problemet summeringen av lång sikt. På grund av användningen av Fourier-summan antar metoden implicit att systemet som studeras är oändligt periodiskt (ett förnuftigt antagande för kristallernas inre). En upprepande enhet i detta hypotetiska periodiska system kallas en enhetscell . En sådan cell väljs som "central cell" för referens och de återstående cellerna kallas bilder .

Långdistansinteraktionsenergin är summan av interaktionsenergier mellan laddningarna i en central enhetscell och alla laddningar i gittret. Därför kan den representeras som en dubbel integral över två laddningstäthetsfält som representerar enhetscellens och kristallgittrets fält

E_{\ell r}=\iint d\mathbf {r} \ ,d\mathbf {r} ^{\prime }\,\rho _{\text{TOT}}(\mathbf {r} )\rho _{uc}(\mathbf {r} ^{\prime })\ \varphi _{\ell r}(\mathbf {r} -\mathbf {r} ^{\prime })

där enhetscellens laddningstäthetsfält $\rho _{uc}(\mathbf {r} )$ är en summa över positionerna $\mathbf {r} _{k }$ av laddningarna $q_{k}$ i centralenhetens cell

\rho _{uc}(\mathbf {r} )\ {\stackrel {\ mathrm {def} }{=}}\ \sum _{\mathrm {avgifter} \ k}q_{k}\delta (\mathbf {r} -\mathbf {r} _{k})

och det totala laddningstäthetsfältet $\rho _{\text{TOT}}(\mathbf {r} )$ är samma summa över enhetscellsladdningarna $q_{k }$ och deras periodiska bilder

\rho _{\text{TOT}}(\mathbf {r} )\ {\stackrel {\mathrm {def} }{=}}\ \sum _{n_{1},n_{2} ,n_{3}}\summa _{\mathrm {avgifter} \ k}q_{k}\delta (\mathbf {r} -\mathbf {r} _{k}-n_{1}\mathbf {a} _{1}-n_{2}\mathbf {a} _{2}-n_{3}\mathbf {a} _{3})

Här är $\delta (\mathbf {x} )$ Dirac deltafunktionen , $\mathbf {a} _{1}$ , $\mathbf {a } _{2}$ och $\mathbf {a} _{3}$ är gittervektorerna och $n_{1}$ , $n_{2}$ och $n_{3}$ intervall över alla heltal. Det totala fältet $\rho _{\text{TOT}}(\mathbf {r} )$ kan representeras som en faltning av $\rho _{uc }(\mathbf {r} )$ med en gitterfunktion $L(\mathbf {r} )$

L(\mathbf {r} )\ {\ stackrel {\mathrm {def} }{=}}\ \sum _{n_{1},n_{2},n_{3}}\delta (\mathbf {r} -n_{1}\mathbf {a} _{1}-n_{2}\mathbf {a} _{2}-n_{3}\mathbf {a} _{3})

Eftersom detta är en faltning är Fourier -transformationen av $\rho _{\text{TOT}}(\mathbf {r} )$ en produkt

{\tilde {\rho }}_{\text{TOT}}(\mathbf {k} )={\tilde { L}}(\mathbf {k} ){\tilde {\rho }}_{uc}(\mathbf {k} )

där Fouriertransformen av gitterfunktionen är en annan summa över deltafunktioner

{\tilde {L}} (\mathbf {k} )={\frac {\left(2\pi \right)^{3}}{\Omega }}\summa _{m_{1},m_{2},m_{3}} \delta (\mathbf {k} -m_{1}\mathbf {b} _{1}-m_{2}\mathbf {b} _{2}-m_{3}\mathbf {b} _{3} )

där de reciproka rymdvektorerna är definierade ${\displaystyle \mathbf {b} _{1}\ {\stackrel {\mathrm {def} }{=}}\ 2 \pi {\frac {\mathbf {a} _{2}\times \mathbf {a} _{3}}{\Omega }}} ($ och cykliska permutationer) där $\Omega \ {\stackrel {\mathrm {def} }{=}}\ \mathbf {a} _{1}\cdot \left(\mathbf {a} _{2}\times \mathbf {a} _{3}\right)$ är volymen av den centrala enhetscellen (om den geometriskt är en parallellepiped , vilket ofta men inte nödvändigtvis är fallet). Observera att både $L(\mathbf {r} )$ och ${\tilde {L}}(\mathbf {k} )$ är verkliga, jämna funktioner.

För korthets skull, definiera en effektiv enpartikelpotential

v(\mathbf {r} )\ {\stackrel {\mathrm {def} } {=}}\ \int d\mathbf {r} ^{\prime }\,\rho _{uc}(\mathbf {r} ^{\prime })\ \varphi _{\ell r}(\mathbf {r} -\mathbf {r} ^{\prime })

Eftersom detta också är en faltning är Fouriertransformationen av samma ekvation en produkt

{\tilde {V}}(\mathbf {k} )\ {\stackrel {\mathrm {def} }{ =}}\ {\tilde {\rho }}_{uc}(\mathbf {k} ){\tilde {\Phi }}(\mathbf {k} )

där Fouriertransformen definieras

{\tilde {V}}(\mathbf {k} )=\int d\mathbf {r} \ v(\mathbf {r} )\ e^{-i\mathbf {k} \cdot \mathbf {r} }

Energin kan nu skrivas som en enstaka fältintegral

E_{\ell r}=\int d\mathbf {r} \ \rho _{\text{TOT}}(\mathbf {r } )\ v(\mathbf {r} )

Med hjälp av Plancherel-satsen kan energin också summeras i Fourierrymden

E_{\ell r}=\int {\frac {d\mathbf {k} }{\left(2\pi \right)^{3}}}\ {\tilde { \rho }}_{\text{TOT}}^{*}(\mathbf {k} ){\tilde {V}}(\mathbf {k} )=\int {\frac {d\mathbf {k} }{\left(2\pi \right)^{3}}}{\tilde {L}}^{*}(\mathbf {k} )\left|{\tilde {\rho }}_{uc} (\mathbf {k} )\right|^{2}{\tilde {\Phi }}(\mathbf {k} )={\frac {1}{\Omega }}\summa _{m_{1}, m_{2},m_{3}}\left|{\tilde {\rho }}_{uc}(\mathbf {k} )\right|^{2}{\tilde {\Phi }}(\mathbf {k} )

där $\mathbf {k} =m_{1}\mathbf {b} _{1}+m_{2}\mathbf {b} _ {2}+m_{3}\mathbf {b} _{3}$ i den slutliga summeringen.

Detta är det väsentliga resultatet. När ${\displaystyle {\tilde {\rho }}_{uc}(\mathbf {k} )} har beräknats,$ är summeringen/integrationen över $\mathbf {k}$ enkel och bör konvergera snabbt. Den vanligaste orsaken till bristande konvergens är en dåligt definierad enhetscell, som måste vara laddningsneutral för att undvika oändliga summor.

Partikelnät Ewald (PME) metod

Ewald summation utvecklades som en metod inom teoretisk fysik , långt innan datorernas tillkomst . Ewald-metoden har dock haft stor användning sedan 1970-talet i datorsimuleringar av partikelsystem, särskilt de vars partiklar interagerar via en omvänd kvadratisk kraftlag som gravitation eller elektrostatik . På senare tid har PME också använts för att beräkna $r^{-6}$ delen av Lennard-Jones potential för att eliminera artefakter på grund av trunkering. Tillämpningar inkluderar simuleringar av plasma , galaxer och molekyler .

I partikelnätmetoden, precis som i standard Ewald-summering, delas den generiska interaktionspotentialen upp i två termer $\varphi (\ mathbf {r} )\ {\stackrel {\mathrm {def} }{=}}\ \varphi _{sr}(\mathbf {r} )+\varphi _{\ell r}(\mathbf {r} )$ . Grundidén med partikelnät Ewald summering är att ersätta den direkta summeringen av interaktionsenergier mellan punktpartiklar

E_{\text{TOT}}=\sum _{i,j}\varphi (\mathbf {r } _{j}-\mathbf {r} _{i})=E_{sr}+E_{\ell r}

med två summeringar, en direkt summa $E_{sr}$ av kortdistanspotentialen i det verkliga rummet

E_{sr}=\sum _{i,j}\varphi _{sr}(\mathbf {r} _{j} -\mathbf {r} _{i})

(detta är partikeldelen av partikelnätet Ewald ) och en summering i Fourierrymden av den långdistansdel

E_{\ell r}=\sum _{\mathbf {k} }{\tilde {\Phi }}_{\ell r}(\mathbf {k} )\left|{\tilde {\ rho }}(\mathbf {k} )\right|^{2}

där ${\tilde {\Phi }}_{\ell r}$ och ${\tilde {\rho }}(\mathbf {k} )$ representerar Fouriertransformer av potentialen och laddningstätheten (detta är Ewald -delen). Eftersom båda summeringarna konvergerar snabbt i sina respektive utrymmen (real och Fourier), kan de trunkeras med liten förlust av noggrannhet och stor förbättring av erforderlig beräkningstid. För att effektivt utvärdera Fouriertransformen ${\tilde {\rho }}(\mathbf {k} )$ av laddningstäthetsfältet, använder man den snabba Fouriertransformen , som kräver att densitetsfältet är utvärderas på ett diskret gitter i rymden (detta är nätdelen ) .

På grund av periodicitetsantagandet som är implicit i Ewald-summering kräver tillämpningar av PME-metoden på fysiska system införande av periodisk symmetri. Metoden är således bäst lämpad för system som kan simuleras som oändliga i rumslig utsträckning. I av molekylär dynamik uppnås detta normalt genom att medvetet konstruera en laddningsneutral enhetscell som kan "beläggas" oändligt för att bilda bilder; Men för att korrekt redogöra för effekterna av denna approximation återförs dessa bilder tillbaka i den ursprungliga simuleringscellen. Den övergripande effekten kallas ett periodiskt randvillkor . För att visualisera detta tydligast, tänk på en enhetskub; det övre ansiktet är effektivt i kontakt med det nedre ansiktet, det högra med det vänstra ansiktet och framsidan med baksidan. Som ett resultat måste enhetscellstorleken väljas noggrant för att vara tillräckligt stor för att undvika felaktiga rörelsekorrelationer mellan två ytor "i kontakt", men fortfarande tillräckligt liten för att vara beräkningsmässigt genomförbar. Definitionen av gränsen mellan kort- och långdistansinteraktioner kan också introducera artefakter.

Begränsningen av densitetsfältet till ett nät gör PME-metoden mer effektiv för system med "jämna" variationer i densitet, eller kontinuerliga potentiella funktioner. Lokaliserade system eller de med stora fluktuationer i densitet kan behandlas mer effektivt med den snabba multipolmetoden från Greengard och Rokhlin.

Dipolterm

Den elektrostatiska energin för en polär kristall (dvs. en kristall med en nettodipol $\mathbf {p} _{uc}$ i enhetscellen) är villkorligt konvergent , dvs. beror på summeringsordningen. Till exempel, om dipol-dipol-interaktionerna för en central enhetscell med enhetsceller placerade på en ständigt ökande kub, konvergerar energin till ett annat värde än om interaktionsenergierna hade summerats sfäriskt. Grovt sett uppstår denna villkorliga konvergens eftersom (1) antalet interagerande dipoler på ett skal med radie $R$ växer som $R^{2}$ ; (2) styrkan hos en enskild dipol-dipol-interaktion faller som ${\frac {1}{R^{3}}}$ ; och (3) den matematiska summeringen $\sum _{n=1}^{\infty }{\frac {1}{n}}$ divergerar.

Detta något överraskande resultat kan förenas med den ändliga energin hos verkliga kristaller eftersom sådana kristaller inte är oändliga, dvs har en speciell gräns. Mer specifikt har gränsen för en polär kristall en effektiv ytladdningstäthet på sin yta $\sigma =\mathbf {P} \cdot \mathbf {n}$ där $\mathbf { n}$ är ytnormalvektorn och $\mathbf {P}$ representerar nettodipolmomentet per volym. Interaktionsenergin $U$ för dipolen i en central enhetscell med den ytladdningstätheten kan skrivas

U={\frac {1}{2V_{uc}}}\int {\frac { \left(\mathbf {p} _{uc}\cdot \mathbf {r} \right)\left(\mathbf {p} _{uc}\cdot \mathbf {n} \right)dS}{r^{ 3}}}

där $\mathbf {p} _{uc}$ och $V_{uc}$ är nettodipolmomentet och volymen för enhetscellen, $dS$ är en infinitesimal area på kristallytan och $\mathbf {r}$ är vektorn från den centrala enhetscellen till den infinitesimala arean. Denna formel är resultatet av att integrera energin $dU=-\mathbf {p} _{uc}\cdot \mathbf {dE}$ där $d\mathbf { E}$ representerar det oändliga elektriska fältet som genereras av en infinitesimal ytladdning ${\displaystyle dq\ {\stackrel {\mathrm {def} }{=}}\ \sigma dS} ($ Coulombs lag )

d\mathbf {E} \ {\stackrel {\mathrm {def } }{=}}\ \left({\frac {-1}{4\pi \epsilon }}\right){\frac {dq\ \mathbf {r} }{r^{3}}}=\ vänster({\frac {-1}{4\pi \epsilon }}\right){\frac {\sigma \,dS\ \mathbf {r} }{r^{3}}}

Det negativa tecknet härrör från definitionen av ${\displaystyle \mathbf {r} } ,$ som pekar mot laddningen, inte bort från den.

Historia

Ewald-summeringen utvecklades av Paul Peter Ewald 1921 (se referenser nedan) för att bestämma den elektrostatiska energin (och därmed Madelung-konstanten ) för jonkristaller.

Skalning

I allmänhet ger olika Ewald-summeringsmetoder olika tidskomplexitet . Direkt beräkning ger $O(N^{2})$ , där $N$ är antalet atomer i systemet. PME-metoden ger $O(N\,\log N)$ .

Se även

Ewald, P (1921). "Die Berechnung optischer und elektrostatischer Gitterpotentiale" . Ann. Phys . 369 (3): 253–287. Bibcode : 1921AnP...369..253E . doi : 10.1002/andp.19213690304 .
Darden, T; Perera, L; Li, L; Pedersen, L (1999). "Nya knep för modellerare från verktygslådan för kristallografi: Ewald-algoritmen för partikelnät och dess användning i nukleinsyrasimuleringar" . Struktur . 7 (3): 55–60 kr. doi : 10.1016/S0969-2126(99)80033-1 . PMID 10368306 . S2CID 40964921 .
Frenkel, D., & Smit, B. (2001). Förstå molekylär simulering: från algoritmer till applikationer, Akademisk press.