Defect-formation energies

Hydrogen atoms were inserted interstitially inside YSZ supercell, at various locations. Minimum-energy equilibrium sites were determined for each charge state of hydrogen by energy minimization. Final formation energies were determined as a function of Fermi energy, E_F, within YSZ band gap, according to the following expression:

(2)

where 𝐸 𝑑𝑒𝑓 𝑡𝑜𝑡 and 𝐸 𝑏𝑢𝑙𝑘 𝑡𝑜𝑡 represent the defect’s total energies (within embedded hydrogen atom) and pristine bulk-lattice supercells, respectively; µ_H denotes the magnitude of hydrogen chemical potential; ∆n accounts for deviations from stoichiometry; q is defect charge state; and E_V stands for the value of valence-band edge, which served as reference energy for Fermi level. Band gap of 3.91 eV was obtained from DFT-PBE, a value that underestimated YSZ experimental gap ⁸.

Final formation energies for different charge states of hydrogen are shown in Fig. 1, as a function of Fermi level in the gap. It can be seen that neutral state is never thermodynamically stable for any position of Fermi level. Instead, monatomic hydrogen is an amphoteric defect in YSZ, with pinning level (+/-) at 1.47 eV below conduction band edge, E_C. In band gap largest part, hydrogen is positively charged; thus, it is thermodynamically stable as a proton defect: H⁺. In Fig. 1, Fermi-level energies ranged from zero (valence-band edge, E_V) to E_C.

Figure 1: Formation energy of interstitial hydrogen in YSZ for hydrogen-rich and hydrogen-poor conditions as a function of Fermi-level position inside band gap. The red vertical line denotes (+/-) pinning level. 

The dependence of formation energies on specific value of hydrogen chemical potential is depicted in Fig. 1. For hydrogen-rich conditions (equilibrium with hydrogen gas, H₂(g)) corresponding formation energy is lower. Nonetheless, hydrogen-poorest conditions are the most realistic conditions in practical applications where YSZ samples are in equilibrium with a surrounding water vapor, H₂O(g). In this case, formation energies were higher, and hydrogen chemical potential is given as:

(3)

where 𝜇 𝐻 2 𝑂 𝑜 are 𝜇 𝑂 𝑜 are chemical potentials of water and oxygen, respectively, in their respective gas states, at standard T and partial-pressure conditions.

Furthermore, assuming an ideal-gas behavior for gas phases, effects of T and partial pressure can also be included to chemical potentials, as follows:

(4)

where µº denotes chemical potential of hydrogen at standard conditions of T = 298.15 K and pº = 1 atm. ∆µ(T, pº) is T-dependent part of chemical potential. These quantities were obtained from thermochemical tables ²².

Structurally, it was observed that H⁺ bound to oxygen ions, forming hydroxyl-type O-H bonds, with bond lengths within a very narrow range, from 0.98 to 1.00 Å. Due to the disordered structure of YSZ lattice, H⁺ binding to various O ions led to not iso-energetic configurations ¹². In fact, lower formation energies were recorded for O-H configurations very near oxygen vacancies (see Fig. 2a).

These results suggest that intrinsic oxygen vacancies of YSZ lattice defined trapping regions within YSZ, where protons preferably bind with oxygen ions which were the nearest neighbors of these vacancies. The structural relaxation pattern shows that O-H hydroxyl bond was accommodated by a breaking of one (or two) Zr-O host bonds, rendering participating O ion under-coordinated (see Fig. 2b).

Figure 2: (a) YSZ lattice near an oxygen vacancy. (b) Proton configuration in YSZ lattice- formation of dative O-H bond and breaking of two Zr-O host bonds. 

After establishing equilibrium sites and formation energies of protons in YSZ lattice, the next step was to determine energy balance of hydration reaction. Accordingly, protons incorporated into bulk YSZ lattice through dissociative dissolution of water molecules from surrounding atmosphere. Hydration reaction can be expressed as:

(5)

A water molecule transferred from gas phase was the source of two protons that entered the oxide lattice, and oxygen atom filled an oxygen vacancy. Hydration enthalpy (T = 0 K) of this reaction can be obtained as follows:

(6)

Considering formation energy of proton defect (see Fig. 1) and of the doubly positively charged oxygen vacancy ¹², net balance of hydration reaction was negative (13 eV). Thus, hydration was an exothermic process and YSZ lattice favorably accommodated proton defects from surrounding water atmosphere, at low T.

Migration pathways and barriers

A number of migration pathways representing plausible possibilities for proton motion in entire YSZ lattice were considered. DFT calculations with NEB method identified two distinct groups of paths with substantially different migration barriers. Pathways spatially confined within the immediate neighborhood of individual oxygen vacancies were found to be short ranged (see Fig. 3), and to have low to moderate barriers ranging from 0.13 to 1.02 eV (Fig. 4).

Figure 3: Migration paths of proton in YSZ lattice, denoted by arrows connecting initial and final proton sites. Oxygen vacancy site is represented by “b”. 

Figure 4: Energy profiles of proton migration paths depicted in Fig. 3. Energy barriers (in eV) for forward/backward motion are shown inside parentheses for specific migration segments. 

These paths comprise both proton-transport and bond reorientation modes of migration. In contrast, vacancy escape paths are exclusively proton-transport type and have a longer range (Fig. 3). These paths should normally be linked to long-range macroscopic diffusion, and are associated with larger barriers, in excess of 1 eV (see Fig. 4). These findings strongly suggested that proton mobility in bulk YSZ should be rather low.

