Journal of Nanomaterials
Volume 2020 (2020), Article ID 8857101, 13 pages
https://doi.org/10.1155/2020/8857101

Research Article

Molecular Dynamics Study on SiO₂ Interfaces of Nonfiring Solids

Tomohiro Sato,¹ Atsuto Kubota,² Ken-ichi Saitoh,¹ Masayoshi Fuji,³ Chika Takai,⁴ Hadi Sena,⁵ Masanori Takuma,¹ and Yoshimasa Takahashi¹

¹Department of Mechanical Engineering, Kansai University, 3-3-35 Yamate-cho, Suita, Osaka 564-8680, Japan
²Graduate School of Science and Engineering, Kansai University, 3-3-35 Yamate-cho, Suita-shi, Osaka 564-8680, Japan
³Advanced Ceramics Research Center, Nagoya Institute of Technology, 10-6-29, Asahigaoka, Tajimi, Gifu 507-0071, Japan
⁴Department of Chemistry and Biomolecular Science, Gifu University, 1-1 Yanagido, Gifu 501-1193, Japan
⁵Institute of Materials and Systems for Sustainability, Nagoya University, Furo-cho, Chikusa, Nagoya 464-860, Japan

Correspondence should be addressed to Tomohiro Sato; tom_sato@kansai-u.ac.jp

Received 26 August 2020; Revised 5 December 2020; Accepted 11 December 2020; Published 23 December 2020

Academic Editor: Filippo Giubileo

Copyright © 2020 Tomohiro Sato et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

As a sustainable ecosystem, the general firing process for ceramics emits large amounts of CO₂ gas; thus in ceramics production, the focus is the nonfiring process; however, the solidification and strengthen mechanism of this nonfiring system, which essentially reacts between surface-activated ceramic particles and a solvent, has not been elucidated to date. The nonfiring process had three steps, i.e., particle surface activate process by grinding process, maintaining the active state until starting nonfiring solidification begins, and nonfiring solidification process. Thus, in this study, the reaction of silica and water was simulated by adapting molecular dynamics based on LAMMPS with ReaxFF potentials. Reproducing the activated silica surface state, three ended models called O model, Si model, and OH model were prepared which indicated ended molecules of each surface. These models and a water molecule as a solvent were bonded in the atomic scale, and the energetic state and mechanical properties were evaluated. A reacted or structured O-H-O bond was reproduced in the nonfiring process in the O-ended model. The bond was a hydrogen bond. A Si-O-Si bond was produced when a Si atom was ended on the interface. The bonded interface was able to tensile. However, the tensile strength was weaker than that of the solid silica model. The nonbonded OH model did not have tensile strength.

1. Introduction

In the powders and powder metallurgy field, particularly for metal or ceramics powders, sintering is a key process; however, industrial sintering furnaces emit CO₂. This problem can be solved using a nonfiring process.

It is a unique surface activation method, and it was proposed as a novel ceramic consolidation method, which does not utilize the firing (sintering) process.

The nonfiring process involves a particle surface active process via mill grinding, maintaining the active state until nonfiring solidification begins, and a nonfiring solidification process. An experimental nonfiring process including these three steps was studied as follows.

The activated particle surface was fabricated using the nonfiring process. A black paper sludge ash surface was activated by a mechanochemical process. Then, these materials could be on a nonfiring process by using KOH solution [1].

The same process was adopted for the combination of alumina particle and absorbed water molecules [2]. Other ceramic materials have been investigated experimentally as candidates for the nonfiring process [3–5]. Generally, the surface energy of the powder is stable because the surfaces are covered oxide. The powder surfaces become active via a grinding process, e.g., milling.

A mechanochemical process is often used for material synthesis because this process can obtain high energy for particle surfaces. For example, barium titanate ceramics was synthesized using this process [6]. It was not a nonfiring process. However, as a mechanochemical process, grinding was useful relative to understanding the solidification process by firing (i.e., the sintering process). By changing its condition, the mechanochemical process could activate only the particle surface, not breaking and not synthesizing the particles.

