Stoichiometry deviation in amorphous zirconium dioxide

Abstract

Amorphous zirconia (a-ZrO₂) has been simulated using a synergistic combination of state-of-the-art methods: employing reverse Monte-Carlo, molecular dynamics and density functional theory together. This combination has enabled the complex chemistry of the amorphous system to be efficiently investigated. Notably, the a-ZrO₂ system was observed to accommodate excess oxygen readily – through the formation of neutral peroxide (O₂²⁻) defects – a result that has implications not only in the a-ZrO₂ system, but also in other systems employing network formers, intermediates and modifiers. The structure of the a-ZrO₂ system was also determined to have edge-sharing characteristics similar to structures reported in the amorphous TeO₂ system and other chalcogenide-containing glasses.

---

1 Introduction

Zirconium dioxide (ZrO₂) is a widely used and important material in engineering applications including thermal barrier coatings within the aerospace industry¹ and transistor components in electronic engineering.² In the nuclear energy sector ZrO₂ has been used as an inert matrix nuclear fuel material³ and has been investigated as a potential nuclear waste-form.⁴ Additionally, ZrO₂ is the protective layer formed when zirconium-based alloys are exposed to air and steam (as in a nuclear reactor environment^5,6). The oxide has three major polymorphs:⁷ the low temperature (T ≲ 1440 K) monoclinic; intermediate temperature (1440 ≲ T ≲ 2640 K) tetragonal; and high temperature (T ≲ 2640 K) cubic. A metastable amorphous phase has also been widely reported^8–12 and studied for potential use in electronics as resistive random-access memory (RRAM) devices¹³ and as a replacement for SiO₂ in complementary metal-oxide-silicon (CMOS) devices.^14,15

The radiation damage and potential amorphization of zirconia is of particular importance when considering the behaviour of zirconium alloys in corrosive, nuclear environments – such as those found in a typical light water reactor (LWR). ZrO₂ is the passivating layer on such alloys and is extremely tolerant to radiation damage,¹⁶ remaining crystalline (in its cubic stabilized state) at low temperatures and high fluence¹⁷ – only amorphizing when the grain size is extremely small (∼50 nm).¹⁸ What has not been fully understood is the potential formation of amorphous films or phases at grain boundaries¹⁹ when exposed to radiation or corrosive environments, and the effect these may have on the behaviour of the bulk ceramic in terms of mechanisms that limit the protective nature of the oxide.

The formal charge of the Zr ion in ZrO₂ is not readily increased from 4+ to 5+ as to do so would require electron removal from a new orbital (an ionization energy of ∼46 eV). In other similar cation systems such as MgO,²⁰ BeO,²¹ (Ba/Sr)ZrO₃,²² and CeO₂,²³ where the cation is at its highest formal charge, it has been predicted that the peroxide ion (O₂²⁻) can accommodate an increase in oxygen content without the need to oxidize a cation in the system. Importantly, the formation of O₂²⁻ ion has been experimentally observed through Raman spectroscopy in both BaZrO₃ and SrZrO₃ as a route to accommodate hyper-stoichiometry.²² Hypo-stoichiometric ZrO_2−x has been widely studied and it has been observed that non-monoclinic phases can be stabilized by the resulting defects and non-stoichiometry.²⁴ Electronic defects may also be a route for charge compensating deviations in oxygen stoichiometry.^25,26

The periodic nature of crystals means that any interruption of their structures is easily described in terms of defects. For instance point-defects such as vacancies and interstitial atoms can be unambiguously defined in terms of their species, effective charge and crystal site, and communicated compactly using Kröger–Vink notation.²⁷ Crystals contain only a limited number of symmetrically distinct atomic sites and are described in terms of a lattice and a generally small repeating atomic motif. These symmetry relationships greatly reduce the number of atomic environments that must be considered to understand the defect chemistry of a crystalline material.

Later in this paper we describe how the structure of amorphous ZrO₂ accommodates changes in the oxygen stoichiometry. In the crystalline material these changes would be easily described in terms of point defects. However, the lack of long-range order that characterises the amorphous state, means that traditional definitions of point-defects break down when there is no lattice to act as a frame of reference. However, certain topological features in the amorphous state often have similarities with point-defects²⁸ and we will draw such comparisons where appropriate.

In this work, the potential for deviations in stoichiometry in amorphous ZrO₂ is investigated using a combination of molecular dynamics and reverse Monte-Carlo calculations using classical potentials and density functional theory (DFT). Both classical and quantum mechanical methods have been used previously to study amorphous systems: effective potential methods having the benefit of being able to simulate large, complex supercells with tens of thousands of atoms,^29–31 while the greater physical rigour of DFT methods and its quantum mechanical description of electronic properties additionally allows subtle chemical effects due to charge transfer and electronic defects to be considered.^32,33

