Um potencial de interação para o estudo de materiais e simulações por dinâmica molecular

An interaction potential for materials and molecular dynamics simulations

interaction potential; molecular-dynamics simulation; structural correlations

Revisão

UM POTENCIAL DE INTERAÇÃO PARA O ESTUDO DE MATERIAIS E SIMULAÇÕES POR DINÂMICA MOLECULAR

José Pedro Rino e Nelson Studart

Departamento de Física, Universidade Federal de São Carlos, 13565-905 São Carlos - SP

Recebido em 4/12/00; aceito em 28/3/01

AN INTERACTION POTENTIAL FOR MATERIALS AND MOLECULAR DYNAMICS SIMULATIONS. The Vashishta-Rahman effective interaction potential, based on the Pauling's concept of "ionic radii", has been successfully employed to investigate structural and dynamical properties of different classes of material. By celebrating Pauling's birth centenary, we review the building up of the Vashishta-Rahman potential and we present molecular-dynamics simulation results for structure and dynamics of superionic materials, chalcogenide glasses and metallic oxides.

Keywords: interaction potential; molecular-dynamics simulation; structural correlations.

INTRODUÇÃO

Simulação computacional tem se constituído numa ferramenta poderosa para a compreensão dos fenômenos físicos e químicos em sistemas da matéria condensada. Podemos apontar alguns motivos que justificam o grande investimento nesta abrangente área de pesquisa:

- Simulações computacionais possibilitam uma estreita ligação entre a teoria e a experiência. É essencial reconhecer que a simulação científica não substitui a teoria e a experiência, mas contribui enormemente para o avanço de ambas. Simulações numéricas bem elaboradas podem conduzir a novos discernimentos teóricos e descobertas experimentais;

- Um grande número de graus de liberdade, quebras de simetrias e efeitos não-lineares são incorporados mais facilmente em experiências computacionais do que através de métodos puramente teóricos. Atualmente, simulações têm sido realizadas para predizer o comportamento em escala atômica de muitos milhões de moléculas ou reações químicas em sistemas complexos;

- O universo acessível à simulação computacional não está limitado a processos que ocorrem na natureza. Neste sentido, problemas cujos princípios científicos básicos não estão ainda bem estabelecidos podem ser abordados por métodos numéricos avançados que podem possibilitar, por exemplo, a melhoria do processamento e desempenho de materiais estruturais bem como a previsão do comportamento de uma ampla gama de sistemas biológicos.

Uma interação profícua entre cientistas das áreas básicas com os matemáticos e os cientistas de computação é essencial para o desenvolvimento de novas arquiteturas computacionais, novos algoritmos e melhores técnicas de visualização que permitirão uma interpretação mais apurada dos dados obtidos na simulação. Diversas análises prospectivas têm sido feitas recentemente tanto no aspecto mais abrangente da pesquisa científica e tecnológica¹ como no caso específico da Física e Ciência dos Materiais^2,3,4.

Dentre as várias áreas de pesquisa onde simulações científicas têm sido largamente utilizadas insere-se, sem dúvida, a investigação das propriedades gerais de materiais e processos artificiais de produção de novos materiais. Uma comparação acurada dos dados obtidos pela simulação com resultados experimentais pode levar ao entendimento das interações microscópicas. Uma simulação realista pode projetar teoricamente as propriedades desejáveis de materiais e somente aqueles potencialmente "promissores" viriam a ser investigados do ponto de vista experimental.

A partir da descoberta da difração de raios-X e da constituição atômica da matéria, uma teoria de materiais começou a ser desenvolvida com base nas interações entre seus constituintes. Modelos teóricos bastante simplificados foram satisfatórios para explicar grande parte das propriedades elétricas e térmicas de materiais mesmo antes do advento da Mecânica Quântica. No entanto, em certos materiais de interesse atual, as interações atômicas são bastante complicadas exigindo modelos mais sofisticados e novas técnicas de simulações atomísticas que requerem computação de alto desempenho.

Os primeiros estudos acerca da matéria condensada através de simulação computacional foram iniciados por Metropolis et al.⁵ que usaram uma primeira versão do atualmente conhecido método de Monte Carlo. Este método permite o acesso aleatório ao espaço de configuração através de uma regra simples, qual seja, o estabelecimento de uma hierarquia de regiões do espaço de configuração determinada pelo fator de Boltzmann. Assim, o método permite calcular médias de quantidades físicas uma vez que estabelece uma estatística bastante eficiente, mas limita-se a propriedades estáticas.

