Li₄P₂Se₆ – structure and properties

Abstract

Non-oxide chalcogenides are among the fastest Li–ion conductors currently available. Here, a new crystalline lithium selenido-phosphate, Li₄P₂Se₆, is reported. Its structure was solved using a combined approach of quantum-chemical structure prediction, powder X–ray diffraction, and solid-state NMR. Li₄P₂Se₆ crystallizes in an orthorhombic unit cell in the space group Pnma and the lattice parameters a = 13.8707(1) Å, b = 11.2115(1) Å, and c = 6.45445(7) Å, representing a novel structure type. ³¹P and ⁷⁷Se magic–angle–spinning NMR spectra are reported and were analyzed for chemical shift and J-couplings. The material's ionic conductivity, characterized by impedance spectroscopy, slightly exceeds that of the non-isostructural crystalline compound Li₄P₂S₆, which contains the homologous complex anion. A computed phase diagram shows that several compounds can be prepared experimentally, which from the computed energy at a temperature of 0 K, feature a positive energy above the convex hull and thus should not exist. This research highlights the influence of configurational and vibrational entropy in stabilizing ionic chalcogenides, emphasizes the complexity of predicting phase stability in related systems by quantum-chemical calculations and contributes to the understanding of non-oxide chalcogenides, including potential fast ionic conductors.

---

1. Introduction

Non-oxide chalcogenides are promising candidates for fast solid ionic conductors.¹ While the Li–P–S system has recently attracted significant attention, its related system Li–P–Se remained relatively unexplored, most likely due to the rarity of Se as comparison to S. An advantageous feature of Se is that its NMR–active isotope, ⁷⁷Se, has a higher natural abundance (7.6%) than ³³S (0.8%). Additionally, ⁷⁷Se has a nuclear spin of I = 1/2, resulting in spectra that are easier to interpret than those of the ³³S quadrupolar spin of I = 3/2. Because the structural characteristics of selenium- and sulfur-containing compounds tend to be similar, meaningful comparisons can be drawn from structural studies, which can probe mobility on the chalcogenide position.

Only a single ternary crystalline phase in the system Li–P–Se is currently known, corresponding to the high–temperature modification of Li₇PSe₆.² This compound crystallizes in the well-known cubic argyrodite structure and exhibits fast Li-ion conductivity, similar to its sulfur-based analogue.³ Efforts to enhance the ionic conductivity of this phase through aliovalent substitution with different elements have been reported.^4,5 Additionally, a previous study provided strong evidence for another crystalline phase, “Li₄P₂Se₆”, based on solid-state NMR and powder X-ray diffraction data, including a suggested unit cell; however, its structure remained unsolved.⁶ A key question to be addressed is whether the crystalline “Li₄P₂Se₆” adopts the same structure as its sulfur analogue Li₄P₂S₆.^7,8

A previous study⁹ on aliovalent doping of Li₄P₂S₆ with magnesium demonstrated that several energetic minima can be accessed by only minor compositional changes, suggesting that the structure of Li₄P₂Se₆ could well differ from the sulfur variant. Furthermore, it has been theoretically proposed that “Li₃PSe₄” may exist in the system, but this has yet to be confirmed experimentally.¹⁰

The energetic hypersurface of existing structures can be explored through a combination of quantum chemical energy calculations and search algorithms like Monte Carlo or genetic algorithm.^11–13 Such structure prediction techniques provide insights into the stability of different crystalline phases by probing the convex hull of all known and hypothesized crystalline structures. Moreover, they allow identification of competing structures for a given composition and estimation of their relative formation energies.

The purpose of this study is to determine the structure and properties of Li₄P₂Se₆ and to relate them to other crystalline phases within the Li–P–Se system.

2. Experimental part

Synthesis

All starting materials were stored and prepared inside a glove box (MBraun, Garching, Germany) under an argon atmosphere. Li₂Se was prepared by reacting selenium vapor and molten lithium metal. Lithium metal was placed on a graphite boat and transferred into a quartz vessel containing selenium powder. The evacuated quartz tube was heated to 600 °C for three days.

Li₄P₂Se₆ was synthesized according to eqn (1). Red phosphorus (4 mmol, 123.9 mg, Alfa Aesar, 99.999%), selenium powder (10 mmol, 790.0 mg, ChemPur, 99.999%) and lithium selenide (4 mmol, 371.0 mg) were thoroughly ground and mixed in an agate mortar and afterwards filled into a graphitized quartz ampule (8 mm outer diameter). The ampule was sealed under vacuum (p = 1.8 × 10⁻² mbar) and heated to 650 °C in a tube furnace for seven days.

NMR spectroscopy

Solid-state NMR spectra were recorded on a Bruker Neo spectrometer, operating at a magnetic flux density B₀ = 14.1T, corresponding to resonance frequencies of 242.96 MHz and 114.47 MHz for ³¹P and ⁷⁷Se, respectively.

Chemical shift values are reported on a deshielding scale. The ¹H resonance of 1% Si(CH₃)₄ in CDCl₃ served as an external secondary reference using the Ξ values for ³¹P and ⁷⁷Se as reported by the IUPAC.¹⁴