2 Methodology

For this work a number of methodologies have been combined to generate representative amorphous structures and then extract relevant property and structural information.

The following section outlines the general workflow used:

(1) Melt-quench: Classical molecular dynamics simulations were used to represent the process of rapidly cooling a system from the molten state to induce an amorphous state. This first step was conducted for a relatively large simulation box which contained 96 000 atoms (in the stoichiometric case). More details of this step are given in Section 2.1.

(2) Reverse Monte-Carlo: In this step much smaller simulation cells (96 atoms when stoichiometric) were generated using the Reverse Monte-Carlo (RMC) method from the description of the amorphous state predicted in the previous step. A particular aim of this study was to compare the electronic structure of crystalline ZrO₂ with that of the amorphous state. This required the use of the DFT method (see the next step). The computational expense of DFT means that, due to their large size, the amorphous structures from step 1 cannot be used directly in DFT calculations. Instead, information obtained from the MD simulations such as density, pair-correlation functions and bond angle distributions, are extracted and used to parametrise RMC simulations. In these, the small systems are refined, through a series of stochastic atom moves, to produce structures that are consistent with the structural information extracted from step 1. A detailed description of this appears in Section 2.2.

(3) DFT: Finally the small cells were structurally relaxed using DFT energy minimisation. This allowed the structure and properties of the system to be better understood and compared with the various polymorphs of crystalline ZrO₂. The DFT method employed is given in Section 2.3.

The method described here has some similarities with previous work by Vanderbilt et al.,³⁴ Zhao¹⁴ and Ceresoli.³⁵ They used a melt-quench technique to study amorphous ZrO₂ and claimed good agreement with experimentally obtained values for the dielectric constant of the material. The melt-quench used by Vanderbilt et al. was performed on small simulation cells (96 atoms) using DFT MD. Here we use much larger simulation cells (96 000 atoms) and classical potentials to perform melt-quenches at much slower quench rates (4 × 10¹¹ K ps⁻¹ versus 3.4 × 10¹⁴ K s⁻¹ of the previous work).

The details of each of the MD, RMC and DFT steps will now be described.

2.1 Molecular dynamics – melt quench

The following procedure was used to simulate rapid cooling from the molten state in order to mimic the type of experimental melt-quench often used to produce amorphous materials. The starting structure for all quenches was based on a 20 × 20 × 20 super-cell of cubic zirconia. As will be described in Section 2.2, hypo- and hyper-stoichiometry was introduced by adding or removing oxygen atoms to this supercell as necessary. The same quench procedure was used for all stoichiometries which is now described. All MD calculations were performed using the 11 Aug 2017 version of the LAMMPS code³⁶ and an integration time-step of 1 fs was used throughout.

The cubic form of ZrO₂ is associated with high temperatures, consequently, the simulation cell was initially equilibrated at a temperature of 2700 K for 100 ps. In order to quickly stabilise the system's temperature and reach the density predicted by the chosen potential model (described below), a Berendsen thermostat and barostat were used (with relaxation times of 0.1 ps and 1 ps respectively).³⁷ This was followed by another 50 ps of MD at 2700 K which now used the Nosé–Hoover thermostat and barostat^38–40 (to allow more representative sampling of the NPT thermodynamic ensemble). The Nosé–Hoover thermostat had a relaxation time of 0.01 ps and the barostat used 0.1 ps. These were maintained for the remainder of the quench. An MD barostat controls pressure by adjusting cell parameters to achieve a target pressure, here this was 0 GPa and, in-line with the isotropic nature of amorphous materials, the barostat was set to maintain a cubic geometry throughout.

To remove any bias due to the crystalline starting structure, the system was heated to a high temperature of 5000 K. This was achieved by ramping the temperature during a 50 ps MD run. Once at temperature, another 50 ps of MD was performed, this can be thought of as the ‘melt’ section of the melt-quench, even though it takes place at a much higher temperature than the experimentally measured melting point of zirconia (2988 K (ref. 41)). Before commencing the quench, the radial distribution function of system was checked to ensure a molten state had been achieved.

The system was then quenched from 5000 K to 300 K by reducing the thermostat's temperature linearly during a 1.88 ns NPT MD run, therefore giving an effective quench rate of 0.4 K ps⁻¹. The melt-quench procedure is summarized graphically in Fig. 1.