Additionally, it is easier to perform a nonfiring solidification process on ceramic powders with active surface energy. However, the strength of ceramics which bonded only surface diffusion in the firing process has not been sufficiently investigated. Generally, during the sintering process, increased grain size results in greater strength due to the diffusion of atoms from the powder’s surface. However, only the surfaces of powders are joined by a nonfiring process. Therefore, the strength of ceramic material that has undergone a nonfiring process may be less than that of sintered ceramics. Ensuring comparable strength is key to making nonfiring processes practical; however, investigating the structural and energetic state of ceramics powder surfaces experimentally at microscale is difficult. Therefore, in this study, interfaces of silica (SiO₂), which is extensively used as a ceramic material, were modeled as nonfiring solids using a molecular dynamics (MD) method based on a Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) with ReaxFF, which is a force field for reactive systems [7].

Many previous MD studies use silica, and most of these studies employ the Tersoff potential [8], which is a popular potential function for Si systems, e.g., SiO₂ and SiC. Recently, reactive MD studies have used ReaxFF, as a reactive force field. For glassy materials, such as glassy silicate, reactive MD has been used to investigate their elastic properties [9].

Another model calculated the acidity of single sites on the humid silica surface represented by a model for the hydroxylated amorphous surface. This recalls the behavior of the out-of-plane silanols on the crystalline (0001) αquartz surface, although the acidity here is even stronger. In this abinitio MD models, details on the solvation of the different surface sites of silica were described [10]. The assessment of the empirical reactive force field ReaxFF to predict the formation of amorphous silica from its crystalline structure and the determination of mechanical properties under tension using MD simulations is presented. Stress relaxation simulations indicate that the transition strain rate occurs when the characteristic time for high-strain rate loading and stress relaxation times are in the same order [11].

In the tribology field, the wear mechanism of silicon at the Si/SiO₂ interface in an aqueous environment has been investigated by MD using ReaxFF [12, 13]. Furthermore, pressure distributions and contact areas were simulated for an octadecylsilane-functionalized silica interface [14]. The mechanical effects of adsorbed films at surface contacts were also investigated, and broadening of the contact area and minimization of direct surface contact were observed. Silica is a very available material; therefore, a self-assembled monolayer (SAM) was stacked on SiO₂ basement material. Then, SiOH was selected as the first reaction component between silica and the SAM, and these reactions were investigated by MD [15]. These studies clarified that MD is a useful method to investigate the energy and structural properties of silica surfaces and interfaces.

In this paper, a nonfiring reaction for some types of interfaces is discussed. Additionally, tensile strength was considered by comparing MD models with and without interfaces. We also considered reactions using ReaxFF where the silica surface is active, and it could be clarified as a part of the mechanism of the nonfiring process of silica.

2. Procedure

2.1. Unit Cell and Potentials

Unit cell properties of β-tridymite [16], the number of atoms, density, cell size, and crystal structure are listed in Table 1. In silica crystal, β-tridymite density is lower, and hence, they are candidate of nonfiring silica experimentally, and the structure was suitable for water adsorption and hydroxylation [17]. A schematic of a β-tridymite unit cell is shown in Figure 1. In this study, ReaxFF is employed as an interatomic potential function that is expressed as follows [18–20]:

where is the bond, is the lone pair, is the over-coordination, is the under-coordination, is the angle, is the torsion, is the Van der Waals, is the Coulomb, and is the hydrogen bond.

Table 1: β-Tridymite unit cell properties and property values.

Number of atoms
12

Density (g/cm³)
2.19

Cell size (Å)

Crystal structure
Hexagonal

Figure 1: β-Tridymite unit cell structure [16].

