Intrinsic Defects and Their Role in the Phase Transition of Na-Ion Anode Na2Ti3O7

The development of high-power anode materials for Na-ion batteries is one of the primary obstacles due to the growing demands for their use in the smart grid. Despite the appealingly low cost and non-toxicity, Na2Ti3O7 suffers from low electrical conductivity and poor structural stability, which restricts its use in high-power applications. Viable approaches for overcoming these drawbacks reported to date are aliovalent doping and hydrogenation/hydrothermal treatments, both of which are closely intertwined with native defects. There is still a lack of knowledge, however, of the intrinsic defect chemistry of Na2Ti3O7, which impairs the rational design of high-power titanate anodes. Here, we report hybrid density functional theory calculations of the native defect chemistry of Na2Ti3O7. The defect calculations show that the insulating properties of Na2Ti3O7 arise from the Na and O Schottky disorder that act as major charge compensators. Under high-temperature hydrogenation treatment, these Schottky pairs of Na and O vacancies become dominant defects in Na2Ti3O7, triggering the spontaneous partial phase transition to Na2Ti6O13 and improving the electrical conductivity of the composite anode. Our findings provide an explanation on the interplay between intrinsic defects, structural phase transitions, and electrical conductivity, which can aid understanding of the properties of composite materials obtained from phase transitions.


Synthesis and physicochemical characterization
prepared as detailed in Ref. [1]. Na 2 Ti 6 O 13 was prepared using an identical synthetic route to Na 2 Ti 3 O 7 but precursors were calcined at 900 °C for 25 h.
Samples were characterised by powder X-ray diffraction (PXRD) at room temperature, using a Smartlab diffractometer (Rigaku Corporation) equipped with a 9 kW Cu rotating anode (λ = 1.54056 Å) operating in reflection mode with Bragg-Brentano geometry. Data were collected in the 5-70 ° 2θ range at a scan speed of 0.02 ° s -1 . Rietveld refinement of the Na 2 Ti 3 O 7 /Na 2 Ti 6 O 13 composite sample was performed from X-ray diffraction data, using the ICSD 250000 ICSD 23877 Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 structures as initial models. The GSAS-EXPGUI software interface 2, 3 was used to perform the Rietveld refinement. Peak shapes were modelled with a Gaussian-Lorentzian function, and the background, lattice parameters, atomic positions and thermal parameters were refined. The occupancy for all atoms was fixed to n = 1. For the structural refinement of the Na 2 Ti 3 O 7 /Na 2 Ti 6 O 13 composite, UISO values for identical atoms were constrained to be equal.
UV-visible spectroscopy was performed to determine the bandgap of the three titanate samples described in this work using a Cary500 spectrometer in the 250-600 nm, coupled with an integrating sphere, to acquire only the diffusive reflectance of the electromagnetic radiation.

