Weak Correlation between the Polyanion Environment and Ionic Conductivity in Amorphous Li − P − S Superionic Conductors

: Amorphous Li − P − S materials have been widely used as solid-state electrolytes for all-solid-state batteries because of their high ionic conductivity (10 − 4 to 10 − 3 S cm − 1 ) as well as good synthetic accessibility and processability. Despite the potential of these materials, their amorphous structures have made it challenging to quantify the relation between the structure and conductivity. In this paper, we use ab initio molecular dynamics simulations to investigate the role of the local structure and density in determining the conductivity of amorphous Li − P − S structures with different polyanion units. We observe similar rates for Li-ion hopping regardless of the local P − S polyanion environment in these amorphous materials, indicating that the path connectivity at a larger length scale may be controlling the overall Li conductivity. This finding will serve as an important guideline in the continued development of amorphous solid electrolytes for advanced all-solid-state batteries.


■ INTRODUCTION
−4 Solid-state batteries have received increased attention as a promising next-generation technology owing to their nonflammable nature and potential for higher energy density, achievable by using a Li−metal anode. 1,2,5The development of solid electrolytes (SEs) with high ionic conductivity and chemical stability is essential for the realization of practical solid-state batteries.
The Li−P−S system is one of the most promising compositional spaces for SEs as numerous crystalline and amorphous phases with high conductivity have been reported in this chemistry.−8 In particular, the 3Li 2 S:1P 2 S 5 (3Li:1P:4S) composition is of interest: the amorphous 3Li 2 S:1P 2 S 5 phase (10 −4 to 10 −3 S cm −1 ) 3,8,9 and nanoporous β-Li 3 PS 4 phase (∼8 × 10 −4 S cm −1 ) 10−12 have both shown high conductivity.Nevertheless, the room-temperature ground-state crystalline phase at this composition, γ-Li 3 PS 4 , exhibits low Li-ion conductivity (10 −8 to 10 −7 S cm −1 ), 11,13 emphasizing the importance of the local structure in determining the ionic conductivity. 14,15It has been suggested that different types of polyanion building blocks (i.e., PS 4 3− , P 2 S 6 4− , and P 2 S 7 4− ) (Figure 1b) exist in amorphous Li−P−S ionic conductors and that they play a key role in determining the Li-ion diffusivity. 8,16However, because of the difficulty in analyzing the structure of amorphous phases, the relationship between the local structure and Li-ion diffusion has not yet been clarified.
In this study, we use ab initio molecular dynamics (AIMD) to investigate the Li-ion diffusion mechanism in amorphous Li−P−S materials with different polyanion building blocks and densities and compare the results to Li transport in four crystalline phases.By analyzing the dynamic and static structures, we find that only a weak correlation between the local structure and conductivity exists in the amorphous phases.

