https://doi.org/10.31527/analesafa.2020.31.1.7
FÍSICA DE LA ATMÓSFERA, LA TIERRA Y EL AGUA. FÍSICA DEL MEDIO AMBIENTE

Movimiento de un borde de grano en láminas delgadas usando método de Monte Carlo

 

Grain boundary migration in thin films using Monte Carlo method

 

C. L. Di Prinzio^1,2
P. I. Achával¹
D. Stoler¹
G. Aguirre Varela^*1,2

¹ FAMAF (Facultad de Matemática, Astronomía, Física y Computación), Universidad Nacional de Córdoba, Medina Allende y Haya de la Torre, (5000) Ciudad Universitaria, Córdoba, Argentina.
² IFEG-CONICET (Instituto de Física "Enrique Gaviola", Universidad Nacional de Córdoba, Medina Allende y Haya de la Torre, (5000) Ciudad Universitaria, Córdoba, Argentina.

Autor para correspondencia: C.L. Di Prinzio, email: carlosdiprinzio@gmail.com; pachaval@famaf.unc.edu.ar

Recibido: 21/11/2019
Aceptado: 29/03/2020

 

---

Resumen:
En este trabajo se estudió la evolución de un borde de grano plano en una muestra delgada mediante un algoritmo basado en el método de Monte Carlo. El borde de grano es impulsado por una fuerza externa y los efectos superficiales sobre su movimiento son estudiados. El movimiento de la parte del borde de grano en la superficie es espasmódico, lo cual significa que tiene periodos alternados de movimiento y de estancamiento. Los periodos de estancamientos son inversamente proporcionales al espesor de la muestra. Los resultados obtenidos computacionalmente se coinciden satisfactoriamente con los resultados teóricos y experimentales obtenidos por diferentes autores.

Palabras clave: láminas delgadas, borde de grano, difusión gaseosa.

Abstract:
This paper presents the evolution of a flat grain boundary in a thin sample, using a numerical algorithm based on the Monte Carlo method. The grain boundary is driven by an external force and the effect of the free surface is studied. The grain boundary migration on the free surface is spasmodic, which means that it has alternating periods of movement and stagnation. Stagnation periods are inversely proportional to the thickness of the sample. The results obtained computationally fitted acceptable with the theoretical results obtained by different authors.

Keywords: thin films, grainboundary, gas difusión.

---

 

I. INTRODUCCIÓN

La variación de tamaño de los cristales en un material es debido mayormente a la migración o movimiento de bordes de grano (BG)⁽¹⁾. La migración de BG depende de la temperatura y de las tensiones en el material entre otros factores, pero además es afectada por las impurezas e inclusiones sólidas⁽²⁾.

FIG. 1: Configuración del BG.

Un BG plano, representado en la Fig. 1, puede ser impulsado por una fuerza constante F. Esta fuerza F puede deberse a una diferencia de tensiones acumuladas⁽³⁾, a un gradiente de temperatura⁽⁴⁾, a diferencias de propiedades magnéticas⁽⁵⁾, etc. En este caso, la coordenada "a" que describe la posición del BG en la superficie (BGS) evoluciona, cuando la superficie no lo afecta, como

donde M es la movilidad del BG.

FIG. 2: Evolución del BG con el surco del BG a velocidad v.

Cuando el BGS se ve afectado por el entorno que rodea a la muestra, se forma un surco superficial y su evolución es afectada profundamente como puede verse en la Fig. 2. En muestras finas^(6,7,8,9)  el BGS se libera del surco y migra durante un período  durante el cual otro surco se va formando y va frenando al BGS hasta estancarlo temporariamente durante un período   [Fig. 2(b) a 2(e)]. La fuerza constante F y la fuerza de capilaridad P, generada en el interior de la muestra, producen una resultante suficientemente grande para desanclar el BG y liberarlo nuevamente. El BGS se mueve por lo tanto en forma espasmódica alternándose períodos de movimiento () y estancamiento () [Fig. 2]. Mullins ⁽¹⁰⁾ sugirió que si el surco del BG se formaba por difusión superficial la coordenada "a", correspondiente a un BG plano como en la Fig. 2, debería moverse durante  como

Donde t es el tiempo y  es una función compleja expresada en serie en . El argumento  es el ángulo medio que forma el BG con la superficie para ambos cristales en donde se ha formado el surco, y  es un ángulo límite dado por  . B esta dado por

donde  es el volumen atómico, σ es la energía superficial, σBG es la energía del BG,  el coeficiente de difusión superficial,  es el número de átomos por unidad de área,  es la constante de Boltzmann y  es la temperatura en grados kelvin.
Gottstein y col⁽¹¹⁾ y Aristov y col.⁽¹²⁾, estudiaron en forma simplificada, el movimiento del BGS de un bicristal plano bajo la acción de una fuerza impulsora constante F, considerando el efecto de la superficie de la muestra durante el tiempo . El surco que se formaba entre el BG y la superficie podía estar producido por efectos de difusión superficial, difusión gaseosa o por el mecanismo de evaporación-condensación.
Si el proceso de formación del surco es por difusión superficial, la coordenada  es:

donde

y  es el espesor de la muestra.
Si el proceso de formación del surco es por evaporación gaseosa, la coordenada  se obtiene a partir de

con β y c2 iguales a

Donde p0 es la presión de vapor de saturación del material y D el coeficiente de difusión gaseosa.
Si el proceso de formación del surco es por evaporación--condensación, la coordenada  se calcula como

con  y

donde Mm es la masa molecular del compuesto.
Se puede demostrar que cuando la migración comienza sin el surco del BG, las Ecs. (4), (8) y (11) pueden reducirse a una única expresión lineal en t dada por la Ec. (1).
Sin embargo, luego de un tiempo suficientemente largo, las Ecs. (4), (8) y (11) se comportan lineales en  , respectivamente,  quedando

Notar que la coordenada "a" evoluciona en el tiempo de una forma relacionada con los mecanismos de formación y de transporte superficial detallados  por Gottstein y col.(11), Aristov y col.(12) y Mullins(13), siendo la evolución más lenta cuando el espesor de la muestra es reducido.
Aristov y col. también dan expresiones del tiempo  que transcurre entre dos sucesivas liberaciones del BGS de su correspondiente surco, para cada proceso físico. Esas expresiones son de la forma

Se ve que para todos los procesos físicos, el tiempo de permanencia del movimiento del BGS con el surco  aumenta cuando el espesor disminuye.

II. MÉTODO COMPUTACIONAL

La migración de BG para el bicristal de la Fig. 1 fue estudiada mediante el método de Monte Carlo (MC) usando los mismos procedimientos que se emplearon en muestras tridimensionales(14,15).

Inicialmente se crea una matriz de N sitios (Nx x Ny x Nz = N) con un bicristal en 3D dentro de la misma. Cada sitio i tiene una orientación única Si > 0. En este modelo, como muestra la Fig. 3, uno de los granos tiene todos los sitios con Si = 40, los sitios en el otro cristal con Si = 20 y el entorno que rodea al bicristal con Si = 30.

FIG. 3. Esquema simplificado usado en la simulación de la cara de arriba del bicristal de la Fig. 1.

El entorno puede ser un gas o el mismo vapor del material del bicristal pero a los fines del cálculo computacional fue considerado como otro cristal.
El algoritmo tiene los siguientes pasos:

a) La energía total W del policristal está dada por