Additional NEB calculations were also performed for proton paths at the core of a high-angle grain boundary ¹². Such extended defects constitute an integral part of Nc-YSZ samples internal structure. Calculations showed that Zr-Zr cation stacking at the core of these interfaces acted as strong obstacle to proton motion, creating very high barriers to overcome. Therefore, it is quite doubtful that grain boundaries can be fast diffusion pathways for proton transport. The only alternative for the reported high protonic conductivity in Nc-YSZ appears to be that protonic defects diffuse unhindered along internal surfaces of Nc-YSZ samples ⁷. This could be accomplished either by executing hopping jumps between oxygen ions as protons, or through a vehicle mechanism where H⁺ ions are part of larger groups (for instance H₃O⁺ ions) that are present from water dissolution to the surface from a wet atmosphere.

MSS

MSS has become a standard technique for studying monatomic hydrogen configurations through modelling with muonium, a light pseudo isotope where the hydrogenic atom possesses a positive muon as nucleus: μ⁺-e^{− (23,24}. Although the muon is only one-ninth the proton mass, muonium reduced mass 99.6% that of hydrogen, so that the respective electronic properties should be nearly identical ^25,26. MSS has also the additional advantage of being restricted to high-dilution limit for muonium impurity, which can thus generally be regarded as isolated, only indirectly affected by other defects or impurities through overall Fermi energy.

MSS experiments took place at EMU instrument of ISIS Facility, Rutherford Appleton Laboratory, UK. A Nc YSZ (Zr_0.92Y_0.08O₂) sample (grain size 13 nm) provided by Innovnano was investigated. In the experiments, 4 MeV positive muons were implanted into the sample. In the first one, a magnetic field (B = 10 mT) was applied perpendicularly to initial MSS polarization, at T of 8.5 K. Fig. 5 shows a Fourier transform of obtained MSS high statistics (using Lomb method ²⁷): the level of probability of false alarm below 0.0001 is shown as a pointed line in Fig. 5.

Figure 5: µSRS Fourier spectrum of YSZ at T = 8.5 K in an external magnetic field of B = 10 mT. The red dashed line is a simulation of a powder spectrum with A_iso = 0 MHz and D = 0.2 MHz. 

The results clearly indicate, as identified previously ¹¹, that MSS signal basically corresponds to positively-charged muons precessing at Larmor frequency. However, the frequency line width clearly exceeds the one expected from nuclear dipolar broadening, which indicates the presence of a small electronic spin density at the muon. Similar results have been obtained in other oxide systems ²⁸, revealing that hyperfine interaction was mainly dipolar. One may immediately estimate its approximate value from the full width at line half maximum, which yields a dipolar term: D ≈ λ/π = 0.11 MHz. A closer inspection of Fig. 5 reveals the existence of a small spectral power above the threshold for false alarm, at frequencies of f1 ≈ 1.3 MHz and f2 ≈ 1.5 MHz. These are likely to correspond to expected steps of a powder-pattern hyperfine spectrum ²⁹. Powder pattern represented as a red dashed line in Fig. 5 is a simulation with the model of ²⁹, using a value of isotropic hyperfine interaction A_iso = 0 and D corresponds approxim. to D ≈ f2 − f1 = 0.2 MHz. Both approaches indicate that hyperfine interaction was of almost pure dipolar origin (A_iso ≈ 0 and with D ≈ 0.1 MHz- at most, 0.2 MHz).

A notable feature is that a significant part of MSS polarization was not seen at low T (missing fraction) ¹¹. The missing fraction is related to muons depolarizing rapidly during thermalization stage, due to the formation of deeply bound muonium. Thus, one assigns the missing fraction to muons thermalizing in interstitial configurations ⁸ (either pure interstitial or oxygen-vacancy site). Nevertheless, hyperfine interaction associated to these configurations can be characterized by means of longitudinal-field repolarization technique. Fig. 6 shows repolarization curve for YSZ, at T = 300 K. This curve is two-stepped, which is a usual sign of anisotropy ³⁰. Therefore, the data were fitted (Fig. 6) with phenomenological repolarization functions proposed in ³⁰. Obtained values were A_iso = 2.1(2) GHz and D = 0.13(2) GHz. The corresponding fit is shown as a red line in Fig. 6.

Figure 6: Repolarization curve for YSZ at T = 300 K. The red line is a fit assuming an axially symmetric anisotropic hyperfine interaction, as discussed in the text. 

However, calculated hyperfine interaction for interstitial hydrogen in cubic bulk YSZ has been calculated in ¹⁰ as being essentially isotropic, with A_iso ≈ 3 GHz. This value is much higher than A_iso = 2.1(1) GHz, which is the present study’s experimental value. However, one must take into account that the sample used in these experiments is Nc, with an average grain size of 13 nm, implying an important role of the grain surface. MSS experiments in Nc II-VI semiconductors present strong evidence of surface states formation ³¹. Segregation of impurities on the nanograins surface is also a well-known effect ³², which may be enhanced with decreasing grain size ³³. Moreover, the first-principles calculations of ¹⁰ showed that hydrogen at the grain boundary core was characterized by different A_iso values from those of hydrogen in bulk regions. In certain cases, the corresponding value was found to be lower by as much as 20%. This finding was attributed to the distinct interface structure of the grain boundary. This could allow (in certain cases) a much larger spilling of valence electron density to neighboring ions, thus leading to a reduction of spin density at hydrogen nucleus and, consequently, to a smaller A_iso value for hydrogen, at the grain boundary core. The corresponding reduced values for muonium, after scaling took into account the 3.183 factor of magnetic moment ratio of the muon and the proton, were around A_iso ≈ 2.3 GHz, which compared well with the experimental value. Therefore, it is probable that anisotropic component seen in Fig.6 corresponds to muonium influenced by the grains surface.