FÍSICA DE LA TIERRA, AGUA Y ATMÓSFERA. FÍSICA DEL MEDIO AMBIENTE

Aplicación de simulación con dinamica molecular al cálculo de energía relativa de borde de grano en hielo

Aplication of molecular dinamics simulation on relative grain boundary energy in ice

 

Di Prinzio C.L.^a,b,*, Pereyra R^a,b., Druetta E^a. y Nasello O.B^a,b.

a-FAMAF (Facultad de Matemática Astronomía y Física. Universidad Nacional de Córdoba) Medina Allende y Haya de la Torre. Ciudad Universitaria (5000)
b-CONICET (Consejo Nacional de Ciencia y Tecnología

 

---

Resumen

En este trabajo se calculó la energía de los bordes de grano en hielo y la energía superficial del hielo mediante simulaciones de dinámica molecular en muestras bicristalinas con bordes de grano simétricos en función de la desorientación cristalina de los mismos. La energía relativa de borde de grano se obtuvo como el cociente entre la energía de los bordes de grano en hielo y la energía superficial del hielo de la misma cara expuesta. Las muestras bicristalinas de hielo creadas para las simulaciones tenían una desorientación alrededor del eje y el ángulo de desorientación 0 se barrió cada 10^o en el intervalo [0^o, 180°]. Los resultados fueron comparados satisfactoriamente con resultados experimentales de energía relativa de borde de grano simétricos en hielo puro obtenidos mediante mediciones de replicas plásticas usando microscopía confocal.

Palabras clave: Hielo; Dinámica Molecular; Gromacs; Borde de grano; Energía superficial.

Abstract

In this work the energies of the ice grain boundary and of the ice surface were calculated by molecular dynamics simulations for symmetric grain boundaries with different crystalline disorientations.. The relative grain boundary energy is obtained as the ratio of the energies of the ice grain boundary and of the ice surface. The ice samples used in the simulations had disorientation around the axis and the disorientation angle 0 in the range [0°, 180°]. The results were compared successfully with experimental relative energies measured in pure ice symmetrical grain boundaries by confocal microscopy, using plastic replicas.

Keywords: Ice; Molecular Dynamics; Gromacs; Grain boundary; Surface energy.

Recibido: 08/10/14;
aceptado: 28/10/15

---

 

I. INTRODUCCIÓN

En el hielo polar se encuentran registros de la temperatura y de los contaminantes atmosféricos de épocas pasadas^[1-2]. Los contaminantes se encuentran en los cristales de hielo y las interfaces entre ellos, llamadas bordes de grano (BG). Los contaminantes pueden moverse en el hielo más fácilmente por los BG o por las uniones triples, denominados caminos de alta difusión (CAD), y en consecuencia su datación podría ser incorrecta^[3].
Los CAD están altamente relacionados con la estructura y energía de los BG. Existen muy pocos trabajos en la literatura que estudien la energía de BG en hielo Ih (hielo de estructura cristalina hexagonal) en forma experimental, teórica o computacional.
En los trabajos de Kechman y Hobbs (1974)^[4] y Suzuki y Koroiwa (1979)^[5] se midió la energía relativa del BG en hielo en función de la desorientación cristalina. En dichos trabajos se comparan los resultados hallados con la teoría de dislocaciones de BG solamente para BG de bajo ángulo.
Recientemente Di Prinzio y col (2014)^[6-7] publicaron datos experimentales de la energía relativa de BG en muestras de hielo puro y contaminado para diferentes concentraciones y temperaturas. Las mediciones de la energía relativa del BG fueron realizadas usando el balance entre la energía superficial del hielo y la energía del BG en la intersección del BG con la superficie libre de la muestra. Para ello se midió el ángulo formado entre el BG y la superficie libre del hielo usando una imagen en tres dimensiones (3D) obtenida mediante un microscopio confocal. Los resultados demostraron que la energía de los BG depende de la orientación cristalina y la inclinación de los mismos. Además, se determinó que los BG simétricos poseen energías relativas más baja que los no simétricos.
La estructura de los BG puede, en un primer paso, ser estudiada mediante la red de sitios de coincidencia (CSL). Para ciertas desorientaciones cristalinas entre cristales existen super redes, formadas por los sitios coincidentes de ambas redes, que, cuando son de alta densidad, pueden indicar inclinaciones donde los BG podrían tener baja energía. En estas inclinaciones los BG tendrían una estructura menos desordenada al tener más sitios coincidentes y por lo tanto su energía sería más baja.
En hielo Ih, las redes de sitios de coincidencias (CSL) fueron obtenidas por Gonzalez et al. (1997)^[8]. En este trabajo se determinaron las orientaciones cristalinas "especiales" con mayor densidad de sitios coincidencias Σ y para estas orientaciones se determinaron las inclinaciones con mayor densidad planar de sitios de coincidencias Γ. Di Prinzio y col^[6-7] encontraron que, en general, había una alta correlación entre 1/Γ y la energía relativa de los BG simétricos con desorientación 0 alrededor del eje . Kobayashi y col^[9-10] encontraron que las interfaces entre cristales de nieve se encuentran mas fácilmente entre algunas orientaciones "especiales" de acuerdo a la teoría de CSL. Además, Hondoh T. (1979)^[11] encontró algunas propiedades "especiales" en la migración de BG en hielo las cuales se producían en BG con alta densidad de sitios de coincidencia.
En otros materiales con estructura hexagonal se ha estudiado la estructura de BG usando dinámica molecular (DM) y los resultados obtenidos muestran una buena concordancia con los resultados experimentales hallados por diferentes técnicas.^[12-18]
Vemos en otros materiales cómo se ha avanzado enormemente en la interpretación de la estructura y energía de los BG y en la correlación entre resultados experimentales y computacionales pero no así en hielo Ih. Es por esa razón que el espíritu de este trabajo es comenzar a calcular mediante DM, las estructuras y las energías de los BG simétricos en hielo Ih y compararlas con las obtenidas experimentalmente en trabajos anteriores. Por lo tanto se darán los primeros pasos en la construcción de estructuras de BG simétricos en hielo con desorientaciones cristalinas obtenidas por rotaciones alrededor del eje prismático .

II. MÉTODOS

Preparación del sistema de hielo Ih:
Inicialmente se construye un cristal de hielo Ih con forma de un paralelepípedo rectangular con las moléculas de agua respetando las reglas de Bernal-Fowler según el algoritmo de Buch^[19] y con condiciones periódicas en todas sus caras.
El plano prismático perpendicular al eje (b) es paralelo a una de las caras (ver figura 1) y coincide con el plano xz. El eje c <0001> forma con el eje x un ángulo θ/2.

Figura 1: Configuración del bicristal estudiado y del BG en relación al sistema de coordenadas xyz.

El número de moléculas del cristal utilizado varió entre aproximadamente 700 y 2200 para las diferentes inclinaciones del eje c. El cristal inicial fue llevado a un estado de equilibrio mediante DM GROMACS [20]. Este proceso llamado "termalización" se hizo a una temperatura de 255 ^oK y a 1 bar de presión durante 1ns. En la próxima sección se dan detalles de la simulación de DM GROMACS^[20] utilizada.
El cristal "termalizado" se replica mediante una reflexión en el plano xy. De esta manera, los dos cristales así formados son puestos en contacto (figura 1) y un bicristal de hielo es formado con un ángulo de desorientación θ alrededor del eje de rotación (b). Elángulo θ se varió entre 10^o y 180^o. Este tipo de bicristal tiene un BG tipo tilt y simétrico denominándose BG tilt θ.

Detalles computacionales de las simulaciones de DM:
Para la relajación del cristal se usa el paquete de simulación de GROMACS v.4.5.5^[20]. Las moléculas de agua interactúan de acuerdo al potencial TIP5P-Ew^[21]. Este modelo de agua funciona muy bien en las simulaciones de hielo, principalmente debido a que su arquitectura presenta dos sitios que representan los pares solitarios del agua, que contribuyen a formar más eficazmente los enlaces de hidrógeno. La simulación fue realizada en una caja con condiciones periódicas de contorno manteniendo los parámetros N, T y P constantes.
El algoritmo de leapfrog fue usado para la integración de las ecuaciones dinámicas con un paso de 0,001ps. Un corte esférico de 0,9 nm fue impuesto al Potencial de Lenard-Jones y a las interacciones electrostáticas de corto alcance. Para las interacciones de largo alcance hemos incluido correcciones usando la aproximación de PME.
La presión y la temperatura del sistema fueron controladas usando el termostato de Nose Hoover y el barómetro de Parrinello Rahman. Las constantes de tiempo en el barómetro y en termostato fueron de 0,5 ps y la compresibilidad fue uniforme e igual a 4,510⁵ bar^-1

Cálculo de la energía interfacial E_I:
La muestra bicristalina obtenida de acuerdo a lo explicado anteriormente fue tratada con DM a fin de llegar a un balance energético. El tiempo de simulación variaba entre 5 a 10 ns de acuerdo a la orientación de los cristales y a su tamaño.
En la figura 2 se presenta un BG tilt 60°en la etapa inicial (a) y después de 2ns (b).

Figura 2: BG tilt 60° / en la etapa inicial (a) y después de 2ns (b)

Inicialmente el BG presenta zonas casi libre de moléculas de agua pero con el correr de la simulación esas zonas son cubiertas por moléculas de agua desordenadas con respecto a las redes cristalinas adyacentes. Hay sin embargo una zona central donde las moléculas se ven siguiendo aproximadamente el orden de las redes adyacentes. Esta disposición de las moléculas podría estar relacionada con alguna estructura periódica de moléculas en el BG pero a esta escala es muy difícil de asegurar.
En la figura 3 se muestra la evolución de la energía total E_T (cinética mas potencial) del sistema en función del tiempo de simulación. Como podemos ver la energía E_T va disminuyendo hasta hacerse en promedio aproximadamente constante. La energía total E_T que se toma es el promedio en el último nanosegundo.

Figura 3: Evolución de la energía total E_T en función del tiempo de simulación.

La energía interfacial E_I se calcula mediante la siguiente ecuación^[22-25]:

Donde ET es la energía total calculada en la simulación en kJ/mol, e_ice e es la energía total por unidad de molécula en el hielo perfecto a la misma temperatura y presión (=45,42 kJ/mol a T= 255 ^oK y p=1 bar), n número de moléculas de cada bloque de hielo, N_A número de Avogadro, d_x y d_y son las dimensiones de la caja paralelas al BG expresadas en nm. El primer factor 2 (en el numerador) es debido a los dos bloques de hielo enfrentados y el segundo factor 2 (en el denominador) es debido a que la muestra estudiada posee dos BG.
En la tabla 1 se presentan los valores de E_T, n, dx, dy y E_I para las 17 orientaciones cristalinas estudiadas.

Tabla 1: Valores de E_T, n, dx,dy y E_I para las 17 orientaciones cristalinas estudiadas

Cálculo de la energía superficie libre ES:
La energía libre de una determinada superficie de hielo en equilibrio con su vapor (ES) se calculó usando una ecuación similar a la ec. (1). Las superficies estudiadas fueron los planos cristalinos que formaban los bordes de grano estudiados anteriormente. En la figura 4 se ilustra la configuración final de una de estas muestras donde se puede observar que la superficie libre presenta un desorden propio de una capa cuasi liquida. En la Tabla 2 se presentan los valores de la energía superficial ES obtenida para las 17 orientaciones cristalinas estudiadas.

Figura 4: Configuración final del plano del BG correspondiente al BG tilt 60°

Tabla 2: Valores de ES y ʎ_IS = E_I/ Es obtenidos para las 17 orientaciones cristalinas estudiadas.

Con los valores de E_I y ES se obtuvo la energía relativa γ_is = E_I/E_s correspondiente a cada una de las 17 orientaciones estudiadas.

III-RESULTADOS OBTENIDOS

Di Prinzio y col (2014)^[6-7] determinaron experimentalmente el valor de γ_is para BG simétricos y no simétricos. En la figura 5 se presentan solo los valores de γ_is para BG simétricos y los obtenidos en este trabajo mediante DM. Se agregaron además los valores 1/Γ para las desorientaciones "especiales" obtenidas en^[8] para la CSL del hielo Ih.

Figura 5: valores de γ_is obtenidos en este trabajo mediante DM y por Di Prinzio et al(2014) [6,7] y valores de 1/Γ.

Podemos ver que los 17 valores de γ_is obtenidos mediante la simulación están distribuidos cada 10° entre 10° y 170°. Los valores experimentales de γ_is son menos y solo han sido medidos alrededor de 40°, 60°, 70°, 140° y 160°. En general se nota que los valores experimentales y simulados tiene el mismo orden de magnitud y son menores que la unidad.
Un hecho notable se da entre 30° y 70° donde los valores experimentales presentan al igual que los simulados un mínimo aproximadamente en 60°. De la misma manera se observan dos máximos relativos en 45° y en 70°. Esto además está en concordancia con el valor de 1/Γ indicando que el BG tilt /60° está ubicado en un plano de mas alta densidad de sitios de coincidencia que el BG tilt / 70° y el BG tilt / 45°. En la simulación con DM también se refleja que el BG tilt / 60° tiene menos desorden que los otros BG vecinos.
Los valores obtenidos experimentalmente para 140° y 160° son comparables con los valores simulados de acuerdo a los errores experimentales obtenidos. Podemos ver que en estos ángulos los valores de energía relativa de los BG no tienen mayores diferencias con los valores de 1/Γ Tambien puede verse una gran correspondencia entre los valores de 1/Γ y los valores ʎ_IS simulados entre 90° y 140° donde no hay mediciones experimentales.
Vemos que es necesario obtener más datos experimentales de ʎ_IS para poder comparar más en detalle los mismos con las simulaciones de DM. Sin embargo la obtención de datos experimentales es muy compleja y requiere una muy buena puesta a punto del método experimental desarrollado en los trabajos previos.

IV-CONCLUSIONES

En este trabajo se usó DM para obtener valores de la energía relativa de los BG (ʎ_IS ) en muestras bicristalina de hielo Ih puro. Los resultados fueron comparados con valores de ʎ_IS obtenidos experimentalmente por Di Prinzio y col (2014)^[6,7].

Podemos en general concluir que:

a) Los valores experimentales y simulados de ʎ_IS obtenidos para todas las desorientaciones son en general menores que 1.
b) Los valores experimentales y simulados de ʎ_IS presentan un mínimo relativo alrededor de 60o y dos máximos locales en 45° y 70°. Este hecho demuestra que los BG con alta densidad de sitios de coincidencia (Γ) son "especiales" y el valor de 1/Γ es directamente proporcional a ʎ_IS. Esta coincidencia entre un modelo geométrico como la CSL y los resultados de DM muestran que las moléculas de hielo tienen muy poca libertad de relajase en las interfaces y casi conservan a baja temperatura sus posiciones originales.
c) Los valores simulados y experimentales son comparables para las desorientaciones de 140° y 160°.
d) Entre 70° y 140° la similitud entre el modelo estático de la CSL y el valor deʎ_IS obtenido por DM es sorprendente lo cual demuestra la poca relajación que sufren las moléculas de agua en las interfaces debido a al carácter del enlace puente hidrogeno.