con S_i y Sj orientaciones de los sitios de la red i y j, respectivamente. J es la energía de interacción entre sitios i y j, V el número de vecinos al sitio i y  la función delta de Kronecker. En este caso V=26, que corresponde a los primeros, segundos y terceros vecinos en una estructura cúbica. Se debe aclarar que la energía J depende de la orientación entre el sitio i y j. La energía entre ambos cristales que formaban el BG fue considerada 0,4kT, pero entre esos cristales y el entorno se consideró 0,5kT. Esto de alguna manera intenta poner de manifiesto que el BG tiene menos energía interfacial que una superficie "externa" en estas simulaciones.

En forma aleatoria \footnote{Los números aleatorios usados en este trabajo fueron generados en cada ciclo por una versión adaptada de la subrutina RANECU escrita por F. James (1990)(16).} se elige un sitio de la red perteneciente a un grano, denominado i, con una orientación Si. Mediante la Ec. (4) se calcula la energía alrededor del sitio i como

donde el supra índice in significa etapa inicial.

b) También es posible, en el caso de BG plano, considerar la presencia de una fuerza impulsora contante F. Eso se modeló calculando la energía inicial

Esta cantidad cuenta cuantos sitios con número (S_i) hay en la vecindad. Con el valor de H se puede aumentar el peso de este término ya que su valor está íntimamente relacionado a la fuerza F.

