Hydrogen Bond Dynamics in Proton-Conducting Lanthanum Arsenate

Article pubs.acs.org/JPCC

Hydrogen Bond Dynamics in Proton-Conducting Lanthanum Arsenate Kazuaki Toyoura*,† and Katsuyuki Matsunaga†,‡ †

Department of Materials Science and Engineering, Nagoya University, Furo, Chikusa, Nagoya 464-8603, Japan Nanostructures Research Laboratory, Japan Fine Ceramics Center, Atsuta, Nagoya 456-8587, Japan

‡

ABSTRACT: The proton incorporation and migration mechanisms in lanthanum ortho-arsenate, LaAsO4, have theoretically been investigated using first-principles calculations in conjunction with the nudged elastic band (NEB) method and the kinetic Monte Carlo (KMC) method. Nine proton sites and twenty migration paths connecting the proton sites were found in the crystal, which result in the anisotropic proton diffusivity and conductivity of LaAsO4. First-principles phonon calculations have also been performed to discuss the relation between the host lattice flexibility and relaxation around a migrating proton. The involved host lattice relaxation during the proton migration mainly consists of rotations and translations of AsO4 units with little distortion of the tetrahedra. This reflects the relatively low eigenfrequencies of vibrational modes attributed to the AsO4 rotations and translations in comparison with the AsO4 bending and stretching modes. Such dynamics of the host lattice and protons should be common in protonconducting oxides with rigid oxygen tetrahedral units (XO4) isolated in the crystals.

1. INTRODUCTION Proton-conducting oxides are of great interest as electrolytes in electrochemical devices involving hydrogen, such as fuel cells, water electrolyzers, and gas sensors. Several classes of proton conductors have already been reported so far, e.g., perovskitetype, pyrochlore-type, monazite-type, scheelite-type, and their related oxides.1−12 Lanthanum ortho-arsenate, LaAsO4, is one of the proton conductors in the monazite-type class, which was reported to exhibit protonic conductivity up to 850 °C in wet atmospheres.12 The hydration behavior was also investigated theoretically from first principles in ref 12, where dominant positively charged defects in LaAsO4 were reported to switch from pyro-arsenate ions, (As2O7)•• 2AsO4, to hydrogen arsenate • ions, (HAsO4)AsO4, depending on dry or wet conditions. They found many proton sites around AsO4 units in the crystal, whose calculated site energies are widely scattering with large dispersion of 0.6 eV. This implies the complicated proton conduction mechanism, which may be one of the major reasons why the conduction paths connecting these proton sites have not been investigated in spite of the highly significant and suggestive information leading to improvement of the proton conductivity. In general, incorporated protons in oxides mainly reside around oxide ions with OH bonds and migrate by repetition of rotations around the oxide ions and hoppings between the rotational orbits. The oxide ions in the crystals form oxygen polyhedral units with cations, e.g., tetrahedral and octahedral units. In some oxides, the oxygen polyhedral units exist in isolation, and in others, they are corner-shared and/or edgeshared to form oxygen polyhedral networks; therefore, oxygen polyhedral types and geometries have a lot of variety reflecting © 2013 American Chemical Society

their crystal structures, which should play a key role for the proton long-range conduction. On the basis of the previous reports using atomic-stale simulations,1,3−6,12 protons in some oxides form hydrogen bonds (OH−O bonds) with the first- and second-nearestneighbor oxide ions (O1NN and O2NN), not only the OH bonds. The hydrogen bond formation involves large relaxation of the host lattices, particularly distortion, rotation, and/or translation of neighboring oxygen polyhedra. Therefore, the host lattice flexibility reflecting the oxygen polyhedral type and geometry should be of importance to determine the proton mobility in the crystal interpreted as hydrogen bond dynamics. Firstprinciples phonon calculations are one of the powerful techniques to quantitatively evaluate the host lattice flexibility. In the present study, we have calculated the vibrational spectrum of the host lattice from first principles, to analyze the components of individual vibrational modes in the host lattice relaxation involved by the proton incorporation and migration. Specifically, the atomic displacements from the original positions were expressed as linear combinations of the calculated vibrational modes. Such analysis using phonon calculations can provide meaningful information on ionic conduction from the structural viewpoint, particularly on proton conductors in which vibrations of protons and the host lattice are considered to be in separation due to their large difference in mass. Received: July 8, 2013 Revised: August 6, 2013 Published: August 8, 2013 18006