Fig. 1 Graphical representation melt-quench procedure.

On reaching 300 K, the data to be used as inputs to the RMC method were collected. This was achieved by performing a final 100 ps MD run during which pair-correlation functions for the Zr–O, Zr–Zr and O–O species pairs were calculated and time averaged over the final 50 ps of the run by collecting data at 5 ps intervals. Atomic densities and relevant bond-angle distributions were also averaged in this way.

The many body potential model derived by Liu, Cooper et al. was used to describe interatomic forces during all MD runs.⁴² Although this potential provides an excellent high temperature description of the tetragonal and cubic forms of ZrO₂ and the transition between them, it must be recognized that the potential model does not reproduce the monoclinic phase.

2.2 Reverse Monte-Carlo

Reverse Monte-Carlo is a structural refinement technique which aims to generate model atomic systems, which are consistent with a set of constraints and structural data. RMC was originally developed to help interpret experimental diffraction data obtained from amorphous systems.^43–45 Of relevance here, Winterer used RMC to examine the structure of amorphous and monoclinic ZrO₂ using experimental data obtained using the EXAFS technique.⁴⁶

Diffraction experiments for crystals can be analysed to give a unit-cell and consequently provide an unambiguous description of the structure. By comparison, the shape-factor and radial distribution functions obtained from the diffraction of amorphous materials provide averaged data about the pair separations of atoms in the material. This means that it is not possible to define a single atomic structure from amorphous diffraction data as multiple structures may be consistent with a particular set of radial distribution functions (RDFs). RMC provides a method to generate atomic configurations that aim to match experimental RDFs. Furthermore, by including additional physical constraints on the material structure such as bond lengths, bond angles and coordination number limits, configurations may be generated that are more physically rigorous than fitting to RDFs alone. A full description of the method is provided elsewhere however a simple description of its operation is now given.^47,48

A set of atoms with the desired composition and density are assigned random positions in a simulation box. The RDF for this model system is calculated and compared with the experimental RDF. The root mean squared (RMS) difference between the two shows the goodness of fit between model and experiment. An atom in the model is then randomly moved and the comparison is made again. If the fit is improved (indicated by a decrease in the RMS difference), the move is accepted – if not, the move may still be accepted based on the Metropolis-Hastings algorithm but is otherwise rejected.^49,50 Similar moves are repeated until a satisfactory match is obtained between the model system and the experimental data. In this way a structure is generated that is consistent with experiment and should therefore be representative of the original material.

Here the RMC method is used, but rather than fitting the model system to an experimentally determined RDF, the pair-correlation functions, for each species pair, obtained from the MD quenches were used as the inputs to RMC, with the aim of generating much smaller cells suitable for DFT.

This method of combining MD derived structural data with RMC, to provide structures for DFT, has been used with success by the authors to study doped ZnO.⁵¹ In comparison to this previous work, several refinements to the method have been made here. First, RMC runs were initialized with cells containing randomised coordinates. For the stoichiometric case these contained 96 atoms (Zr₃₂O₆₄). Secondly, after every 500 accepted Monte-Carlo moves, an energy minimisation was performed on the RMC structure using the same potential model as for MD. This was then used as the input for the next round of RMC tests. In this way, the RMC algorithm was made to converge more quickly. Optimized structures were only accepted when the energy difference between consecutive energy minimised structures was less than 0.01 eV. Finally, the O–Zr–O and Zr–O–Zr bond angle distributions were included in RMC as additional constraints to the system.

For each composition, twenty distinct RMC cells were generated for use in DFT by initiating each run with a different random seed. The rmc++ code (version 1.4) was used for all RMC runs.⁵² During RMC single atom displacements of up to 2 Å were made during each Monte-Carlo loop.

2.3 Density functional theory

The DFT method was used to structurally minimise the unit cells obtained using the MD + RMC process. The Vienna Ab-initio Simulation Package (VASP)⁵³ was used with the projector augmented wave (PAW)⁵⁴ library supplied with the software. The GGA-PBE exchange correlation functional⁵⁵ was used with the highest number of electrons treated as valence supplied for each species.

The cut-off energy was set at 450 eV for all calculations with a smearing width of 0.1 eV (using Gaussian smearing). A 4 × 4 × 4 k-point grid was set automatically using the Monkhorst–Pack scheme. The self-consistent field (SCF) stopping criterion was set to 1 × 10⁻⁴ eV and the geometry optimization stopping criterion was set at 1 × 10⁻³ eV. Atomic positions and supercell size and shape were all free to fully relax.