Agradecimientos:

Este trabajo fue subsidiado por CONICET (PIP 11220090100537) y SECYT-UNC ("Estudio de las interfaces en hielo: aplicación a hielo natural" 2014). Los autores agradecen el apoyo técnico del señor José Barcelona.

V. REFERENCIAS

1. Alley, R. B., Perepezko J. H., and. Bentley C. R. Grain Growth in Polar Ice: I. Theory, J. Glaciol., 32(112), 415-424, (1986a).         [ Links ]

2. Alley, R. B.,. Perepezko J. H, and Bentley C. R. Grain growth in polar ice: II. application, J. Glaciol., 32(112), 425- 433, (1986b).         [ Links ]

3. Barnes, P.R.F., Wolff E.W., Mader H.M.,. Udisti R, Castellano E. and Rothlisberger R. The evolution of chemical peak shapes in the Dome C, Antarctica, ice core, J. Geophys. Res., 108(D3), 4126, (2003)         [ Links ]

4. Kelcham, W. M.; Hobbs, P. V. An experimental determination of the surface energies of ice, Philosophical Magazine 19,162, 1161-73. (1969),         [ Links ]

5. Suzuki S. Y.; Kuroiwa D. Grain-Boundary Energy and Grain-Boundary Groove Angles in Ice, Journal of Glaciology, 11, 265-277. (1972),         [ Links ]