Um método diferente, denominado Dinâmica Molecular (DM), consiste em determinar explicitamente as trajetórias de pontos representativos do espaço de fase através da solução numérica das equações do movimento. As primeiras simulações de DM foram realizadas por Alder e Wainwrigth⁶ com o propósito de estudar o conhecido paradoxo da reversibilidade: um sistema clássico de muitas partículas é governado pelas equações temporais reversíveis enquanto a descrição macroscópica (termodinâmica) do mesmo sistema está baseada por leis irreversíveis. Eles mostraram que a distribuição de velocidades do sistema de 100 esferas impenetráveis convergia rapidamente ao equilíbrio. A primeira aplicação do método de DM ao estudo de materiais foi feita por Vineyard et al.⁷ através da investigação do processo de dano no material por radiação usando um potencial repulsivo de curto alcance e um potencial responsável pela coesão do cristal. Em seguida, Rahman⁸ foi o primeiro a investigar sistemas descritos por potenciais contínuos simulando o argônio líquido. Foi surpreendente observar que um sistema de 864 partículas, com condições periódicas de contorno, poderia reproduzir satisfatoriamente as propriedades termodinâmicas de sistemas reais. Mais informações sobre a história e uma coletânea de artigos sobre a aplicação de DM em sólidos podem ser encontradas na Ref. 9.

As propriedades de materiais têm sido investigadas usando métodos que podem ser divididos em duas classes: métodos ab-initio (ou primeiros princípios) que visam obter as propriedades a partir da solução aproximada da equação de Schrödinger de muitos corpos que descreve o material; métodos empíricos ou semi-empíricos em que valores experimentais de quantidades físicas são incluídas nos cálculos.

Os métodos ab-initio são úteis para a determinação de propriedades de novos materiais com estruturas bastante complexas, mas demandam um altíssimo tempo computacional para simulações realísticas. Dentre os vários métodos ab-initio, destaca-se a teoria do funcional de densidade (TFD), desenvolvido por Kohn e colaboradores¹⁰, que consiste na solução numérica aproximada da equação de Schrödinger para N elétrons fazendo uso do fato de que o estado fundamental do sistema é uma função apenas da densidade eletrônica¹¹. Esta teoria apresenta bons resultados para constantes de rede, posições atômicas e propriedades elásticas, mas falha na previsão de largura de gaps de semicondutores e exige considerável esforço computacional. Outro método consiste em proceder iterativamente através dos seguintes passos: a partir de uma configuração de posições e velocidades atômicas, calcula-se a estrutura eletrônica para esta configuração, e a seguir, calcula-se as forças interatômicas usando o teorema de Hellmann-Feynman, move-se os átomos de acordo com as equações clássicas do movimento e as novas coordenadas atômicas são consideradas como input para um novo ciclo¹².Embora conceitualmente factível, o esquema é difícil de ser implementado devido a problemas de instabilidade dinâmica, dentre outros, embora alguns avanços tenham sido obtidos recentemente. Uma técnica mais efetiva foi desenvolvida por Car e Parinello¹³ que engenhosamente combinaram o problema da estrutura eletrônica com a dinâmica molecular dos átomos através da solução de equações dinâmicas derivadas não diretamente das equações de Newton mas a partir de equações de movimento geradas por uma Lagrangeana fictícia do sistema e através de processo de otimização dinâmica¹⁴.

Toda esta área de Mecânica Quântica aplicada à Ciência dos Materiais exige métodos numéricos sofisticados, um esforço computacional incomensurável para simulações realísticas, paralelismo em larga escala, etc¹⁵. Porém são imprescindíveis quando tratamos de materiais onde a distribuição eletrônica ou a ligação química em torno do átomo muda significativamente durante os movimentos dos átomos que são de interesse na simulação. Por outro lado, uma grande quantidade de materiais pode ser simulada e muitas de suas propriedades dinâmicas interpretadas em termos dos movimentos de seus constituintes que interagem entre si através de potenciais efetivos dependentes explicitamente de suas posições nucleares. Estes potenciais podem ser construídos com uma boa intuição física em termos das distribuições de carga e características químicas dos átomos envolvidos. Neste caso, o velho e tradicional método da Dinâmica Molecular Clássica^16,17 ainda se constitui numa poderosa ferramenta na investigação das propriedades de materiais numa extensa classe que inclui, dentre outros, os materiais estruturais cerâmicos (sólidos iônicos, óxidos, silicatos) e mesmo supercondutores de alta temperatura crítica. Em particular, a classe de cerâmicas é industrialmente importante como materiais funcionais que vão de condutores super-iônicos, eletrólitos em células combustíveis e baterias, a materiais ferroelétricos e piezoelétricos.

