Andreas R. Stamminger*ab,
Benedikt Ziebartha,
Matous Mrovecb,
Thomas Hammerschmidtb and
Ralf Drautzb
aRobert Bosch GmbH, Corporate Research, 71272, Renningen, Germany. E-mail: andreas.stamminger@de.bosch.com
bRuhr-Universität Bochum, Interdisciplinary Centre for Advanced Materials Simulation, Universitätsstraβe 150, 44801, Bochum, Germany
First published on 13th March 2020
The synthesis of Li superionic conductor Li7P3S11 may be accompanied by the formation of a detrimental Li4P2S6 phase due to a high mixing sensitivity of precursor materials. This phase exhibits a poor ionic conductivity whose origins are not fully understood. Recently Dietrich et al. investigated the energetics of Li ion migration in Li4P2S6 with nudged elastic band (NEB) calculations. The observed large migration barrier of 0.51 eV for purely interstitial diffusion leads to an interpretation of the low ionic conductivity by kinetic limitations. Based on ab initio molecular dynamics simulations (AIMD) we propose a new and energetically much more favorable diffusion path available to interstitial Li ion charge carriers that has not been considered so far. It consists of a concerted process in which a second lithium atom is pushed out from its equilibrium lattice position by the diffusing lithium ion. A detailed analysis with NEB calculations shows that the energy barrier for this concerted diffusion is only 0.08 eV, i.e. an order of magnitude lower than the previously reported value for purely interstitial diffusion. Therefore, the observed low ionic conductivity of Li4P2S6 is likely not originating from kinetic limitations due to high diffusion barriers but rather from thermodynamic reasons associated with a low concentration of free charge carriers. We therefore expect that increasing the charge carrier concentration by doping is a viable design route to optimize the ionic conductivity of this material.
Recently, chalcogenide-based compounds,7–9 such as Li7P3S11, have been identified as the most promising SSE candidate materials.10 A glass-ceramic phase of Li7P3S11 was found to exhibit a remarkably high ionic conductivity of up to 10−3 S cm−1 at room temperature.4,11–14 This glass-ceramic phase can be prepared by a careful heat treatment from a mixture of Li2S and P2S5 precursor materials with a ratio of 70:30.15–17 Moreover, small deviations from this optimal composition lead to a dramatic decrease of ionic conductivity in the resulting materials. For instance, a slightly altered ratio of 67:33 leads to the formation of glass-type Li4P2S7. Although synthesized from the same precursors this material shows ionic conductivities of only up to 10−5 S cm−1, i.e., two orders of magnitude below the conductivities reported for Li7P3S11.1,4 An even lower ionic conductivity of about 10−7 S cm−1 at room temperature has been found for the thermodynamically stable Li4P2S6 phase that crystallizes readily from the Li4P2S7 glass phase at elevated temperatures and is associated with a decrease of the sulfur content.4,14,17–21
As the formation of these detrimental phases is unavoidable during large-scale production (e.g., due to local inhomogeneities in the composition or temperature fluctuations), it is important to understand the origins of their low ionic conductivities. Recent studies showed that the ionic conductivity in several electrolyte materials is primarily determined by the correlated motion of charge carriers. For instance, He et al.22 showed that the diffusion of lithium ions in the high-conducting Li7La3Zr2O12 (LLZO), Li1.3Al0.3Ti1.7(PO4)3 (LATP) and Li10GeP2S12 (LGPS) is mitigated by a complex migration process. They showed that the occupation of high-energy sites by mobile charge carriers induces concerted migration processes with significantly lower migration barriers than direct diffusion between interstitial sites. Similar results for Li2CO3 (ref. 23) and β-Li3PS4 (ref. 24–26) showed that the observed high lithium ion conductivity in these materials can be explained by the simultaneous movement of lithium ions. Likewise, Zhang et al.27 reported a cooperative diffusion mechanism for Na+ ions in Na3Zr2Si2PO12 (NaSICON). Despite these findings of the paramount role of concerted migration in several electrolyte materials, surprisingly, up to now only purely interstitial diffusion processes have been considered in the estimation of migration barriers for Li4P2S6. Only Hood et al.18 already suggested the possible importance of more complex migration paths in this material, involving both interstitial and host lattice ions.
Inspired by this gap in literature and the recent reports on low migration barriers for concerted diffusion in other superionic conductor materials, we investigate the possibility of concerted diffusion of Li ions in the most detrimental, crystalline Li4P2S6 phase. Our theoretical investigation employs extensive first-principles calculations based on density functional theory. We combine ab initio molecular dynamics (AIMD) and the nudged elastic band (NEB) method to obtain an in-depth insight into Li diffusion. While AIMD simulations enable the sampling of a large number of Li migration trajectories in an unbiased way, the NEB calculations provide a detailed and quantitative analysis of minimum energy paths. This combined approach allowed us to identify an additional migration process in Li4P2S6 that cannot be described by independent jumps of a single interstitial or vacancy but rather by concerted interaction of two Li ions.22,28,29 Moreover, the identified process is found to create a two dimensional diffusion network with a low energy barrier, available for all present charge carriers. This finding reveals a high mobility for interstitial lithium ions, in addition to the high mobility found for vacancies by Dietrich et al.,1 ruling out kinetic limitations as reason for the observed low ionic conductivity of Li4P2S6.
The paper is organized as follows. A detailed discussion of the computational methods and parameters is provided in Section 2. The results obtained from our simulations are presented in Section 3 followed by a detailed discussion of the results in Section 4. Finally, we conclude in the last Section 5.
To ensure the validity of our results for the investigated Li4P2S6 structure, we first performed static calculations to converge the total energy of the structure with respect to the plane-wave energy cutoff and the k-point sampling of the Brillouin zone. For a primitive unit cell containing a total of 12 atoms, an energy cutoff of Ecut = 420 eV and a (3 × 3 × 2) k-point mesh are sufficient to converge the total energy to within 0.001 eV per atom. Using the determined cutoff and k-point parameters both the volume and atomic positions of the initial structure were relaxed to obtain an equilibrated configuration for the subsequent NEB and AIMD simulations. To achieve a well relaxed unit cell, tight convergence thresholds of 10−9 eV and 10−4 eV Å−1 have been used for the self-consistency cycle and the forces, respectively.
For the AIMD simulations, the convergence criterion for the total energy was set to 0.01 eV per atom for the sake of computational efficiency. From convergence tests we verified that for a supercell containing 96 atoms, used in our AIMD simulations, a Γ-point only calculation with a plane-wave cutoff of Ecut = 300 eV suffices to fulfill this criterion.
Apart from the validation of the plane-wave cutoff and the k-point sampling, we also used our NEB calculations to test the performance of the available Li/Li_sv pseudo potentials treating 1s-electrons either as core (Li) or valence (Li_sv) electrons. From this benchmark we found that both potentials yield the same results and hence, the computationally less demanding Li pseudo potential was employed in all our calculations.
Species | Wyckoff | x | y | z |
---|---|---|---|---|
Li1 | 2c | 0.3333 | 0.6667 | 0.0 |
Li2 | 2d | 0.6667 | 0.3333 | 0.5 |
P | 2e | 0.0 | 0.0 | 0.1715 |
S | 6k | 0.3237 | 0.0 | 0.25 |
Fig. 1 Unit cell structure of Li4P2S6 with constituent atoms placed according to the Wyckoff positions defined in Table 1 (P: light purple, S: yellow, Li: green). Additionally, a single Li atom (marked in red) was placed at the interstitial ‘1b’ Wyckoff position used as initial position for AIMD and NEB calculations of Li migration. |
(1) |
The mean squared displacements obtained from 100 ps AIMD simulations at temperatures between 400 K and 500 K, shown in Fig. 2, exhibit several discontinuities corresponding to Li migration events. Such a step-like behavior adds a large uncertainty to the evaluation of the average MSD(Li)(t) values. To obtain more reliable results we increased the number of migration events recorded in the AIMD trajectories. To this end, we performed AIMD simulations at elevated temperatures ranging from 800 K to 1000 K to sample a larger number of migration events in our trajectories. To reduce statistical noise, the total simulation time was increased to 300 ps for the calculations. However, despite the larger temperatures and the simulation time being tripled, the recorded MSD(Li)(t) dependencies still show sizable fluctuations (see ESI Fig. S1†).
Fig. 2 Mean squared displacements obtained from AIMD simulations between 400 K and 500 K. Only few migration events corresponding to the discontinuities are observed. |
The reason for this behavior is found from analyzing the MSD components along the directions of the crystal axes. From this analysis, we obtain information about the anisotropy of Li migration, which is shown in Fig. 3 for the AIMD simulation at 1000 K. This decomposition reveals that there is no diffusion component in the [001]-direction, i.e. the Li diffusion is purely two-dimensional and confined to the {001} planes only. This strong anisotropy greatly reduces the number of migrating lithium ions only to those located in the {002} planes, where also the additional Li atom was added (cf. Fig. 1). As a consequence, only 9 out of the total of 33 Li ions in the supercell contribute to the MSD, enhancing statistical noise.
Fig. 3 MSD components along the [100], [010] and [001] lattice directions (cf. Fig. 1) obtained from an AIMD simulation at 1000 K. |
Fig. 4 Temperature dependence of Li diffusivity in Li4P2S6 obtained from AIMD simulation data following eqn (1). Linear regression gives a low activation energy of about 200 meV for lithium migration. |
To get a rough estimate for the activation barrier from the MSD data we assume an equilibration time of 10 ps. For a better statistical significance we split the remaining total MSD evaluation into 10 sub-samples MSD(Li)i(τ) with a length of τ = 29 ps for each sample. Diffusivities for the mobile lithium ions are then obtained from linear regression to the averaged MSD dependence
(2) |
(3) |
Based on our simulation results, we estimate the corresponding activation barrier for Li interstitial diffusion in Li4P2S6 to be EAIMDa = (200 ± 50) meV, with prefactor D0 = 6.1 × 10−5 cm2 s−1.
Although our AIMD simulations at lower temperatures are not reliable for extraction of diffusion data, the small number of migration events makes these AIMD trajectories valuable for tracking and in-depth analysis of the underlying migration paths. To this end, we further analyzed the MSD dependence recorded at 310 K, shown in Fig. 5. In this AIMD simulation, only a single jump event occurred during the 100 ps run, greatly simplifying the task of identifying the underlying processes and tracking the corresponding migration path. Our analysis reveals a migration of the interstitial Li ion between two ‘1b’ Wyckoff positions, that cannot be described entirely by independent migration events of a single Li ion. Instead, the complete process is composed of several subsequent steps. In the first step, an interstitial Li ion migrates from its initial ‘1b’ Wyckoff position to an intermediate ‘6k’ Wyckoff position. This is followed by the knock-out of a bound lithium atom from its ‘2d’ Wyckoff position (Li2, cf. Table 1) to another ‘6k’ site (shown schematically in Fig. 6), and the incorporation of the migrating ion into the vacated ‘2d’ site. In the last step, the released lithium ion migrates to a nearest interstitial ‘1b’ Wyckoff position that is energetically most favorable.
To explore the energetics of this migration mechanism, we analyzed the whole pathway using the NEB method. For the NEB calculation, the complete migration process was divided into three segments. The first segment corresponds to the migration of the interstitial Li ion from its initial ‘1b’ Wyckoff site to a neighboring ‘6k’ Wyckoff site as shown in Fig. 6a and b. The second segment deals with the actual concerted process displayed in Fig. 6b–d. The final migration of the displaced lattice ion to the next ‘1b’ Wyckoff position corresponds to Fig. 6d and e. The computed minimum energy path (MEP) is shown in Fig. 7. It is symmetric, as expected from the geometry of the Li4P2S6 crystal structure, and contains two local minima corresponding to configurations where the diffusing Li ion occupies the ‘6k’ Wyckoff sites. The migration barrier for the forward and backward jumps between the ‘1b’ Wyckoff site and the metastable ‘6k’ Wyckoff site is about 80 meV and 40 meV, respectively. The barrier associated with the concerted process is symmetric with a height of about ENEBconcerted = 80 meV. The overall barrier for the complete migration process equals to ENEBa = 120 meV and is significantly lower compared to previously reported barriers for purely interstitial diffusion (see Fig. 8 for a direct comparison).
Fig. 7 Energy profile of an interstitial lithium ion migrating between the most stable interstitial positions at the ‘1b’ Wyckoff sites via the ‘6k’ Wyckoff sites. The energy barrier necessary to trigger the concerted process, allowing for migration between different P2S6 structures, is shown within the gray box (data points annotated by (a), (b), (c), (d) and (e) correspond to the intermediate structures shown in Fig. 6). |
Fig. 8 Energy barriers for the diffusion path ‘1b’–‘6k’–‘6k’ with and without the concerted process calculated in this work compared to available literature results for purely interstitial diffusion. |
By analyzing the process using NEB calculations, we found similar activation barriers of about ENEB1b–6k = ENEBconcerted = 80 meV for the forward ‘1b’ → ‘6k’ interstitial jump and the ‘6k’ → ‘2d’ concerted process. However, since the ‘6k’ interstitial site is energetically less favorable than the ‘1b’ site, the migration barrier for the first jump is not symmetric (cf. Fig. 7). As the backward ‘6k’ → ‘1b’ jump is associated with a barrier of only 40 meV, i.e., about half of that for the forward jump, the probability for backward jumps significantly exceeds that of the forward jumps. It is therefore likely that a single concerted migration event will include a large number of preceding back and fourth jumps between the ‘1b’ and ‘6k’ sites. This presumption is fully consistent with our AIMD simulations. In Fig. 5 four smaller peaks, corresponding to the jumps between ‘1b’ and ‘6k’ Wyckoff sites, are visible between 0.0 and 60.0 ps before a concerted process happens at about 62.0 ps.
Our predicted activation barrier of ENEB1b–6k = 80 meV for the ‘1b’ → ‘6k’ migration is in excellent agreement with the result found by Dietrich et al.1 However, the energy barrier for the complete process of concerted diffusion (see Fig. 6) of about ENEBa = 120 meV is much lower than the energy barrier of about 500 meV for the in-plane Li migration found in their work.1 This difference can be attributed to the two distinct migration mechanisms. The barrier investigated by Dietrich et al. corresponds to a direct interstitial jump between two ‘6k’ sites, where the migrating Li ion has to pass between two closely-spaced Li atoms occupying the ‘2d’ lattice sites. Due to strong steric and Coulombic interactions the corresponding barrier is significantly higher than the barrier associated with the concerted mechanism investigated in this work. The low energy barrier of the concerted process is consistent with the low estimated activation energy of EAIMDa = (200 ± 50) meV derived directly from our AIMD simulations (cf. Fig. 4).
Our study illustrates limitations in predicting diffusion properties based solely on NEB calculations since these are usually restricted to anticipated migration mechanisms. Based on AIMD simulations, we were able to show that Li4P2S6 features concerted interstitial migration leading to a low energy barrier compared to purely interstitial migration.
From the energy barriers obtained by our NEB calculations, we assume that migration paths enabling long range 2D diffusion with low activation barriers exist in crystalline Li4P2S6. To support this assumption we analyze the lithium probability density calculated from the AIMD simulation at 1000 K. The calculated distribution of lithium ions shows a 2D diffusion network of adjacent Wyckoff ‘1b’ sites connected by the concerted migration process. Furthermore, a negligible probability is found for the passage of an interstitial lithium ion through two neighboring Wyckoff ‘2d’ positions (see diagonal path in Fig. 9). This confirms that the correlated jumps, described in this work, are the dominant migration process while purely interstitial migration remains inactive even at elevated temperatures due to their associated large activation barriers.
Assuming an Arrhenius behavior over the whole temperature range, we estimate the room temperature diffusivity to D(300 K) ≈ 3.2 × 10−8 cm2 s−1. If we consider the interstitial Li ions to be the only active charge carriers, the theoretical room temperature conductivity amounts to σ(300 K) ≈ 1.2 × 10−4 S cm−1. This value is only one order of magnitude less than the experimentally determined room temperature conductivity of the super-ionic conductor Li7P3S11 (σexp0(300 K) ≈ 3.2 × 10−3 S cm−1 corresponding to an activation barrier of about Eexpa ≈ 0.19 eV (ref. 12 and 17)). However, this is in sharp contrast to the low conductivities measured for crystalline Li4P2S6.1,4,18 This large discrepancy may have several reasons. As mentioned above, the fast diffusion channels are limited to the {001} planes only forming a 2D diffusion network. In addition, the estimated theoretical conductivity implies a large charge carrier concentration of about 9.7 × 10−4 mol cm−3. However, this required charge carrier concentration cannot be reached naturally due to the large formation energy of about 1.0 eV for Frenkel defects in Li4P2S6 (ref. 1 and 44) suppressing an efficient formation of free charge carriers. Hence, even though migration paths with low activation barrier exist in crystalline Li4P2S6, its ionic conductivity is ultimately controlled by thermodynamics leading to an intrinsically small charge carrier concentration. A possible remedy to this problem could be to increase the amount of interstitial Li ions artificially by aliovalent doping of Li4P2S6, which may be subject of a future work. This would bypass the problem of the intrinsically low defect concentration and would allow to exploit the discovered low energy migration path to increase the ionic conductivity of the material.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/d0ra00932f |
This journal is © The Royal Society of Chemistry 2020 |