6. Di Prinzio C. L.; Druetta, E; Nasello, O. B. Experimental Determination of... Tilt Grain Boundary Energies in Ice, Journal of Materials Science Research vol. 3 p. 69 -77, (2014)         [ Links ]

7. Di Prinzio C. L.; Druetta, E; Nasello, O. B., Journal of Physical Chemistry B, in press. (2014)         [ Links ]

8. B. J. Gonzalez Kriegel, C. L. Di Prinzio and O. B. Nasello. Exact Coincidence Site Lattice in Ice Ih, J. Phys. Chem. B, 101 (32), pp 6243-6246, (1997),         [ Links ]

9. Kobashashi T.; Furukawa Y. On twelve-branched snow crystals, Journal of Crystal Growth, 28, 21-28, (1975)         [ Links ]

10. Kobashashi T.; Furukawa Y. Epitaxial relationships during the formation of three-dimensional snow dendrites, Journal of Crystal Growth, 45, 48-56, (1978)         [ Links ]

11. Hondoh T.; Higashi A. Anisotopy of migration and faceting of large-angle grain boundaries in ice bicrystals, Phil. Magazine A 39, 137-149. (1979),         [ Links ]

12. Bullough R. The Dislocation Content of a Large Angle Tilt Boundary, Phil. Mag. 12 1139, (1965)         [ Links ]