A partir de seu trabalho pioneiro sobre o condutor super-iônico a-AgI^18,19, Vashishta e colaboradores iniciaram um exaustivo estudo das propriedades estruturais das fases cristalinas e da fase amorfa de outros sólidos super-iônicos²⁰ e materiais óxidos binários AX₂^21,22. Nestes trabalhos o mesmo potencial de interação proposto na Ref. 18 foi usado, com algumas variações, nas simulações de DM. O modelo de potencial interatômico foi inspirado nas idéias e conceitos introduzidos por Pauling nos primórdios da Mecânica Quântica ao final da década de 20 e consiste de um potencial de dois corpos que surge devido à repulsão estereométrica dos íons, interações Coulombianas, interação de carga-dipolo e um potencial de três corpos responsável pelas interações covalentes. Este modelo foi usado ainda com sucesso na simulação de outros materiais como BaBiO₃ e seu parceiro supercondutor BaKBiO₃^23,24.Mais recentemente, as propriedades estruturais de óxidos metálicos foram estudadas pelos autores usando a mesma abordagem²⁵.

Neste artigo, dedicado ao centenário de nascimento de Pauling, objetivamos rever alguns tópicos do nosso trabalho em que a combinação de um simples e bastante efetivo potencial interatômico, construído a partir do conceito de "raios iônicos" de Pauling, com a técnica de DM produziu excelentes resultados na investigação de vários materiais e contribuiu para o avanço de nosso entendimento dos processos microscópicos que ocorrem nestes materiais.

O trabalho está dividido em 5 partes. Na primeira, discutimos como derivar um potencial entre pares baseado em idéias fenomenológicas e interações simples entre átomos. Uma descrição geral do método de DM é apresentada na 2ª parte e não pretende ser exaustiva, mas concentrada nos pontos-chave do formalismo, remetendo o leitor interessado à extensa bibliografia sobre o assunto. Na terceira, introduzimos as expressões para as diversas grandezas físicas que descrevem as propriedades estruturais e dinâmicas. Na quarta, apresentamos os resultados da simulação computacional de diversos materiais estudados e na quinta encerramos com um sumário.

CONSTRUINDO A FORMA DO POTENCIAL DE INTERAÇÃO DE ACORDO COM AS IDÉIAS DE PAULING

Na construção de um esquema para obter uma forma explícita de um potencial de interação devemos incluir as várias contribuições como as interações Coulombianas resultantes da transferência de carga, a repulsão estereométrica devido aos tamanhos atômicos, e interações de carga-dipolo e dipolo-dipolo para incluir os efeitos da polarizabilidade eletrônica dos átomos bem como as interações covalentes de três corpos.

A abordagem adotada aqui será fenomenológica no sentido de que serão omitidas considerações detalhadas de Química Quântica acerca das contribuições eletrônicas em torno dos íons em sua fase condensada e apenas conceitos físicos simples serão empregados.

O modelo de potencial de Vashishta-Rahman foi construído usando as idéias de Pauling acerca dos raios iônicos¹⁸.Na sua derivação, vamos seguir os mesmos passos empregados por Pauling na sua descrição da estrutura dos cristais iônicos²⁶.

O nosso sistema é formado por uma coleção de íons não superpostos que são mantidos juntos devido às interações de Coulomb. A palavra íon deve ser interpretada aqui no sentido de que as ligações são majoritariamente iônicas mas não necessariamente correspondem ao caso iônico extremo.

A grandes distâncias, os íons são atraídos ou repelidos pela interação de Coulomb entre as cargas Z_ie e Z_je dada por

em que i e j denotam as espécies de íons e a carga efetiva Ze é definida em termos da densidade eletrônica integrada dentro de um volume do íon. Como a escolha do volume é de algum modo arbitrária, o conceito de carga iônica não é precisamente determinado. Por exemplo, podemos com boa precisão atribuir Z=+1 para Na⁺ e Z=-1 para Cl^-, mas não devemos estabelecer Z=-2 para o íon O^-2. Este íon não se forma no espaço livre, porém em um cristal, é estabilizado devido ao potencial eletrostático dos íons positivos vizinhos. Portanto os valores das cargas dos íons são parâmetros ajustáveis no modelo.

Em adição ao potencial de Coulomb, devemos considerar uma contribuição atrativa devido à polarização de cada íon devido ao campo elétrico do outro. Este termo, que chamaremos de carga-dipolo, é bastante pequeno em comparação à repulsão ou atração Coulombiana, e em geral tem sido desprezado em vários modelos de potencial de interação. Entretanto, vários íons negativos, presentes nos sistemas que estudamos estão entre os maiores da Tabela Periódica e portanto altamente polarizáveis. Devido à presença do campo elétrico produzido por um átomo, a distribuição de carga do outro é distorcida resultando em um momento de dipolo induzido dado por onde a é a polarizabilidade eletrônica. Como conseqüência, a sua energia diminui, em primeira aproximação, de -(1/2)aE^{2 27}. A interação carga-dipolo resultante pode então ser escrita como