■ RESULTS
Preparation of Model Structures.To investigate the effect of the P−S local environment on the Li-ion diffusion, we selected two descriptors: (1) the amount of different polyanion units (PS 4 ) and (2) the density.The composition is fixed to 3Li:1P:4S to unambiguously elucidate the effect of the local structure of polyanions and to remove indirect effects of the composition.
The nature of the P−S polyanions is believed to have a significant effect on the Li-ion diffusion in amorphous phases 8,16 because some crystalline superionic conductors such as β-Li 3 PS 4 , 12  , and P 2 S 7 4− polyanions. 8,16Li−P−S SEs are generally pelletized under a high pressure (∼400 MPa) to minimize voids between SE particles.−25 In total, 16 different amorphous structures were generated to examine the role of the polyanion environment and density on the Li ionic conductivity.The structure that only consists of PS 4 3− polyanions (128 atoms, 48 Li + ions, and 16 PS 4 3− ions) was selected as the reference amorphous structure.(See the Methods section for the detailed explanation.)Note that the known crystalline phases at this composition (β-, γ-, and α-Li 3 PS 4 ) only contain PS 4 3− polyanions (Figure 1a).In experiments, amorphous Li 3 PS 4 has been reported to consist of 70−90% PS   4− + S 2 2− .Moreover, off-stoichiometry amorphous structures without the free sulfur ions (Li 3 PS 4-δ , δ ≤ 0.25) were also studied because the existence of free sulfur ions or their dimers (S 2− or S 2 2− ) has not been experimentally clarified.For these structures, background charges were applied to preserve charge neutrality.Consequently, eight amorphous structures (Figures 1a and S1) with different P−S local environments were prepared, and the polyanion configurations for each are summarized in Figure 1b.After annealing at 1000 K, all the structures were found to be successfully amorphized without destroying the polyanion units.Using these structures, AIMD calculations were performed using the NVT ensemble.Two densities (1.8 and 2.0 g cm −3 ) were considered for all eight amorphous structures to investigate the effect of density.The densities of fully densified samples have been reported to be in the range 1.8 g cm −3 < ρ < 2.0 g cm −3 , and a 1 order of magnitude conductivity difference was observed between higher-density samples (∼2.0 g cm −3 ) and lower-density samples (∼1.8 g cm −3 ), even when the porosity was minimized. 8,21,22igure 1c shows the neutron-weighted pair distribution functions (PDFs, G′(r)) of the eight amorphous structures with a density of 1.8 g cm −3 .Only small differences in the shape of the PDFs can be observed among the amorphous structures with different polyanion environments (gray shade in Figure 1c).The experimental PDF of amorphous 3Li 2 S:1P 2 S 5 measured by Ohara et al. 16 is plotted with green stars.The PDFs of the generated model structures match the experimental curve well, except for the intensity of the first peak, which corresponds to the P−S bond distance.This discrepancy can be attributed to the significant sulfur deficiency in the experimental sample.While the P:S ratio should be 4 in the ideal structure, inductively coupled plasma (ICP) analysis 16 indicated a ratio of 3.6 in the synthesized material.
Li-Ion Conductivity.The Li conductivity was extracted from the AIMD Li + trajectories using the standard methodology as described in the Methods section.For the 16 amorphous and 4 crystalline structures, the activation energy of  16 diffusion and conductivity at 300K (Table 1 and Figure 2) were calculated from Arrhenius fits to the temperaturedependent conductivity between 600 and 1000 K (Figures S2−S4).
The predicted conductivity of all the structures corresponds well to the conductivity range observed in experiments 8,[10][11][12]21,22,26,27 and other theoretical studies.18,26,28−32 As observed in Figure 2, the conductivity of amorphous LPS with a density of 1.8 g cm −3 is in the range 10 −3 to 10 −2 S cm −1 (red squares), whereas that of amorphous LPS with a density of 2.0 g cm −3 is 1 order of magnitude lower, between 10 −4 and 10 −3 S cm −1 (blue squares). Experimental vlues have been reported to be 1 × 10 −3 (ref 21 ) and 2.8 × 10 −4 (ref 8 ) for 1.80 and 1.88 g cm −3 samples, respectively (gray stars). Considering that the experimentally measured conductivity is usually somewhat lower than the intrinsic conductivity due to interface effects, the results from the AIMD simulations on the amorphous structures agree reasonably well with the experimental findings. Throom-temperature conductivities of crystalline β-Li 3 PS 4 and Li 7 P 3 S 11 are calculated to be 2 × 10 −3 and 6 × 10 −2 S cm −1 , respectively, which are comparable to the measured values for nanoporous β-Li 3 PS 4 (∼8 × 10 −4 S cm −1 ) 10−12 and Li 7 P 3 S 11 (2 × 10 −2 S cm −1 ).26,27 These conductivity values are also similar to the results from other AIMD studies.18,26,28−31 The room-temperature conductivity of α-Li 3 PS 4 has not been experimentally measured because this phase is only stable above 740 K; however, the value in our calculation (1 × 10 −2 S cm −1 ) is similar to the result from another AIMD work of 8 × 10 −2 S cm −1 .,33 Because of the low conductivity of this phase, the number of Li hops in our simulation was not large enough to obtain a statistically meaningful mean squared displacement (MSD) within 500 ps.Thus, only structural analysis was performed for this compound.
Consistent with the experimental results, 8,21,22 Figure 2 shows that density has a pronounced influence on the conductivity of the amorphous structures, with the structures with a density of 2.0 g cm −3 exhibiting up to an order of magnitude lower room-temperature conductivity compared to that of the 1.8 g cm −3 structures.
Although about 1 order of magnitude variation (Figure 2 and Table 1) is observed among amorphous structures with the same density, it is difficult to identify a correlation between the conductivity and polyanion type.For example, 1P 2 S 6 (4.04 × 10 −2 S cm −1 ) and 2P 2 S 6 −4S (1.14 × 10 −3 S cm −1 ) are, respectively, the best and worst Li ionic conductors among the simulated amorphous structures with a density of 1.8 g cm −3 , despite both having P 2 S 6 4− anions.Similar examples can be found for structures containing P 2 S 7 4− anions: 2P 2 S 7 −2S (5.68 × 10 −3 S cm −1 ) and 1P 2 S 7 −1S (1.11 × 10 −4 S cm −1 ) are the best and worst Li ionic conductors, respectively, among the model structures with a density of 2.0 g cm −3 .Furthermore, the existence of sulfur deficiency does not significantly affect the Li-ion conductivity (e.g., ).Therefore, our results do not support the argument made in previous studies that the nature of the polyanion building blocks (PS 4 ) mostly affects the ionic conductivity. 8,16DISCUSSION