Magic angle spinning (MAS) was performed using a commercial 3.2 mm triple resonance Bruker probehead within zirconia rotors (3.2 mm, ZrO₂). Rotors were packed under an argon atmosphere in a glove box. All spectra were recorded at a sample spinning frequency of ν_rot = 20 kHz, except for the ⁷⁷Se and ³¹P MAS spectra used for the determination of chemical shift tensors, which were recorded at ν_rot = 15 and 5 kHz, respectively.

For the ³¹P and ⁷⁷Se MAS NMR measurements, repetition delays of 600 s and 820 s were used respectively.

For the homonuclear ³¹P–³¹P 2D zero-quantum (ZQ) exchange spectroscopy measurement, the super-cycled R6₆² sequence¹⁵ was used with a repetition delay of 32 s, 72 R-elements, an increment of 50 μs and rotor-synchronized data sampling of the indirect dimension accumulating 8 transients per FID. The R-element consisted of a 90° pulse and a 270° pulse.

X-ray powder diffraction

X-ray powder diffraction (XRD) patterns were recorded using a STOE STADI P powder diffractometer (Stoe, Darmstadt, Germany) with Ge(111)–monochromated CuK_α1 radiation (λ = 1.54056 Å) in the Debye–Scherrer geometry (capillary inner diameter: 0.28 mm). A 2Θ range between 5° and 90° with a step size of 0.01° was measured, using a counting time of 90 s per step.

Rietveld refinements were carried out using TOPAS academic V7.¹⁶ The structure model with the space group Pnma (no. 62) was used as the starting model. In the given order, the following parameters were refined: scale factor and background coefficients using a Chebyshev function with 12 free parameters, the peak shape using the fundamental parameter approach, the zero-shift error, the lattice constants, the atom positions and the isotropic atomic displacement parameters. Fit indicators of R_wp, R_exp, and GOF were used to assess the quality of the refined structural models.¹⁷

Quantum-chemical calculations

Density functional theory (DFT) calculations and the computation of NMR parameters have been performed using the Quantum ESPRESSO¹⁸ package and GBRV pseudopotentials¹⁹ with the Perdew–Burke–Ernzerhof (PBE) functional as the exchange–correlation functional. A 7 × 7 × 7 Monkhorst–Pack grid of k points was used. The force convergence threshold was set to 10⁻¹⁰ Ry per a.u., the convergence threshold for self-consistency was set to 10⁻¹⁴ Ry and the total energy convergence threshold was set to 10⁻¹³ Ry. The kinetic energy cutoff for wavefunctions was set to 100 Ry. The calculated results were converted via XCrysDen ²⁰ and VESTA 3 ²¹ into a SHELX file in order to analyze the symmetry of the calculated structure with KPLOT.²²

Structure prediction was done with USPEX version 9.4.4, including predefined molecular building elements.^23–26 The formula Li₄P₂Se₆ was decomposed into Li and P₂Se₆ units. Calculations with different numbers of formula units Z ∈ {1, 2, 3, 4} and different population sizes N ∈ {20, 30, 50, 100} were started individually assuming a temperature of 0 K. For each individual structure suggestion, three consecutive calculations were carried out with increasingly stricter energy convergence criteria (0.001, 0.0001 and 0.00001 a.u.) and a finer k-point density (0.14, 0.10 and 0.04 2π per Å), respectively. Quantum chemical calculations within USPEX were carried out with the plane-wave self-consistent field package from Quantum ESPRESSO (version 6.7)¹⁸ using Gaussian smearing and a degauss parameter of 0.02 Ry. Calculations used PBE-pseudopotentials from the GBRV library.¹⁹

Electrochemical impedance spectroscopy

The Li₄P₂Se₆ powder sample was thoroughly ground in an agate mortar and pressed into a cylindrical pellet (diameter d = 13 mm) applying a load of 4 tons (equivalent to p = 295 MPa) for 20 minutes. The thickness of the pellet was measured with an electronic micrometer screw gauge. Afterwards, the pellet was brought into contact with silver foil (thickness: 0.1 mm, R. Götze GmbH & Co. KG, silver content 935) by pressing silver sheets (diameter d = 13 mm) onto both sides of the pellet (load of 4 tons, equivalent to p = 295 MPa for 20 minutes). All preparatory steps were carried out inside a glove box (MBraun, Garching, Germany) under an argon atmosphere.

Electrochemical impedance spectroscopy (EIS) measurements were recorded using a NEISYS electrochemical impedance analyzer (Novocontrol Technologies, Montabaur, Germany) in a home-built cell, which was calibrated before actual measurements based on short/load calibration standards with a 100 Ω resistor as the load. The impedance measurements were recorded in potentiostatic mode, with an amplitude of 7.1 mV_rms, in a frequency range from 50 mHz to 1 MHz. The temperature was controlled using a variable temperature and flow controller (NMR Service GmbH) with a constant nitrogen gas flow. Each temperature point was held for 20 min to ensure thermal equilibrium with an accuracy of ±0.1 K throughout the EIS measurement. The data analysis was performed using a home-written Python script.²⁷

3. Results and discussion

In the initial phase of the study, Li₄P₂Se₆ was characterized by solid-state NMR, which provided key information about the molecular fragments. These fragments were subsequently used as constraints for the structure solution. The structure was solved by comparison with theoretically predicted structures. Finally the refined structure of Li₄P₂Se₆ is related to other structures within the Li–P–Se phase triangle. Finally, the ionic conductivity of Li₄P₂Se₆ was investigated through variable temperature impedance spectroscopy.