13. Hirth J P, Pond R C and Lothe. Diconnections in Tilt Walls, J Acta Mater. 54 4237, (2006)         [ Links ]

14. Sutton A P and Balluffi R W 1995 Interface in Crystalline Materials (Oxford: Oxford University Press)         [ Links ]

15. Hirth J P and Pond R C. Contability and accommodation in displacive phase transformations, Progress. Mater. Sci. 56 586, (2011)         [ Links ]

16. Wang J, Woo C H and Huang H. Destabilization of dislocation dipole at high velocity, Appl. Phys. Lett. 79 3621, (2001)         [ Links ]

17. Wang J and Huang H C Novel deformation mechanism of twinned nanowires, Appl. Phys. Lett. 88 203112, (2006)         [ Links ]

18. Liu X Y, Adams J B, Ercolessi F and Moriarty J A. EAM potential for magnesium from quantum mechanical forces, Modelling Simul. Mater. Sci. Eng. 4 293, (1996).         [ Links ]

19. Wang Y, Chen L Q, Liu Z K and Mathaudhu S N , Scr. Mater. 62 646, (2010)         [ Links ]

20. Zope R R and Mishin Y. Interatomic potentials for atomistic simulations of the Ti-Al system, Phys. Rev. B 68 024102, (2003)         [ Links ]

21. Buch V.; Sandler P.; Sadlej J. Simulations of H₂O Solid, Liquid, and Clusters, with an Emphasis on Ferroelectric Ordering Transition in Hexagonal Ice, J. Phys. Chem. B, 102, 8641-8653, (1998)         [ Links ]

22. Chen, F., & Smith, P.E. Simulated surface tensions of common water models, Journal of chemical physics, 126(22), 3. (2007).         [ Links ]

23. Ikeda, K, Yamada, K, Takata, N, Yoshida, Nakashima,H, Tsuji,N. Grain Boundary Structure of Ultrafine Grained Pure Copper Fabricated by Accumulative Roll Bonding, Materials Transactions l49, 1,24-30. (2008)         [ Links ]

24. González-Romero R,L., Meléndez J. J., Gómez-García D, Cumbrera F. L., Domínguez-Rodríguez A. A molecular dynamics study of grain boundareis in YSZ: structure, energetics and diffusion of oxygen, Solid State Ionics 219, 1-10. (2012)         [ Links ]

25. Zhang J, Zhao J, Lu J. Intrinsic strength and failure behaviors of graphene grain boundaries. ACS Nano 6, 3, 2704-11. (2012)         [ Links ]