ReaxFF, a calculation program in the reactive force field molecular dynamics method (ReaxFF method), was originally developed by van Duin, Goddard, and coworkers. The LAMMPS introduces the many-body potential required for the ReaxFF method, i.e., the ReaxFF potential. Here, the potential is calculated in the same way as conventional classical MD calculations. In addition, the ReaxFF potential is based on the bond order concept and is calculated using the instantaneous bond order from the interatomic distance. However, by itself, the ReaxFF potential cannot accurately represent the bond order for each atom, and a partial over-coordination state appears. Thus, the correct bond order can be expressed by following the bond atomization theory. Then, the bond breakage and generation between atoms can be expressed. This is one of the major features of ReaxFF, and as a result, it is classified as a reaction force field and used in systems where chemical reactions occur. In addition, the obtained bond order is incorporated in the bond length term and angle term between three atoms, the torsion term between four atoms, and all atoms regardless of the existence of bonds between atoms; i.e., nonbonded interactions (e.g., Van der Waals and Coulomb interactions) are calculated between them, and nonbonded interactions at short distances are avoided by shielding. Nonfiring solidification involves a chemical reaction between SiO₂ and H₂O molecules that have undergone surface activation. In this process, H₂O molecules dissociate and bonds form between dissociated O and H atoms and SiO₂. Therefore, in this study, ReaxFF was adopted as the interatomic potential.

ReaxFF force field parameters have been developed for various compounds, such as silicon oxide and nitramines, starting with those made for hydrocarbons, and are roughly categorized into three systems: (1) a system that considers O and H atoms as gases only (because the target system temperature exceeds the boiling point of water); (2) a combustion system force field; and (3) an aqueous system force field that can treat O and H atoms as liquids. It is an independent force field that does not use atoms. This branching occurs because O and H atoms parameterized as gases cannot be directly handled as liquids. Consequently, the same parameters cannot always be used for the same atomic species, and the versatility per parameter is not high. Therefore, many parameters have been developed, and parameters that can be used for ~60% of atomic species have been determined. The system used in this study involved three types of atoms, i.e., Si, O, and H. van Duin et al. adapted them to silicon and silicon oxide as parameters including them. In this study, H₂O was used in the liquid state; therefore, the abovementioned aqueous force field was required, and the parameter was expanded by Fogarty et al. In this study, ReaxFF parameters used by Kulkarni et al. are adopted and extended to gas-surface interactions.

2.2. Simulation Models

Three models were developed to reproduce the nonfiring reaction between silica interfaces. The ended atoms of each surface had another kind of atoms. The O (oxygen) model is shown in Figure 2. The ended atom of the silica surfaces was O (oxygen). A schematic of the interface is shown in Figure 3. The O atom-ended interfaces included water molecules.

Figure 2: Structure of the O model.

Figure 3: Schematic of the interface structure of the O model.

In the same manner as the O model, a Si (silicon) model was prepared (Figure 4). A schematic of this model is shown in Figure 5. Si-ended atoms were neighbored water molecules as an initial configuration.

Figure 4: Structure of the Si model.

Figure 5: Schematic of the interface structure of the Si model.

The third silica interface model was an OH model shown in Figure 6. This configuration was modeled as surface grind and an adsorbed model experimentally [5]. As discussed in the literature, the surface of the silica was adsorbed by OH (Figure 7). The distance between the ended-atom of each interface atom and the oxygen of a water molecule was set to 1 Å. Experimentally, these models were able to exist because of surface activation [5].

Figure 6: Structure of the OH model.

Figure 7: Schematic of the interface structure of the OH model.

OH groups are added to the SiO₂ surface to maintain the increased activity of particles due to the grinding, and this model was reproduced in a simulation. Here, it is assumed that the surface properties remain the same as when the activated particles are held after compression molding.

For all the three models, the surface structure of the particles was such that the interface was in the [0001] direction of the β-tridymite structure. In this direction, the surface energy (reactivity) of β-tridymite is high, and the surface activity due to grinding is indicated using this direction as the interface. In addition, the periodic boundaries are taken in all directions; thus, the upper and lower ends of the calculation model do not take discontinuous structures. The size of the supercell is ~25 mm in the [001 0] direction, ~13 mm in the [1 100] direction, and~45–50 mm in the [0001] direction.