Effect of Polyanion Groups on Li Hopping Behavior.
To unambiguously examine the atomic scale effect that polyanions have on Li-ion diffusion, the movement of Li ions near different polyanion environments was statistically analyzed.The preference of Li binding to certain polyanions can be evaluated by measuring the time required for Li to escape from each specific P environment.Using the AIMD

D
trajectories of the amorphous structures at 600 K, we recorded the times for Li ions to migrate in and out of a 3 Å radius of specific P atoms.The difference between the two times is defined as the time required for Li to escape the environment (inset of Figure 3a).Note that the 3 Å radius was selected because it is the average distance of the first-neighbor Li−P (Figures 1c and 4d).The cumulative fraction of Li ions that have jumped away from the environment within a certain time range is plotted in Figures 3 and S5.We separately compared the Li escapes from PS 4 3− versus those from P 2 S 7 4− (Figures 3a and S5a−c) or PS 4 3− versus P 2 S 6 4− (Figures 3b and S5d,e) in the same structure to confirm whether a specific polyanion works as an obstacle for Li-ion diffusion.
The data in Figures 3 and S5 can be modeled by an exponential decay function (eq 1), indicating a statistical escape mechanism consistent with hopping motion where F, a, and t are the Li fraction that has escaped, the decay constant, and time, respectively.Here, the decay constants a can be interpreted as a Li hopping frequency Γ, E a , k, and T being the attempt frequency, activation barrier of hopping, Boltzmann constant, and temperature, respectively.The Li escape as a function of time can be successfully fitted to the decay function (eq 1) for all the polyanion types (R 2 > 0.95, Table 2).Under the assumption that Γ does not change in the same structure, one can evaluate the difference between the activation barrier for Li escape from two polyanion types A and B We find that the difference between the activation barrier for Li escape from PS 4 3− and from P 2 S 6 4− /P 2 S 7 4− is less than 6 meV for all the amorphous model structures (Table 2).This statistical analysis confirms that the Li-ion hopping behavior is not locally modified in any significant way by the nature of the polyanion group.
Local Environment Analysis.Since we observe that the polyanion environment does not significantly modify the local Li hopping, we examined in more detail the local structures by constructing the partial radial distribution functions (pRDFs) of amorphous and crystalline Li−P−S structures from the AIMD trajectories.Figures 4a and S6 show the P−P pRDFs of amorphous structures with different polyanion compositions.Consistent with the presence of different P−S groups, we find distinct characteristic peaks representing P−P bonding in different polyanion groups, for example, P−P in P 2 S 6 4− at 2.3 Å or P−P in P 2 S 7 4− at 3.5 Å, as well as distinct interpolyanion peaks in the 5−8 Å region.
The Li-related pRDFs, that is, Li−Li (Figure 4b), Li−S (Figure 4c), and Li−P (Figure 4d), on the other hand show virtually no variation among the different amorphous phases with the same density (1.8 g cm −3 ), indicating that Li interacts with them in the same way.The black line gives the average RDF with the gray shading indicating maximum and minimum values.
The lack of any variation among the amorphous structures is in direct contrast to the distinct and well-defined peak positions and intensities of Li-related pRDFs in crystalline phases (Figure 4b−d).A similar observation is made for the amorphous structures with a higher density (2.0 g cm −3 ), which also show little variation among them, although the position of the first-neighbor peaks shifts to the left compared with those of the lower-density structures because of the shortened atomic distances (Figure S6).Such analogous local Li environments in amorphous Li−P−S SEs have also been reported in another theoretical study. 34The paper showed similar shapes of Li-related pRDFs among amorphous phases, even for different stoichiometries, for example, a-Li 3 PS 4 , a-Li 7 P 3 S 11 , a-Li 4 P 2 S 7 , a-Li 4 P 2 S 6 , and so forth. 34elocity autocorrelation functions, calculated from the Li trajectories, are also similar among the amorphous structures (Figure S7).Together, these results qualitatively suggest that the Li local environments are not the origin of the variations of diffusivity among amorphous structures of the same density.
Self-Correlation Effect.Since Li-ion diffusion in amorphous structures occurs through a set of pseudorandom paths, both local mobility changes and longer-range variations in path connectivity can modify the diffusivity.The relation between Li diffusion, which characterizes the long-range displacement of ions, and local Li hopping can be described using the selfcorrelation factor 35 where D diffusion is the tracer diffusivity calculated from the MSD and D hop is an equivalent random walk diffusivity with the same number of Li hops.The self-correlation represents both topological and physical correlation factors.The shape and connectivity of the Li pathways create a topological correlation, while correlation in the hopping of different lithium ions can create physical contributions.For example, if a lithium ion strongly prefers to return to its original position immediately after hopping, the f value is expected to be low, as is, for example, the case in systems with strong Li-vacancy ordering, or highly varying site energies.In this study, we evaluate D hop from the sum of the travel distances of all Li hops 30 Here, the number of Li hops (N hop ) is counted by finding the number of times Li ions become displaced by some minimal distance (R H ) (see the Methods section for details), while n Li and t are the number of Li ions in the structure and simulation time, respectively.
Because D hop does not consider the shape of the Li path, it is related to local environmental descriptors that control the hopping rate.An activation barrier for hoppings (E hop ) can be calculated from an Arrhenius fit of D hop obtained at various temperatures (Figures S8−S10).Comparable values of E hop are observed for amorphous phases with the same density, as seen in Table 1.All the values of E hop lie within the narrow range of 160−200 meV, whereas a comparatively larger range of E diffusion from 180 to 340 meV can be observed among amorphous structures with a density of 1.8 g cm −3 (Table 1).This result again confirms our finding that the local environment of Li in amorphous structures is similar, which is in good agreement with the results of Figure 4.It is also noteworthy that the activation barriers for hopping and diffusion both increase when the density increases to 2 g cm −3 , indicating that density