A pequenas distâncias, as camadas eletrônicas dos átomos começam a se sobrepor, e como decorrência do princípio de exclusão de Pauli, surge uma forte interação repulsiva que cresce drasticamente quando a distância entre íons decresce e diminui rapidamente quando r_ij aumenta. Born foi o primeiro a sugerir que este potencial repulsivo deveria decair algebricamente com a potência n (na verdade os argumentos de Born eram baseados na estabilidade do cristal). Pauling incorporou esta sugestão assumindo a forma b_ij/r_ijⁿ para este potencial, onde os parâmetros b_ij e n seriam determinados a partir dos valores experimentais da distância mínima entre os íons no equilíbrio (parâmetro de rede) e da compressibilidade do cristal. Pauling introduziu o conceito de "raios iônicos" para estimar os tamanhos dos íons e formulou um conjunto semi-empírico de raios usando como ponto de partida os dados experimentais dos parâmetros de rede de uns poucos cristais conhecidos²⁸. Uma discussão detalhada do procedimento adotado por Pauling pode ser encontrada em seu livro²⁶.A forma de potencial adotada por Pauling para descrever a repulsão estereométrica é dada por

onde s representa o raio iônico. Note que cada "contato" iônico contribui com uma energia A pois a intensidade da interação repulsiva é escalada pela soma dos raios iônicos. Portanto, como afirmou Pauling, esta é uma forma conveniente porque as forças repulsivas aumentam de intensidade com o aumento dos tamanhos dos íons. Assim, o conceito de Pauling de raios iônicos constitui um modo de expressar o comprimento da ligação química. Este conceito é usado literalmente quando se define (s_i+s_j) como a menor distância entre íons vizinhos na estrutura cristalina.

Uma outra expressão para a repulsão estereométrica bastante utilizada é devido a Born e Mayer que, a partir de cálculos das forças entre os íons baseados na Mecânica Quântica, sugeriram que a simples lei de potência r^-n deveria ser substituída por uma exponencial na forma V_R=C_ijexp(-r_ij/r _ij) ²⁹.

Finalmente, devemos considerar as interações de van der Waals (forças causadas por deformações mútuas de átomos e moléculas), introduzidas por London e Margeneau³⁰, no estudo dos cristais de gases nobres. Apesar de sua origem dever-se essencialmente a efeitos quânticos (a existência de uma energia de nível zero finita) são comumente descritas em termos da interação clássica dipolo-dipolo induzida pelas polarizações produzidas pelos íons³¹.Usando um modelo simplificado de dois osciladores (dipolos representando os átomos) acoplados pela força de Coulomb, mostra-se facilmente que a interação de van der Waals V_W é atrativa e varia com a sexta potência da distância entre eles³².É evidente que a intensidade da interação é diretamente proporcional às polpolarizabilidades. Deste modo, podemos escrever

Reunindo todas as contribuições, obtemos a forma da interação de pares

Uma vez que os s's, a's e W's são conhecidos da literatura, as demais constantes são os parâmetros ajustáveis do modelo. Para determiná-los, consideramos dados experimentais acerca da estrutura cristalina, energia de coesão e compressibilidade disponíveis em tabelas-padrão.

Na Figura 1, apresentamos o potencial de interação usado em cálculos de DM para a determinação das propriedades estruturais do TiO₂ com os parâmetros dados na Tabela 2 da Ref. 25.

A expressão mais geral para a energia de coesão de um cristal formado de N íons com parâmetro de rede R pode ser dada como uma soma de contribuições de um corpo, dois corpos, três corpos e termos de mais alta ordem:

As interações de três corpos são essenciais quando se considera as ligações covalentes em semicondutores elementais ou semicondutores e isolantes do tipo AX₂ (A = Si ou Ge e X = O, S, Se e Te)²¹. Uma forma para a interação de três corpos é

onde B_jiké a intensidade da interação de três corpos e as funções f e Q representam os efeitos de estiramento (stretching) e entortamento (bending) da ligação respectivamente. Definimos e q_jikcomo o ângulo formado por e com vértice em i. A contribuição de três corpos anula-se quando o ângulo q_jik = , sendo positivo nos outros casos. A parte angular da Eq. (7) favorece uma distribuição angular com uma dada geometria em cujo vértice encontra-se a partícula i. Uma das expressões usadas para f e Q podem ser escritas como³³

e

sendo r₀ a distância de corte para a interação de três corpos. Observe que, na Eq. (8), não há descontinuidade nas derivadas com respeito a r quando r_ij ou r_ik=r₀. Nas interações de três corpos X-A-X e A-X-A, a parte radial f é a mesma, mas a parte angular Q difere nos valores dos ângulos e . As intensidades da interação de três corpos são obviamente diferentes.

O MÉTODO DE DINÂMICA MOLECULAR

O objetivo básico da técnica de DM, assim como a de Monte Carlo, é observar a evolução do sistema dado através da determinação do movimento das partículas individuais. Devido às interações entre partículas, o sistema é capaz de manter tanto o equilíbrio mecânico quanto térmico, e no caso de perturbações externas o sistema pode atingir uma nova configuração de equilíbrio. A diferença essencial entre DM e MC está na natureza da evolução do sistema.