In the actual state, the activated particles are compression-molded and placed in a solvent to obtain a compact. In this process, a solvent enters in each particle and promotes a chemical reaction between the particles. Then, solidification occurs due to this reaction. The computational model used in this study focuses on a set of particles, and the model has a structure in which a solvent was placed between the particles.

2.3. Simulation Conditions

Table 2 shows the reaction (relaxation) simulation conditions for each model. Here, using the NVT ensemble, the temperature was maintained at 300 K, which is consistent with a nonfiring state. Note that periodic boundary conditions were also adopted in all models. Table 3 shows the tensile simulation conditions of each model. To facilitate comparison, a solid β-tridymite model was simulated. The constant strain rate was set to 2.5 ps after NPT simulation in 2.5 ps as relaxation after the reaction simulation. This is a very high value compared to the strain rate of the actual system; however, the tensile properties of silica have ranged from to in previous simulation studies, and the strain rate in this study is within that range.

Table 2: Calculation condition on relaxation calculation.

Calculation model
O model
Si model
OH model

Numbers of atoms
1125
1035
1215

Number of H₂O
30

Potential
ReaxFF

Boundary condition
Periodic boundary

Ensemble
NVT

(Å)
25.235
25.235
25.235

(Å)
13.1125
13.1125
13.1125

(Å)
47.572
46.4149
49.3578

Temperature (K)
300

Calculation time (ps)
10

Table 3: Calculation condition on tensile simulation.

Calculation model
O model
Si model
OH model
β-Tridymite

Numbers of atoms
1125
1035
1215
1080

Potential
ReaxFF

Boundary condition
Periodic boundary

Ensemble
NPT → NVT

(Å)
25.235
25.235
25.235
25.235

(Å)
13.1125
13.1125
13.1125
13.1125

(Å)
47.572
46.4149
49.3578
49.572

Pressure (Pa)
0

Temperature (K)
300

Tensile direction
direction

Strain rate (fs⁻¹)
⁻⁵

Calculation time (ps)
2.5→5.0

3. Results and Discussion

3.1. Reaction in O Model

Figure 8 shows the configurations of the interface atoms of the O model after the simulation. Some atoms in the interface had reacted water molecules, which put the interface of the silica. One of the reacted atoms is illustrated in Figure 9. As can be seen, the H atom that bonded with the silica surface was part of the water molecule prior to the relaxation calculation. However, it was completely dissociated from the water molecule after relaxation; i.e., it was not part of the water molecule. For comparison, interface atoms that did not react and bond between interfaces and internal atoms are also shown in Figure 9. As shown in Figure 9, the interface bonded through the H (hydrogen) atom from the water molecule. In this model, five types of these bonds were observed at the interfaces. Consequently, the O model exhibited O-H-O bonds between the interfaces.

Figure 8: Structure of interface on the O model.

Figure 9: Characteristic structure of interface on the O model.

Figure 10 shows the changes of the potential energy of the atoms shown in Figure 9. Here, two types of potential energy are shown. As can be seen, bonded atoms indicated lower energy, which indicates more stable potential energy than that of nonbonded atoms. Internal atoms were the most stable, and a reference value [21] had a similar bond value. The potential energy of β-tridymite and a reference value [21] are compared in Figure 11. As can be seen, they were nearly the same value.

Figure 10: Comparison of potential energy in the O model.

Figure 11: Potential energy of β-tridymite.

The bonded length was evaluated for five bonds that reacted as nonfiring state (Table 4). Five water molecules underwent this process, and the entire system contains 30 water molecules; thus, the ratio is ~17%. Here, each length was <2.4 Å. In a previous study, the length of the hydrogen bond (H-OH bond) was <2.2–2.5 Å [22]. From the configuration, potential energy, and bond length, the bonds of the nonfiring O model were hydrogen bonds.

Table 4: Distance of bonded parts on the O model.

Bonding part
1
2
3
4
5

O-H…O distance (Å)
2.303
2.267
2.324
2.258
2.378

