Density Functional Tight Binding Theory Approach for the CO2 Reduction Reaction Paths on Anatase TiO2 Surfaces

Herein, we have investigated the CO2 reduction paths on the (101) anatase TiO2 surface using an approach based on the density functional tight binding (DFTB) theory. We analyzed the reaction paths for the conversion of carbon dioxide to methane by performing a large number of calculations with intermediates placed in various orientations and locations at the surface. Our results show that the least stable intermediate is CO2H and therefore a key bottleneck is the reduction of CO2 to formic acid. Hydrogen adsorption is also weak and would also be a limiting factor, unless very high pressures of hydrogen are used. The results from our DFTB approach are in good agreement with the hybrid functional based density functional theory calculations presented in the literature.


INTRODUCTION
The phenomenon of "global warming" caused by emission of CO 2 to the earth atmosphere puts a demand to achieve a fossilfree economy and make use of sustainable energy sources such as sun, wind, geothermal, and hydrothermal. 1 The importance of such strive has been stressed in the recently announced Paris agreement in which the United Nations Framework Convention on Climate Change (UNFCCC) is dealing with greenhouse gas emission mitigation. CO 2 reduction by semiconductor catalyzers 2,3 and/or electrochemical 4,5 and photoelectrochemical methods 6,7 along with direct capture of CO 2 from air 8 are all well-known technologies that can both reduce the CO 2 amount in the atmosphere and produce valuable chemicals or fuels. 9,10 They can be classified into the so-called "artificial photosynthesis" because the processes mimic plants, which convert CO 2 into O 2 and store energy in useful hydrocarbons. Here, the overall reaction for conversion of carbon dioxide to methane can be represented by the following reaction (1) However, the exact mechanism of the conversion of CO 2 to methane is not yet understood and is also dependent on the reaction medium, the catalyst, and its interfaces, as well as other experimental conditions. 11−14 The reported kinetics so far are slow and the efficiencies are rather poor even for electrocatalytically boosted reactions. Therefore, there is a clear need for new experimental work and theoretical simulations that allow us to understand and predict new roots and ways to optimize CO 2 reduction devices.
TiO 2 is one of the most attractive and widely investigated materials in this field 13,15,16 because of its chemical stability, biocompatibility, 17,18 photocatalytic ability, and abundance. 19 The ongoing research for the CO 2 reduction on anatase titanium dioxide is exploring different reaction paths, 16 as well as impacts of catalyst, 20 semiconductor crystal planes and surface chemistry, 21 rate-limiting steps, 16 and particle size 22 in order to quantify, understand, and develop the CO 2 reduction process. Here, theoretical simulations can provide valuable and detailed information about the mechanisms of reduction processes. In a recent publication, such mechanisms, when taking place at the pristine anatase (101), were explored by means of density functional theory (DFT). 16 Although accurate and widely transferable, DFT is typically too computationally demanding to widely screen a large number of adsorption geometries of structures such as nanoparticles or complex surfaces with multiple defects, surface kinks, impurities, and so forth. Alternative, approximate methods such as density functional tight binding (DFTB) are capable of targeting such systems, but their accuracy and limitations need to be established for the current application. In this paper, we address this issue and validate the quality of DFTB for the methanization reaction on the anatase (101) surface. Herein, we have employed a rather unique strategy to combine lowcost computational methods with standard DFT simulations to study the details of the above reaction. Our approach is especially useful because we can optimize the "computational costs". This will allow us to expand the "implemented framework" and tackle systems of much larger complexity. As an example, we have analyzed the most critical steps in the reaction in an anatase nanoparticle model, which cannot be easily pursued by a similar DFT approach.

METHODS
2.1. DFTB Calculations. The DFTB+ code 23 in its implementation with self-consistent charges (SCC) 24 has been used for all DFTB calculations. The implemented anatase (101) surface slab consisted of 96 atoms and a 10 Å vacuum gap. The following molecules have been considered on the anatase surfaces according to the possible reaction paths for the conversion of carbon dioxide to methane: CO, CH 3 , CH 3 OH, CH 4 , COH, COOH, CH 2 , COHOH, CHOH 2 , COH 3 , COH 2 H 2 , COHH, COOH 2 , H, and H 2 . Ti 3d-orbitals, O 2p-orbitals, C 2p-orbitals, and H 1s-orbitals have been considered for the tight binding basis sets. A Hubbard model with U = 3.5 eV was applied in order to correctly describe the localized d-electrons of Ti. The force and self-consistent charge conversion thresholds for relaxing the atomic coordinates and  the self-consistent cycles were 0.0002 H/Bohr and 0.001 H, respectively. Spin of the electrons has been taken into account by comparing the results using different numbers of unpaired electrons in the system. For each relevant molecule, many different possible adsorption configurations on the anatase surface were generated by varying the angle and distance of the molecules on top of the various sites at the surface. We also considered dissociative adsorption configurations. In total, over 500 different configurations have been optimized.
The final DFTB results have been compared to a previously reported DFT study using the Heyd−Scuseria−Ernzerhof (HSE) hybrid functional. We also performed our own DFT calculations using the GGA+U functional. These latter simulations were done with the VASP software 25−28 using the Perdew−Burke-Ernzerhof (PBE) functional 29,30 with an energy cutoff of 400 eV. We used the standard projector augmented wave 31,32 potential provided with the VASP software.
In order to allow a more feasible exploration of the various adsorption configurations, we developed a small octave software that assembles initial structures for calculations of adsorbates on a surface. The primary input to the code is the structure of the molecule to be adsorbed and the structure of the surface slab in the VASP POSCAR format. An origin and two vectors describing the direction of the molecule have to be provided. The code will then automatically place the origin on top of the designated points at the surface at a specified distance and align the first vector along the z-direction and the second along the x-axis of the slab system. Next, the code will rotate the molecule along the z-axis in steps specified by the user. At each step, a subfolder containing a POSCAR-file of the slab with the adsorbed molecule will be created. The octave code, including a simple example, is provided in the Supporting Information.