Em DM, calculamos a trajetória de fase do sistema que obedece à dinâmica de Newton-Hamilton, enquanto que em MC o sistema evolui de acordo com uma dinâmica estocástica, mesmo que a energia do sistema seja especificada pela mesma Hamiltoniana. A partir de uma configuração inicial , i=1...N, para um dado tempo t₀ (6N condições iniciais), as soluções das equações clássicas do movimento das N (10² - 10⁹) partículas interagindo através de potencial conhecido, tornam possível conhecer todas as configurações sucessivas t_j=jDt, com j=1,... para uma seqüência de tempos posteriores t_j. Se a Hamiltoniana do sistema for invariante por translação e rotação e for explicitamente independente do tempo então as correspondentes componentes dos momentos linear e angular, bem como a energia total, são conservadas. Esta lei de conservação de energia se aplica independentemente da existência ou não de um campo externo. A condição essencial é que a força que atua no sistema não dependa explicitamente do tempo ou da velocidade.

O procedimento tradicional na resolução das equações de Newton consiste em discretizar as equações diferenciais acopladas, ou seja, transformá-las em diferenças finitas. Partindo do pressuposto de que o potencial de interação, e portanto as forças entre partículas, são funções contínuas e diferenciáveis, dadas as condições iniciais em um instante t₀, a posição, a velocidade e qualquer outra variável dinâmica, pode ser obtida em um instante de tempo posterior t+Dt com a precisão adequada. A escolha de Dt é fundamental: não deve ser tão pequeno que o sistema não consiga evoluir (ou demore muito tempo para evoluir), e também não tão grande de modo que as constantes de movimento não se mantenham invariantes. Na maioria das vezes, Dt está relacionado com alguma freqüência característica do sistema (freqüência de Einstein, por exemplo). Existe uma variedade imensa de algoritmos, todos baseados em expansões do tipo da série de Taylor. Os mais empregados são os de Beeman, Runge-Kutta, Verlet, Gear (predictor-corrector). O compromisso entre a precisão e o armazenamento de dados na memória do computador, na maioria das vezes, é determinante na escolha do algoritmo. De modo a eliminar efeitos de superfície torna-se necessário o uso de condições periódicas de contorno. Um dos passos mais importantes e o que consume mais tempo computacional é o cálculo das forças entre as partículas. O leitor interessado deve consultar o livro de Allen e Tildesley¹⁶ para se familiarizar com procedimentos úteis (truques) na simulação.

Uma vez conhecidas as trajetórias das partículas no espaço de fase, o próximo passo consiste em obter as propriedades macroscópicas do sistema. Para isto nos baseamos na hipótese ergódica. Ou seja, a média de uma grandeza f(r^N,p^N) sobre um ensemble estatístico determinado pela função densidade r(r^N,p^N)

é equivalente à média sobre o tempo

Dado que, com a solução das equações do movimento, podemos realizar médias temporais, as grandezas macroscópicas podem assim ser calculadas. Deste modo, podemos determinar, a temperaturas finitas, a estrutura microscópica do sistema bem como as propriedades termodinâmicas. Com relação à sua estrutura, podemos obter, por exemplo, a ordem de curto alcance, a distância entre primeiros vizinhos, comprimento de ligação, número de coordenação e distribuição de ângulos de ligação. Quanto à sua dinâmica, podemos determinar o coeficiente de difusão, o espectro de estados vibracionais, o fator de estrutura dinâmico. Com respeito à sua termodinâmica, a capacidade térmica, a energia livre, a entropia e pressão, entre outras grandezas, podem ser explicitamente calculadas.

FUNÇÕES DE CORRELAÇÃO

Estruturais

As correlações de dois corpos são determinadas em termos das funções de correlações de pares parciais g_ab(r) definidas através de

onde é o número de partículas da espécie b dentro de uma casca esférica entre r e r+Dr em torno de uma partícula da espécie representa tanto a média de ensemble como a média sobre todas as partículas da espécie b. Deste modo, g_ab(r) representa a probabilidade de se encontrar uma partícula do tipo b a uma distância r de uma dada partícula do tipo a. A densidade total é r=N/W com N=N_a+N_b, a concentração c_b=N_b/N, e W é o volume do sistema. A função correlação dos pares é definida como

O número de coordenação C_ab(r) definido como o número de partículas b, vizinhas mais próximas de uma partícula a, é obtido calculando-se a área sob a função de distribuição de pares parciais correspondente, e é dado por

As funções fator de estrutura densidade-densidade e de raios-X podem ser facilmente obtidas a partir de fatores de estrutura parciais, obtidos a partir das transformadas de Fourier das correspondentes funções de correlações de pares parciais, ou seja,

