Somatório de Ewald

O somatório de Ewald (ou, às vezes, a soma de Ewald ) é um método de calcular as energias de interação de sistemas periódicos (e, em particular, cristais ), e especialmente energias eletrostáticas. O somatório de Ewald é um caso especial da fórmula do somatório de Poisson, com a substituição do somatório das energias de interação no espaço real por um somatório equivalente em um espaço de Fourier. A vantagem dessa abordagem é a rápida convergência da soma no espaço de Fourier, em comparação com seu equivalente no espaço real, quando as interações ocorrem a longa distância. Como as energias eletrostáticas incluem termos de interação de curto e longo alcance, é muito interessante decompor o potencial de interação em termos de curto alcance - que é somado no espaço real - e longo range - que é somado no espaço de Fourier. O método foi desenvolvido por Paul Peter Ewald, em 1921 (ver referências) a fim de determinar a energia eletrostática (e, portanto, a constante de Madelung) dos cristais iônicos.

Princípio

O somatório de Ewald reescreve o potencial de interação como a soma de dois termos  :

$\phi (\mathbf {r} )\ {\stackrel {\mathrm {def} }{=}}\ \phi _{sr}(\mathbf {r} )+\phi _{lr}(\mathbf {r} )$

onde $\phi _{sr}(\mathbf {r} )$ representa o termo de curto alcance, cujo somatório é rápido no espaço real, enquanto o $\phi _{lr}(\mathbf {r} )$ é o termo de longo alcance cujo somatório é rápido no espaço de Fourier. A parte de longo alcance deve ser finalizada para todos os argumentos (particularmente para $r=0$ ), mas pode ter qualquer forma matemática adequada, mais comumente uma distribuição gaussiana. O método postula que a parte de curto alcance pode ser somada facilmente, portanto, o problema é lidar com a soma do longo prazo. Em razão do uso de uma integral de Fourier, o método postula implicitamente que o sistema estudado é infinitamente periódico (postulado relevante para um cristal perfeito - isto é, para o interior de um cristal real). A unidade básica reproduzida pela periodicidade é chamada célula unitária. Essa malha é escolhida como a “malha central” de referência e as células reproduzidas periodicamente são chamadas de imagens .

A energia de interação de longo alcance é a soma das energias de interação entre as cargas de uma célula da unidade central e todas as cargas da rede. Dessa forma, pode ser descrita como uma dupla densidade de carga dois integrais dos campos correspondentes aos campos na célula unitária e que da rede cristalina.

$E_{lr}=\int \int d\mathbf {r} d\mathbf {r} ^{\prime }\rho _{\text{tot}}(\mathbf {r} )\rho _{uc}(\mathbf {r} ^{\prime })\ \phi _{lr}(\mathbf {r} -\mathbf {r} ^{\prime })$

em que o campo de densidade de carga da célula unitária $\rho _{uc}(\mathbf {r} )$ é uma integração em posições $\mathbf {r} _{k}$ das cargas $q_{k}$ na célula da unidade central  :

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

e o campo de densidade total de carga $\rho _{\text{tot}}(\mathbf {r} )$ é a mesma integração nas cargas da malha unitária $q_{k}$ e suas imagens periódicas  :

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

Aqui, $\delta (\mathbf {x} )$ é a função Dirac δ, $\mathbf {a} _{1}$, $\mathbf {a} _{2}$ e $\mathbf {a} _{3}$ são os vetores de malha e $n_{1}$, $n_{2}$ e $n_{3}$ são inteiros. O campo total $\rho _{\text{tot}}(\mathbf {r} )$ pode ser descrito como uma convolução de $\rho _{uc}(\mathbf {r} )$ com uma função de rede $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})$

Como temos uma convolução, a transformação de Fourier de $\rho _{\text{tot}}(\mathbf {r} )$ é um produto:

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

na qual a transformada de Fourier da função de rede é outra integração nas funções δ:

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

onde os vetores espaciais recíprocos são definidos da seguinte maneira: $\mathbf {b} _{1}\ {\stackrel {\mathrm {def} }{=}}\ \mathbf {a} _{2}\times \mathbf {a} _{3}/\Omega$ (exceto permutações cíclicas) onde $\Omega \ {\stackrel {\mathrm {def} }{=}}\ \mathbf {a} _{1}\cdot \left(\mathbf {a} _{2}\times \mathbf {a} _{3}\right)$ é o volume da célula da unidade central (se for geometricamente um paralelepípedo, o que geralmente é o caso, mas não necessariamente). $L(\mathbf {r} )$ e ${\tilde {L}}(\mathbf {k} )$ são funções reais pares.

