Phase separation in amorphous tantalum oxide from first principles

The transition between Ta2O5 and TaO2 governs resistive switching in tantalum oxide-based resistive random access memory. Despite its importance, the Ta2O5–TaO2 transition is scarcely described in the literature, in part because the tantalum oxide layer in devices is amorphous, which makes it difficult to characterize. In this paper, we use first-principles calculations to construct the convex hull of the amorphous Ta2O5−x system for 0 ≤ x ≤ 1 and show that oxygen deficiency in tantalum oxide leads to phase-separation into Ta2O5 and TaO2. In addition, our work challenges the conventional interpretation of X-ray Photoelectron Spectroscopy (XPS) spectra of the Ta 4f orbitals. Specifically, we find that TaO2 exhibits both the Ta4+ peak associated with TaO2 and the Ta5+ peak normally associated with Ta2O5. While our simulated Ta2O5 peak originates from a narrow range of oxidation states, the TaO2 peak comes from disproportionated Ta atoms with Bader charges ranging from +3 to +1, the lowest of which are well below Ta atoms in crystalline TaO. Finally, we demonstrate that the XPS blueshift of around 1 eV observed experimentally in amorphous Ta2O5 with respect to crystalline Ta2O5 comes from both the presence of under-coordinated Ta atoms and longer Ta–O bond distances in the amorphous system. Our simulated XPS analysis shows that amorphous XPS spectra may be more complex than previously thought, and hence, caution should be applied when assigning XPS peaks to oxidation states.