c) A continuación se remplaza la orientación del sitio i (S_i) por una orientación (Sj) del sitio j obtenida aleatoriamente de sus vecinos que pertenecen a un grano.

d) Se calcula nuevamente la energía  del sitio i y también  la energía  donde el supra índice fi significa final. La cantidad  calcula nuevamente cuantos sitios de orientación (S_i) quedaron después del cambio.

e) Luego se calcula la diferencia de energías  como

El primer término está relacionado con la capilaridad y el segundo término con la fuerza externa. Si la cantidad de sitios con orientación (S_i) aumentó y (S_i) está asociada al cristal que uno quiere hacer avanzar, entonces la cantidad del segundo término será negativa. Por lo tanto, el término asociado a la fuerza externa produce una mayor disminución de la energía total del sistema y aumentando la probabilidad de que el sitio con orientación (S_i) no cambie.

f) Si el valor obtenido con la Ec. (22) resulta negativo o nulo, el cambio se produce permanentemente, y si es positivo, se calcula una probabilidad P dada por

Para permitir que el sistema produzca cambios por activación térmica, se elige un número aleatorio Z, entre 0 y 1, y se compara con P. Si P es más grande que Z, entonces se hace el cambio de S_i por Sj, en caso contrario no.

III. RESULTADOS COMPUTACIONALES

En este trabajo se simularon muestras 3D bicristalinas como las mostradas en la Fig. 1, donde las dimensiones de la tapa son 100 x 100 pixeles y se han considerado los siguientes espesores: 20, 40, 60, 100 y 200 pixeles.