A função janela, W(r/L)= sen(2p r/L)/(2pr/L), onde L é o comprimento da aresta do cubo de simulação, é usualmente introduzida para reduzir o efeito de truncar a integral com um limite superior finito. Uma vez obtidos S_ab(q), o fator de estrutura é dado por

e o fator de estrutura estático de raios-X por

sendo b_a os fatores de forma de raios-X. A função de estrutura estático para espalhamento por nêutrons pode ser obtida de forma análoga a partir dos fatores de estrutura parciais ponderados por comprimentos do espalhamento coerente de nêutrons. Com o conhecimento das funções de correlação parciais, outros fatores de estrutura podem ser definidos e facilmente obtidos²¹.

Dinâmicas

O conhecimento do espaço de fase permite-nos obter além das propriedades estruturais do sistema, como descrito acima, todas as propriedades dinâmicas que se queira. Como exemplo discutiremos agora algumas destas funções dinâmicas.

A função de auto-correlação de velocidades (VAF) normalizada é definida como¹⁶

sendoa velocidade da partícula i do tipo a no tempo t. A densidade de estados vibracionais de fônons, G(w), pode ser obtida através da transformada de Fourier da função de auto-correlação de velocidades, ou seja,

O valor finito de G(w), para w =0, está relacionado com a constante de auto-difusão dada por

As trajetórias de fase permitem-nos calcular também a função de auto-correlação corrente-corrente, cuja transformada de Fourier fornece a condutividade iônica dependente da freqüência s(w), ou seja,

sendo a densidade de corrente definida por .

A relação de Einstein permite, então, calcular o coeficiente de difusão através de

Além das grandezas dinâmicas acima definidas, outras funções dependentes do tempo podem ser obtidas, tais como a condutividade térmica e a viscosidade de cisalhamento, a partir do nosso conhecimento das trajetórias no espaço de fase.

MATERIAIS ESTUDADOS

Nesta seção, apresentamos sucintamente alguns resultados de nosso estudo de materiais usando simulação DM.

Super-Iônicos

Materiais com alta condutividade iônica, chamados condutores de íons rápidos ou super-iônicos, consistem de duas "sub-redes" interpenetrantes. Sais de prata, fluoretos e b-aluminas fazem parte de uma grande classe de materiais super-iônicos, sendo AgI, Ag₃SI, Ag₂S e Ag₂Se os mais estudados. Nestes materiais, os anions na fase super-iônica formam uma rede cúbica de corpo centrado (ccc) estável mas altamente anarmônica, por onde os cátions se difundem.

As fases super-iônica e líquida do Ag₂Se, usando a técnica de DM, foram investigadas na Ref. 20. Na Figura 2, o fator de estrutura estático, calculado via DM, é mostrado para o Ag₂Se na fase líquida. Observe que os resultados de DM estão em boa concordância com as experiências de espalhamento por nêutrons. A partir dos fatores de estrutura parciais, foi possível determinar o papel das correlações inter-íons que dão origem aos vários picos que aparecem nas experiências.

A Figura 3 mostra a dependência da constante de auto-difusão, dada pela Eq. (20), para o íon Ag⁺ na fase super-iônica do Ag₂Se. O quociente de Haven, definido como H(T)=D(T)/D_a(w=0,T), pode ser determinado e encontramos, para T=804 K, o valor 0,58 enquanto que o resultado experimental encontra-se na faixa de [0,42 - 0,50].

Vidros Calcogênicos

Dos vidros calcogênicos, sem dúvida, o dióxido de silício tem sido o material mais estudado pelo físicos da matéria condensada, químicos, cientistas de materiais e engenheiros. O SiO₂ é encontrado com cerca de 40 diferentes estruturas cristalinas, embora as mais estáveis, para o cristal quimicamente puro, sejam cristobalita, quartzo, coesita, stishovita e quiatita. É natural que tal riqueza de estruturas permaneça nas fases líquidas e vítreas do SiO₂ a diferentes temperaturas e pressões. As estruturas destas fases têm sido investigadas usando as técnicas de difração de raios-X e de nêutrons (veja Ref. 21 para uma relação da literatura pertinente).

O potencial de interação, dado pela Eq. (5), juntamente com o método de DM foi utilizado também para estudar as correlações estruturais do SiO₂, bem como de outros vidros semicondutores e isolantes do tipo AX₂, onde A denota Si ou Ge e X representa O, S, Se, ou Te.

A Figura 4 mostra um diagrama esquemático dos processos de resfriamento (quenching) e tratamento térmico usados nas simulações de DM para preparação do estado vítreo a partir do estado líquido.

A Figura 5 representa os resultados de DM do fator de estrutura estático e sua comparação com resultados experimentais de difração de nêutrons. Como vemos, a concordância das alturas e larguras dos picos com a experiência é excelente. Os resultados para os fatores de estrutura parciais tornam evidente que o segundo pico envolve contribuições das correlações da ligação Si-Si e O-O com um cancelamento parcial das anti-correlações de Si-O, enquanto o terceiro e quarto surgem de contribuições de Si-Si, Si-O e O-O.