I. INTRODUCTION
Resistive switching devices constitute an important research topic within the general area of random access memory (RAM) technology. 1,2 Since Chua proposed the "memristor"-memory resistor-in 1971, 3 researchers have demonstrated memory resistive properties in metal/metal-oxide/metal devices 4,5 using semiconducting and insulating transition metal oxides such as TiO 2 , 4 HfO 2 , 6,7 and Ta 2 O5 8 as the active layer. The switching mechanism in these layers is initiated by applying a large electric field across the device. This leads to the creation and subsequent migration of oxygen vacancy (vO) defects and eventually to the growth of nanoscale filaments of oxygen-deficient and conductive regions through the device. [9][10][11] Tantalum pentoxide (Ta 2 O5) is one of the key candidate materials for the switching layer of memristive switching devices. In a direct comparison to TiO 2 , which was the first material to be linked to memristive switching, 4 Ta 2 O5 shows five times greater ionic mobility 4,12 and, by extension, may yield greater switching speed and lower power consumption. 13 Additionally, a switching endurance exceeding 10 12 cycles has been demonstrated without device breakdown. 14 Structurally, Ta 2 O5 thin films fabricated at low to room temperature are amorphous, [15][16][17] which makes them difficult to study on a theoretical level. Nevertheless, previous ab initio Molecular Dynamics (AIMD) studies based on Density Functional Theory (DFT) have provided in-depth investigations of amorphous tantalum oxide. Utilizing the atomic-level information offered by DFT, researchers have described atomic arrangements surrounding singular vO's, 18 density-of-states spectra of several TaOx stoichiometries, 19 electrical conductivity of TaOx, 20 and formation energies of vO's in Ta 2 O5 close to a Ta-Ta 2 O5 interface. 21 To our knowledge, however, neither experiment nor simulation has thus far attempted to map the Ta 2 O5-TaO 2 transition, which governs resistive switching in TaOx devices.
One reason for this is that determining the exact composition of amorphous materials is a non-trivial task, which is highly dependent on fabrication parameters. 17 Several studies have applied X-ray Photoelectron Spectroscopy (XPS) to obtain binding energies (BEs) especially for the Ta 4f orbitals. 17,[22][23][24][25][26][27] Usually, the resultant spectra are deconvoluted using Gaussian/Lorentzian functions 17,24,26,27 assigned either to different oxidation states of Ta or to different stoichiometries. The two assignments can be used interchangeably since the stoichiometries Ta 2 O5, TaO 2 , Ta 2 O 3 , TaO, Ta 2 O, and Ta correspond to the oxidation states +5, +4, +3, +2, +1, and 0, respectively. A compositional analysis can be performed based on the relative areas underneath the fitting functions. However, the XPS spectra of amorphous materials are more complex compared to crystalline materials, which exhibit well-defined patterns.
In this study, we zoom in on the transition between Ta 2 O5 and TaO 2 , which is the range of compositions where resistive switching takes place. We use DFT to investigate amorphous Ta 2 O5−x (a-Ta 2 O5−x), sampling 14 different snapshots for each of six different compositions with 0 ≤ x ≤ 1. We draw the convex hull for the Ta 2 O5-TaO 2 transition and confirm that the hull consists of these two stoichiometries, with Ta 6 O 13 being another possible meta-stable intermediate phase. Furthermore, we simulate XPS spectra for each stoichiometry and show that the XPS spectrum of TaO 2 shows both Ta 5+ and Ta 4+ peaks, the latter of which comes from disproportionated atoms with a wide range of oxidation states. Finally, by fitting Gaussian components to our oxidation state spectra and our XPS spectra, we establish a linear correlation between the Gaussian peak locations. These conclusions aid the understanding of XPS results of a-Ta 2 O5 and amorphous materials in general.
Amorphous systems of Ta 2 O5−x, x ∈ {0.00, 0.33, 0.50, 0.67, 0.80, 1.00}, were generated using AIMD. Beginning with the primitive cell of λ-Ta 2 O5 suggested by Lee et al., 34 we created a 2 × 2 × 3 orthogonal supercell containing 48 Ta atoms and 120 O atoms, the same size as those used by Bondi et al. 18,20 and Guo and Robertson. 35 Oxygen atoms were subsequently deleted randomly to obtain 8, 12, 16, 20, and 24 vacancies. We expanded each lattice vector by a factor of 1.06, which was chosen to match the experimental density for a-Ta 2 O5 of 6.88 × 10 3 kg/m 3 . 36 The expanded cell vectors were kept fixed throughout the duration of the MD simulations.
The amorphous structure generation used the following meltand-quench procedure: (1) equilibration at 5000 K for 10 ps, (2) further equilibration and extraction of in total 14 structure snapshots per stoichiometry in intervals of 2 ps following the initial equilibration, (3) quenching of each snapshot from 5000 K to 300 K, and (4) equilibration for 2 ps at 300 K. The MD simulations used a low-accuracy/high-speed setup with a cutoff energy of 300 eV and sampling only the Γ point. The time step was 1 fs for equilibrations and 2 fs for quenching.
Following the AIMD procedure, we increased the cutoff energy to 520 eV and k point sampling to (2 × 2 × 2) and performed full DFT relaxations of ionic positions and lattice vectors. Despite the relaxation of the lattice vectors, the average density of the final Ta 2 O5 geometry was still similar to the experimental density (see Fig. S1 of the supplementary material). Atomic forces were converged to within 0.025 eV/Å, which resulted in little change to the cell volume and thus density. We performed a high-precision calculation of the optimized structures using the hybrid functional of Heyd-Ernzerhof-Scuseria (HSE06). [37][38][39] Due to their high computational demands, these calculations sampled only the Γ point. Following these calculations, the BEs of the Ta 4f orbitals were computed using the initial state approximation, which involves recalculating the Kohn-Sham eigenvalues of the core states after a self-consistent calculation of the valence charge density. The HSE06-calculated wave functions served as the basis for a Bader charge analysis 40 using the program developed by Henkelmann et al. 41 In addition to our amorphous systems, we performed PBE relaxation followed by HSE06 calculations of the following crystalline structures: λ-Ta 2 O5, rutile TaO 2 , corundum Ta 2 O 3 , rocksalt TaO, and BCC Ta. These calculations used the same setup as the DFT calculations of the amorphous systems with a k point sampling density of at least 3.60 k points per Å −1 .