Stoichiometry was varied in the supercells initially by varying the Zr : O ratio in the empirical calculations. This was carried out to provide a hyper-stoichiometric DFT modelled system containing 32 Zr atoms and 66 O atoms (ZrO_2.0625) while the hypo-stoichiometric DFT model contained 32 Zr atoms and 62 O atoms (ZrO_1.9375). Oxygen gas was modelled as a dimer in a large box in order to accurately assess the bond strength enabling calculations investigating the drive to accommodate both hyper-stoichiometry and hypo-stoichiometry in the amorphous ZrO₂ system.

Monoclinic ZrO₂ was modelled to compare the energies of amorphous ZrO₂ with the crystalline form – as well as understand the difference in energy and mechanism to accommodate deviations in stoichiometry compared to the amorphous ZrO₂ system. Deviations in stoichiometry in crystalline ZrO₂ were performed in a supercell containing 96 lattice sites (32 ZrO₂ units). Deviations in stoichiometry were considered by accounting for vacancies on all symmetrically distinct sites as well as the incorporation of interstitial species onto ten symmetrically distinct sites.

3 Results

3.1 Structure

A simulated XRD pattern was created for the amorphous ZrO₂ structure obtained from the molecular dynamics melt-quench routine (Fig. 2(ii) – grey pattern). When compared to experimental data by Sugiyama et al. (Fig. 2(iii) – red pattern) the patterns can be seen to be similar – showing a broad, amorphous hump with at peak at ∼30° and a second, lower intensity hump for 2θ in the 40–60° range. This is a strong indication that the MD model is a good representation of the experimentally observed amorphous ZrO₂ structure. The density of the amorphous structure generated by MD was 6.32 g cm⁻³.

Fig. 2 X-ray diffraction patterns obtained from simulated and experimental structures: (i) shows calculated XRD trace averaged over the 20 DFT structures, (ii) is the calculated XRD pattern for MD structure following melt-quench, (iii) shows the experimental XRD pattern for amorphous ZrO₂ obtained by Sugiyama et al. at T = 423 K,⁵⁶ (iv, v and vi) are the patterns for the cubic,⁴⁴ tetragonal⁴⁵ and monoclinic⁴³ polymorphs of ZrO₂ respectively (Cu K-α X-rays used across all datasets). The simulated diffraction patterns (i and ii) were obtained using the user/diffraction package in LAMMPS.^36,57

The RMC derived supercells were relaxed using DFT. Simulated XRD patterns were then calculated for all twenty structures and then averaged to give the trace seen in Fig. 2(i). It can be seen that unlike the pattern generated from the supercell after the molecular dynamics routine, the features of the XRD patterns are sharper but still vary significantly from the fully crystalline ZrO₂ polymorphs Fig. 2(iv–vi). It is expected that this difference in pattern profile can be attributed to a number of factors including: the limited size of the supercell causing imaging effects near 10 Å distances; the use of a relatively small total number of ZrO₂ units (640 – limiting the ability to completely simulate a bulk amorphous solid); and the static nature of the simulations that results in a lack of thermal effects. The average density for the supercells after relaxation was 5.95 g cm⁻³.

Despite the differences in the patterns, the local structures in the DFT supercells compared to experimental are similar and therefore the local chemical effects can be reasonably predicted. The radial distribution function of the DFT data, the molecular dynamics data and analysis performed on experimental data are shown in Fig. 3 highlighting the similar local structures. Fig. 4 reports the similarity between the local morphology in the amorphous phase to the crystalline phase.

Fig. 3 Radial distribution functions, g(r), of simulated amorphous zirconia (where r is atomic separation). Grey data: full amorphous structure obtained through molecular dynamics. Pink data: average data from 20 × 96 atom cells produced through RMC technique. Blue data: average RDF for the 15 lowest energy relaxed supercells obtained through DFT.

Fig. 4 An example of a DFT relaxed amorphous ZrO₂ cell. A section of the amorphous network corresponding to the outlined atoms from (a) is shown in sub-figure (b). This shows that the network is primarily composed of oxygen tetrahedra in edge sharing configurations. The dashed lines indicate the position of these shared edges. Sub-figure (c) shows that this edge sharing topology is also evident in the crystal polymorphs of zirconia.

In monoclinic ZrO₂ the Zr is coordinated by 7 oxygen atoms whilst in both tetragonal and cubic ZrO₂ the each Zr is coordinated by 8 oxygen atoms. Oxygen in the monoclinic polymorph is coordinated by either 3 or 4 Zr atoms whilst in the tetragonal and cubic polymorphs each O is coordinated by 4 Zr atoms. The coordination of both the Zr atoms in the amorphous ZrO₂ structures generated in this work was compared to the previous work reported by Vanderbilt et al.³⁴ in Fig. 5. Bonds were considered up to a cut-off distance of 3 Å to aid comparison with the previous work.