Na Tabela 1, apresentamos os resultados de DM para os comprimentos e ângulos de ligação bem como as respectivas larguras totais a meia altura (LTMA) para as fases líquidas e amorfas do SiO₂.

Outros aspectos da estrutura dos estados vítreo e líquido do SiO₂, a partir da estatística e energética dos anéis de conectividade, foram ainda detalhadamente discutidos na Ref. 34. As propriedades dinâmicas e o processo de densificação foram analisados nas Refs. 35 e 36.

Óxidos Metálicos

O dióxido de titânio é outro material exaustivamente estudado ao longo dos tempos devido principalmente a suas aplicações tecnológicas, tendo sido empregado em catálise, eletrocromagem e sensores. Existe na forma cristalina em três fases estruturais, rutilo, anatásio e broquita. A sua estrutura mais estável é o rutilo, que exibe propriedades de isolante (gap 3,0 - 3,3 eV), mas pode se tornar um semicondutor quando existe deficiência de oxigênio. O cristal de rutilo é o protótipo de muitos óxidos metálicos com coordenação octaédrica. Embora não exista na fase vítrea por si próprio³⁷, TiO₂ torna-se um vidro quando combinado com alguns outros óxidos (por exemplo SiO₂). Uma forma amorfa instável pode ser obtida por métodos de sol-gel.

A partir de simulações de DM com o potencial acima proposto, investigamos as propriedades estruturais do TiO₂. Na Figura 6, mostramos o comportamento da energia total por número de partículas do sistema, evidenciando na região hachurada, a transição cristal-líquido. Experimentalmente, esta transição foi observada a T=2128 K.

Através da análise dos fatores de estrutura, estudamos o papel das correlações estruturais de curto e intermediário alcance que são responsáveis pelo primeiro pico de difração em experiência de espalhamento de raios-X e de nêutrons. As distribuições angulares das ligações estão mostradas na Figura 7. Vemos que a fase amorfa surge a partir das distorções dos octaedros Ti(O_1/3) tendo em vista o aparecimento de picos bem definidos em 90º e 120º na distribuição angular Ti-O-Ti, 60º na distribuição O-O-O, e também porque cada íon Ti possui, em média, 5,4 átomos de O vizinhos próximos. A conectividade entre octaedros, descritos pelos ângulos da ligação Ti-O-Ti, parecem com a da estrutura do rutilo.

CONCLUSÃO

Neste artigo, usando um potencial de interação construído a partir das idéias seminais de Pauling, apresentamos resultados acerca das propriedades estruturais e dinâmicas de materiais de diferentes classes obtidos através de simulações de DM. Apesar dos resultados aqui exibidos terem sido limitados àqueles obtidos pelos autores, existe uma vasta literatura concernente a estudos em outros materiais realizados por outros grupos.

Nossa mensagem final é que embora as recentes abordagens baseadas em métodos da Mecânica Quântica Aplicada e avançadas técnicas computacionais tenham alcançado grande sucesso, existe ainda espaço para estudos que se baseiam em potenciais simples e que permitem obter resultados qualitativos mais rapidamente e obter tendências para novas simulações em maior escala, bem como para a investigação de propriedades gerais de sistemas complexos que são intratáveis por métodos mais rigorosos e que não dependem de detalhes finos da energética do sistema.

AGRADECIMENTOS

Os autores agradecem aos Profs. P. Vashishta e R. Kalia por tê-los ensinado a técnica de DM e pela sua colaboração na maioria dos trabalhos aqui apresentados. Somos ainda gratos pela hospitalidade e camaradagem durante a nossa permanência tanto no Laboratório Nacional de Argonne quanto na Universidade de Louisiana. Agradecemos à FAPESP pelo financiamento da infra-estrutura computacional que permitiu a realização dos cálculos descritos na parte de Óxidos Metálicos e ao CNPq pelas bolsas de Produtividade em Pesquisa. O nosso muito obrigado a André R. Studart e ao Prof. Antonio José Roque da Silva pela leitura crítica deste ensaio.

REFERÊNCIAS

1. Veja o Relatório Final do DOE-NSF National Workshop on Advanced Scientific Computing em http://www.er.doe.gov/ssi/LangerReport.pdf

2. Ceperley, D. M.; Rev. Mod. Phys. 1999, 71, S438.

3. Langer, J. S.; Phys. Today 1999, July, 11.

4. Bernholc, J.; Phys. Today 1999, September, 30.

5. Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E.; J. Chem Phys. 1953, 21, 1087.

6. Alder, B. J.; Wainwrigth, T. E.; J. Chem. Phys. 1957, 27, 1207.

7. Gibson, J. B.; Goland, A. N.; Milgram, M.; Vineyard, G. H.; Phys. Rev. 1960, 120, 1229.