dx.doi.org/10.1021/jp406701w | J. Phys. Chem. C 2013, 117, 18006−18012

The Journal of Physical Chemistry C

Article

more precisely optimized until the residual forces became less than 1 × 10−5 eV/Å, since phonon calculations are sensitive to the residual forces of the initial structure. Each atom in the supercell was displaced by ±0.01 Å in each of x, y, and z directions to obtain all the interatomic force constants. The eigenvectors, ei, and eigenfrequencies, ωi, were obtained by diagonalizing the dynamical matrix consisting of the interatomic force constants and the atomic masses.27,28 The atomic displacements in the relaxed structures around a migrating proton can be expressed using the eigenvectors at the Γ-point of the supercell consisting of 2 × 2 × 2 unit cells due to the orthonormal basis set with the same degrees of freedom reflecting the periodic boundary condition. Specifically, the displacement column vector, x, with the 3n dimensions (n = 192: number of atoms in the supercell) and the 3n × 3n diagonal matrix of mass, M, are related to the linear combination of the eigenvectors as follows:

In this paper, we have demonstrated the above analysis of ionic conductions using phonon calculations, taking LaAsO4 with isolated AsO4 units as a model system. First, the proton incorporation and migration behaviors in the crystal were analyzed in a first-principles manner using the nudged elastic band (NEB) method13 and the kinetic Monte Carlo (KMC) method,11,14,15 respectively. Then, the involved host lattice relaxation around a migrating proton was evaluated by firstprinciples phonon calculations based on the supercell approach with the finite displacement method,16 to discuss the characteristics of the proton conduction in the related oxides with isolated oxygen tetrahedral units.

2. COMPUTATIONAL PROCEDURES The present study was based on first-principles calculations with the projector augmented wave (PAW) method17 implemented in the VASP code.18−22 The generalized gradient approximation (GGA) parametrized by Perdew, Burke, and Ernzerhof23 was used for the exchange-correlation term. The 5s, 5p, 6s, and 5d orbitals for lanthanum, 4s, 4p, and 5s for strontium, 4s and 4p for arsenic, 2s and 2p for oxygen, and 1s for hydrogen were treated as valence states. The plane wave cutoff energy was set to be 400 eV. A supercell consisting of 2 × 2 × 2 unit cells of LaAsO4 was used with a single k-point sampling at the Γ-point. The analysis of the proton conduction in LaAsO4 was performed in the same first-principles manner as we previously reported.11,24,25 First, the local energy minima of a proton in the crystal, i.e., proton sites, were determined by construction of the potential energy surface with the fixed atomic positions and the subsequent structural optimizations. The evaluated potential energy surface of a proton was here limited to spheres with a radius of 1 Å around each oxide ion, due to the strong OH bonds of typically 1 Å in oxides. The spherical grid includes 184 points per oxide ion. The structural optimizations were subsequently performed using the local energy minima in the PES as the initial structures. The atomic positions were fully optimized until the residual forces became less than 0.02 eV/Å. The proton conduction paths and their energy profiles were evaluated using the NEB method for all the possible paths connecting the proton sites within 3.5 Å distance. The proton diffusion simulations were finally performed on the basis of the KMC methods using the migration paths with low potential barriers below 1 eV with reference to that at the most stable site. The migrating protons in the crystal were assumed to be independent particles without their interactions due to the low solubility limit of dopants, ∼1%.12 The KMC steps were carefully determined in the range from 1.5 × 104 to 3 × 108 by checking convergence of the diffusion coefficient at each temperature. Association energies between protons and dopants must be evaluated for comparison with experimentally reported conductivity data, because the association can significantly reduce the concentration of mobile protons in the crystal. The association energy with strontium dopants was here evaluated from the energy difference between the two supercells including both a proton and a strontium ion with the nearest and furthest H+−Sr2+ distances. Phonon calculations for the perfect crystal of LaAsO4 were performed using the phonopy code26 based on the supercell approach with the finite displacement method, to evaluate the host lattice relaxation by proton incorporation and migration in the crystal. The atomic positions before displacements were