A. Energetic stability and structural analysis
We first examine the phase change between Ta 2 O5 and TaO 2 as well as differences in the Ta-O bonding. Figure 1 presents the convex hull, cumulative distribution function (CDF) for Ta-O bonds limited to the sixth nearest neighbor, and a histogram of Ta coordination numbers for all simulated systems (radial distribution functions can be found in Fig. S2 of the supplementary material). The total energies making up the hull are extracted from the HSE06 calculations and subsequently normalized to 24 formula units (f.u.). Our high sampling of 14 structures per stoichiometry allows us to visualize the energy variance, providing reassurance that our convex hull is accurate.
The convex hull in Fig. 1(a) shows a preference for phase separation into Ta 2 O5 and TaO 2 . However, while Ta 3 O7 (x = 0.33), Ta 4 O 9 (x = 0.50), and Ta 12 O 25 (x = 0.80) lie on average more than 0.6 eV above the hull, Ta 6 O 13 (x = 0.67) is considerably lower at only 0.11 eV above the hull. Our calculations show a clear division into high-energy and low-energy stoichiometries, and we expect some degree of similarity between the three low-energy stoichiometries.
To investigate what the three low-energy stoichiometries (Ta 2 O5, Ta 6 O 13 , and TaO 2 ) have in common, we consider the CDF for the Ta-O bonds in the atomic structures of Fig. 1(a). Since Ta is sixfold coordinated in most, if not all, stable and metastable crystalline TaOx phases, we limit our investigation to exactly the sixthnearest neighbor. The resultant spectra are shown in Fig. 1(b). All systems are characterized by an ascent beginning at ∼2.0 Å and ending at ∼2.5 Å. Therefore, we use 2.5 Å as the cutoff radius when an O atom counts toward the coordination number of a Ta atom. The height of the initial ascent thus describes the ratio of fully sixfold coordinated Ta atoms, and the low-energy systems are characterized by relatively larger quantities of fully coordinated Ta atoms than are seen for the high-energy systems. Within this picture, it is also clear that over 50% of Ta atoms are under-coordinated in all stoichiometries, including a-Ta 2 O5. We further illustrate the trends in coordination numbers in Fig. 1(c) by plotting coordination numbers 6, 5, 4, and <4 in a cumulative histogram. In particular, Ta 6 O 13 and TaO 2 show higher sixfold coordination and lower fourfold coordination than stoichiometries with slightly greater oxygen content. It is worth mentioning that two of the three low-energy stoichiometries, TaO 2 and Ta 2 O5, exist in crystalline form. 42 The crystalline form of the other low-energy stoichiometry, Ta 6 O 13 , has not been reported in the literature. However, tantalum shares periodic table Group 5 with vanadium, and V 6 O 13 has been found in a crystalline form where vanadium is also sixfold coordinated. 43 In contrast, the high-energy stoichiometries indicated in Fig. 1 either have no corresponding crystalline Group 5 metal-oxide phase (Ta 12 O 25 ), 46 or have corresponding crystalline phases with lower cation coordination numbers (Ta 3 O7 and Ta 4 O 9 ). 47,48 This analogy with crystalline phases underlines the importance of conserving the Ta sixfold coordination to keep the system stable. This explains why amorphous Ta 2 O5, Ta 6 O 13 , and TaO 2 exhibit a relatively high degree of Ta sixfold coordination and are thus particularly stable, while Ta 3 O7, Ta 4 O 9 and Ta 12 O 25 do not.

B. XPS spectra
In this section, we present simulated XPS spectra for our systems and compare them to the experimental XPS spectra reported by Benito and Palacio, 26 Simpson et al., 27 and Li et al. 17 To this end, we calculate the energies of the core electrons of the Ta 4f orbital, which correspond to measuring their BEs, thus simulating an XPS spectrum.
We "calibrated" our XPS results in two ways prior to making any comparisons. We benchmarked our 4f binding energy spectra for the crystalline λ phase 34 of Ta 2 O5 (19.1 eV) against the Ta 7/2 5+ peak reported by Ho, Contarini, and Rabalais 24 for crystalline Ta 2 O5. The reported XPS peak value is 26.2 eV; hence, we shifted all our XPS spectra uniformly by 26.2 eV − 19.1 eV = 7.1 eV. Such a shift is customary, and even necessary, in order to compare DFT-calculated XPS results to experiment. 44 Additionally, we remark that realistic XPS spectra are usually obtained at room temperature, whereas DFT calculations are performed at 0 K. We applied Gaussian broadening to our data to overcome this difference using a broadening width obtained from the energy spectrum of 2000 AIMD iterations at 300 K using a 520 eV cutoff energy and (2 × 2 × 2) k point grid. With these two alterations, we proceed to discuss our XPS spectra and compare them to the aforementioned experimental results. Figure 2(a) presents our calculated XPS spectrum for a-Ta 2 O5 with a single Gaussian function fitted to it. The fit is near perfect, signifying that the XPS spectrum comes from atoms with uniform oxidation states. The spectrum of a-TaO 2 in Fig. 2(b), by contrast, shows two peaks and has been deconvoluted using two Gaussian functions. We find that our calculated peaks' positions lie slightly above those reported by Benito and Palacio, Simpson et al., and Li et al. The calculated highest energy peak (i.e., the one usually assigned to Ta 5+ in the literature) and the second-highest energy peak (i.e., the one usually assigned to Ta 4+ in the literature) are  1.6 eV). The slight blueshift of our results compared to the experimental values could arise from using the initial state approximation to calculate the core-level energies. Overall, we believe that our peak locations and spacing are in good agreement with experiment.
Apart from matching experiments, we highlight the XPS spectrum of TaO 2 as a particularly important result. Even in a sample with a nominal Ta oxidation state of +4, the Ta 5+ peak is visibly present. Additionally, the Ta 4+ curve is significantly wider than the Ta 5+ curve, which contradicts the deconvolutions performed in the previous experimental studies. 17,26,27 For this reason, we recommend caution when using deconvoluted XPS spectra to perform compositional analyses on amorphous a-TaOx samples and likely for other amorphous metal-oxide samples.