8. Rahman, A.; Phys. Rev. 1964, 136, A405.

9. Ronchetti, M. e Jacucci, G. in Simulation Approach to Solids, Kluwer; Dordrecht, 1990.

10. Hohenberg, P. ; Kohn, W.; Phys. Rev. 1964, 136, B864; Kohn, W. ; Sham, L.; Phys. Rev. 1965, 140, A1133.

11. Dreizler, R. M. ; Gross, E. K. U., in Density Functional Theory; Springer-Verlag, Berlim 1990.

12. Payne, M. C.; Allan, D. C.; Teter, M. P.; Arias, T. A.; Joannopoulos, J. D.; Rev. Mod. Phys. 1992, 64, 1045.

13. Carr, R. ; Parinello, M.; Phys. Rev. Lett. 1985, 55, 2471.

14. Remler, D. K. ; Madden, P. A.; Molec. Phys. 1990, 70, 921.

15. Para uma recente revisão da aplicação de métodos computacionais à Ciência dos Materiais, veja Ohno, K.; Esfarjani, K.; Kawazoe, Y.; Computational Materials Science; Springer; Berlim 1999.

16. Allen, M. P. ; Tildesley, D. J.; Computer Simulations of Liquids; Oxford Univ. Press; London 1989.

17. Ciccotti, G.; Frenkel, D.; McDonald, I. R.; Simulations of Liquids and Solids; North-Holland; Amsterdam 1987.

18. Vashishta, P. ; Rahman, A.; Phys. Rev. Lett. 1978, 40, 1337.

19. Rahman, A. ; Vashishta, P. in The Physics of Superionic Conductors and Electrode Materials; Perram, W., Ed.; Plenum; New York, 1983; p83.

20. Rino, J. P.; Hornos, Y. M. M.; Antonio, G. A.; Ebbsjö, I.; Kalia, R. K. ; Vashishta, P.; J. Chem. Phys. 1988, 89, 754.

21. Vashishta, P.; Kalia, R. K.; Rino, J. P.; Ebbsjö, I.; Phys. Rev. B 1990, 41, 12197.

22. Nakano, A.; Kalia, R. K.; Vashishta, P.; J. Non-Cryst. Solids 1994, 17, 157.

23. Loong, C. K.; Hinks, D. G.; Vashishta, P.; Jin, W.; Kalia, R. K.; Degani, M. H.; Price, D. L.; Jorgensen, J. D.; Dabrowski, B.; Mitchel, A. W.; Richards, D. R.; Zheng, Y.; Phys. Rev. Lett. 1991, 66, 3217.

24. Degani, M. H.; Braz. J. Phys. 1994, 24, 956.

25. Rino, J. P.; Studart, N.; Phys. Rev. B 1999, 59, 6643.

26. Pauling, L.; The Nature of the Chemical Bond; 3^rd Edition; Cornell University Press; New York, 1960.

27. Veja, por exemplo, Morrison, M. A.; Estle, T. L.; Lane, N. F.; Quantum States of Atoms, Molecules, and Solids; Prentice-Hall; New Jersey, 1976.

28. Pauling, L.; J. Am. Chem. Soc. 1927, 49, 765.

29. Born, M. ; Mayer, J. E.; Z. Physik 1932, 75, 1; Mayer, J. E.; J. Chem. Phys. 1933, 1, 270.

30. London, F.; Z. Physik 1930, 63, 245; Margeneau, H.; Phys. Rev. 1931, 38, 747; Margeneau, H.; Rev. Mod. Phys. 1939, 11, 1.

31. Kittel, C.; Introdução à Física do Estado Sólido; 5a. Edição; Guanabara Dois; Rio de Janeiro, 1978.

32. Born, M.; Física Atômica; 4a. Edição; Fundação Calouste Gulbenkian; Lisboa, 1986; Apêndice XL.

33. Stillinger, F. H. ; Weber, T. A.; Phys. Rev. B 1985, 31, 5262.

34. Rino, J. P.; Ebbsjö, I.; Kalia, R. K.; Nakano, A.; Vashishta, P.; Phys. Rev. B 1993, 47, 3053.

35. Jin, W.; Kalia, R. K.; Vashishta, P.; Rino, J. P.; Phys. Rev. Lett. 1993, 71, 3146; ibid., Phys. Rev. B 1993, 48, 9359.

36. Jin, W.; Kalia, R. K.; Vashishta, P.; Rino, J. P.; Phys. Rev. B 1994, 50, 118.

37. Zachariasen, W. H.; J. Am. Chem. Soc. 1932, 54, 3841.

Datas de Publicação

• Publicação nesta coleção
  08 Fev 2002
• Data do Fascículo
  Dez 2001

Histórico

• Recebido
  04 Dez 2000
• Aceito
  28 Mar 2001