Chemistry of Materials
affects the Li local environments, as shown in Figures 4, S6, and S7.The density effect can be explained by the previous observation that the Li-ion diffusion barrier generally increases when the volume per S atom decreases regardless of the sulfur framework. 14In contrast to our observation in the amorphous phases, we find that E hop varies a lot among the crystalline phases.For example, E hop equals 251 meV for β-Li 3 PS 4 but only 142 meV for Li 7 P 3 S 11 , consistent with their diverse Li local environments.
We believe that the large difference between D hop and D diffusion results from the tortuous Li path by which Li diffuses in amorphous materials.Tables 3 and S1−S3 show the selfcorrelation factor calculated at 600 K and 800 K for amorphous and crystalline Li−P−S structures.The selfcorrelation factor increases with increasing temperature in all the structures, which implies that higher energy channels become accessible at high temperatures, consistent with there being a distribution of energy barriers in an amorphous material.Therefore, E diffusion − E hop (the gap between the activation barrier of diffusion and that of hopping) can be interpreted as the energy required to shorten the diffusion paths.Indeed, as observed in Table 3, a larger difference between E diffusion and E hop results in a larger self-correlation factor difference between the 800 and 600 K values.

■ CONCLUSIONS
In this work, the relationship between Li ionic conductivity and local structures in amorphous Li−P−S was studied using AIMD simulations.The dynamic and static analysis reveals that the Li-ion diffusivity is similar among various amorphous structures with a given density regardless of the P−S local environment, disproving the conventional argument designating the polyanion environments as a key parameter determining the ionic conductivity.Instead, in our simulated structures, the diffusivity of the amorphous structures appears to be controlled by the complexity of the Li path as measured by the self-correlation factor.While, in principle, one could design an "optimal" local environment in amorphous compounds that facilitates Li hopping and try to make it percolate throughout the amorphous structure, it is likely to be challenging due to the lack of specific structure control when making amorphous materials.
The fact that it is the connectivity of pathways rather than the presence of specific local environments that controls the macroscopic Li conductivity may explain why the performance of amorphous Li−P−S conductors appears to be strongly processing-dependent.It seems likely that different degrees of randomness are induced when changing processing methods, and this will certainly affect how hopping paths connect.
The diffusion behavior of amorphous materials is in apparent contrast to that of crystalline phases, which have distinct Li local environments and diffusion behaviors based on strictly defined crystallographic sites.Unlike amorphous phases, crystalline phases show at best a few characteristic hops originating from their well-defined Li sites.Thus, an improvement of the conductivity in crystalline phases has been achieved by controlling the Li site energies [23][24][25]29 or exploring new phases with small Li site energy differences 14,36 to realize a flat Li local energy landscape (small E hop ) as well as a short diffusion path (small E diffusion ).
■ METHODS AIMD Simulations.The Li dynamics were studied using the Vienna ab initio simulation package (VASP) 37 based on the non-spinpolarized density functional theory (DFT) calculations with the projector augmented-wave method. 38The generalized gradient approximation 39 was used as the exchange−correlation functional.The kinetic energy cutoff of the plane-wave basis set was set to 300 eV.The Γ-only k-point AIMD calculations were performed with the NVT ensemble using a Nose−Hoover thermostat 40 and a 2 fs time step.For amorphous materials, the r-PS 4 (random) structure was first generated through the following steps: (1) selecting the lowest electrostatic energy configuration between 100 structures with randomly placed 16 P 5+ atoms; (2) generating PS 4  3− anions based on positions of P 5+ ions and placing 48 Li cations in the void space; and then, (3) relaxing the structure using a 20 ps AIMD at 1000 K.The volume and lattice vectors of the simulation cell were fixed to retain the cubic symmetry and reproduce the experimental density, 1.8 and 2.0 g cm −3 .Other amorphous structures were constructed based on the r-PS 4 structure by replacing two PS 4 3− anions with a single P 2 S 6 4− or P 2 S 7 4− ion and corresponding S 2− or S 2 2− ions.The structures were heated from 100 to 1000 K at a heating rate of 300 K/ ps and equilibrated for 20 ps.We confirmed that the polyanions were not destroyed during the equilibration process.Files with the atomic positions (POSCAR files) are provided in the Supporting Information.The equilibrated structures were quenched to the target temperatures, 600−1000 K.For the crystalline materials, the initial structures were imported from the Materials Project 41 and experimental reports. 33,421 × 2×2, 2 × 2×2, 1 × 2×2, and 2 × 1×1 supercells were used for β-, γ-, α-Li 3 PS 4 , and Li 7 P 3 S 11 phases, respectively.After 0 K relaxation, the structures were heated up to the desired temperatures in 3 ps.In all the AIMD simulations, an equilibration time of 3 ps was excluded from the MSD calculation, and the simulation was performed for more than 200 ps at each temperature to obtain statistically converged MSDs.The diffusivity D diffusion was calculated as where R i stands for the displacement of the i-th Li ion and the bracket denotes the time average.To evaluate the statistical variances of the calculated diffusivity, the AIMD simulations were divided into several pieces, from which the average and variances of diffusivity were calculated.This process was performed using the aimd.diffusionmodule developed by He et al. 43 The Li-ion conductivity was calculated using the Nernst−Einstein relationship where n is the number of mobile ions per unit volume, q is the charge of the moving ion, and D σ is the diffusion coefficient, which includes the interparticle correlation effect. 35However, a converged evaluation of D σ requires much longer Li trajectories than those needed for D diffusion and is thus computationally more demanding.In this study, D diffusion was used instead of D σ , ignoring the interparticle correlation effect.Atomic Structure Analysis.The pymagen 44 package was used for structural manipulations, generation of inputs for DFT calculations, and Li trajectory analysis.Neutron-weighted PDF and pRDF curves were calculated using the Diffpy package 45 and the pymatgen diffusion module, 46 respectively.The curves were generated by averaging every snapshot of the 600 K AIMD trajectory.
Jump and Self-Correlation Factor Analysis.The number of Li hops is counted by tracking individual Li ions.The displacement between the initial position and the position at a certain time step , the event is counted as a jump.Then, the origin of displacement is initialized to the current position, = = R R t t t 0 , trying to evaluate a new displacement value to find the next jump of the Li atom.The cutoff radius R H is set to 3 Å in this study.By counting the number of jumps for all Li ions during the simulation (N hop ), D hop was calculated from eq 5.The self-correlation factor was calculated from the ratio between D diffusion and D hop , as expressed in eq 4. The activation barrier for diffusion (E diffusion ) and that of hopping (E hop ) was calculated from an Arrhenius fit of D hop and D diffusion obtained at various temperatures where D diffusion,0 and D hop,0 denote pre-exponential factors for Li diffusion and hopping, respectively.In evaluating D diffusion,0 , we used the aimd.diffusionmodule developed by He et al., 43 which can calculate an uncertainty value based on the average and variances of diffusivities at different temperatures.

