The Study of Dynamics Heterogeneity in SiO2 Liquid

A molecular dynamics simulation has been carried out to investigate the dynamics heterogeneity SiO2 liquid at 2600 K 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 happened, but 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. We found that the clustering become more pronounced at ambient pressure. This evidences the dynamics heterogeneity in the SiO2 liquid.


Introduction
Dynamical heterogeneity (DH) which has been mentioned in many studies is spatiotemporal fluctuations in local dynamical behavior [1][2][3][4][5][6][7][8][9][10][11][12]. Experimental studies found the phenomenon in polymers [3] and organic compounds [4]. Sen (2008) [5] confirmed the existence of DH in K 2 SiO 3 by the nuclear magnetic resonance experiment. However, these experiments have not provided information about the spatial arrangement of mobile or immobile atoms. More detailed about features of these atoms can be given by computer simulations, because it can observe the motion of individual atoms. For example, Hoang et al. (2007) [6] investigated liquid SiO 2 at some pressures, they indicated that the displacement distribution of atoms at a density about 5.35 g/cm 3 is non-Gaussian that is completely different from dynamical homogeneous systems. They believed that it is evidence of heterogeneous dynamic. Molecular dynamic simulation for Lennard-Jones systems found the evidence of DH based on the multi-point correlation equation [7,8]. This equation only indirectly detects the DH as well as does not mention the cause of the existence of faster and slower regions. However, experimental results have found the evidence of the breakdown and recombination of Si-O bonds [9] in liquid SiO 2 . Hung et al. (2019 and 2020) [10][11][12] continuously informed that O atoms in the coordination number unit SiO x can move to the next unit (Si for OSi y ) correponding when the transition occurs between the units as follow: 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 least. When atom is classified in the group of atoms that are stable, unstable, most mobile or immobile, respectively.

2
The atoms in each of these groups are distributed heterogeneously in the system and form cluster in space which is 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 set of the fastest and slowest atoms [8,13]. When, the structure of this material comprises 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 SiO 2 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][17][18]. In general, the DH is always expected to relate to the structural heterogeneity [8,17,[19][20][21][22]. However, the origin of DH is still not properly understood yet. Therefore, this problem needs more systematic and extensive studies. So, in present paper, we focus on the reactions of DH. And then we determine the clustering from specified sets of atoms to clarify the relating of the nature of DH as well as the structural heterogeneity in the liquid SiO 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 SiO x SiO x1 and OSi y OSi y1 . Other types also occur, for instance SiO x SiO x+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.

(b)
We first determine sets M i consisting of 50 silicon and 100 oxygen (5 % total of atoms); i = 1, 2…5. The atoms of M 1 are chosen randomly from the system. Other M i is specified by the mobility of atoms or the frequency of reactions. The M 2 , M 3 composes of most mobile and immobile atoms. We regard the most mobile and immobile atoms to ones that they displace for the time t D over a distance larger or smaller than remaining atoms in the system, respectively. Thereby the set M 2 , M 3 is determined from the positions of atom in the starting configuration and the configuration at time t D . The M 4 , M 5 consists of atoms with which the reactions happen most frequently or rarely, respectively. The set M 4 , M 5 is determined from the number of reactions which happen with the atoms of these sets within the time t D .

Results and Discussion
We first examine pair radial distribution function (PRDF) as seen in Figure 2 (left) and compared to previous works [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 SiO 2 liquid. This tetrahedral network structure can be seen in the Figure 2 (right). One can see mainly the units SiO 4 and linkages OSi 2 and only has some the units SiO 5 as well as the linkages OSi 3 .

Figure 2. The PRDF (left) and snapshots of arrangement of atoms (right) in the SiO 2 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 G s (r,t) of the van Hove correlation function for the particles [8,24]. The deviation of G s (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: Where 2 () rt  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 SiO 2 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 having the same n. As seen from Figure 3, MSD for central atom with small n is larger than one with big n. Therefore,

Figure 3. The MSD of SiO x and OSi y units in SiO 2 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 M 2 and M 3 indicated DH in the liquid. To analyze the trajectory of atoms of M i (i=2, 3) we first determine M i in the configuration at time t D = 10 ps and at ambient pressure. Then we determine the clusters and mc Mi in the configurations at different time. The data is shown in Figure 4. It can be seen that the mc Mi found varies from 55 to 68 that significantly smaller than mc M1 which is in the interval from 99 to 114. This can be interpreted by that the atoms of M 2 and M 3 reside in the separate regions. Within the time considered the atoms of M 2 and M 3 do not move in random directions, but they displace so that those regions are internally rearranged. On the other hand, the clustering of atoms in sets of M 4 and M 5 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 M 1 , M 4 and M 5 in the configuration at ambient pressure. One can see that the atoms in set of M 1 move in average over a distance larger than M 5 , but significantly smaller than M 4 . This implies that the mobility of atoms of M 1 , M 4 and M 5 increases in the order from M 5 , M 1 to M 4 . The clustering and significant difference in the mobility of atoms from M 4 , M 5 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 6 shows the snapshot of 5 % the most mobile (a), 5 % the most immobile (b) atoms and 5% the randomly atoms (c) in SiO 2 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.

Conclusion
We have shown 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 fix positions. The 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 M i 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 reaction is responsible for DH.