Synthesis and NMR of Li₄P₂Se₆

The synthesis of Li₄P₂Se₆ was straightforward, using P, Se and Li₂Se as starting agents, and resulted in a fine brown-red powder. Li₄P₂Se₆ was analyzed by solid-state ³¹P and ⁷⁷Se MAS NMR in order to obtain information on its structural composition and confirm the suggested building units.⁶ The ³¹P MAS NMR spectrum of Li₄P₂Se₆ (Fig. 1) displays two sharp signals at δ_iso = 46.6 ppm and 55.1 ppm with a peak area ratio of exactly 1 : 1. Additionally, both signals exhibit splittings due to homonuclear J-coupling. The isotropic chemical shift values are consistent with the previously reported values of 46.5 ppm and 55.1 ppm.⁶ A slower spinning spectrum was recorded to determine the ³¹P chemical shift tensors and the homonuclear J-coupling constants (Fig. 2). The simulation was performed using a ³¹P–³¹P direct dipole–dipole coupling constant D = −2038 Hz, corresponding to a P–P distance of 2.13 Å determined from the crystal structure.

Fig. 1 Experimental ³¹P MAS NMR spectrum of Li₄P₂Se₆ recorded at a spinning frequency ν_rot = 20 kHz. Isotropic chemical shifts δ_iso are marked with asterisks. The peak corresponding to a minor impurity of Li₃PO₄ is marked with a cross.

Fig. 2 Experimental and simulated ³¹P magic angle spinning NMR spectrum of Li₄P₂Se₆, simulated spectrum and difference. The simulated spectrum is based on a model including chemical shift tensors and homonuclear direct and indirect dipolar ³¹P–³¹P interactions using a two-spin model. The spectrum was recorded at a spinning frequency of 5 kHz. Isotropic spinning side-bands are marked with asterisks. The peak corresponding to a minor impurity of Li₃PO₄ is marked with a cross.

The ³¹P chemical shift tensors were assumed to be collinear with the internuclear vector of the dipolar interaction. This assumption is based on the idealized D_3d molecular geometry of [P₂Se₆]⁴⁻, in which the principal C₃ axis is aligned along the P–P bond. According to pseudosymmetry arguments, the presence of a C_n axis (with n ≥ 3) suggests that the principal axis of the chemical shift tensor should be parallel to the rotation axis. Although the molecular symmetry in the crystal structure deviates from the idealized case, a significant misalignment between the principal axis of the chemical shift tensor and the dipolar vector is not expected.

The results are sensitive to the sign of the J-coupling, as seen in the lineshape and relative intensities of the spinning sidebands. Simulations suggest that a negative ¹J(³¹P,³¹P) coupling provides a slightly better match to the experimental spectrum, consistent with the literature values for similar systems.²⁸ A comparison of simulations using positive and negative J-coupling values is provided in Fig. S1.† The fitted parameters are summarized in Table 1.

Table 1 Fitted chemical shift tensor parameters and eigenvalues for ³¹P and ⁷⁷Se in Li₄P₂Se₆, determined according to the Haeberlen–Mehring–Spiess convention as implemented in SIMPSON.^29,30 The parameters δ_iso, δ_aniso, and η were determined from simulations as shown in Fig. 2 (³¹P) and Fig. 4 (⁷⁷Se). Additional simulation details are provided in the ESI.† Errors of the parameters δ_iso(³¹P), δ_aniso(³¹P), η(³¹P), δ₁₁(³¹P), δ₂₂(³¹P), δ₃₃(³¹P), δ_iso(³¹P), δ_aniso(⁷⁷Se), η(⁷⁷Se), δ₁₁(⁷⁷Se), δ₂₂(⁷⁷Se), δ₃₃(⁷⁷Se), ¹J_iso(³¹P,⁷⁷Se) and ¹J_iso(³¹P,³¹P) are estimated to be 0.1, 4, 0.1, 4, 4, 4, 1, 10, 0.1, 10, 10, 10, 10, and 15, respectively

	^31P NMR		^77Se NMR
	Spin 1	Spin 2	Spin 1	Spin 2	Spin 3	Spin 4
δiso/ppm	46.4	55.2	83	103	166	105
δaniso/ppm	−109	–109	−436	−446	557	478
η	0.3	0	0.7	1.0	0.7	1.0
δ11/ppm	115	110	455	550	724	585
δ22/ppm	85	110	149	104	80	100
δ33/ppm	−62	–54	−352	−342	−305	−367
^1Jiso(^31P,^77Se) coupling/Hz	—	—	488	532	531	496
^1Jiso(^31P,^31P) coupling/Hz	−227	−227	—	—	—	—

---

The homonuclear ³¹P–³¹P 2D ZQ exchange MAS NMR spectrum of Li₄P₂Se₆ (Fig. 3) exhibits two cross peaks for the ³¹P doublets, consistent with a J-coupled two-spin system originating from the same crystalline phase. This indicates that the crystal structure of Li₄P₂Se₆ must contain at least two crystallographic orbits for P with equal multiplicities. The presence of homonuclear J-coupling indicates a P–P bond, suggesting the presence of a [P₂Se₆]⁴⁻ unit featuring two different P sites. These findings confirm the results of the previous study.⁶