The reaction process of the O model is shown in Figure 12. When hydrogen atom parts from one water molecule, the remaining OH atoms parted from interface. In the O model, some nonfiring hydrogen bonds were observed. However, the reason for such bonding was not clarified. As will be described in detail in the tensile test section, the interface created by the non-firing process has lower strength than the other sections. Certainly, this strength is inferior to that of the original powder; however, any shape can be made by a nonfiring process.

Figure 12: Reaction between SiO₂ and H₂O on the O model.

3.2. Reaction in Si Model

Figure 13 shows the configurations of interface atoms of the Si model after the relaxation simulations. Certain Si atoms in the upper and lower parts of the interface were bonded through the oxygen atom. Here, 13 bonds were confirmed after the relaxation simulations. Some of the bonded part, nonbonded part, and internal structures are shown in Figure 14. Initially, in the bonded model, the O atoms that tied the upper and lower Si atoms were a water molecule. Moreover, OH and H atoms were observed, as shown in Figures 13 and 14. Some of these atoms were bonded to a silica atom. Other atoms appeared to become hydrogen molecules. When a Si-H bond was constructed on the surface, it might be difficult to bond to each other, as shown in Figure 14. The change of potential energies through the relaxation simulation is shown in Figure 15. Similar to the O model shown in Figure 10, here, the bonded atoms rather than nonbonded atoms are in a stable state. The reaction path of the Si model is shown in Figure 16. First, an O atom parted from a water molecule, and then, a Si-O-Si bond was formed between silica and hydrogen molecules. Then, one of the hydrogen atoms that was initially a water molecule was bonded with a Si atom. Another hydrogen atom becomes a hydrogen molecule by reacting with another hydrogen atom. Table 5 shows the measured length of bonded atoms was for 13 bonded combinations. As can be seen, siloxane bonds with [23] were formed because its distance is ~1.5–1.6 Å. The configurations and changes of potential energy and distance were evidence of siloxane bonds in this nonfiring reaction of the Si model.

Figure 13: Structure of interface on the Si model.

Figure 14: Characteristic structure of interface on the Si model.

Figure 15: Comparison of potential energy in the Si model.

Figure 16: Reaction between SiO₂ and H₂O on the Si model.

Table 5: Distance of bonded parts on the Si model.

Bonding part
1
2
3
4
5
6
7
8
9
10

Si-O-Si distance (Å)
1.578
1.543
1.594
1.606
1.562
1.570
1.553
1.593
1.528
1.549

3.3. Reaction in OH Model

Figure 17 shows the configurations of the OH model after relaxation simulations. Most water molecules seemed to remain its state. By comparing the potential energy of bonded atoms in the O and Si models and the nonbonded atoms of the OH model, the nonbonded atoms of the OH model indicated lower energy (Figure 18). In addition, the lower energy of the OH model was the same as the energy of the bonded atoms. Thus, the surfaces of the OH model were in a stable energy state from the beginning of the relaxation simulations. As a result, a reaction between water molecules and the interfaces did not occur. Figure 19 shows the configurations of interface atoms. Here, the surface structure that had OH groups was indicated. Table 6 shows the electrical charge of some surface atoms of the OH model. Table 7 also shows the electrical charge of some of the inner atoms of the OH model.

Figure 17: Structure of interface on the OH model.

Figure 18: Comparison of potential energy.

Figure 19: Surface structure in the OH model.

Table 6: Charge of surface structure.

Si₁
O₁
O₂
O₃
O₄
H₁
H₂
H₃
Total

Charge (e)
0.527
-0.262
-0.137
-0.148
-0.133
-0.019
-0.011
-0.017
0.14

Table 7: Charge of internal structure.

Si₂
O₁
O₅
O₆
O₇
Total

Charge (e)
0.613
-0.262
-0.294
-0.316
-0.322
0.016

By comparing these tables, we find that the surface atoms had a positive charge despite also having stable potential energy. Experimentally, powders before nonfiring maintained this surface state which was the OH model to prevent oxidization of the surface as raw materials. After that, compressing and water adsorption were performed. Then a part of surface state might become the O model and/or Si model.