Ab initio calculations
All the calculations performed in this work employed density functional theory (DFT) as implemented in the Vienna Ab initio Simulation Package code. 4,5 Interactions between core and valence electrons were described using the projector augmented wave (PAW) method. 6 The electron configurations Na (3s 1 ), Ti (3d 3 4s 1 ), and O (2s 2 2p 4 ) were treated as the valence electrons.
Convergence to plane wave energy was checked, with a cut-off of 500 eV found to be sufficient to converge the total energy to within 0.01 eV atom -1 . Brillouin zones for all compounds were sampled such that the k-points were converged in an accuracy of the total energy in 0.001 eV atom -1 . All calculations were deemed to be converged when the forces on all atoms were less than 0.01 eV Å -1 .
To accurately predict the phase stability of Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 , we used hybrid functional to recalculate enthalpies of phases of the Na-Ti-O system that are previously investigated with PBEsol functional. 1 In this calculation, the screened hybrid functional (HSE06), 7 in which 25% of exact non-local Fock exchange is added to the PBE 8 functional, is selected to investigate phase stability of 20 different phases of the Na-Ti-O system that are reported to be stable in the Materials Project database ( Figure S1a). 9 To take into account the temperature effect on the thermodynamic stability of phases, we employed vibrational entropies of nine phases (Na, Ti, Na 2 Ti 3 O 7 , Na 2 Ti 6 O 13 , Na 2 O, Na 4 TiO 4 , Na 8 Ti 5 O 14 , Na 4 Ti 5 O 12 , and TiO 2 ) that were calculated in our previous study. 1 For the gas phase of O 2 , the empirical gas phase thermodynamic data from the NIST data 10 are utilized instead to evaluate the temperature dependence of the gas free energies. In addition to these phases, we also calculated vibrational entropies of two phases of NaTi 2 O 4 , NaTi 5 O 10 that have chemical compositions close to Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 to better predict the chemical potential region of stability for Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 . This was achieved by generating 3×6×2 and 2×2×2 supercells of NaTi 2 O 4 and NaTi 5 O 10 , respectively, followed by calculating the harmonic force constants from atomic displacements of supercells using quasi-harmonic equations as implemented in Phonopy package. 11 However, the calculation of vibrational entropies for these large supercells is challenging with computationally intensive HSE06 functional. Furthermore, preliminary calculations revealed that HSE06 functional can largely overestimate the vibrational entropy of phases ( Figure S2). For this reason, we selected PBEsol functional to calculate vibrational entropies, for it has been proven to accurately reproduce lattice parameters and lattice dynamics in solid systems while maintaining a relatively low computational cost. 12,13 Based on the calculated enthalpy (H) and entropy (S) of phases, Gibbs free energies were then calculated as: where T is absolute temperature. The free energy of formation for Na x Ti y O z can then be estimated according to: Calculated free energies of formation were used to predict the chemical potential spaces that can be accessible during synthesis condition. This was done by comparing free energies of formation at a synthesis temperature of 1070 K, 1  using Python packages of pymatgen 15 and bsym 16 (Table 1). All constructed defective cells were geometrically optimized using HSE06 functional, where integrations over the Brillouin zone
Based on the predicted chemical potential regions and defect formation energies, the formation energy of a defect D in the charge state q at 0 K can be calculated according to: 17 (3) where is the energy of the defective supercell and is the energy of the host supercell. 17, , 18 The second term represents the energy change during removal/addition of defect atoms, i, where n i corresponds to the number of atoms exchanged, is the reference energy of a defect atom, and Δµ i indicates the chemical potential of the defect atom. In the third term, E f is the Fermi level and represents the energy needed to add or remove an electron from the VBM to the Fermi level. E corr is a correction term to account for the artificial electrostatic interaction due to periodic boundary conditions. 19,20 The calculation of defect formation energy at elevated temperature Elemental chemical potentials also differs from ones evaluated at 0 K ( ) because of the ∆ ∆ changes in the chemical potential stability region at elevated temperature. Among the free energy contributions in eq. (4), we assumed the difference in vibrational entropy arising from a point defect in a supercell structure, i.e., , to be zero, considering that the entropy ,changes for adsorption and desorption of a single atom is negligible. 22,23 . The other free energy terms from elemental and competing phases were calculated using quasi-harmonic approximation as implemented in Phonopy, 11 except O 2 gas of which entropy was obtained from NIST thermodynamic database. 10 The thermodynamic transition level, , which represents the defect formation ( , / ') energies at which a defect with differing charge states, (q, q'), have the same formation energy was predicted according to: Plotting thermodynamic transition level of all defects and their formation energy along Fermi energy, gives transition level diagram, which conveys overall defect chemistry of materials. For a given temperature, defect and charge carrier concentrations are given as a function of Fermi energy and defect formation energy and thus, varies towards the x-axis of transition level diagram. Thus, by calculating charge carrier concentrations for all defects along Fermi energy, one can predict the SC Fermi level in the middle of band gap where overall charge carriers satisfy the charge S5 neutrality condition, which was done in this study using the SC-FERMI code. 24 In this calculation, defect concentrations at the SC Fermi level were calculated for two selected chemical potential sets of A and B in Figure 1a at synthesis temperature of 1070 K.
The temperature dependent electrical conductivity of Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 were calculated using the AMSET package. 25 Table S5. Figure S1. (a) Na-Ti-O phase diagram defined by Na, Ti and O obtained from structures with the lowest formation enthalpy calculated using HSE06 functional. Target phases of Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 are highlighted by red rhombus, whereas competing phases with compositions close to target phases are selected for high temperature stability and denoted by orange rhombus. Among 17 phases denoted in (a), Na 2 Ti 3 O 7 is found to be unstable at athermal limit.           Tables   Table S1. List of lattice parameters (a, b, c) and angles (α, β, γ) of Na 2 Ti 3 O 7 and Na 2 Ti 6 O 13 optimized using PBEsol and HSE06 functionals. The lattice parameters a, b, and c are given in Å, whereas angles α, β, and γ are given in degrees. All Table S2. Thermodynamic transition levels (ε) of native defects in Na 2 Ti 6 O 13 obtained with respect to VBM. Note that most defects show deep ε(0/-) or ε(0/+) transition levels similar to those predicted from Na 2 Ti 3 O 7 , preventing Na 2 Ti 6 O 13 to show polaronic conduction at room temperature. One exception is Ti Na with its ε(0/+) transition level being positioned above CBM, but it is less likely to contribute to polaronic conduction as Ti Na has the highest formation energy among all other donor defects (see Figure S5).