Fig. 5 Zirconium coordination number distribution in amorphous ZrO₂. The top two panes (shown in grey) reveal the distributions obtained by Vanderbilt et al.³⁴ in previous work using the activation relaxation method and a 96 atom DFT melt-quench. The lower three panes (in blue) show those predicted by the current work. The MD results are for the 96 000 atom melt-quench with classical potentials. The MD + RMC distribution is averaged across the twenty RMC optimized cells (each containing 96 atoms). Finally, the MD + RMC + DFT results show the RMC cells following relaxation using DFT. Note, the considerable differences in the magnitude on the frequency axes is due to the varying total number of atoms considered for each method, not necessarily the supercell size (see Section 2).

The average coordination of both Zr and O in the current work is higher than the reported coordination values in the work of Vanderbilt et al. who considered the same number of atoms in each DFT supercell as were used in this work³⁴ (see Fig. 5 and 6). The MD + RMC + DFT data from this work possesses a lower average coordination relative to the MD and MD + RMC generated amorphous ZrO₂ structures in this work indicating that the smaller system size and/or the change in atomic description may reduce the average Zr coordination. Similarly, the oxygen coordination in the work of Vanderbilt et al.³⁴ (a combination of 3 and 4) is more akin to the coordination environment in the monoclinic phase, whilst the oxygen coordination in the structures from this work are predominantly 4 coordinated, akin to the tetragonal and cubic polymorphs of ZrO₂.

Fig. 6 Oxygen coordination number distribution in amorphous ZrO_2. The top two panes (shown in grey) reveal the distributions obtained by Vanderbilt et al.³⁴ in previous work using the activation relaxation method and a 96 atom DFT melt-quench. The lower three panes (in blue) show those predicted by the current work. The MD results are for the 96 000 atom melt-quench with classical potentials. The MD + RMC distribution is averaged across the twenty RMC optimized cells (each containing 96 atoms). Finally, the MD + RMC + DFT results show the RMC cells following relaxation using DFT. Note, the considerable differences in the magnitude on the frequency axes is due to the varying total number of atoms considered for each method, not necessarily the supercells size (see Section 2).

3.2 Deviations in stoichiometry

The ability for amorphous ZrO₂ to accommodate deviations in stoichiometry is now investigated and compared to the crystalline ZrO₂ systems. Deviations in stoichiometry in ionic materials impact a number of key properties including intrinsic diffusion mechanisms, mechanical properties and thermal transport properties. It is not expected that Zr will oxidize to a higher charge state than its formal 4+ charge in the stoichiometric ZrO₂ system.

For excess oxygen, one may expect the formation of neutral oxygen defects such as the peroxide ion. The peroxide defect is an oxygen interstitial that combines with an oxygen on a lattice site to form an O₂²⁻ ionic species where the two constituent oxygen atoms are covalently bonded at a distance of 1.49 Å from each other. The strong covalent bond is distinct from the other ionic bonding where it has been previously reported.²²

Twenty hyper-stoichiometric (ZrO_2.0625) amorphous supercells were relaxed using the same parameters as the stoichiometric supercells allowing the following simple reaction to be considered:

O₂ + Zr₃₂O₆₄ → Zr₃₂O₆₆ (1)

The reaction was computed between each of the 15 lowest energy stoichiometric supercells and the 15 lowest energy hyper-stoichiometric supercells: totalling 225 reactions. The average reaction energy was 0.02 eV (0.01 eV per extra oxygen) indicating that excess oxygen can be readily accommodated into the structure via this interstitial-like topological feature. The range of energies is large (see Fig. 7) indicating that portions of amorphous ZrO₂ will accommodate excess oxygen exothermically.

Fig. 7 Histogram illustrating the reaction energy for O₂ to enter the amorphous ZrO₂ supercells.

The fifteen most stable hyper-stoichiometric amorphous structures were all found to accommodate either one or both of the excess oxygen species through the formation of a peroxide ion, characterized by the 1.49 Å bond. This bond is identifiable on the average RDF pattern for the hyper-stoichiometric system (Fig. 8).

Fig. 8 Average radial distribution function (RDF) for the hyper-stoichiometric ZrO₂ system after DFT structural minimization. *Indicates the peak at 1.49 Å typical of the peroxide ion.