Fig. 3 Homonuclear ³¹P–³¹P 2D zero-quantum exchange MAS NMR spectrum of Li₄P₂Se₆ at ν_rot = 20 kHz. The 1D projection from the one-pulse experiment is shown on top of the 2D spectrum. The diagonal line (autocorrelation diagonal) refers to the hypothetical peak positions of two isochronous spins.

The ⁷⁷Se MAS NMR spectrum reveals four ⁷⁷Se resonances at δ_iso = 83.8, 103.8, 105.8 and 166.0 ppm, with intensity ratios of 1 : 1.1 : 1.9 : 1.8 (Fig. 4). These observations provide valuable constraints for the minimum number and multiplicities of the Se sites in the crystal structure of Li₄P₂Se₆. The chemical shift parameters for all signals of the experimental ⁷⁷Se MAS NMR spectrum are listed in Table 1.

Fig. 4 Experimental and fitted ⁷⁷Se MAS NMR spectrum of Li₄P₂Se₆. The simulated spectrum accounts for ⁷⁷Se chemical shift tensors, heteronuclear ³¹P–⁷⁷Se direct dipolar interactions, and ¹J_iso(³¹P,⁷⁷Se) couplings using a two-spin model. The spectra of the individual spins are shown separately for clarity. Isotropic chemical shifts are marked with asterisks. A small unidentified impurity at approximately 800 ppm is marked with a cross. Recorded at a spinning frequency ν_rot = 15 kHz.

Structure solution

NMR data indicate the presence of two distinct P-sites within the [P₂Se₆]⁴⁻ unit, suggesting a relatively small unit cell. This made the compound an ideal candidate for structure determination using structure prediction employing a genetic algorithm in combination with quantum-chemical calculations.

The algorithm was provided with the composition Li₄P₂Se₆ and the number of formula units Z ∈ {1, 2, 3, 4} as the input. During the search, the candidate structures proposed by the algorithm were optimized using density functional theory (DFT) under periodic boundary conditions. The search space was drastically reduced by using the fragment [P₂Se₆] and Li as the basic building elements.

The X-ray diffraction patterns for the 20 best structures generated for each value of Z were obtained and visually compared to the experimental powder X-ray diffraction pattern. Among these, only the orthorhombic structure described below provided a good match with the experimental pattern. The structure that matches the experimental data is among the three structures, which are lowest in energy (Table 4). It is less than 2 meV per atom higher in energy than the most stable structure.

Additionally, the structural model was required to reflect the number of sites observed by solid-state NMR; at least two and four sites should be present for P and Se, respectively. The proposed structural model satisfied this condition, with two P-sites both with Wyckoff labels 4c and four Se-sites with Wyckoff labels 8d, 8d, 4c and 4c, consistent with the experimental NMR data.

The main reflections from the powder XRD pattern and the ³¹P NMR peaks reported in a previous study are in agreement with the proposed structure in this work.⁶ This indicates that the presented structure corresponds to the same phase, although the previously suggested unit cell was different.

The structural parameters were refined using the Rietveld method (Fig. 5). The refined parameters are listed in Table 2. The refined structure is novel and has not been classified previously in the ICSD database. It consists of hexaselenidohypodiphosphate units, with Li atoms occupying octahedrally coordinated positions (Fig. 6). The bond lengths and angles are within typical ranges (Table 3).

Fig. 5 Rietveld refinement of Li₄P₂Se₆: experimental powder pattern (black points), fitted curve (red line), difference (blue line), and Bragg reflections of Li₄P₂Se₆ (green lines) and of a small impurity of Li₂Se (dark green lines).

Fig. 6 The structure of Li₄P₂Se₆. Lithium (gray), phosphorus (black) and selenium (red) positions are shown. The octahedral Li polyhedra are shown in purple shades.

Table 2 Comparison of the crystallographic unit cell information (cell dimensions, crystal system, and space group) of the predicted structural model and the structural model of Li₄P₂Se₆ at room temperature obtained from Rietveld refinement of X-ray powder data (λ = 1.54056 Å), as well as the crystallographic data (atomic coordinates, B_eq, and occupancy) of the refined model. Fit residuals (R_wp, R_exp, and GOF): 4.21, 2.20, and 1.92

			Predicted structural model		Refined structural model
a/Å			14.081		13.8707(1)
b/Å			11.256		11.2115(1)
c/Å			6.494		6.45445(7)
Crystal system			Orthorhombic		Orthorhombic
Space group			Pnma (No. 62)		Pnma (No. 62)
Z			4		4
Atom	Wyckoff site	x	y	z	Occ.	Beq/Å^2
Li1	8d	0.625(2)	0.397(2)	0.0932(49)	1	1
Li2	8d	0.127(2)	0.893(2)	0.945(4)	1	1
P1	4c	0.7867(4)	3/4	0.0912(12)	1	0.72(9)
P2	4c	0.9401(4)	3/4	0.1101(10)	1	0.72(9)
Se1	8d	0.9829(1)	0.9137(2)	0.2716(3)	1	0.29(2)
Se2	8d	0.7490(1)	0.5862(4)	0.9243(4)	1	0.29(2)
Se3	4c	0.9957(2)	3/4	0.7882(4)	1	0.29(2)
Se4	4c	0.7387(2)	3/4	0.4197(7)	1	0.29(2)