Additional notes
In addition to its capability to predict the defect concentrations at synthesis temperature, defect formation energies can also provide insights on the number of mobile Na ions in sodium titanate anodes. 30 In case of Na 2 Ti 3 O 7 , the Na stoichiometries calculated under the "frozen-in" assumption are x Na,n = 1.92 and x Na,p = 1.93 for oxygen-poor and oxygen-rich conditions (Figure 4), respectively, the number of which are often considered to be the mobile charge carrier populations. However, the mobility of Na + ions can be affected by Coulombic attractions and configurational entropies of other defect species, which may alter the concentration of mobile Na + ions.

Effect of Na Ti antisites
Among intrinsic defects calculated in this study, one possible case that can lower the

Effect of clustering of Na i and V Na
Another factor that influences the mobility of Na is the clustering of Na interstitial and vacancies owing to their Coulombic attraction to each other. If this Coulombic interaction surpass the tendency for defects to be disordered state (manifested by configurational entropy), the Na i and V Na will bind together to form defect pairs and lowers the mobility of Na.
The binding energy of the two defects can be estimated by the difference in energy before and after dissociation, which can be done by either comparing the energies between the closest and farthest Na i -V Na pairs or calculating dissociation energy of Na i -V Na pairs  Na i + V Na . The preliminary calculations using PBEsol functional showed that the shortest Na i -V Na distance without recombination during geometry optimization was 4.0 Å whereas the largest separation of Na i -V Na pair possible in a 72-atom-supercell was 7.3 Å. Subsequent structural optimizations with HSE06 functionals revealed the difference in energy between these two supercells, S18 dissociation energy of Na i -V Na pair, was 1.01 eV. On the other hand, the dissociation energy calculated from Na i -V Na pairs  Na i + V Na gives 1.78, 2.09 and 3.09 eV for the dissociated Na interstitial in site 1, 2 and 3, respectively ( Figure 2). Overall, the calculated dissociation energies are lower than the formation energies of individual Na i and V Na (2.10 -2.75 eV) when Na interstitials are in between Ti-O layers (Na i,1 and Na i,2 in Figure 2), suggesting that Na i and V Na located in between Ti-O layers of Na 2 Ti 3 O 7 act as independent defects with disordered configurations. However, when Na interstitial is trapped inside Ti-O layers (Na i,3 in Figure 2), these interstitials are difficult to be removed and thus, prefers the formation of Na i -V Na pair and impedes the Na migrations.
Above analyses on the Na mobility suggests that the clustering of Na i -V Na pair in Na 2 Ti 3 O 7 may reduce the number of mobile Na ions. The reduction in the Na mobility can lead to the fewer generation of V Na,2 -V O,5 Schottky pairs, which results in the less phase transition from Na 2 Ti 3 O 7 to Na 2 Ti 6 O 13 . From this perspective, we believe that the amount of precipitated secondary Na 2 Ti 6 O 13 phase during synthesis will be lower than those (7-8 %) predicted from Na stoichiometries in Figure 4.