Accommodation of excess oxygen in crystalline, monoclinic ZrO₂ was investigated to compare to the amorphous system. It was found that an excess oxygen atom preferentially sits at a (0.21, 0.44, 0.16) interstitial site which bonds to the (0.05, 0.314, 0.354) oxygen lattice site to produce a peroxide ion (as reported by Lyons et al.⁵⁸). The enthalpy required to take ½O₂ into solution was computed to be 1.09 eV – significantly larger than the value for solution into the amorphous structure.

Hypo-stoichiometry in amorphous ZrO₂ was investigated in a similar manner to hyper-stoichiometry. Twenty supercells with 32 Zr atoms and 62 O atoms were created using the same routine using data from a molecular dynamics simulation. The following reaction was then considered to understand the drive to release oxygen:

Zr₃₂O₆₄ → Zr₃₂O₆₄ + O₂ (2)

The average reaction energy for this (considering the 15 lowest energy structures) was 10.46 eV or 5.23 eV per oxygen removed. The energy to produce an oxygen vacancy in crystalline monoclinic ZrO₂ was computed to be 11 eV, meaning that hypo-stoichiometry in amorphous ZrO₂ is preferred, despite the significant energy penalty for hypo-stoichiometry in either system.

The predicted energy penalty to go from monoclinic ZrO₂ to amorphous ZrO₂ is 0.54 eV per ZrO₂ unit. This energy is lower when considering deviations in stoichiometry: the energy penalty per ZrO_2.0625 unit to go from crystalline to amorphous is 0.48 eV implying that the presence of any excess oxygen in crystalline ZrO₂ will aid the formation of an amorphous phase when exposed to radiation damage. The same was found when considering hypo-stoichiometry – the energy penalty for a ZrO_1.9375 unit to go from crystalline to amorphous was computed to be 0.45 eV. Larger deviations from stoichiometry are expected to reduce the energy penalty for amorphization, accordingly.

4 Summary

The amorphous structure of ZrO₂ has been modelled using a combination of empirical potentials, reverse Monte-Carlo refinement and density functional theory. This produces results that compare well with previous atomic scale descriptions of the system. The coordination environment in the present work is more closely related to the tetragonal and cubic polymorphs of ZrO₂, whilst the structures reported in the work of Vanderbilt et al.³⁴ can be regarded as more similar to the coordination environment in the low temperature monoclinic polymorph. It should be noted that Vanderbilt et al. conducted their quenches exclusively using DFT which would have the advantage of being able to capture the low temperature monoclinic form. The use of DFT imposes some computational limitations, especially related to system size and the need for very rapid quench rates. By comparison our method allows large system sizes and considerably slower, more reasonable, quench rates (although still far higher than experiment). These methodological differences will lead to slight differences in the final structures but are unlikely to have a major consequence on the conclusions related to stoichiometry.

Unlike silica based glasses which consist of corner-sharing SiO₄ tetrahedra, amorphous zirconia can be defined by tetrahedra in edge sharing configurations (akin to chalcogenide glasses⁵⁹) as illustrated in Fig. 4. Previous investigations have focused on how the edge-sharing topology of chalcogenide glasses may change their behaviour in comparison to corner sharing glasses,^60,61 these have included molecular dynamics studies.⁶² Amorphous TeO₂ has been reported to be edge sharing by Brady et al.⁶³ who also discuss the glass forming ability of corner-sharing glasses with edge-sharing and face-sharing glass and the increasing difficulty in glass formability with increasing structural restriction. The edge sharing nature of amorphous ZrO₂ is consistent with the low amorphous to crystalline transition temperature observed experimentally.⁹

It is predicted that deviations in the stoichiometry are more readily accommodated in amorphous zirconia compared to the crystalline structure. The results suggest that the same mechanism accommodates hyper-stoichiometry in both the crystalline and amorphous systems but the morphological differences in the a-ZrO₂ system have enabled the O₂²⁻ dumbbell to be accommodated more readily in the non-crystalline form. It is expected that electronic disorder will also play a role in the accommodation of non-stoichiometry (especially with charged extrinsic defects)²⁵ and may be a fruitful avenue for further work.

In a situation where an amorphous region of material is in contact with crystalline material, it is expected that the amorphous region will getter deviations in stoichiometry assuming the migration of the defects is kinetically viable. This is particularly relevant when considering amorphous phases that may exist at grain boundaries induced by radiation damage – for example in the semi-protective layer that forms as a result of oxidation of zirconium alloys used as cladding materials in light water reactors (LWRs). Migration of neutral peroxide species through grain boundaries may be the route for Zr metal oxidation as oxygen mobility through the bulk is slow and deviations from stoichiometry are very unfavourable. The presence and evolution of glassy phases along grain boundaries may therefore impact the oxidation behaviour of alloys with zirconia semi-protective layers. The same behaviour is expected for other HCP metals such as titanium alloys and warrants further work.