Los resultados encontrados en esta sección podrían explicar ciertos comportamientos observados en la migración de muestras de hielo puro obtenidas por Di Prinzio y col. (1997)(9). Para obtener la equivalencia entre el MCS y un valor temporal real se debe conocer el valor de  del BG y la longitud real que tiene el pixel (unidad de longitud en el programa de Monte Carlo), como se había comentado en otro artículo publicado en esta revista (Di Prinzio y col. (2019)⁽¹⁷⁾. En una muestra real estudiada por Di Prinzio y col. (1997) se tenía un espesor entre 3 mm y 5 mm aproximadamente. Por lo tanto podemos decir, por ejemplo, que la muestra computacional con un espesor de 60 pixel es equivalente a una muestra real de 3 mm.

En el trabajo de 2019 de Di Prinzio y col.⁽¹⁷⁾ se encontró que el valor de , mediante la simulación del movimiento de un BG esférico. Para la temperatura de T = -5ºC, el valor promedio entre todas las muestras bicristalinas analizadas por Di Prinzio y col (1997) fue . Por lo tanto, resulta entonces que

A partir de la Fig. 8 se puede afirmar que la configuración final a 10000 MCS representa, aproximadamente, una migración de 2100 horas. Di Prinzio y col. (1997) observaron la migración de los BG y hielo puro por 1500 horas aproximadamente, por lo cual, el valor máximo de MCS cubre correctamente el periodo del experimento real.

En la Fig. 4 se presenta la evolución de la muestra bicristalina en forma 3D (100 x 100 x 40 pixeles) entre los 1600 MCS y los 1900 MCS (MCS: paso de Monte Carlo) con H = 4kT.
En dicha figura, el cristal 1 está representado de color blanco y el cristal 2 de color gris claro. El medio ambiente que rodea al bicristal está representado por un color gris oscuro (ver también Fig. 3). La tapa y la cara del costado izquierdo del bicristal muestran la presencia del surco. El surco tiene un espesor de 1 a 2 pixeles y no entra al interior del bicristal, como se mostró en el esquema de la Fig. 2. En la Fig. 5 se muestra una sección correspondiente a la parte media del bicristal de la Fig. 3(a).

FIG. 4: Evolución del bicristal 3D de 100 x 100 x 20 pixeles para los 1600 MCS (a), 1700 MCS (b), 1800 MCS (c)  y 1900 MCS (d).

Se puede observar que todo el bicristal está rodeado del medio ambiente, que en el interior del bicristal no hay presencia del surco y que el surco tiene escasa profundidad.

FIG. 5: Corte en el medio de la muestra presentada en la Fig. 4(a).

En esta simulación, el cristal  que avanza va perdiendo puntos que son ganados por el medio ambiente. Este efecto se puede considerar como un proceso de difusión gaseosa o evaporación. En la Fig. 6(a) se presentan la forma del BG en el corte de la Fig. 5 y su digitalización. También en la Fig. 7 se presentan las formas del BG, a diferentes MCS, para las muestras delgadas con los diferentes espesores En la Fig. 6(b) se puede ver que la forma coincide aproximadamente con la de una función catenaria como lo reportado por Mullins⁽¹⁰⁾ y Frost y col.⁽⁷⁾. La catenaria es la forma adoptada por una cuerda sostenida fijamente de sus extremos, con una determinada tensión y bajo el efecto de una fuerza constante. El BG adopta en este caso esta forma porque tiene las mismas condiciones que dicha cuerda.

FIG. 6: (a) Forma del BG a los 1600 MCS para la muestra 100 x 100 x 20 pixeles. (b) Coordenadas (x,y) de los puntos situados sobre el BG de la Fig. 6(a) y ajuste por cuadrados mínimos de los mismos con una función catenaria.

La catenaria para un BG, por lo tanto, se puede encontrar haciendo la analogía con la forma catenaria de una cuerda. De esta manera, la catenaria de un BG resulta

donde (x,y(x)) representa las coordenadas de los puntos sobre el BG y . Dicha ecuación fue ajustada por cuadrados mínimos a los puntos (x,y) del corte  de la Fig. 5 y representados en la Fig. 6(b), resultando  pixel,  pixel. Se puede ver que el BG claramente presenta una forma catenaria.

En la Fig. 7, el BG en muestras con menor espesor adopta una forma similar a una catenaria debido a la influencia de la superficie. A medida que el espesor de la muestra aumenta, la influencia de la superficie es menos evidente y eso puede notarse por la forma plana del BG al avanzar. En la Fig. 8 se presenta la coordenada "a" en función del MCS para la muestra 100 x 100 x 50 pixeles con H = 4kT y H = 3,7kT. También se presenta la evolución temporal de "a" de

FIG. 7: Forma del BG en la zonas mostrada en la Fig. 5, a diferentes MCS, para las muestras con los espesores (a) 20 pixeles (2300 MCS), (b) 25 pixeles (2500 MCS), (c) 30 pixeles (1600 MCS), (d) 35 pixeles (900 MCS), (e) 40 pixeles (800 MCS) y (f) 100 pixeles (700 MCS).

Si se multiplica G y MF, se obtiene  . En un trabajo anterior de Di Prinzio y col.⁽¹⁷⁾ se encontró que el valor de , mediante la simulación del movimiento de un BG esférico. Con ambos valores se puede calcular  que es comparable con el valor de  de la simulación en la Fig. 6(a). Por este resultado se puede decir que la simulación describe con su algoritmo la forma que va adoptando el BG en su movimiento bajo efecto del surco superior.

FIG. 8: Coordenada "$a$" en función del MCS para la muestra 100 x 100 x 50 pixeles sin efecto superficial (sin) y con H = 4kT y H = 3,7kT.

Otro importante resultado hallado con este algoritmo computacional es la migración espasmódica del BGS. La migración espasmódica se caracteriza por el movimiento y frenado sucesivos del BGS, caracterizado por la coordenada "a". En la Fig. 8, el movimiento espasmódico se incrementa al reducir la fuerza externa, como lo demuestra el tiempo de estancamiento  encontrado por Aristov y col.⁽¹²⁾ [Ecs. (16), (17) y (18)]. En el caso de la muestra 100 x 100 x 50  pixeles con H = 3,7kT, la coordenada "a" registra un movimiento y frenado sucesivo a medida que trascurre el tiempo. Este efecto es debido, como se explicó en la Fig. 2, a efectos del BG al encontrase con la superficie libre de la muestra.

Para poder ponderar el efecto de estancamiento en la migración de BGS, se estudiaron muestras bicristalinas de diferentes espesores con H = 3,7kT. En la Fig. 9 se presenta la coordenada "a" de bicristales, como los presentados en la Fig. 3, en función del MCS para los siguientes espesores 20, 25, 30, 35, 40, 60 y 100 pixeles.

FIG. 9: Coordenada "a" de bicristales del tipo presentado en la Fig. 3, en función del MCS para diferentes espesores 20, 25, 30, 35, 40, 60 y 100 pixeles.

Aquí se tomaron los resultados promedios de 3 simulaciones para cada espesor. Se observa que:

    a) El desplazamiento neto del BGS es menor cuando menor es el espesor.
   
    b) En todos los casos el movimiento del BGS es espasmódico.
   
    c) El BGS se va frenando en forma suave en algunos de los casos. Sin embargo, el movimiento puede interrumpirse también abruptamente.
   
    d) Cuando el movimiento del BGS se va frenando suavemente, el mismo podría ser descripto mediante una función del tipo , como lo expresado en las Ecs. (14)-(16).
   
    e)  varía con el espesor, notándose que este tiempo se reduce cuando el espesor disminuye. Sin embargo, no fue posible encontrar una prueba explícita de las ecuaciones para  reportadas por Aristov y col.