---

Table 3 Bond-lengths and angles found in the refined structural model of Li₄P₂Se₆

Atom 1	Atom 2	Atom3	Bond length/Å	Bond angle/°
Se1	P2	—	2.193(4)	—
Se2	P1	—	2.192(6)	—
Se3	P2	—	2.216(7)	—
Se4	P1	—	2.222(9)	—
P1	P2	—	2.131(8)	—
Se4	P1	P2	—	104.2(4)
P2	P1	Se2	—	105.4(2)
Se1	P2	P1	—	107.3(2)
P1	P2	Se3	—	107.1(3)
Se4	P1	Se2	—	113.4(2)
P2	P1	Se2	—	105.4(2)
Se1	P2	Se3	—	110.6(2)
P1	P2	Se1	—	107.3(2)
Se2	P1	Se2	—	113.8(3)
Se1	P2	Se1	—	113.7(3)

---

Comparison and phase diagram

The three most likely candidate structures suggested by the genetic algorithm for Li₄P₂Se₆ differed by less than 200 J per formula unit or 2 meV per atom. Interestingly, the energetically most favorable structure is isotypic to that of Li₄P₂S₆,^7,8 but it does not correspond to the experimentally determined structure for Li₄P₂Se₆. Such small differences in energy are challenging to accurately resolve using quantum chemical calculations at the density–functional theory level. Nevertheless, it is promising that the experimentally favored structure was included among the top candidates suggested by the genetic algorithm.

The new structure enables the construction of a preliminary phase diagram (Fig. 7), which considers only ordered structures, whether experimentally confirmed or theoretically predicted in the absence of experimental data (Table 4). The set of thermodynamically stable structures, those that form the convex hull based on their quantum-chemical formation energy at 0 K, are separated from metastable phases. Phases above the convex hull decompose into adjacent phases on the diagram and are thus not stable.

Fig. 7 Phase diagram Li–P–Se based on calculated energies (Table 4).

Table 4 Computed energies per formula unit E_f of ordered crystalline compounds in the composition space Li–P–Se and their energy ΔE_ah above the convex hull. The structures of the predicted structures are available in the ESI† as cif files

Phase	Pearson symbol	Source	Ef (kJ mol^−1)	ΔEah/atom/meV	Source	Collection code
Pblack	oS8	XRD	0.0	0	ICSD	25253
Se8	mP32	XRD	0.0	0	ICSD	2718
Li	cI2	XRD	0.0	0	ICSD	44759
Li2Se	cF12	XRD	−370.9	0	ICSD	60433
Li3P	hP8	XRD	−280.3	0	ICSD	26880
Li3P7	oP40	XRD	−383.8	0	ICSD	60774
LiP	mP16	XRD	−107.7	0	ICSD	100465
LiP5	oP24	XRD	−136.0	10	ICSD	88710
LiP7	tI128	Unclear	−152.4	0	ICSD	23621
P2Se5	mP28	XRD	−20.3	23	ICSD	74546
P4Se3	oP112	XRD	−48.6	7	ICSD	26483
P4Se4	mP32	XRD	−71.1	0	ICSD	74878
P4Se5	oP36	XRD	−65.5	6	ICSD	16140
Li3PSe4	cP8	Predicted	−626.7	0	OQMD	7772844
Li7PSe6	cP56	Ag7PSe6	−1309.8	43	ICSD	54055
Li4P2Se6	oP48	XRD, predicted	−874.0	0.15	This work
Li4P2Se6	hP12	Predicted	−866.2	6.92	This work
Li4P2Se6	oP24	Predicted	−873.8	0.30	This work
Li4P2Se6	mP24	Predicted	−873.8	0.33	This work
Li4P2Se6	hP36	Predicted	−874.2	0	This work
Li4P2Se6	oP48′	Predicted	−874.1	0.06	This work

---

At a first glance, it is surprising that certain compounds, such as LiP₅, which can be synthesized at high temperatures, are not predicted to be stable according to the calculated energies. This is unexpected, as it would be assumed that thermodynamic minima are easily achieved for such phases. Upon closer inspection, the small excess energies above the convex hull energy stand out (Table 4). These energy differences are minor when compared to the energy contributions from vibrational or configurational entropy.

For instance, at the synthesis temperature of 600 °C, the contribution from configurational entropy for Li₇PSe₆ amounts to approximately 5 kJ mol⁻¹. This perspective suggests that using energy above the convex hull as a strict criterion for determining new stable structures needs to be applied with caution.³¹ Thus, structure prediction is not only limited by the limited search space that can be explored in practice but also by uncertainties in free energy. Kinetic control is still another point which may be difficult to consider in computational chemistry but may be less relevant for synthesis at high temperatures.

Ionic conductivity

Variable temperature impedance spectroscopy measurements of Li₄P₂Se₆ were performed to determine its macroscopic lithium ion conductivity. A representative Nyquist diagram (Fig. 8) shows a single half circle assumed as the bulk region and a double-layer formation at the sample/electrode interface.

Fig. 8 Nyquist (imaginary part of the impedance against real part) plot of impedance measurements of Li₄P₂Se₆. Pellet (diameter = 13 mm, thickness = 1.3 mm) at a temperature of 393 K. The gray cycles represent the experimental data while red stars represent the result of a fit using the shown model. R_b is the bulk resistance while Q_b and Q_el represent the constant phase elements for the bulk and the electrode polarization, respectively.