C. Bader analysis
Our CDF spectra showed how a significant fraction of Ta atoms in all samples are under-coordinated, i.e., less than 6 oxygen neighbors within 2.5 Å. We now examine the connection between the Ta coordination numbers, BEs, and Bader charges (which provide an indication of the oxidation state). Figures 3(a)-3(c) shows the ascalculated Ta 4f binding energies plotted against their corresponding Bader charges for the three low-energy stoichiometries (Ta 2 O5, Ta 6 O 13 , and TaO 2 ). Atoms are color-and shape-coded according to the number of oxygen neighbors within 2.5 Å. Also shown are five points corresponding to crystalline Ta-O phases and elemental Ta along with kernel density estimates of the Bader charges (horizontal axis) and BEs (vertical axis). Discussions henceforth are limited to the low-energy stoichiometries.
For stoichiometric Ta 2 O5, the observed XPS peak comes from a narrow range of Bader charges centered at roughly +3, which is consistent with the results of Xiao and Watanabe. 19 Any slight variation in Bader charge is clearly related to the coordination number of Ta, visible as discernible regions of different colors. As O deficiency increases, the Bader charge spectrum becomes smeared out in direction toward lower positive charge, i.e., Ta atoms retain more of their valence electrons. However, a significant part of the Ta 2 O5 spectrum-highlighted by the dashed ellipsis-remains even for TaO 2 . Evidently, for x > 0, a-Ta 2 O5−x becomes disproportionated, which is the origin of the dual peak spectrum we observe for TaO 2 .
In comparison with crystalline phases, the Ta atoms of the amorphous structures exhibit much higher BEs even at lower oxidation states. For Ta 2 O5, the calculated difference between the BE of the Ta 4f orbitals in the crystalline phase and the center of the Gaussian peak of the amorphous phase is around 1 eV, in excellent agreement with the experimental difference (compare the results of the amorphous phase by Benito and Palacio, 26 Simpson et al., 27 and Li et al. 17 with those of the crystalline phase by Ho, Contarini, and Rabalais 24 ). For lower Ta oxidation states, the differences in BEs between amorphous samples and their crystalline counterparts become even more pronounced. This phenomenon can be understood by comparing the coordination number of Ta atoms between the crystalline and amorphous phases. While Ta is always sixfold coordinated for crystalline phases, every considered amorphous The change in coordination number is not the only cause of the difference in BEs between amorphous and crystalline phases. In Fig. 3(d), we show the variation of the BE peak for crystalline Ta 2 O5 when varying the lattice parameter of the crystal (which is equivalent to varying the average Ta-O distances). We scan the variation of the Ta 4f binding energy in the range from ε = −0.03 to ε = +0.03, where ε is the relative expansion factor applied to the lattice parameters. We observe a quasi-linear correlation in the investigated region between the BE of Ta 4f orbitals and average Ta-O bond length, which is consistent with the expected reduction of the Coulombic repulsion between the Ta 4f electron and the oxygen ligands when the Ta-O distance increases. We observe that the average Ta-O distances for the sixfold coordinated Ta atoms in Ta 2 O5 are ∼2 pm larger than those in crystalline λ-Ta 2 O5 (at ε = 0). Thus, we do not only attribute the blueshift in BEs in the amorphous samples to under-coordinated Ta atoms but also to the increase in the Ta-O distances for sixfold coordinated Ta atoms in amorphous samples relative to crystalline ones.
We now examine the relation between the Bader charge and BE hinted in Fig. 3 where subscripts A and i represent the Ta ions and their 6 nearest O neighbors in our implementation, respectively, q is the Bader charge, r is the distance between i and A, and β 0,1,2 are fitting parameters to be determined. This model improves the correlation coefficient for a-Ta 2 O5 to 0.492 but has no impact on a-TaO 2 (r = 0.881). We conclude that the correlation between BE and Bader charge in both a-Ta 2 O5 and a-TaO 2 is more complex than can be accounted for in the above models (in contrast, the model of Ebadi et al. 45 performs well for the crystalline phases; see Fig. S4 of the supplementary material). Instead, we decompose the Bader charge and XPS spectra of Fig. 3 (black dashed curves) into Gaussian components in order to compare the peak locations. We find that the best fit is achieved using a single Gaussian function for a-Ta 2 O5 and three Gaussian functions for a-Ta 6 O 13 and a-TaO 2 . We plot the decomposed BE peak values vs the Bader charge peak values in Fig. 5(g) and achieve a clear linear correlation with r = 0.993. We conclude that the XPS/Bader charge spectra of amorphous tantalum oxide are best understood as statistical entities, even when complete knowledge of the atomic arrangements is available.