Por questões de brevidade, definimos um potencial efetivo com uma única partícula  :

$v(\mathbf {r} )\ {\stackrel {\mathrm {def} }{=}}\ \int d\mathbf {r} ^{\prime }\rho _{uc}(\mathbf {r} ^{\prime })\ \phi _{lr}(\mathbf {r} -\mathbf {r} ^{\prime })$

Sendo essa equação também uma convolução, a transformação de Fourier dessa equação também é um produto:

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

na qual a transformação de Fourier é definida da seguinte maneira:

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

A energia pode, então, ser escrita como parte integrante de um único campo:

$E_{lr}=\int d\mathbf {r} \ \rho _{\text{tot}}(\mathbf {r} )\ v(\mathbf {r} )$

Usando o teorema de Parseval, a energia também pode ser somada no espaço de Fourier:

$E_{lk}=\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 }}\sum _{m_{1},m_{2},m_{3}}\left|{\tilde {\rho }}_{uc}(\mathbf {k} )\right|^{2}{\tilde {\Phi }}(\mathbf {k} )$

onde $\mathbf {k} =m_{1}\mathbf {b} _{1}+m_{2}\mathbf {b} _{2}+m_{3}\mathbf {b} _{3}$ na última integral.

Este é o resultado principal. Uma vez que ${\tilde {\rho }}_{uc}(\mathbf {k} )$ é calculado, a soma/integração em $\mathbf {k}$ é simples e deve convergir rapidamente. A razão mais comum para a ausência de convergência é uma definição deficiente da célula unitária, que deve ser eletricamente neutra para evitar somas infinitas.

Método de malha especial de Ewald

O somatório de Ewald foi desenvolvido como um método da física teórica, muito antes do advento dos computadores e da ciência da computação. No entanto, é muito difundido desde a década de 1970 nas simulações numéricas de sistemas de partículas, particularmente aqueles que interagem com as leis de força em inverso do quadrado, como gravidade ou electrostática. As aplicações dessas simulações incluem objetos tão diversos quanto plasmas, galáxias ou biomoléculas .

Como em um somatório "normal" de Ewald, um potencial genérico de interação é separado em dois termos  : $\phi (\mathbf {r} )\ {\stackrel {\mathrm {def} }{=}}\ \phi _{sr}(\mathbf {r} )+\phi _{lr}(\mathbf {r} )$ - um termo de interação de curto alcance $\phi _{sr}(\mathbf {r} )$ à integração rápida no espaço real é um termo de interação de longo alcance $\phi _{lr}(\mathbf {r} )$ que se integra rapidamente ao espaço de Fourier. A ideia básica do somatório de uma malha de Ewald específica é substituir o somatório direto das energias de interação entre partículas pontuais:

$E_{TOT}=\sum _{i,j}\phi (\mathbf {r} _{j}-\mathbf {r} _{i})=E_{sr}+E_{lr}$

com dois somatórios, integração direta $E_{sr}$ no potencial de curto alcance no espaço real:

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

(Porção particular da malha de Ewald) e integração no espaço de Fourier da porção de longo alcance:

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

onde ${\tilde {\Phi }}_{lr}$ e ${\tilde {\rho }}(\mathbf {k} )$ representam as transformadas de Fourier do potencial elétrico e da densidade de carga (parte de Ewald). Devido à sua rápida convergência em seus respectivos espaços (real e de Fourier), eles podem ser truncados sem nenhuma perda real de precisão e um ganho significativo no tempo de computação necessário. Para avaliar efetivamente a transformação de Fourier ${\tilde {\rho }}(\mathbf {k} )$ no campo de densidade de carga, pode-se usar a rápida transformação de Fourier, que requer que a densidade de campo seja avaliada sobre uma grade discreta do espaço (malha à parte). Por causa do postulado de periodicidade implícita do somatório de Ewald, as aplicações do método de malha particular de Ewald (particular mesh Ewald em inglês, abreviada pela sigla PME) aos sistemas físicos exigem a imposição de simetria por periodicidade. Por isso, o método é muito mais eficiente para sistemas que podem ser simulados como infinitos no espaço. Nas simulações de dinâmica molecular, essa configuração é normalmente obtida construindo deliberadamente uma célula unitária de carga neutra que pode ser infinitamente reproduzida por ladrilhos afim de formar imagens; no entanto, para levar em conta os efeitos dessa aproximação, essas imagens são reincorporadas na caixa de simulação original. O efeito geral é semelhante a uma versão tridimensional do jogo Asteroids, na qual cada dimensão "se dobra" em si mesma. Essa propriedade da malha é chamada de condição de limite periódica. Para melhor visualizar isso, podemos imaginar um cubo unitário; o lado superior está em contato com o lado inferior, o lado direito com o lado esquerdo e a frente com o lado traseiro. Isso resulta em que o tamanho da caixa deve ser escolhido com cuidado para evitar gerar correlações de movimento impróprias entre duas faces "em contato", mas também suficientemente pequenas para serem processadas numericamente. A definição do raio de corte, a distância na qual se passa da parte de interação de curta distância da parte de interação de longa distância, também pode introduzir artefatos.