Figure 1 .
Figure 1.(a) Amorphous Li 3 PS 4 structure generated by randomly distributing Li + /PS 4 3− ions and melt quenching the structure at 1000 K using AIMD.(b) Polyanion composition of amorphous Li 3 PS 4−δ and known crystalline phases simulated in this study.(c) Averaged neutron-weighted PDFs (G′(r)) of the model amorphous structures (black line, gray shade covers between maximum and minimum values) compared with that of amorphous Li 3 PS 4 measured in the experiments (green stars).16

Figure 2 .
Figure 2. Conductivity at 300 K calculated via extrapolation of the Arrhenius plot of conductivity evaluated at five−seven temperatures.The experimental values of amorphous materials (gray stars) are taken from [a] ref 21 and [b] ref.8

a
The activation barrier difference (| | E a ) of Li escape between PS 4 3− and P 2 S 7 4− or PS 4 3− and P 2 S 6 4− is evaluated using eq 3.

Figure 4 .
Figure 4. (a) P−P, (b) Li−Li, (c) Li−S, and (d) Li−P pRDFs of amorphous structures with a density of 1.8 g cm −3 compared with those of crystalline Li−P−S phases.Black lines stand for the average of pRDFs of all amorphous model structures, while the gray shades cover the region between maximum and minimum pRDF values. 4

Table 1 .
Estimated Li-Ion Conductivity at 300 K and Activation Barriers of Li-Ion Diffusion and Local Hop Calculated for Amorphous Li 3 PS 4−δ and Crystalline Phases

Table 2 .
Decay Constants (a) Obtained by Fitting the Li Fraction vs Required Escape Time Curves to eq 1 ,a

Table 3 .
Difference between the Activation Barriers of Diffusion and That of the Hops and Self-Correlation Factor Calculated at 600 and 800 K for 1.8 g cm −3 Amorphous Li 3 PS 4−δ and Crystalline Phases https://doi.org/10.1021/acs.chemmater.2c02458Chem.Mater.XXXX, XXX, XXX−XXX is recorded for each Li ion.If a displacement becomes larger than the cutoff hopping distance (R H ) at a certain time step t″,