Available online at www.HighTechJournal.org HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 1 The Study of Dynamics Heterogeneity in SiO2 Liquid G. T. T. Trang a*, N. H. Linh a, N. T. T. Linh a, P. H. Kien a* a Thai Nguyen University of Education, 20 Luong Ngoc Quyen, Thai Nguyen, Viet Nam. Received 20 January 2020; Revised 17 February 2020; Accepted 19 February 2020; Published 01 March 2020 Abstract A molecular dynamics simulation has been carried out to investigate the dynamics heterogeneity of SiO2 liquid at 2600 Kelvin and ambient pressure. We indicate that the diffusion in the liquid is realized by the rate of effective reaction, SiOxSiOx’ and OSiyOSiy’. Moreover, the reactions are non-uniform: they are spatially clustered. In addition, we found the clustering from different sets of atoms specified by the mobility of atom or frequency of reactions. Also, results show that the clustering becomes more pronounced at ambient pressure. This evidences the dynamic heterogeneity in the SiO2 liquid. Keywords: Wolf Pack Algorithm; Improvement; Adaptive; Levy Flight; Structural Optimization. 1. Introduction Dynamical heterogeneity (DH), which has been mentioned in many studies, is spatiotemporal fluctuations in local dynamical behavior [1, 2]. Experimental studies found the phenomenon in polymers [3] and organic compounds [4]. The nuclear magnetic resonance experiment confirmed the existence of DH in K2SiO3 by Sen (2008) [5]. However, these experiments have not provided information about the spatial arrangement of mobile or immobile atoms. More detailed information about the features of these atoms can be given by computer simulations because they can observe the motion of individual atoms. For example, Hoang et al. (2007) [6] investigated liquid SiO2 at some pressures and indicated that the displacement distribution of atoms at a density of about 5.35 g/cm3 is non-Gaussian, which is completely different from dynamical homogeneous systems. They believed that it was evidence of heterogeneous dynamics. Molecular dynamic simulation for Lennard-Jones systems found evidence of DH based on the multi-point correlation equation [7, 8]. This equation only indirectly detects the DH and does not mention the cause of the existence of faster and slower regions. However, experimental results have found evidence of the breakdown and recombination of Si-O bonds [9] in liquid SiO2. Hung et al. (2019 and 2020) [10-12] continuously informed us that O atoms in the coordination number unit SiOx can move to the next unit (Si for OSiy) when the transition occurs between the units as follows: SiOx↔SiOx’ và OSiy ↔SiOy’ (1) The heterogeneous distribution of transformations (1) leads to the existence of some faster and slower regions whose conversion frequency is larger or smaller. Also, due to the breakdown and recombination of these Si-O bonds, the list of coordinated atoms of each atom in the system over time may not change, change, change the most or the * Corresponding author: tranggt@tnue.edu.vn; phamhuukien@dhsptn.edu.vn http://dx.doi.org/10.28991/HIJ-2020-01-01-01  This is an open access article under the CC-BY license (https://creativecommons.org/licenses/by/4.0/). © Authors retain all copyrights. https://creativecommons.org/licenses/by/4.0/ https://orcid.org/0000-0003-1002-9678 HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 2 least. When an atom is classified into the group of atoms that are stable, unstable, most mobile, or immobile, respectively, the atoms in each of these groups are distributed heterogeneously in the system and form clusters in space that are stable, unstable, most mobile or immobile regions, respectively. In addition, by tracking the atoms in a liquid, it is possible to group the atoms in the system into a set of the fastest and slowest atoms [8, 13]. When the structure of this material comprises the fastest and slowest regions, the results in Xu et al. (2012) and Vargheese et al. (2010) studies [14, 15] have visualized the fastest and slowest atoms in liquid models. The formation of heterogeneous dynamic regions in liquid SiO2 mentioned above is similar to the heterogeneous regions which are proposed in cooperatively rearranging regions theory by Adam and Gibbs (1965), mode coupling theory by Gotze and Sjogren (1992), or the theory of two-order-parameter by Tanaka (2005) [16-18]. In general, the DH is always expected to relate to structural heterogeneity [8, 17, 19-22]. However, the origin of DH is still not properly understood yet. Therefore, this problem needs more systematic and extensive studies. So, in the present paper, we focus on the reactions of DH. And then we determine the clustering from specified sets of atoms to clarify the relationship between the nature of DH and the structural heterogeneity in the liquid SiO2. 2. Computational Procedure The models composing 1000 silicon and 2000 oxygen atoms have been generated by means of MD simulation. We used the van Beest–Kramer–van Santen (BKS) potential. The MD step is equal to 0.478 fs. Initial configuration is generated by randomly placing all particles in a simulation box. This configuration first is heated to 5000 K and cooled down to 2600 K. Then, the sample has been compressed to specified density. Next the sample is relaxed in N- T-V ensemble (the constant volume and temperature). Then to collect dynamical data we have relaxed the obtained model in N-E-V ensemble (the constant volume and energy) for long times. The temperature and pressure is calculated by averaging over 1000 configurations separated by 10 steps. The volume V and temperature T are maintained several times to reach the desired temperature and pressure. More details about the preparing BKS model can be found elsewhere [23, 24]. Next, the coordination cell which consists of a central atom and neighbors is considered as follows. The distance between central atom and neighbor is less than cutoff distance which is equal to 2.40 Å. The bond is formed by a pair of O and Si which locate in the distance less than the cutoff distance. As a reaction happens, the current bond breaks or new bond is created. Most reactions are SiOxSiOx1 and OSiyOSiy1. Other types also occur, for instance SiOxSiOx+2. However, they occur extremely rarely. Figures 1a and 1b shows how atoms rearrange when two reactions occur. For the first case, a bond is replaced by new bond. In second case, a bond is broken and then restored. The difference between two cases concerns that the atoms in the second case almost vibrate around fix positions. The reaction in first case is called effective reaction. Figure 1. Schematic illustration of two reactions: a) the reactions SiO4SiO5, SiO5SiO4 happen and the bond Si12- O2 is replaced by Si12- O6; b) the reactions SiO5SiO4, SiO4SiO5 happen and the bond Si10- O5 is broken, then restored (a) 1 2 3 4 15 1 2 3 4 15 6 1 3 4 15 6 1 2 3 4 11 5 1 2 3 4 11 1 2 3 4 11 5 (b) HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 3 We first determine sets Mi consisting of 50 silicon and 100 oxygen (5 % total of atoms); i = 1, 2…5. The atoms of M1 are chosen randomly from the system. Other Mi is specified by the mobility of atoms or the frequency of reactions. The M2, M3 composes of most mobile and immobile atoms. We regard the most mobile and immobile atoms to ones that they displace for the time tD over a distance larger or smaller than remaining atoms in the system, respectively. Thereby the set M2, M3 is determined from the positions of atom in the starting configuration and the configuration at time tD. The M4, M5 consists of atoms with which the reactions happen most frequently or rarely, respectively. The set M4, M5 is determined from the number of reactions which happen with the atoms of these sets within the time tD. 3. Results and Discussion We first examine pair radial distribution function (PRDF) as seen in Figure 2 (left) and compared to previous studies [23, 24]. The result is consistent with the experiment in the positions and heights of PRDF peaks. From PRDF, we calculated the averaged coordination numbers for pair Si-O, O-Si for the configuration at ambient pressure, they is closed to 4 and 2, respectively. It shows the tetrahedral network structure of SiO2 liquid. This tetrahedral network structure can be seen in the Figure 2 (right). One can see mainly the units SiO4 and linkages OSi2 and only has some the units SiO5 as well as the linkages OSi3. Figure 2. The PRDF (left) and snapshots of arrangement of atoms (right) in the SiO2 liquid at 2600 K and ambient pressure; blue (large) and red (small) spheres represent Si and O, respectively As we known, the DH in the liquid system is often detected by the time dependence of the self-part Gs(r,t) of the van Hove correlation function for the particles [8, 24]. The deviation of Gs(r,t) from a Gaussian form at intermediate times reflected the existence of DH in the system, where r is the distance traveled by a particle in a time t. Hence, deviations can be determined by the non-Gaussian parameter which has the form as below: 4 2 2 2 3 ( ) ( ) 1, 5 ( ) r t t r t        (2) Where 2 ( )r t  is the mean squared displacement, if the system is dynamically homogeneous α2(t)=0 at high temperatures due to the dynamical homogeneous and it has a maximum at lower temperatures due to the dynamical heterogeneities. In this paper, we detect the DH by mobility of atoms. Namely, we examine the types of atomic motion and DH in SiO2 liquid. In Figure 3 we show the mean square displacement (MSD) of central atom at the ambient pressure. We denote n is the number of Init-bonds which are linked to the central atom. Because central atom and n neighbors move together during the time considered, hence the set of central atom and n neighbors behaves like a super-molecule (SM) which flows in the liquid. The SM consists of n+1 atom. The MSD is obtained by averaging over all central atoms 0 3 6 9 12 0 3 6 9 12 T h e P R D F , g (r ) r (angstrom) Si-O Si-Si O-O HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 4 having the same n. As seen from Figure 3, MSD for central atom with small n is larger than one with big n. Therefore, the mobility of central atom is correlated with the size of SM. Moreover, the MSD of SiO4 units is unchanged with the time. This is caused by the stable of SiO4 network-structure at ambient pressure. Figure 3. The MSD of SiOx and OSiy units in SiO2 liquid According to Hung et al. (2016) [25], the DH connected to the clustering from most mobile and immobile atoms in the system. Therefore, we also use the approximation to observe clustering from sets of M2 and M3 indicated DH in the liquid. To analyze the trajectory of atoms of Mi (i=2, 3) we first determine Mi in the configuration at time tD = 10 ps and at ambient pressure. Then we determine the clusters and mcMi in the configurations at different time. The data is shown in Figure 4. It can be seen that the mcMi found varies from 55 to 68 that significantly smaller than mcM1 which is in the interval from 99 to 114. This can be interpreted by that the atoms of M2 and M3 reside in the separate regions. Within the time considered the atoms of M2 and M3 do not move in random directions, but they displace so that those regions are internally rearranged. Figure 4. The MSD (left) and number of clusters (right) for atoms of M2 and M3. The set M2 and M3 is determined for the configuration at tD = 10 ps and at ambient pressure 2 4 6 8 10 0 5 10 15 20 25 The time, ps T h e m e a n s q u a re d is p la c e m e n t, Å 2 T h e m e a n s q u a re d is p la c e m e n t, Å 2 SiO0; SiO1 SiO2; SiO3 SiO4 2 4 6 8 10 0 5 10 15 20 25 OSi1 OSi2 OSi3 The time, ps 2 4 6 8 10 55 60 65 70 75 T h e M S D o f a to m s T h e n u m b e r o f c lu s te rs The time, ps M2 M3 0 2 4 6 8 10 0 3 6 9 12 15 M2 M3 The time, ps HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 5 On the other hand, the clustering of atoms in sets of M4 and M5 shows that the reactions are not uniformly distributed in the liquid, but they happen with quite different frequency in separate regions. Expectedly, Figure 5 shows the MSD and the mean number of reaction per atom for M1, M4 and M5 in the configuration at ambient pressure. One can see that the atoms in set of M1 move in average over a distance larger than M5, but significantly smaller than M4. This implies that the mobility of atoms of M1, M4 and M5 increases in the order from M5, M1 to M4. The clustering and significant difference in the mobility of atoms from M4, M5 indicates that there are mobile regions where reactions frequently happen and immobile ones where reactions occur rarely. This result also indicates the DH which is caused by the non-uniform spatial distribution of reactions. Figure 5. The mean number of reactions per atom (right) and MSD for the configuration (left) at ambient pressure Figure 6 shows the snapshot of 5 % the most mobile (a), 5 % the most immobile (b) atoms and 5% the randomly atoms (c) in SiO2 liquid. Compared to the randomly atoms, the most mobile and immobile are non-uniform distributed in space, but instead they lead to cluster with each other into the separate groups/clusters. This trend also shows from snapshot of the positions where the reactions happen as seen Figure 7. It can be seen that the reactions happen rarely as well as frequently non-uniform distributed in the space. Namely, they are distributed corresponding with the immobile regions and with the mobile regions, respectively. The Figures 6 and 7 informs that atoms of similar mobility are spatially correlated and that atoms with different mobility tend to be anti-correlated. It is noted that the correlation is quantitatively considered by calculating static pair correlation functions between atoms belonging to the mobile/immobile regions [8]. These results again demonstrated that the DH which is caused by the non-uniform spatial distribution of the most mobile or immobile. Figure 6. Snapshot of the positions of 5 % the most mobile (M3), 5 % the most immobile (M2) atoms and 5 % the randomly atoms (M1) in SiO2 liquid. Here the blue sphere is Si or O atom 2 4 6 8 10 2 4 6 8 T h e n u m b e r o f re a c ti o n s T h e M S D o f a to m s The time, ps M1 M4 M5 2 4 6 8 10 0 50 100 150 200 250 M1 M4 M5 The time, ps (M3) (M2) (M1) HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 6 Figure 7. Snapshot of the positions where the reactions happens rarely (M4) and frequently (M5) in the space; here the purple sphere presents the positions of reaction 4. Conclusion The results show that the reactions are not uniformly distributed in the liquid. We reveal that within the time studied, the liquid comprises separate mobile and immobile regions. Large rigid clusters consist of immobile atoms which vibrate around fixed positions. Mobile atoms move by free or correlation motion type. The reactions happen rarely in the immobile regions, but frequently in the mobile regions. We also show that the non-uniform spatial distribution of reactions causes the DH in the liquid. Our simulation reveals the clustering from sets Mi, which are specified by the atom mobility or frequency of reaction. The clustering is most pronounced for the configuration at ambient pressure. This shows DH in the liquid. We also found that the non-uniform spatial distribution of reactions is responsible for DH. 5. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. 6. References [1] Mizuno, H., & Yamamoto, R. (2011). Dynamical heterogeneity in a highly supercooled liquid: Consistent calculations of correlation length, intensity, and lifetime. Physical Review E, 84(1). doi:10.1103/physreve.84.011506. [2] Franz, S., Parisi, G., Ricci-Tersenghi, F., & Rizzo, T. (2011). Field theory of fluctuations in glasses. The European Physical Journal E, 34(9). doi:10.1140/epje/i2011-11102-0. [3] Tracht, U., Wilhelm, M., Heuer, A., Feng, H., Schmidt-Rohr, K., & Spiess, H. W. (1998). Length Scale of Dynamic Heterogeneities at the Glass Transition Determined by Multidimensional Nuclear Magnetic Resonance. Physical Review Letters, 81(13), 2727–2730. doi:10.1103/physrevlett.81.2727. [4] Qiu, X. and Ediger, M. D. (2003). Length scale of dynamic Heterogeneity in supercooledd-Sorbitol: comparison to Model predictions. The Journal of Physical Chemistry B, 107(2), 459–464. doi:10.1021/jp021888b. [5] Sen, S. (2008). Differential mobility and spatially heterogeneous dynamics of oxygen atoms in a supercooled glass-forming network liquid. Physical Review B, 78(10). doi:10.1103/physrevb.78.100201. [6] Hoang, V. V., Zung, H., & Hai, N. T. (2007). Diffusion and dynamical heterogeneity in simulated liquid SiO2under high pressure. Journal of Physics: Condensed Matter, 19(11), 116104. doi:10.1088/0953-8984/19/11/116104. [7] Flenner, E., & Szamel, G. (2007). Anisotropic spatially heterogeneous dynamics in a model glass-forming binary mixture. Journal of Physics: Condensed Matter, 19(20), 205125. doi:10.1088/0953-8984/19/20/205125. [8] Vogel, M., & Glotzer, S. C. (2004). Temperature dependence of spatially heterogeneous dynamics in a model of viscous silica. Physical Review E, 70(6). doi:10.1103/physreve.70.061504. (M5) (M4) HighTech and Innovation Journal Vol. 1, No. 1, March, 2020 7 [9] Mei, Q., Benmore, C. J., & Weber, J. K. R. (2007). Structure of Liquid SiO2: A Measurement by High-Energy X-Ray Diffraction. Physical Review Letters, 98(5). doi:10.1103/physrevlett.98.057802. [10] Hung, P. K., Hong, N. V., Trang, G. T. T., & Iitaka, T. (2019). Topological analysis on structure and dynamics of SiO2 liquid with the help of Si-particle and O-particle statistics. Materials Research Express, 6(8), 085201. doi:10.1088/2053- 1591/ab1c73. [11] Hung, P. K., Vinh, L. T., Hong, N. V., Trang, G. T. T., & Nhan, N. T. (2019). Insight into microstructure and dynamics of network forming liquid from the analysis based on shell–core particles. The European Physical Journal B, 92(8). doi:10.1140/epjb/e2019-100255-8. [12] Hung, P. K., Vinh, L. T., Ha, N. T., Trang, G. T. T., & Hong, N. V. (2020). Domain structure and oxygen-pockets in the silica melt under pressure. Journal of Non-Crystalline Solids, 530, 119780. doi:10.1016/j.jnoncrysol.2019.119780. [13] Vogel, M., & Glotzer, S. C. (2004). Spatially Heterogeneous Dynamics and Dynamic Facilitation in a Model of Viscous Silica. Physical Review Letters, 92(25). doi:10.1103/physrevlett.92.255901. [14] Xu, W.-S., Sun, Z.-Y., & An, L.-J. (2012). Effect of attractions on correlation length scales in a glass-forming liquid. Physical Review E, 86(4). doi:10.1103/physreve.86.041506. [15] Vargheese, K. D., Tandia, A., & Mauro, J. C. (2010). Origin of dynamical heterogeneities in calcium aluminosilicate liquids. The Journal of Chemical Physics, 132(19), 194501. doi:10.1063/1.3429880. [16] Adam, G., & Gibbs, J. H. (1965). On the Temperature Dependence of Cooperative Relaxation Properties in Glass‐Forming Liquids. The Journal of Chemical Physics, 43(1), 139-146. doi:10.1063/1.1696442. [17] Gotze, W., & Sjogren, L. (1992). Relaxation processes in supercooled liquids. Reports on Progress in Physics, 55(3), 241–376. doi:10.1088/0034-4885/55/3/001. [18] Tanaka, H. (2005). Two-order-parameter model of the liquid–glass transition. II. Structural relaxation and dynamic heterogeneity. Journal of Non-Crystalline Solids, 351(43-45), 3385–3395. doi:10.1016/j.jnoncrysol.2005.09.009. [19] Weeks, E. R., J. C. Crocker, A. C. Levitt, A. Schofield, and D. A. Wei (2000). Three-Dimensional Direct Imaging of Structural Relaxation Near the Colloidal Glass Transition. Science, 287(5453), 627–631. doi:10.1126/science.287.5453.627. [20] Weeks, E. R., Crocker, J. C., & Weitz, D. A. (2007). Short- and long-range correlated motion observed in colloidal glasses and liquids. Journal of Physics: Condensed Matter, 19(20), 205131. doi:10.1088/0953-8984/19/20/205131. [21] Tanaka, H., Kawasaki, T., Shintani, H., & Watanabe, K. (2010). Critical-like behaviour of glass-forming liquids. Nature Materials, 9(4), 324–331. doi:10.1038/nmat2634. [22] Kawasaki, T., & Tanaka, H. (2011). Structural signature of slow dynamics and dynamic heterogeneity in two-dimensional colloidal liquids: glassy structural order. Journal of Physics: Condensed Matter, 23(19), 194121. doi:10.1088/0953- 8984/23/19/194121. [23] Hung, P. K., Hong, N. V., & Vinh, L. T. (2007). Diffusion and structure in silica liquid: a molecular dynamics simulation. Journal of Physics: Condensed Matter, 19(46), 466103. doi:10.1088/0953-8984/19/46/466103. [24] Tatarinova, L. I. (1983). Structure of solid amorphous and liquid substances. From: Moscow, Nauka, 619097. [25] Hung, P. K., Kien, P. H., San, L. T., & Hong, N. V. (2016). The study of diffusion in network-forming liquids under pressure and temperature. Physica B: Condensed Matter, 501, 18–25. doi:10.1016/j.physb.2016.07.033.