The measured impedance spectra were fitted with an equivalent circuit (Fig. 8, inset) consisting of a resistor R_b and a constant phase element Q_b, depicting the conductive and non-ideal capacitive characteristics of the bulk region, respectively. The Q_el depicts the non-ideal double layer formation at the sample/electrode interface. The impedance of a constant-phase element is given by Z_Q = 1/Q(iω)^α, with pre-factor Q and α < 1 exponent.

The effective capacitance C_b of the bulk region is calculated from the expression C_b = (QR^(α−1))^1/α,³² which is 55(9) pF at 393 K (Fig. 8).

The conductivity σ of the bulk region was determined from the fitted bulk resistance values R_b at corresponding temperatures using the expression , where d = 0.13 cm and A = 1.33 cm² denote the sample thickness and the electrode area, respectively.

The activation energy of the macroscopic ionic conductivity was obtained using the given expression σ·T = a₀·exp(−E_A/(k_BT)), where a₀ is the pre-exponential factor, E_A is the activation energy, k_B is Boltzmann's constant, and T is the temperature (Fig. 9).

Fig. 9 Plot of the product of electrical conductivity σ and temperature T versus inverse temperature for Li₄P₂Se₆. The results were obtained from impedance spectra in the temperature range from 333 K to 393 K. Linear fits are indicated by dashed lines. The data points are indicated by grey circles of the bulk conductivities.

The activation energy and the corresponding pre-exponential factors are shown in Table 5. The difference of 1% between the conductivities and activation energies of the repeated measurements (ESI Fig. S2†) confirms the reproducibility of the results.

Table 5 The macroscopic conductivity at room temperature σ_{298 K}, pre-exponential factor a₀, and activation energy E_A of Li₄P₂Se₆ as compared to those of Li₄P₂S₆ ⁹

	Li4P2S6	Li4P2Se6
a These values have been recalculated from ref. 9 using σ·T = a0·exp(−EA/(kBT)).
σ298 K/(S cm^−1)	1.28(12) × 10^−12 ^a	1.85(44) × 10^−12
a0/(S K cm^−1)	3.51(2) × 10^5 ^a	6.98(1) × 10^3
EA/eV	0.88(2)^a	0.78(1)

---

The contribution of the grain boundary could not be distinctly separated from that of the grain itself. Consequently, the measured activation energy of 0.78(1) eV for the ionic motion of Li⁺ in Li₄P₂Se₆ likely represents an upper limit for the bulk activation energy. This value remains slightly lower than the activation energy determined for Li₄P₂S₆, which is 0.88(2) eV.⁹ This reduction in activation energy is expected because Se²⁻ has a higher size and polarizability than S²⁻, and the larger ionic radii of Se²⁻ weaken the attractive force between the ligand and the Li-cation, which enhances the ion mobility.³³ Moreover, the lattice softness of the anion sub-lattice has been shown to decrease the activation energy barrier of the cation transport in the crystal system.³⁴ To enhance the ion transport mechanism, vacancies and interstitial lattice sites could be created by aliovalent doping in anionic as well as cationic lattice sites by substituting Se and Li atoms, but it is rather unlikely that conductivities in the mS cm⁻¹ range can be achieved this way.⁹

In the context of the paddle-wheel mechanism,³⁵ it is interesting to determine the activation energy for the rotation of a PSe₃ fragment around the P–P bond. The calculated activation energy for the rotational motion in Li₄P₂Se₆ (Fig. 6) was estimated by nudged elastic band calculations to be 1.29 eV. These results imply that the Li-ion motion is decoupled from the rotational motion of the complex anion. Hence the movement of ions can be treated independently from the rotation of PSe₃ fragments.

4. Conclusions

The structure of Li₄P₂Se₆ was successfully determined using a combination of powder X-ray diffraction data, quantum-chemical structure prediction, and solid-state NMR. This combined approach proved particularly valuable in addressing challenges associated with the low scattering contrast between heavy atoms like Se and light atoms like Li. In addition to structural determination, the ionic conductivity of Li₄P₂Se₆ and the ³¹P and ⁷⁷Se NMR tensors were characterized, which can serve as reference data for future dynamic studies. Exploration of the Li–P–Se phase diagram indicates the existence of phases above the calculated convex hull at a temperature of 0 K, which indicates that configurational and vibrational entropies play an important role in the prediction of phases in similar material systems.

Author contributions

Conceptualization: JSADG; investigation: SN, NM, SCA, and HJCS; writing – original draft: SN, NM, SCA, HJCS, and JSADG; writing – review and editing: JSADG; visualization: SN, NM, SCA, and HJCS; funding acquisition: JSADG; and supervision: JSADG.

Conflicts of interest

There are no conflicts to declare.

Data availability

Crystal structure: deposited at CCDC, deposition number 2419766.† Predicted structures: see the ESI.† All experimental data will be made available if required.

References

