Microscopic theory of ionic motion in crystals


 The drive towards an electricity based economy and emerging green technologies has resulted in a tremendous push to create safer and more efficient energy storage devices.1–3 The development of solid-state batteries is a major effort in this direction4,5. Unlike the case of a traditional electrochemical apparatus, in solid-state batteries ions move through a solid crystalline electrolyte. Ionic motion is thus intimately linked to the condensed matter description of the system – that is, the periodic electronic and ionic properties of the crystal - and is not adequately described by the existing electrochemical tenet. In the present article, we propose a microscopic, first-principles, description of the ionic conduction in crystals. This allows us to understand the ideal characteristics of materials for ionic conduction in general, and for solid-electrolyte applications in particular. Using ab initio calculations, we show that our formalism results in ionic mobilities consistent with experiments for several materials. Our work opens the possibility for the development of solid electrolytes based on fundamental physical principles rather than empirical descriptions of the underlying processes.

The semiconductor technology revolution that occurred in the middle of the 20 th century has its origin in the development of the microscopic theory of electron motion in crystals in the early days of quantum mechanics. Trailblazers such as F. Bloch showed the importance of taking into account the periodicity of the crystal structure in order to understand electronic behaviour and the properties of many different types of materials, from metals to semiconductors 6 . Not only did this theoretical framework become important for the understanding of naturally occurring materials, it also allowed for the development of new materials with tailored properties that did not exist before. Electronics turned from a heuristic discipline, based on trial and error, into a predictive science and technology. With the development of superior semiconductor-based technologies, such as the one based on complementary metal-oxide-semiconductors (CMOS), these basic concepts became a common language between scientists and engineers. Such extraordinary developments were driven by the necessity of replacing the obsolete vacuum tube technology that was ubiquitous to the electronic devices of that era.
Arguably, we find ourselves in a similar situation in the area of energy storage. The electrochemical science that has been the basis of battery technology for the last 200 years is faced with the necessity of reinventing itself in order to fulfil the societal needs of today, namely, safe, efficient, reliable, fast, and long-lasting energy storage devices that can be seamlessly incorporated into a modern environmentally conscious lifestyle. The pressure felt by important industrial sectors, such as the automotive, has inspired scientists and engineers to look for alternative solutions to traditional approaches. As such, the field of solid-state batteries has emerged as a possible solution to the conundrum of developing commercially viable electric vehicles 4,5,7 .
In a solid-state battery, the movable ions, such as lithium (Li + ), traverse a crystal, that is, a periodic structure consisting of fixed ions (that we call the framework) along with their electrons, and interact with these elements via strong Coulomb forces. Unlike the traditional electrochemical problem in a liquid medium, the ionic motion in crystals depends strongly on their periodicity and symmetries, as in the case of Bloch's theorem. The problem at hand is akin to the famous many-body problem found in strongly interacting electron materials (where phenomena such as magnetism and superconductivity occur) with some fundamental differences, namely, the mass of an ion is approximately 10,000 times that of an electron and while electrons have a strong wave-like characteristics, ions are essentially particles with atomic size. The intricate dance between particles and waves in a crystalline environment is what determines the ionic conductivity and the ultimate efficiency of a solid electrolyte in a battery.
The objective of this work is to develop the basic principles and microscopic formalism that describe ionic motion in crystalline electrolytes. Instead of looking at the problem from a traditional electrochemical perspective, we take a modern condensed matter approach and include the basic elements (mobile and fixed ions and their electrons, in addition to the periodicity of the crystal) from the very beginning. We obtain the steady state equation for the ion motion in a crystal and show that it obeys a Langevin dynamics, consistent with the fluctuation-dissipation theorem, where the ion mobility is determined by the curvature profile of the potential in the atomic lattice and the low frequency phonons of the framework.
In order to substantiate our results, we make use of state-of-the-art density functional theory (DFT) and ab initio molecular dynamics (AIMD) simulations in order to calculate the ion mobility for some crystalline electrolytes and show that our results are consistent with experimental values in polycrystalline samples.
The microscopic Hamiltonian for the problem is given by: where K refers to the kinetic energy and V to the Coulomb potential of interactions for electrons (subscript e) and ions (subscript i). Given that the mass of the ions is much larger than that of the electrons, we can consider mobile and framework ions to be static in the time scale of motion of electrons (the Born-Oppenheimer approximation). Furthermore, a solid-state electrolyte is characterized to be an insulator for electrons, which implies that the material has a large band gap in the electronic structure so that electron-hole excitations are not created during the ionic motion. This implies that one can use the adiabatic theorem and assume that the ionic motion only leads to smooth modifications of the electron-ion interaction. Within these two standard approximations one is left with an effective ion-ion interaction, U(r,u), which depends on the position of the mobile ions, r, and the framework ions, u (mathematical details will be presented elsewhere).
By definition, framework ions are the ones that remain in their crystal lattice position during the flow of the mobile ions. Hence, at any temperature T below the melting point of the crystal, the framework ions undergo oscillatory motion around their equilibrium positions, as illustrated in Fig. 1. Once again, we can use standard solid-state language to describe the framework ions in terms of their phonons, which are characterised by their frequencies Ωs, where the subscript s labels the phonon modes. Of particular importance are the acoustic phonons, or sound modes, with speed of propagation vL,T (where L and T labels the longitudinal and transverse modes). One should stress that these modes do not exist in a liquid or a traditional electrolyte. They only exist in a crystal because a periodicity is induced by the presence of the lattice.
We use the non-equilibrium Keldysh formalism in the path integral representation in order to get the equation of motion for the mobile ions in the form of Newton's equation: where t is time, u0 is the equilibrium position of the framework ions (∇ is the gradient with respect to r) M is the ion mass, F is an applied external force (i.e., electric field), × is the framework mass matrix for ions with mass and D is the system dimensionality, εs is the phonon polarization vector). In Equation (2), the first term on the r.h.s. is a classical term that describes the periodic potential of the static lattice with the framework ions fixed to their equilibrium positions. The second term on the r.h.s. describes the dissipation of energy due to the ion motion when the mobile ion interacts with the phonons of the framework ions. We note that unlike the dynamics in liquids, the dissipative term is position and direction dependent and described in terms of a tensor. Hence, dissipation in a crystal is not isotropic as in a fluid and has quantum nature (as can be seen from the presence of Planck's constant in its definition). This term is rather non-trivial, as we will see below. ̃( ) is the fluctuation force due to the vibrations of the framework. Its origin is the same as the dissipative term. One can show that its correlation function is given by: where KB is the Boltzmann constant. At high temperatures, ≫ ℏ Ω , we can readily see that: which is a generalized fluctuation-theorem for ionic motion in crystals. Hence, dissipation and fluctuation are intimately related. Once again, unlike in liquids, the fluctuation forces are position and direction dependent.
Having established the existence of a Langevin dynamics of mobile ions in a crystal lattice we can simplify the problem further by considering the steady-state motion in the presence of a constant applied electric field, E. In this case, the acceleration and fluctuation terms vanish, and Eq.(2) can be solved for the ion velocity, ( ) = ̇ ( ): (7) where q is the ion charge. We now can use the fact that the lattice potential and the dissipation have the periodicity of the lattice and, hence, can be expanded in a Fourier series in terms of reciprocal lattice vectors, K, in order to get: which can be thought of as the equivalent of Bloch's theorem for ionic motion in a crystal. If we are interested only in the average drift velocity we can take → 0 in the above equation and find: where 〈 −1 〉 = (2 ) 3/2 ∑ −1 (10) is the lattice-averaged dissipation coefficient. From the above expression we can readily obtain the ion mobility: where 〈 −1 〉 can be computed from first principles for any crystal lattice.
In order to gain more insight into the problem we can simplify the problem considerably by assuming that in Eq. (5) only the acoustic phonon modes contribute to the problem. In this case one can show that: where ρ is the mass density of the crystal and Hr is the Hessian operator. We now can see very clearly the dissipative mechanism of ionic motion in a crystal, namely, the softer the crystal (smaller sound velocity) the more dissipative the ionic motion. Furthermore, the movement of the ion is dissipationless in regions of the potential where the Hessian vanishes, which are the saddle point regions of the periodic potential in the unit cell of the crystal. Furthermore, this expression gives us clues regarding what kinds of crystals would be good electrolytes, namely, hard crystals with smooth potential configurations.
The variations of the potential energy surface can be quantified by computing U(r) from first principles, which we do, as illustration, for the metal-halide electrolytes AgCl, LiCl, LiI, α-AgI, and α-CuBr, as shown in Fig. 2 for AgCl and α-AgI. From here, we can obtain μ, per mobile ion, for each compound via the calculation of γ(r), as defined in Eq. (12). These ion mobilities, μcalc, assuming q=e, are listed in Table 1 . For comparison, we also extract the ionic mobilities per mobile ion, μexpt, from available experimental data. The mobility values obtained from fitting the experimental data are lower than the theoretical values, as expected, due to the non-ideal conditions in the experiments. Experimental samples are polycrystalline with reported size and space-charge effects; samples may have more than one mobile defect or ionic species; and the conductivity is usually measured using AC measurements, with sample impedance depending on the frequency. The higher theoretical values show that there is still ample room for improvement in the experimental mobility. Despite these many complications, we note that µcalc, computed using only Eq. (11) and simple first-principles inputs, shows a remarkable consistency with µexpt. In summary, we have developed a microscopic theory for ionic motion in crystals and derived the equivalent of Bloch's theorem for ions. We found that the ionic mobility depends essentially on the lattice softness (via the third power of the sound velocity) and the curvature of the atomic potential felt by the ions; namely, hard materials with smooth atomic potentials are the best candidates for high ionic mobility. Our calculations for the ionic mobility of several ionic conductors using ab initio methods provide an upper bound for the mobility in single crystals and indicate that the currently measured ionic motilities found in the literature are dominated by extrinsic effects such as impurities, grain boundaries, and other types of defects that are already well-known in solid-state physics. Although in the last century we have developed a powerful theoretical framework to study the effect of defects and interfaces in the motion of electrons in solids, the same cannot be said for the case of ions. The understanding of how ions interact with defects and interfaces in solids is an unexplored land. Any further progress in the development of solid-state electrolytes, which are the key elements of solid-state batteries, depends fundamentally on progress in this area of research.