3.4. Tensile Simulation

Figure 20 shows the stress-strain curve of each model. Here, the reaction state of the three models and the solid β-tridymite model, which has no interface, was compared. The results indicate that solid β-tridymite was the strongest. Then, the Si model, which had 13 siloxane bonds, exhibits the next best strength. The O model, which had five hydrogen bonds, showed weaker tensile strength. However, the OH model had no tensile strength, which means that the interface of the OH model was not bonded. It is a brittle material; thus, it is difficult to compare it in tensile test experiments.

Figure 20: Stress-strain curves for each model.

However, it can be compared to the results of previous simulations. Using the same potential, the literature is amorphous (number of atoms: 32928), and this study differs from β-tridymite (number of atoms: 1080). However, the tensile strength is greater by approximately 15 GPa.

Figure 21 shows the relationship between strain and atomic force in the direction for fractured and nonfractured atoms on the Si model. The fractured atoms were bonded by the nonfiring process. When the strain was >0.175, the bonded atoms demonstrate a large force value; i.e., fracture occurred. However, nonbonded atoms did not exhibit these phenomena. Figure 22 shows the compares the Si model just prior to breaking and the β-tridymite fracture surface. Here, the fracture surface of the Si model was the interface formed by the Si-O-Si bond, and the fracture surface of β-tridymite was a Si-O-Si bond parallel to the z direction. As shown in Table 3, the area of the plane of β-tridymite was equal for each calculation model. Since the applied strain rate was equal for each model, the stresses acting on the plane were also equal. The β-tridymite was replicated in five directions in the direction and three directions in the direction; thus, there were 15 Si-O-Si bonds parallel to the direction in the entire system. However, because the number of bonds generated in the Si model is 13; thus, the intensity of β-tridymite is considered to be greater than the number of bonds.

Figure 21: Force of bonded parts in the Si model.

Figure 22: Comparison of fracture surface.

4. Conclusion

In this study, nonfiring models of ceramics were constructed using MD simulations with ReaxFF. Three types of silica interface models with water molecules were simulated in a reaction process and tensile. The simulation results clarify the following findings:

(1)An activated silica interface and water molecule reacted, and some bonds between interfaces were observed(2)When O atom was ended on the interface, an O-H-O bond was reproduced in the nonfiring process. The bond was a hydrogen bond. However, when an Si atom was ended on the interface, a Si-O-Si bond was produced in the nonfiring process, and this bond was siloxane bond(3)Bonded interface was able to have tensile strength. However, the tensile strength was weaker than that of the solid silica model. The non-bonded OH model did not have tensile strength

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This study was supported by Japan Science and Technology Adaptable and Seamless Technology Transfer Program through target-driven R&D (JST A-STEP).

References

A. Eiad-Ua, T. Shirai, T. Kato et al., “Novel fabrication route for porous ceramics using waste materials by non-firing process,” Journal of the Ceramic Society of Japan, vol. 118, no. 1380, pp. 745–748, 2010.

T. Shirai, A. Eiad-Ua, T. T. T. Hien, and M. Fuji, “Surface characterization of activated alumina powder through the mechano-chemical treatment for fabrication of non-fired ceramics,” Journal of the Ceramic Society of Japan, vol. 120, no. 1406, pp. 438–441, 2012.

T. T. T. Hien, T. Shirai, and M. Fuji, “Mechanochemical treatment of amorphous silica powder for glasses without firing,” Journal of Aamrican Ceramics Society, vol. 96, no. 12, pp. 3708–3711, 2013.

Y. Nakashima, H. Razavi-Khosroshahi, C. Takai, and M. Fuji, “Non-firing ceramics: Activation of silica powder surface for achieving high- density solidified bodies,” Advanced Powder Technology, vol. 29, no. 8, pp. 1900–1903, 2018.

Y. Nakashima, H. Razavi-Khosroshahi, H. Ishida, C. Takai, and M. Fuji, “Non-firing ceramics: activation of silica powder surface by a planetary ball milling,” Advanced Powder Technology, vol. 30, no. 2, pp. 461–465, 2019.