FIG. 10: Coordenada "a" del BGS para la muestra 100 x 100 x 60 pixeles de la Fig. 7 entre 2500 y 7000 MCS.

En la Fig. 10 se presenta el tramo de movimiento del BG para la muestra 100 x 100 x 60 pixeles de la Fig. 8 entre 2500 y 7000 MCS. Se puede ver que la relación entre "a" y MCS es lineal al comienzo, siguiendo la Ec. (1) con  Después puede ser descripta adecuadamente por una relación como la Ec. (14), propuesta por Aristov y col.⁽¹²⁾, relacionada con un proceso de difusión gaseosa. El factor . Es probable que esta buena coincidencia entre la teoría y la simulación esté relacionada con el hecho que los puntos del cristal que avanzan van pasando a ser parte del entorno del bicristal. Como se mencionó anteriormente, esto puede estar vinculado a un proceso de tipo evaporativo.

IV. CONCLUSIONES

El algoritmo utilizado para describir el movimiento del BGS de una muestra bicristalina delgada coincide en gran medida con la descripción teórica descripta por Gottstein y col.⁽¹⁰⁾, Aristov y col.⁽¹¹⁾ y  Mullins⁽¹²⁾  y con los datos experimentales aportados por Nasello y col.  Los BGS simulados presentaron los siguientes aspectos relevantes: 

1- La velocidad del BGS en muestras delgadas disminuye con el espesor; siendo casi nula cuando el espesor es muy pequeño.

2- La migración del BGS es espasmódica cuando el mismo es afectado por el surco del BG. Existen en general períodos  de estancamiento total y períodos  de migración.