M1/2x =

∑ Ciei i

(1)

Hereafter, the coefficient Ci is defined as the component of the ith vibrational mode to reproduce the host lattice relaxation, given by Ci = eiTM1/2x

(2)

where the superscript T means transposition from column vectors to row vectors.

3. RESULTS AND DISCUSSION 3.1. Proton Sites. Figure 1a shows the crystal structure of LaAsO4 with the monoclinic monazite structure (space group: P21/n). The calculated lattice parameters, a, b, c, and β, are 7.11 Å, 7.30 Å, 6.85 Å, and 104.95°, respectively, which coincide with the experimental values29 within 1.5%. Arsenic and oxygen atoms form AsO4 tetrahedral units, which are isolated in the crystal and crystallographically equivalent. Four oxygen sites at the corners of AsO4 tetrahedra are inequivalent, which are classified into four, O1, O2, O3, and O4, divided by color in the figure. Figure 1b shows the determined local energy minima of a proton (proton sites) in LaAsO4. The AsO4 tetrahedral unit is shown by the original shape before the structural optimizations due to the different structural relaxations for the individual proton sites. The energies in parentheses are the potential

Figure 1. (a) The crystal structure of LaAsO4 with the monazite structure. The green balls and gray tetrahedra denote La ions and AsO4 tetrahedral units, respectively. The blue, pink, black, and red small balls at the corners of the tetrahedra are four oxygen sites, O1, O2, O3, and O4. (b) The calculated proton sites around each oxide ion in LaAsO4. The energies in parentheses are the potential energies with reference to the most stable. 18007

dx.doi.org/10.1021/jp406701w | J. Phys. Chem. C 2013, 117, 18006−18012

The Journal of Physical Chemistry C

Article

structures around the most stable site, O4-3, before and after the structural optimization. The O1NN−O2NN distance changes from 2.99 to 2.54 Å by the structural relaxation. The local structures around protons at the other proton sites except O3-1 also change to make the O1NN−O2NN distances decrease, 2.4− 2.7 Å in the relaxed structures. This means protons are more stabilized in the crystal by forming OH−O bonds, i.e., hydrogen bonds, which was also reported to be found in several perovskites (ABO3).1,3−5 The unique natures in oxides with isolated oxygen tetrahedra are that the hydrogen bond formation involves the two neighboring oxygen tetrahedral units, and that the changes in local structures are mainly rotations and translations of the two tetrahedra with their little distortions, in contrast to the intra-polyhedral hydrogen bonds with large distortions of BO6 octahedra in the perovskites. The characteristic structural relaxations around the proton sites are originated from the isolation of the AsO4 tetrahedral units in the crystal. The upper chart in Figure 2b shows the calculated vibrational spectrum of the supercell of 2 × 2 × 2 LaAsO4 unit cells without a proton at the Γ-point, which can be divided into three regions, as also seen in LaPO4.30 The first region has low eigenfrequencies up to 6 THz, which is attributed to lattice vibrations consisting of rotations and translations of AsO4 units in the crystal. The second region with middle frequencies in the range 7−14 THz corresponds to AsO4 bending modes with the O−As−O angles in the tetrahedra changed. The third is the high-frequency region of 22−26 THz corresponding to AsO4 stretching modes with the As−O bond lengths changed. The lower chart in Figure 2b shows the components of the vibrational modes to reproduce the atomic displacements in the relaxed structure around the O4-3 site. The vibrational components are relatively high in the first low-frequency region up to 6 THz corresponding to the rotation and translation of AsO4 units in comparison with those in the middle- and high-frequency regions. The similarity in the vibrational components can be seen for the other proton sites in LaAsO4, which is reasonable in that the local structures around the proton sites are relaxed with the smallest energy loss for the energy gain by hydrogen bond formation. Such the local structural relaxations around proton sites by oxygen tetrahedral rotations and translations should be the unique nature of oxides with relatively rigid oxygen polyhedra. Thus, proton sites are determined predominantly by hydrogen bond formation involving local structure relaxations. However, other minor factors also influence the positions of proton sites in oxides, because protons are not exactly located on the lines connecting two neighboring oxide ions, as shown