C. Miclea, C. Tanasoiu, I. Spanulescu et al., “Microstructure and properties of barium titanate ceramics prepared by mechanochemical synthesis,” Romanian Journal Of Information Science and Technology, vol. 10, no. 4, pp. 335–345, 2007.

MOLSIS, 2019, https://www.molsis.co.jp/materialscience/adf/reaxff/.

J. Tersoff, “New empirical approach for the structure and energy of covalent systems,” Physical Review B, vol. 37, no. 12, pp. 6991–7000, 1988.

Y. Yu, B. Wang, M. Wang, G. Sant, and M. Bauchy, “Revisiting silica with ReaxFF: towards improved predictions of glass structure and properties via reactive molecular dynamics,” Journal of Non-Crystalline Solids, vol. 443, pp. 148–154, 2016.

M. Pfeiffer-Laplaud, D. Costa, F. Tielens, M. P. Gaigeot, and M. Sulpizi, “Bimodal acidity at the amorphous silica/water interface,” The Journal of Physical Chemistry C, vol. 119, no. 49, pp. 27354–27362, 2015.

S. C. Chowdhury, B. Z. Haque, and J. W. Gillespie Jr., “Molecular dynamics simulations of the structure and mechanical properties of silica glass using ReaxFF,” Journal of Material Science, vol. 51, no. 22, pp. 10139–10159, 2016.

T. P. Senftle, S. Hong, M. M. Islam et al., “The ReaxFF reactive force-field: development, applications and future directions,” Computational Materials, vol. 2, no. 1, p. 15011, 2016.

J. Wen, T. Ma, W. Zhang et al., “Atomic insight into tribochemical wear mechanism of silicon at the Si/SiO₂ interface in aqueous environment: molecular dynamics simulations using ReaxFF reactive force field,” Applied Surface Science, vol. 390, pp. 216–223, 2016.

B. W. Ewers and J. D. Batteas, “Utilizing atomistic simulations to map pressure distributions and contact areas in molecular adlayers within nanoscale surface-asperity junctions: a demonstration with octadecylsilane-functionalized silica interfaces,” Langmuir, vol. 30, no. 40, pp. 11897–11905, 2014.

O. M. Roscioni, L. Muccioli, A. Mityashin, J. Cornil, and C. Zannoni, “Structural characterization of alkylsilane and fluoroalkylsilane self-assembled monolayers on SiO₂ by molecular dynamics simulations,” The Journal of Physical Chemistry C, vol. 120, no. 27, pp. 14652–14662, 2016.

National Institute for Materials Science Atom Works, 2017, http://crystdb.nims.go.jp/.

A. G. Grivtsov, L. T. Zhuravlev, G. A. Gerasimova, and L. G. Khazin, “Molecular dynamics of water: adsorption of water on β-tridymite,” Journal of Colloid and Interface Science, vol. 126, no. 2, pp. 397–407, 1988.

“ReaxFF User Manual,” 2019, https://www.scm.com/wp-content/uploads/ReaxFF-users-manual-2002.pdf.

A. C. T. van Duin, S. Dasgupta, F. Lorant, and W. A. Goddard, “ReaxFF: a reactive force field for hydrocarbons,” Journal of Physical Chemistry A, vol. 105, no. 41, pp. 9396–9409, 2001.

A. C. T. van Duin, A. Strachan, S. Stewman, Q. Zhang, X. Xu, and W. A. Goddard, “ReaxFFSiO Reactive force field for silicon and silicon oxide systems,” Journal of Physical Chemistry A, vol. 107, no. 19, pp. 3803–3811, 2003.

National Institute for Materials Science, “CompES-X,” 2018, http://compes-x.nims.go.jp/.

The chemical society of Japan, Kagakubinran Kisohen Kaitei 5 Han(Chemical Handbook, Basic, 5th), Maruzen, 2004.

International Baccalaureate, Chemistry book Diploma Programme(DP)Version2, International Baccalaureate Organization, 2014.