IV. CONCLUSIONS
We have used AIMD and DFT to sample amorphous and oxygen-deficient Ta 2 O5. We have analyzed our a-Ta 2 O5−x systems with respect to relative formation energy, CDF, XPS, Bader charge, and coordination numbers. By drawing a convex hull over all 84 included structures, we show that oxygen-deficient a-Ta 2 O5 will preferentially separate into a-Ta 2 O5 and a-TaO 2 . Additionally, the hitherto unknown stoichiometry of a-Ta 6 O 13 showed energies significantly lower than other intermediate stoichiometries, which is linked to the existence of a vanadium-oxide phase with stoichiometry V 6 O 13 . Our CDF spectra showed that stoichiometries with a relatively low energy contain higher proportions of Ta atoms with six O neighbors within a radius of 2.5 Å.
By using Bader charges to model oxidation states explicitly, we found that Ta atoms in amorphous TaO 2 were divided into two different Bader charge spectra: one representing oxidation state +5 contributing to the Ta 5+ peak and the other with a wide range of oxidation states contributing to the Ta 4+ peak. Unlike the crystalline TaO 2 with a nominal Ta oxidation state +4, amorphous TaO 2 possesses disproportionated Ta atoms that make such a description less accurate. Our finding that the XPS spectrum of amorphous TaO 2 contains both Ta 4+ and Ta 5+ peaks challenges the conventional view of a-TaOx spectra. Our results suggest that only the Ta 2 O5 stoichiometry may be reduced to a single oxidation state (+5). This finding is important for researchers wishing to use XPS spectra to estimate the composition of amorphous TaOx samples and conceivably samples of other amorphous metal oxides as well. Finally, we determined the origin of the experimentally observed blueshift of the Ta 4f XPS peaks in amorphous Ta 2 O5 relative to its crystalline counterpart. The presence of under-coordinated Ta atoms and the longer Ta-O distances in amorphous TaOx leads to lower Coulomb repulsion and hence higher binding energies.
As a final remark, we decompose our calculated XPS and Bader charge spectra in Gaussian components and find a linear correlation between peak locations. The correlation strength of this Gaussian model exceeds the strength of a linear model accounting for Coulomb interactions with the sixth nearest neighbors, suggesting that the more accurate way to analyze amorphous materials is by using statistical models. For future analysis of the Ta 6 O 13 system, we suggest studying the kinetic routes to reach this stoichiometry and other possible disproportionated states of the TaOx system. To supplement this, it is also necessary to investigate the sensitivity of the sixfold Ta coordination with respect to addition and removal of oxygen, particularly for the Ta 6 O 13 phase. Finally, while we have studied tantalum oxide exclusively, our approach can be expanded to any amorphous metal-oxide and shed light on the XPS spectra of these disordered structures.

SUPPLEMENTARY MATERIAL
See the supplementary material for more detailed information on amorphous tantalum oxides such as their densities, radial distribution function, XPS spectra of solid solution and phase separate phases, and a linear relation between binding energy and Bader charge for crystalline tantalum oxides.