A restrição do campo de densidade a uma malha torna o método PME ainda mais eficiente para sistemas com variações de densidade “suaves” ou funções de potenciais contínuos. Os sistemas localizados (por exemplo, sistemas radicais) ou com fortes variações de densidade podem ser tratados com mais eficácia pelo método rápido multipolar de Greengard e Rokhlin.

Termo dipolar

A energia eletrostática de um cristal polar (isto é, um cristal com um dipolo identificado $\mathbf {p} _{uc}$ na célula unitária) é semi-convergente, isto é, depende da ordem da soma. Por exemplo, se as interações dipolo-dipolo de uma célula unitária central com células unitárias estiverem localizadas em um cubo sempre crescente, a energia converge para um valor diferente daquele obtido para uma integração esférica. Em outros termos, essa semi-convergência ocorre porque: em primeiro lugar, o número de dipolos que interagem numa concha de raio R aumenta em R² ; segundo, a força de uma interação dipolo-dipolo simples diminui conforme ${\frac {1}{R^{3}}}$ ; e terceiro, o somatório matemático $\sum _{n=1}^{\infty }{\frac {1}{n}}$ diverge. Este resultado relativamente surpreendente pode ser reconciliado com o fato de que a energia dos cristais reais é finita, porque estes não são infinitos e, portanto, têm limites particulares. Mais especificamente, a fronteira de um cristal polar mostra uma densidade de carga superficial efetiva em sua superfície de $\sigma =\mathbf {P} \cdot \mathbf {n}$, onde n é o vetor normal na superfície e P representa o momento dipolar líquido volumétrico. A energia de interação do dipolo em uma célula unitária central com essa densidade de carga superficial pode ser escrita como:

$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}}}$

onde $\mathbf {p} _{uc}$ e $V_{uc}$ são os momentos de dipolo líquidos e o volume da célula unitária, $dS$ é uma área infinitesimal na superfície do cristal e r é o vetor da célula unitária com essa área infinitesimal. Esta fórmula vem da integração de energia $dU=-\mathbf {p} _{uc}\cdot d\mathbf {E}$ onde $d\mathbf {E}$ representa o campo elétrico infinitesimal gerado por uma carga superficial infinitesimal $dq\ {\stackrel {\mathrm {def} }{=}}\ \sigma \,dS$ (Lei de Coulomb):

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

O sinal negativo vem da definição de r, que aponta para a carga, não da carga.

Escalabilidade do método

Em geral, diferentes métodos de somatório de Ewald geram diferentes escalabilidades. Assim, os cálculos diretos fornecem uma escalabilidade em O(N²), onde N é o número de átomos no sistema. O PME fornece escalabilidade em O(N log N).

Ver também

• Paul Peter Ewald
• Constante de Madelung
• Fórmula de somatório de Poisson

Referências

• Ewald P. (1921) " Die Berechnung optischer und elektrostatischer Gitterpotentiale ", Ann. Phys. 64, 253-287.
• Darden T, Perera L., Li L. e Pedersen L. (1999) " Novos truques para modeladores do kit de ferramentas de cristalografia: o algoritmo de Ewald da malha de partículas e seu uso em simulações de ácidos nucleicos ", Estrutura 7, R55-R60.
• Schlick T. (2002). Modelagem Molecular e Simulação: Um Guia Interdisciplinar Matemática Aplicada Interdisciplinar de Springer-Verlag, Biologia Matemática, vol. 21 Nova Iorque, NY.
• Este artigo foi inicialmente traduzido, total ou parcialmente, do artigo da Wikipédia em inglês cujo título é «Ewald summation», especificamente desta versão.