Development of a ReaxFF Reactive Force Field for Titanium Dioxide / Water Systems
Sung-Yup Kim1, Nitin Kumar2, Petter Persson3, Jorge Sofo2, Adri C. T. van Duin1,(*) and James D. Kubicki4 1
Department of Mechanical and Nuclear Engineering, The Pennsylvania State University, University Park, PA 16802 2
Department of Physics, The Pennsylvania State University, University Park, PA 16802 3
Department of Theoretical Chemistry, Lund University, Box 124, S-221 00 Lund, Sweden
4
Department of Geosciences and the Earth & Environmental Systems Institute, The Pennsylvania State University, University Park, PA 16802
*Corresponding author. Email:
[email protected] Bulk modulus The bulk modulus of a substance measures its resistance to uniform compression. Based on the EOS, bulk moduli were fitted against these results by using the Rose equation-of-state4 for the three most important phases. As is shown in Figure S1, ReaxFF results reproduce the experimental values within 17% error. Depending on the DFT or experimental methods, the bulk moduli vary within those values. Considering that, a 17% discrepancy is reasonable. Also, the trend from DFT results is also 1
reproduced by ReaxFF as rutile has highest value followed by brookite and anatase. This can be understood by considering the atomic structure of these phases. In other words, the atoms for rutile are packed in unit cell more densely compared to anatase and brookite, which is responsible for the higher zero-pressure density of rutile.
Bulk Modulus 400 350
ReaxFF QM (Muscat et al.,2004, Beltran et al.,2006)
Bulk Modulus (Gpa)
Expt.(Gerward et al.,1997, Swamy et al.,2001, Luo et al.,2005)
300 250 200 150 100 50 0
1
Rutile
2
Anatase
3
Brookite
Figure S1: The bulk moduli of rutile (Expt = 230 GPa3), anatase (Expt = 178±1 GPa2) and brookite (Expt = 255 GPa3) as calculated with DFT (Muscat et al.5, 2004 for rutile and anatase); Beltran et al.6, 2006 (for brookite ) and ReaxFF (this work).
Energy differences between the structures
2
Energy difference (kcal/mol)
35 ReaxFF QM (Muscat et al.,2004) (Beltran et al.,2006)
30 25 20 15 10 5 0
Anatase 1
Brookite 2
Columbite 3
Baddeleyite 4
Pyrite 5
6 Fluorite
7 Cotunnite
Figure S2: Energy differences between the structures (reference: rutile) Surface formation energies Figure S3 shows surface energies of 5 different faces, 3 for anatase and 2 for rutile. ReaxFF agreement with DFT is as good as the agreement between DFT calculations from different authors. As is shown above, anatase (101) has the lowest value among the anatase surfaces and rutile (110) has the lower value compared to rutile (101), which makes sense because the most common faces are (101) for anatase and (110) for rutile. The discrepancy between ReaxFF and Hummer et al. reaches up to 57%. Interestingly, ReaxFF overestimates the surface formation energy for anatase, while underestimating it for rutile.
3
Surface Energy 1.4 ReaxFF Barnard et al.,2004 Hummer et al.,2009
2
Surface energy (J/m )
1.2 1 0.8 0.6 0.4 0.2 0
1 (001) Anatase
2 (100) Anatase
Anatase (001)
3 (101) Anatase
Rutile4(110)
Anatase (100)
Rutile (110)
Rutile5(101)
Anatase (101)
Rutile (101)
Figure S3: Surface formation energies and structures which were used in the DFT references for the 3 anatase and 2 rutile surfaces
4
Water dissociation preference dependence on slab thickness
Liu et al.1 demonstrated that the calculated dissociation preference of water adsorbed on the rutile (110) surface can depend on the thickness of the slab used in the simulation. Figure S4 shows the energy difference between fully associated and partially dissociated configuration for a 2x1 surface unit cell with 1 monolayer (ML) of water as the thickness of the slab changes from 3 to 15 Ti-O layers. A positive or negative energy difference corresponds to fully associated or partially dissociated configuration being preferred, respectively. Figure S4 shows that the DFT1 predicts fully associated configuration as the preferred method for water adsorption as the thickness of the slab is increased. Our ReaxFF results also predicts fully associated configuration for the adsorbed water molecules on the rutile (110) surface. ReaxFF predicts much smaller fluctuation, than the DFT, in the energy difference which is around 1.1~2.2 kcal/mol per H2O. This can be due to electronic structure effect in the DFT calculations that is not present in the ReaxFF simulations7.
4 Liu et al.,2010 ReaxFF
Energy (kcal/mol)
3
2
1
0
-1
-2 2
4
6
8 10 Number of Layers
12
14
16
Figure S4: Energy difference between fully associated and partially dissociated configuration, for a 2x1 surface unit cell simulation box with 1 ML water, with respect to number of Ti-O layers calculated from DFT1 and from ReaxFF. DFT and our ReaxFF both predict fully associated configuration as the preferred method of water adsorption on the rutile (110) surface for thicker layers. 5
TiO2 clusters A total of 150 bare and hydrated TiO2 cluster structures, mainly based on stoichiometric TiO2 cluster cores of the type (TiO2)x with x in the range 1 – 28, were included in the optimization of the force field. Suitable cluster structures were identified from the literature on small TiO2 clusters8, 9 and model TiO2 nanocrystals10. Structural isomers for some of the larger clusters, e.g. (TiO2)16 and (TiO2)28, were initially generated during the force field development during optimizations, low-temperature MD simulations, or simulated annealing with preliminary force fields, and subsequently calculated by quantum chemical calculations and incorporated into the force field fitting procedure. Quantum chemical cluster data was generally either taken from the literature, or obtained using density functional theory cluster calculations performed with the Jaguar 7.0 program [Jaguar 7.0; Schrodinger Inc.: Portland, OR, 2007]. using the B3LYP hybrid functional11, 12 and a LACVP** basis set and effective core potential (ECP) combination that provides an all electron 6-31G** basis set for the H and O atoms13, 14, while the Ti atoms are treated using a Hay-Wadt ECP15 for the inner core (Ne) electrons together with a double zeta description of the 12 remaining outer core and valence electrons. This level of quantum chemical theory is consistent with that recently employed for the development of ReaxFF force fields for closely related metal oxide systems such as vanadium oxides16. Two types of 3(TiO2) and 4(TiO2) were picked and bond lengths and angles are presented for comparison between ReaxFF and QM. Also, potential energy differences between pairs of clusters and water binding energies for the selected structures are shown in the Table 10~12. Figure S5 (d) to (g) shows the 4 configurations of TiO2 clusters which were chosen for the comparison between ReaxFF and QM results. Measured bond lengths and angles of selected clusters are shown and comparison is made to the quantum-minimized clusters for each case. Overall, agreement is good (± 0.1 Å and ± 1.2°) for cluster structures, potential energies and water binding energies except in a few cases.
6
(a) 2(TiO2)-type I
(b) 2(TiO2)-type II
(c) 3(TiO2)-type IV
[Config.2 in Figure S6]
A
A C B
(d) 3(TiO2)-type(I)
D
D C
B
C B
D
A
(e) 3(TiO2)-type(II)
(f) 4(TiO2)-type(I)
(g) 4(TiO2)-type(II)
(h) 10(TiO2)-type I
(i) 10(TiO2)-type II
[Config.11 in Figure S6]
[Config.14 in Figure S6]
[Config.5 in Figure S6]
C
B D
A
(j) 12(TiO2)-type II [Config.24 in Figure S6]
(k) 12(TiO2)-type V
(l) 15(TiO2)-type II [Config.28 in Figure S6] 7
(m) 16(TiO2)-type I [Reference Structure in Figure S6] Figure S5: Configurations of the TiO2 clusters Figure S6 shows that the potential energy differences between pairs of structures for Ti-O clusters. The reference structure is (m)16(TiO2)-type(I) in Figure S5 and the differences between the energies of 29 different structures divided by numbers of TiO2 units and that of reference structure divided by 16 are shown. Agreement is good overall. The potential energy difference between QM and ReaxFF decreases as the clusters increase in size, which means that large-scale simulations of TiO2-H2O using ReaxFF should be accurate.
Potential energy difference (kcal/mol)
120 ReaxFF
100
QM 80 60 40 20 0 -20
0
5
10 15 20 Configurations of the TiO2 clusters
25
30
8
Figure S6: Potential energy difference between the selected configurations and the reference structure[(m), in Figure S5]
(a) (TiO2) (H2O)
(d) 4(TiO2)5(H2O)
(b) 2(TiO2)2(H2O)
(c) 4(TiO2)2(H2O)
(e) 2(TiO2)3(H2O)
Figure S7: Configurations of the hydroxide TiO2 clusters for the water binding energy
9
P o te n tia l e n e rg y (k c a l/m o l)
-1.08
x 10
5
-1.082 -1.084 -1.086 -1.088 -1.09 -1.092 -1.094 -1.096 -1.098 0
5
10
15
20 25 Time (ps)
30
35
40
(a) 1ML
P o te n tia l e n e rg y (k c a l/m o l)
-1.41
x 10
5
-1.415 -1.42
-1.425 -1.43
-1.435 0
2
4
6 8 Time (ps)
10
12
14
(b) 3ML Figure S8: Energy conservations with respect to time for 1ML and 3ML of rutile (110)
10
References 1.
Liu, L. M.; Zhang, C. J.; Thornton, G.; Michaelides, A., Structure and dynamics of liquid water
on rutile TiO(2)(110). Phys. Rev. B 2010, 82, (16), 4. 2.
Swamy, V.; Gale, J. D.; Dubrovinsky, L. S., Atomistic simulation of the crystal structures and
bulk moduli of TiO2 polymorphs. J. Phys. Chem. Solids 2001, 62, (5), 887-895. 3.
Luo, W.; Yang, S. F.; Wang, Z. C.; Wang, Y.; Ahuja, R.; Johansson, B.; Liu, J.; Zou, G. T.,
Structural phase transitions in brookite-type TiO2 under high pressure. Solid State Commun. 2005, 133, (1), 49-53. 4.
Li, J. H.; Liang, S. H.; Guo, H. B.; Liu, B. X., Four-parameter equation of state and
determination of the thermal and mechanical properties of metals. J. Alloy Compd. 2007, 431, (1-2), 2331. 5.
Muscat, J.; Swamy, V.; Harrison, N. M., First-principles calculations of the phase stability of
TiO2. Phys. Rev. B 2002, 65, (22), 15. 6.
Beltran, A.; Gracia, L.; Andres, J., Density functional theory study of the brookite surfaces and
phase transitions between natural titania polymorphs. J. Phys. Chem. B 2006, 110, (46), 23417-23423. 7.
Kowalski, P. M.; Meyer, B.; Marx, D., Composition, structure, and stability of the rutile
TiO2(110) surface: Oxygen depletion, hydroxylation, hydrogen migration, and water adsorption. Phys. Rev. B 2009, 79, (11). 8.
Albaret, T.; Finocchi, F.; Noguera, C., Density functional study of stoichiometric and O-rich
titanium oxygen clusters. J. Chem. Phys. 2000, 113, (6), 2238-2249. 9.
Jeong, K. S.; Chang, C.; Sedlmayr, E.; Sulzle, D., Electronic structure investigation of neutral
titanium oxide molecules TixOy. J. Phys. B-At. Mol. Opt. 2000, 33, (17), 3417-3430. 10.
Lundqvist, M. J.; Nilsing, M.; Persson, P.; Lunell, S., DFT study of bare and dye-sensitized
TiO2 clusters and nanocrystals. Int. J. Quantum Chem. 2006, 106, (15), 3214-3234.
11
11.
Becke, A. D., DENSITY-FUNCTIONAL THERMOCHEMISTRY .3. THE ROLE OF EXACT
EXCHANGE. J. Chem. Phys. 1993, 98, (7), 5648-5652. 12.
Lee, C. T.; Yang, W. T.; Parr, R. G., DEVELOPMENT OF THE COLLE-SALVETTI
CORRELATION-ENERGY FORMULA INTO A FUNCTIONAL OF THE ELECTRON-DENSITY. Phys. Rev. B 1988, 37, (2), 785-789. 13.
Radom, L.; Hariharan, P. C.; Pople, J. A.; Schleyer, P. V., MOLECULAR-ORBITAL THEORY
OF ELECTRONIC-STRUCTURE OF ORGANIC COMPOUNDS .22. STRUCTURES AND STABILITIES OF C3H3+ AND C3H+ CATIONS. J. Am. Chem. Soc. 1976, 98, (1), 10-14. 14.
Francl, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; Gordon, M. S.; Defrees, D. J.; Pople, J.
A., SELF-CONSISTENT MOLECULAR-ORBITAL METHODS .23. A POLARIZATION-TYPE BASIS SET FOR 2ND-ROW ELEMENTS. J. Chem. Phys. 1982, 77, (7), 3654-3665. 15.
Hay, P. J.; Wadt, W. R., ABINITIO EFFECTIVE CORE POTENTIALS FOR MOLECULAR
CALCULATIONS - POTENTIALS FOR K TO AU INCLUDING THE OUTERMOST CORE ORBITALS. J. Chem. Phys. 1985, 82, (1), 299-310. 16.
Chenoweth, K.; van Duin, A. C. T.; Persson, P.; Cheng, M. J.; Oxgaard, J.; Goddard, W. A.,
Development and application of a ReaxFF reactive force field for oxidative dehydrogenation on vanadium oxide catalysts. J. Phys. Chem. C 2008, 112, (37), 14645-14654.
12