Methods
Density functional theory simulations at 0 K DFT calculations are performed using the Quantum ESPRESSO 12,13 code. Structural relaxations and total energy calculations are performed using a PAW basis 14,15 and the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional 16 . The kinetic energy cutoffs of the charge density and wavefunctions are set to at least the minimum recommended values of the PAW pseudopotential 15 . The Brillouin zones for all materials are sampled using a uniform 6×6×6 grid of K-points.
For the calculation of U(r) we allow one ion of the mobile metal species to move while keeping all other ions fixed. The mobile ion is moved within the cubic unit cell by intervals of 1/64 of the lattice parameter, a. Only configurations in which the distance from the mobile ion to any fixed ion is greater than (5/12)a (α -AgI and α -CuBr) or (1/3)a (AgCl, LiCl, LiI) are permitted. In the cases of α -AgI and α-CuBr we note that the metal ions have partial occupation and thus multiple possible positions. We therefore compute the total energies of all possible permutations of the positions of the mobile ions within the unit cell. The resulting minimum energy configurations are those in which the metal ions are located at the tetrahedral positions on adjacent faces, in agreement with AIMD calculations for α -AgI 17 . Accordingly, to compute U(r) for these materials we fix one of the tetrahedral metal ions and allow the other to move to all other permitted positions, as described above.
Phonon calculations are performed using SG15 optimized norm-conserving Vanderbilt (ONCV) pseudopotentials 18,19 and a PBE exchange-correlation functional, with a 60 Ry kinetic energy cutoff for wavefunctions. Sound velocities are derived from the phonon dispersions along the Γ-X path for cubic cells, as defined in Ref. 20 .
The volumetric images shown in Fig. 2 were generated in VESTA. 21 Ab initio molecular dynamics simulations AIMD simulations are carried out using the SIESTA code. 22 The forces are calculated using the local density approximation (LDA) of density functional theory, 23 and a Harris functional is used for the first step of the self-consistency cycle.
The core electrons are represented by pseudopotentials of the Troullier-Martins scheme. 24 The basis sets for the Kohn-Sham states are linear combinations of numerical atomic orbitals, of the polarized double-zeta type. 25,26 The Γ -point is used for Brillouin zone sampling.
α-AgI AIMD calculations were are performed in 256 atom supercells. AIMD calculations for rocksalt structures are performed in 216 atom supercells.
The temperature is controlled by means of a Nosé thermostat 27 . The integration time step used is 1 fs and the total integration time is 26 ps. The equilibration time varies between different temperatures and systems and is determined by examining the mean square displacement.
Calculations of the ionic mobility from existing experimental data Experimental mobilities are found by fitting existing experimental data in the literature. We have assumed that there is only one mobile ion species. The conductivity is fitted using σ = qnµ, where n = Nexp(−Ea /kT ) is the number of mobile ions. Here, N is the total number of atoms of the mobile species in the crystal, Ea is the activation energy necessary to make the ion mobile, k is the Boltzmann constant and T is the temperature.