energies of proton sites with reference to that of the most stable O4-3 site. What determines the positions of proton sites in the arsenate? OH bond formation is one of the major factors, which restricts proton sites on the several spheres with radius ∼1 Å around oxide ions. Another factor is the location of neighboring oxide ions in other tetrahedra because the OH bonds at all the proton sites head toward a neighboring oxide ion. It is interesting to note that the local structure around a proton changes to shorten the distance between the 1NN and 2NN oxide ions. As an example, Figure 2a shows the local

Figure 2. (a) The change in local structure around the most stable proton site, O4-3. The tetrahedra shown by dashed lines with white balls and the gray solid tetrahedra with colored balls are the AsO4 tetrahedra before and after the structural optimization, respectively. (b) (upper) The calculated vibrational spectra of the supercell of 2 × 2 × 2 LaAsO4 unit cells without a proton at the Γ-point. The contributions of La, As, and O atoms are also shown by green, blue, and red regions, respectively. (lower) The components of the vibrational modes to create the atomic displacements in the relaxed structure around the O4-3 site.

Figure 3. (a) The proton rotational paths around oxide ions in LaAsO4. The colored and white small balls denote the proton sites and the rotational orbits connecting the proton sites around the same oxide ions. (b) The inter-tetrahedral hopping paths in the crystal shown by the solid lines. The number with each path corresponds to the path number in Table 1. 18008

dx.doi.org/10.1021/jp406701w | J. Phys. Chem. C 2013, 117, 18006−18012

The Journal of Physical Chemistry C

Article

Table 1. The Calculated Migration Paths with Low Potential Barriers below 1 eV with Reference to the Potential Energy at the Most Stable Site (O4-3) potential barrier height (eV)a

path no.

initial site

final site

distance (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

O4-3 O4-3 O1-1 O3-2 O4-2 O4-1 O4-3 O3-2 O3-3 O1-1 O4-3 O4-1 O1-1 O4-1 O3-2 O2-1 O2-1 O4-2 O4-1 O4-1

O1-1 O2-1 O1-1 O3-1 O3-1 O4-2 O4-2 O3-3 O3-3 O3-1 O2-2 O4-2 O4-1 O3-3 O2-2 O2-2 O2-1 O4-2 O4-1 O3-3

0.57 1.45 2.08 0.71 1.02 1.48 1.86 1.37 2.47 2.32 1.87 2.92 3.18 3.23 2.48 2.26 2.97 3.23 3.32 2.31

classification inter-tetrahedral inter-tetrahedral inter-tetrahedral rotation inter-tetrahedral rotation rotation rotation inter-tetrahedral inter-tetrahedral inter-tetrahedral inter-tetrahedral inter-tetrahedral inter-tetrahedral inter-tetrahedral rotation inter-tetrahedral inter-tetrahedral inter-tetrahedral inter-tetrahedral

hoppingb hopping hopping hoppingb

hopping hopping hopping hopping hopping hopping hopping hopping hopping hopping hopping

vs O4-3

vs initial site

vs final site

0.08 0.35 0.36 0.54 0.55 0.56 0.59 0.67 0.68 0.71 0.71 0.76 0.77 0.80 0.81 0.82 0.84 0.84 0.84 0.88

0.08 0.35 0.28 0.08 0.08 0.32 0.59 0.21 0.12 0.64 0.71 0.52 0.70 0.56 0.34 0.55 0.56 0.37 0.60 0.64