1. Z. Zhang, Y. Shao, B. Lotsch, Y.-S. Hu, H. Li, J. Janek, L. F. Nazar, C.-W. Nan, J. Maier, M. Armand and L. Chen, New horizons for inorganic solid state ion conductors, Energy Environ. Sci., 2018, 11, 1945–1976,  10.1039/C8EE01053F.
2. S. T. Kong, Ö. Gün, B. Koch, H. J. Deiseroth, H. Eckert and C. Reiner, Structural Characterisation of the Li Argyrodites Li₇ PS₆ and Li₇ PSe₆ and their Solid Solutions: Quantification of Site Preferences by MAS–NMR Spectroscopy, Chem. – Eur. J., 2010, 16, 5138–5147,  DOI:10.1002/chem.200903023.
3. J. Hartel, A. Banik, J. M. Gerdes, B. Wankmiller, B. Helm, C. Li, M. A. Kraft, M. R. Hansen and W. G. Zeier, Understanding Lithium-Ion Transport in Selenophosphate-Based Lithium Argyrodites and Their Limitations in Solid-State Batteries, Chem. Mater., 2023, 35, 4798–4809,  DOI:10.1021/acs.chemmater.3c00658.
4. B. Koch, S. T. Kong, Ö. Gün, H.-J. Deiseroth and H. Eckert, Site preferences and ion dynamics in lithium chalcohalide solid solutions with argyrodite structure: I. A multinuclear solid state NMR study of the system Li₆ PS_5−xSe_xI and of Li₆AsS₅I, Z. Phys. Chem., 2022, 236, 853–874,  DOI:10.1515/zpch-2021-3135.
5. B. Koch, S. T. Kong, Ö. Gün, H.-J. Deiseroth and H. Eckert, Site preferences and ion dynamics in lithium chalcohalide solid solutions with argyrodite structure: II. Multinuclear solid state NMR of the systems Li₆PS_5−xSe_x Cl and Li₆PS_5−xSe_xBr, Z. Phys. Chem., 2022, 236, 875–898,  DOI:10.1515/zpch-2021-3139.
6. R. H. P. Francisco, T. Tepe and H. Eckert, A Study of the System Li–P–Se, J. Solid State Chem., 1993, 107, 452–459,  DOI:10.1006/jssc.1993.1369.
7. S. Neuberger, S. P. Culver, H. Eckert, W. G. Zeier and J. Schmedt auf der Günne, Refinement of the crystal structure of Li₄P₂S₆ using NMR crystallography, Dalton Trans., 2018, 47, 11691–11695,  10.1039/C8DT02619J.
8. H. B. Yahia, K. Motohashi, S. Mori, A. Sakuda and A. Hayashi, Twinned single crystal structure of Li₄P₂S₆, Z. Kristallogr. – Cryst. Mater., 2023, 238, 209–216,  DOI:10.1515/zkri-2023-0013.
9. S. Neuberger, N. Mathew, S. C. Adediwura and J. Schmedt auf der Günne, Influence of Mg on the Li ion mobility in Li_4–2xMg_x,P₂S₆, Dalton Trans., 2023, 52, 16894–16902,  10.1039/D3DT02624H.
10. H.-C. Wang, S. Botti and M. A. L. Marques, Predicting stable crystalline compounds using chemical similarity, npj Comput. Mater., 2021, 7, 12,  DOI:10.1038/s41524-020-00481-6.
11. J. Pannetier, J. Bassas-Alsina, J. Rodriguez-Carvajal and V. Caignaert, Prediction of crystal structures from crystal chemistry rules by simulated annealing, Nature, 1990, 346, 343–345,  DOI:10.1038/346343a0.
12. J. C. Schön and M. Jansen, First Step Towards Planning of Syntheses in Solid–State Chemistry: Determination of Promising Structure Candidates by Global Optimization, Angew. Chem., Int. Ed. Engl., 1996, 35, 1286–1304,  DOI:10.1002/anie.199612861.
13. S. M. Woodley, P. D. Battle, J. D. Gale, C. Richard and A. Catlow, The prediction of inorganic crystal structures using a genetic algorithm and energy minimisation, Phys. Chem. Chem. Phys., 1999, 1, 2535–2542,  10.1039/a901227c.
14. R. K. Harris, E. D. Becker, d. M. S. M. Cabral, P. Granger, R. E. Hoffman and K. W. Zilm, Further conventions for NMR shielding and chemical shifts (IUPAC Recommendations 2008), Pure Appl. Chem., 2008, 80, 59–84,  DOI:10.1351/pac200880010059.
15. A. Brinkmann, J. Schmedt auf der Günne and M. H. Levitt, Homonuclear zero-quantum recoupling in fast magic-angle spinning nuclear magnetic resonance, J. Magn. Reson., 2002, 156, 79–96,  DOI:10.1006/jmre.2002.2525.
16. A. A. Coelho, TOPAS Academic V7, 2020.
17. The Rietveld Method, ed. R. A. Young, Oxford University Press, Oxford, New York, 1995.
18. P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo, A. D. Corso, S. Fabris, G. Fratesi, S. de Gironcoli, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. Seitsonen, A. Smogunov, P. Umari and R. M. Wentzcovitch, Quantum ESPRESSO: a modular and open-source software project for quantum simulations of materials, J. Phys.: Condens. Matter, 2009, 21, 395502,  DOI:10.1088/0953-8984/21/39/395502.
19. K. F. Garrity, J. W. Bennett, K. M. Rabe and D. Vanderbilt, Pseudopotentials for high-throughput DFT calculations, Comput. Mater. Sci., 2014, 81, 446–452,  DOI:10.1016/j.commatsci.2013.08.053.
20. A. Kokalj, XCrySDen—a new program for displaying crystalline structures and electron densities, J. Mol. Graphics Modell., 1999, 17, 176–179,  DOI:10.1016/S1093-3263(99)00028-5.
21. K. Momma and F. Izumi, VESTA 3 for three-dimensional visualization of crystal, volumetric and morphology data, J. Appl. Crystallogr., 2011, 44, 1272–1276,  DOI:10.1107/S0021889811038970.
22. R. Hundt, C. J. Schön and D. Zagorac, KPLOT A Program for Plotting and Analysing Crystal Structures, Version 9, Technicum Scientific Publishing, Stuttgart, Germany, 2016. https://opentechnicum.com/?mbt_book=kplot-a-program-for-plotting-and-analysing-crystal-structures-version-9 (accessed February 5, 2018).
23. A. R. Oganov, A. O. Lyakhov and M. Valle, How Evolutionary Crystal Structure Prediction Works—and Why, Acc. Chem. Res., 2011, 44, 227–237,  DOI:10.1021/ar1001318.
24. A. R. Oganov and C. W. Glass, Crystal structure prediction using ab initio evolutionary techniques: Principles and applications, J. Chem. Phys., 2006, 124, 244704,  DOI:10.1063/1.2210932.
25. A. O. Lyakhov, A. R. Oganov, H. T. Stokes and Q. Zhu, New developments in evolutionary structure prediction algorithm USPEX, Comput. Phys. Commun., 2013, 184, 1172–1182,  DOI:10.1016/j.cpc.2012.12.009.
26. Q. Zhu, A. R. Oganov, C. W. Glass and H. T. Stokes, Constrained evolutionary algorithm for structure prediction of molecular crystals: methodology and applications, Acta Crystallogr., Sect. B: Struct. Sci., 2012, 68, 215–226,  DOI:10.1107/S0108768112017466.
27. S. C. Adediwura, N. Mathew and J. Schmedt auf der Günne, Combining NMR and impedance spectroscopy in situ to study the dynamics of solid ion conductors, J. Mater. Chem. A, 2024, 12, 15847–15857,  10.1039/D3TA06237F.
28. M. A. M. Forgeron, M. Gee and R. E. Wasylishen, A Theoretical Investigation of One-Bond Phosphorus–Phosphorus Indirect Nuclear Spin–Spin Coupling Tensors,¹J (³¹P,³¹P), Using Density Functional Theory, J. Phys. Chem. A, 2004, 108, 4895–4908,  DOI:10.1021/jp040034j.
29. M. Bak, J. T. Rasmussen and N. C. Nielsen, SIMPSON: A General Simulation Program for Solid-State NMR Spectroscopy, J. Magn. Reson., 2000, 147, 296–330,  DOI:10.1006/jmre.2000.2179.
30. Z. Tošner, R. Andersen, B. Stevensson, M. Edén, N. C. Nielsen and T. Vosegaard, Computer-intensive simulation of solid-state NMR experiments using SIMPSON, J. Magn. Reson., 2014, 246, 79–93,  DOI:10.1016/j.jmr.2014.07.002.
31. Y. Gao, A. M. Nolan, P. Du, Y. Wu, C. Yang, Q. Chen, Y. Mo and S.-H. Bo, Classical and Emerging Characterization Techniques for Investigation of Ion Transport Mechanisms in Crystalline Fast Ionic Conductors, Chem. Rev., 2020, 120, 5954–6008,  DOI:10.1021/acs.chemrev.9b00747.
32. G. J. Brug, A. L. G. Van Den Eeden, M. Sluyters-Rehbach and J. H. Sluyters, The analysis of electrode impedances complicated by the presence of a constant phase element, J. Electroanal. Chem. Interfacial Electrochem., 1984, 176, 275–295,  DOI:10.1016/S0022-0728(84)80324-1.
33. J. C. Bachman, S. Muy, A. Grimaud, H.-H. Chang, N. Pour, S. F. Lux, O. Paschos, F. Maglia, S. Lupart, P. Lamp, L. Giordano and Y. Shao-Horn, Inorganic Solid-State Electrolytes for Lithium Batteries: Mechanisms and Properties Governing Ion Conduction, Chem. Rev., 2016, 116, 140–162,  DOI:10.1021/acs.chemrev.5b00563.
34. M. A. Kraft, S. P. Culver, M. Calderon, F. Böcher, T. Krauskopf, A. Senyshyn, C. Dietrich, A. Zevalkink, J. Janek and W. G. Zeier, Influence of Lattice Polarizability on the Ionic Conductivity in the Lithium Superionic Argyrodites Li₆ PS₅ X (X = Cl, Br, I), J. Am. Chem. Soc., 2017, 139, 10909–10918,  DOI:10.1021/jacs.7b06327.
35. Z. Zhang and L. F. Nazar, Exploiting the paddle-wheel mechanism for the design of fast ion conductors, Nat. Rev. Mater., 2022, 7, 389–405,  DOI:10.1038/s41578-021-00401-0.

---

Footnote

† Electronic supplementary information (ESI) available. CCDC 2419766. For ESI and crystallographic data in CIF or other electronic format see DOI: https://doi.org/10.1039/d5dt00227c

---