3- Durante el período  el BGS se mueve abruptamente o suavemente. Cuando el movimiento del BGS se va frenado suavemente, el mismo se puede describir mediante una función del tipo   En las simulaciones presentadas, el surco fue formándose por reconversión de los puntos pertenecientes al cristal que avanzaba. Este proceso fue asociado a una difusión gaseosa y eso se reflejó en el exponente encontrado en la función  donde n fue cercano a 3.

4- El BG adopta una forma catenaria como lo indica la teoría.

Agradecimientos

Este trabajo fue posible gracias a la colaboración de José Barcelona (FAMAF, Univ. Nacional de Córdoba) y del apoyo económico de SeCyT (Secretaría de Ciencia y técnica de la UNC).

Referencias

1. A. Sutton y R. W. Balluffi. Interfaces in Crystalline Materials ISBN: 0-19-851385-2 (Oxford: Clarendon Presss, Oxford, Reino Unido, 1995).         [ Links ]

2. G. G. Gosttein y L. S. Shvindlerman. Grain boundary migration in metals ISBN: 084938222X (CRC Press, Boca Raton; London, 1999).         [ Links ]

3. E. A. Ceppi y O. B. Nasello. Computer simulation of bidimensional grain growth. Scr. Metall. 18, 1221-1225 (1984).         [ Links ]

4. B. Light, R. E. Brandt y S. G. Warren. Hydrohalite in cold sea ice: Laboratory observations of single crystals, surface accumulations, and migration rates under a temperature gradient, with application to "Snowball Earth". J. Geophys. Res. 114, C07018 (2009).         [ Links ]

5. D. A. Molodov y L. S. Shvindlerman. Interface Migration in Metals (IMM): "Vingt Ans Après" (Twenty Years Later). Int. J. Mater. Res. 100, 461-482 (2009).         [ Links ]

6. J. S. Lee, B. J. Lee, Y. M. Koo y S. Kim. 3D simulations of grain growth in polycrystalline sheets and rods: Effects of free surface and shape of specimen. Scr. Mater. 110, 113-116 (2016).         [ Links ]

7. H. Frost, C. Thompson y D. Walton. Simulation of thin film grain structures-I. Grain growth stagnation. Acta Metall. Mater. 38, 1455-1462 (1990).         [ Links ]

8. A. Novick-Cohen, O. Zelekman-Smirin y A. Vilenkin. The effects of grain grooves on grain boundary migration in nanofilms. Acta Mater. 58, 813-822 (2010).         [ Links ]

9. C. L. Di Prinzio y O. B. Nasello. Study of Grain Boundary Motion in Ice Bicrystals. J. Phys. Chem. B 101, 7687-7690 (1997).         [ Links ]

10. W. W. Mullins. The effect of thermal grooving on grain boundary motion. Acta Metall. 6, 414-427 (1958).         [ Links ]

11. G. Gosttein y L. S. Shvindlerman. Grain boundary migration in metals: Thermodynamics, Kinetics, Applications 2.a ed. pag. 157. ISBN: 9780429147388 (CRC Press, Boca Raton, 2009).         [ Links ]

12. V. Y. Aristov, V. E. Fradkov y L. S. Shvindlerman. Phys. Met. Metall. 45, 83 (1978).         [ Links ]

13. W. W. Mullins. Flattening of a Nearly Plane Solid Surface due to Capillarity. J. Appl. Phys. 30, 77-83 (1959).         [ Links ]

14. C. L. Di Prinzio, E. Druetta y O. B. Nasello. More about Zener drag studies with Monte Carlo simulations. Modelling Simul. Mater. Sci. Eng. 21, 025007 (2013).         [ Links ]

15. C. L. Di Prinzio y O. B. Nasello. Development of a model for ice core dating based on grain elongation. Polar Science 5, 319-326 (2011).         [ Links ]

16. F. James. A review of pseudorandom number generators. Comput. Phys. Commun. 60, 329-344 (1990).         [ Links ]

17. C. L. Di Prinzio, P. I. Achaval y C. I. Rodriguez. Evolución cinética de un cristal esférico 3D con partículas móviles usando Monte Carlo. Anales AFA 30, 25-30 (2019).         [ Links ]