It is unclear whether a general rule can be formulated to understand whether or not an amorphous oxide system is able to deviate in stoichiometry. Other ionic glass network components such as Al, Ti, Be and Zn may be expected to behave similarly to amorphous ZrO₂ compared with the more covalent Si and Ge glasses. The influence of multiple network formers, modifiers and intermediates is also likely to impact the overall propensity for a glassy system to vary its oxygen stoichiometry.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was carried out as part of the Sêr Cymru II programme funded through the Welsh European Funding Office (WEFO) under the European Development Fund (ERDF). Computing resources were made available by HPC Wales and Supercomputing Wales and special thanks are due to Aaron Owen, Benjamin Nash and Adrian Fewings. SCM would like to acknowledge the support of the members of the MUZIC-3 programme, particularly Magnus Limbäck, Antoine Claisse, Javier Romero and Jonna Partezana.

References

1. D. Clarke and S. Phillpot, Mater. Today, 2005, 8, 22–29.
2. G. Adamopoulos, S. Thomas, P. Wöbkenberg, D. Bradley, M. McLachlan and T. Anthopoulos, Adv. Mater., 2011, 23, 1894–1898.
3. L. Wang, S. Wang and R. Ewing, Philos. Mag. Lett., 2000, 80, 341–347.
4. W. Weber, R. Ewing, C. Catlow, T. de la Rubia, L. Hobbs, C. Kinoshita, H. Matzke, A. Motta, M. Nastasi, E. Salje, E. Vance and S. Zinkle, J. Mater. Res., 1998, 13, 1434–1484.
5. A. Motta, A. Couet and R. Comstock, Annu. Rev. Mater. Res., 2015, 45, 311–343.
6. M. Lindgren and I. Panas, RSC Adv., 2013, 3, 21613–21619.
7. C. Howard, B. Hill and B. Reichert, Acta Crystallogr., Sect. B: Struct. Sci., 1988, 44, 116–120.
8. T. Kokubo, Y. Teranishi and T. Maki, J. Non-Cryst. Solids, 1983, 56, 411–416.
9. H. Naguib and R. Kelly, J. Nucl. Mater., 1970, 35, 293–305.
10. D. Kundu, P. Biswas and D. Ganguli, Thin Solid Films, 1988, 163, 273–278.
11. Z. Yanwie, G. Fagherazzi and S. Polizzi, J. Mater. Sci., 1995, 30, 2153–2158.
12. J. Chang and T.-S. Lin, Appl. Phys. Lett., 2001, 79, 3666.
13. P. Parreira, G. Paterson, S. McVitie and D. MacLaren, J. Phys. D: Appl. Phys., 2016, 49, 0951111.
14. X. Zhao, D. Ceresoli and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 2005, 71, 085107.
15. D. Ceresoli and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 2006, 74, 125108.
16. K. Sickafus, H. Matzke, T. Hartmann, K. Yasuda, J. Valdez, P. Chodak, M. Nastasi and R. Verrall, J. Nucl. Mater., 1999, 274, 66–77.
17. C. Degueldre, P. Heimgartner, G. Ledergerber, N. Sasajima, K. Hoju, T. Muromura, L. Wang, W. Gong and R. Ewing, MRS Proceedings, 1996, 439, 625.
18. F. Lu, J. Wang, M. Lang, M. Toulemonde, F. Namavar, C. Trautmann, J. Zhang, R. Ewing and J. Lian, Phys. Chem. Chem. Phys., 2012, 14, 12295–12300.
19. P. Cantwell, M. Tang, S. Dillon, J. R. G. Luo and M. Harmer, Acta Mater., 2014, 62, 1–48.
20. S. Middleburgh, K. Lagerlof and R. Grimes, J. Am. Ceram. Soc., 2013, 96, 308–311.
21. S. Zhang, F. Li, H. Xu and G. Yang, Inorg. Chem., 2017, 56, 5233–5238.
22. S. Middleburgh, I. Karachevtseva, B. Kennedy, P. Burr, Z. Zhang, E. Reynolds, R. Grimes and G. R. Lumpkin, J. Mater. Chem. A, 2014, 2, 15883–15888.
23. H. Näfe and D. Guatam, J. Phys. Chem. C, 2016, 120, 10523–10529.
24. S. Middleburgh, G. Lumpkin and R. Grimes, Solid State Ionics, 2013, 253, 119–122.
25. K. Sasaki and J. Maier, Solid State Ionics, 2000, 134, 303–321.
26. B. D. C. Bell, S. T. Murphy, P. A. Burr, R. W. Grimes and M. R. Wenman, J. Appl. Phys., 2015, 117, 084901.
27. F. A. Kröger and H. J. Vink, Solid State Phys., 1956, 3, 307.
28. M. Ojovan and W. Lee, J. Phys.: Condens. Matter, 2006, 18, 11507–11520.
29. T. Stechert, M. Rushton, R. Grimes and A. Dillon, J. Non-Cryst. Solids, 2012, 358, 1917–1923.
30. T. Stechert, M. Rushton and R. Grimes, J. Am. Ceram. Soc., 2013, 96, 1450–1455.
31. J. Fogarty, H. Aktulga, A. Grama, A. van Duin and S. Pandit, J. Chem. Phys., 2010, 132, 174704.
32. D. King, S. Middleburgh, A. Liu, H. Tahini, G. Lumpkin and M. Cortie, Acta Mater., 2015, 83, 269–275.
33. U. Schwingenschlögl, A. Chroneos, C. Schuster and R. W. Grimes, Appl. Phys. Lett., 2010, 96, 242107.
34. D. Vanderbilt, X. Zhao and D. Ceresoli, Thin Solid Films, 2005, 486, 125.
35. D. Ceresoli and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 2006, 74, 125108.
36. S. Plimpton, J. Comp. Physiol., 1995, 117, 1–19.
37. H. Berendsen, J. Postma, W. van Gunsteren, A. DiNola and J. Haak, J. Chem. Phys., 1984, 81, 3684.
38. S. Nosé, J. Chem. Phys., 1984, 81, 511.
39. S. Nosé, Mol. Phys., 1984, 52, 255.
40. W. Hoover, Phys. Rev. B: Condens. Matter Mater. Phys., 1985, 31, 1695.
41. H. Scott, J. Mater. Sci., 1975, 10, 1527.
42. X. Liu, M. Cooper, K. McClellan, J. Lashley, D. Byler, B. Bell, R. Grimes, C. Stanek and D. Andersson, Phys. Rev. Appl., 2016, 6, 1.
43. C. Howard, E. Kisi, R. Roberts and R. Hill, J. Am. Ceram. Soc., 1990, 73, 2828–2833.
44. G. Katz, J. Am. Ceram. Soc., 1971, 54, 531.
45. N. Igawa and Y. Ishii, J. Am. Ceram. Soc., 2001, 84, 1169.
46. M. Winterer, J. Appl. Phys., 2000, 88, 5635–5644.
47. R. L. McGreevy and L. Pusztai, Mol. Simul., 1988, 1, 359.
48. R. L. McGreevy, J. Phys.: Condens. Matter, 2001, 13, R877.
49. N. Metropolis, A. Rosenbluth, M. Rosenbluth, A. Teller and E. Teller, J. Chem. Phys., 1953, 21, 1087.
50. W. Hastings, Biometrika, 1970, 57, 97.
51. D. Muñoz Ramo, A. Chroneos, M. Rushton and P. Bristowe, Thin Solid Films, 2014, 555, 117.
52. O. Gereben, P. Jóvári, L. Temleitner and L. Pusztai, J. Optoelectron. Adv. Mater., 2007, 9, 3021.
53. G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558.
54. P. Blochl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953.
55. J. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865.
56. K. Sugiyama, Y. Waseda and S. Kudo, ISIJ Int., 1991, 31, 1363.
57. S. P. Coleman, D. E. Spearot and L. Capolungo, Modell. Simul. Mater. Sci. Eng., 2013, 21, 055020.
58. J. Lyons, A. Jonotti and C. Van de Walle, Microelectron. Eng., 2011, 88, 1452–1456.
59. H. Eckert, J. Kennedy, A. Pradel and M. Ribes, J. Non-Cryst. Solids, 1989, 113, 287–293.
60. M. Micoulaut, Eur. Phys. J. B, 1998, 1, 277–294.
61. S. Sugai, Phys. Rev. B: Condens. Matter Mater. Phys., 1987, 35, 1345–1361.
62. P. Vashishta, R. Kalia and I. Ebbsjo, J. Non-Cryst. Solids, 1988, 106, 301–304.
63. G. Brady, J. Chem. Phys., 1957, 27, 300–303.

---