RESULTS AND DISCUSSION
The anatase TiO 2 nanoparticles, which are typically implemented in energy and catalytic applications, are predominantly exposing the most thermodynamically stable (101) facets along with a small fraction of the (001) 33,34 facet. Here, the periodic slab models are adopted to simulate the adsorption at the majority (101) anatase facets. The full reaction network including all considered intermediates in the conversion of CO 2 to methane is presented in Figure 1 and the most favorable pathway is given in Figure 2. Figure 3 shows a comparison of the SCC−DFTB data to the HSE−DFT data of ref 16. These recent HSE−DFT simulations focused on the photocatalytic CO 2 reduction where the source of hydrogen is not H 2 but rather H + in the solution.
In order to make a comparison, we have therefore used a common reference point for hydrogen in our figures/ comparison, namely, that of an isolated hydrogen atom bound to the surface. In the photocatalytic reaction, this specie is formed from a proton (H + ) in the solution and a photogenerated electron, while in our scheme, it originates from a dissociated H 2 molecule. For typical intermediates, our DFTB data show a very good agreement with HSE−DFT data, which is a proof of concept for our DFTB approach (see Figure 3).
Moreover, another possible reaction path was derived from our calculations in which carbon is an intermediate product.
Such a pathway can be relevant for experiments in the gas phase such as the work of Dimitrijevic et al. 35 Our data also show that the conversion of COH 2 to C (reaction 3) is associated with prohibitively high energies in some of the intermediates.

Adsorption Configurations and Adsorption
Energies. The optimized geometries of intermediates in the most favorable pathway are given in Figure 2. DFTB and DFT both agree that the first bottleneck in the methanization reaction is likely to be the weak adsorption of H 2 on the surface. This means that, even though the dissociation of H 2 is highly exothermic, we are typically more likely to desorb the H 2 molecule than to dissociate it. The problem is likely alleviated when using a metal component in conjunction with TiO 2 as is typically done in the experiments reported thus far or otherwise by using a high hydrogen gas pressure.
The accuracy of the DFTB calculated adsorption energies are high and are in agreement with the previous reports. For example, the adsorption energies of CO 2 , CO, and H 2 were estimated to be 0.439, 0.200, and 2.01 eV by DFTB in agreement with the reported values of 0.481, 0.262, and 2.3 eV, respectively, in refs, 15,18 (see Table 1). 3.3. CO 2 Reduction on Anatase Nanoparticles. Finally based on our DFTB approach, one can consider the reduction of CO 2 to the formic acid, which is a key bottleneck on the (101) facet, on more realistic nanoparticle models. Such simulations are feasible with our SCC−DFTB approach. Our calculation shows that the energy of the CO 2 H intermediate can be reduced ever so slightly when taking place at the tip of a small nanoparticle (see Figure 4). However, this is merely one of many possible "nanomotifs", and the methodology presented here could be used to explore many more configurations in the future. As discussed in the last two sections, we believe that conclusions regarding the stability of intermediates on such nanoparticle models should be reliable and that they could shed more light into the CO 2 reduction pathways on titania surfaces.

CONCLUSIONS
In this paper, we have investigated the reaction path for CO 2 reduction on the TiO 2 anatase (101) surface using DFTB method. The results show that DFTB can be utilized as a "computationally cheap" and "scientifically insightful" method to explore the reaction paths of CO 2 reduction on anatase TiO 2 surfaces. We have studied and optimized about 500 configurations of different relevant adsorbed hydrocarbons on the anatase (101) surface using the DFTB method in which the computational time is two orders of magnitude smaller than that with standard DFT. The calculated adsorption energies and predicted reaction path are in good agreement with the results of hybrid functional-based DFT presented in the literature. Based on our calculations, for the gas-phase reaction path, two bottlenecks were clearly identified: the week adsorption of H 2 on the anatase surface and the high energy barrier of CO 2 conversion to formic acid. The DFTB approach shows a good reliability when compared to conventional computationally expensive approaches for carbon dioxide chemical conversion paths on anatase titania. We emphasize here that the calculations presented in this work are most relevant for the gas-phase reactions where there is no applied voltage, incoming light, solvent molecule, or additives present.
OCTAVE_CODE_WITH_EXAMPLE.zip: Octave code to assemble the initial structures for calculations of adsorbates on a surface with an example (ZIP)