Computational investigation of the Mg-ion conductivity and phase stability of MgZr4(PO4)6

Solid electrolyte materials exhibiting high Mg-ion conductivity are required to develop Mg-ion batteries. In this study, we focused on a Mg-ion-conducting solid phosphate based electrolyte, MgZr4(PO4)6 (MZP), and evaluated the ionic conductivity of NASICON-type and β-iron sulfate-type MgZr4(PO4)6 structures via density functional theory calculations. The calculations suggest that the migration energy of Mg is 0.63 eV for the NASICON-type structure and 0.71 eV for the β-iron sulfate-type one, and the NASICON-type structure has higher ion conductivity. Although the NASICON-type MZP structure has not been experimentally realised, there is only an energy difference of 14 meV per atom with respect to that of the β-iron sulfate-type structure. Therefore, in order to develop a synthesis method for the NASICON-type structure, we investigated pressure- and temperature-dependent variations in the free energy of formation using density functional perturbation theory calculations. The results suggest that the formation of the NASICON-type structure is disfavoured under the 0–2000 K and 0–20 GPa conditions.


Introduction
Lithium ion batteries (LIBs) are widely used in portable electronic devices, such as mobile phones and laptop computers, owing to their high energy density and long cycle life. In recent years, the use of (hybrid) electric vehicles (EVs) has gradually increased to address the environmental and energy issues. 1 Therefore, LIBs are also attracting attention as secondary batteries for on-vehicle use. However, the following two major technical drawbacks of current LIBs need to be overcome: (i) insufficient driving range of the vehicle compared to the conventional petro-powered vehicle due to the small energy density and (ii) safety concerns arising from inammable organic electrolytes. In this respect, developing all solid-state Mg-ion batteries (MIBs) is believed to be one of the ultimate solutions to resolve the above two issues. 2 The Mg ion is divalent and can carry two electrons, so that the energy density is expected to double. All solid-state battery technology also contributes to the fabrication of a bipolar-type battery, thus enabling increased energy density, and replacement of organic electrolyte by an inammable solid electrolyte can solve the safety concerns.
However, the generally poor Mg diffusivity in the solid state prevents the realisation of MIB-based devices, since divalent Mg ions strongly interact with counter anions. Therefore, it is necessary to understand the factors affecting the Mg-ion conductivity and nd fast Mg-ion conductors to develop future all solid-state MIBs.
As solid electrolyte materials, oxide-based Na super ionic conductor (NASICON)-type structures are attractive for LIBs because NASICON-type ion conductors show high ion conductivity in both Na 3-9 and Li 10-17 systems. So far, one of the best Mgion conductors is rhombohedral NASICON-type (Mg 0.1 Hf 0.9 ) 4/ 3.8 Nb(PO 4 ) 3 ( Fig. 1(a and b)), in which the ionic conductivity of Mg 2+ is 2.1 Â 10 À6 S cm À1 at 573 K and its activation energy is 0.68 eV. 18 This compound consists of (Nb,Hf)O 6 octahedra and PO 4 tetrahedra alternately bridged at the two-coordinated edge O atoms, and Mg ions migrate through the interstitial sites of the corner-shared network. This local bonding rule leads to another type of structure with a different network topology, known as the monoclinic b-iron sulfate-type structure, and MgZr 4 (PO 4 ) 6 belonging to this structure type ( Fig. 1(c and d)) also exhibits relatively high Mg-ion conductivity. 19,20 An example of a material with this type of structure is Mg 0.7 (Zr 0.85 Nb 0.15 ) 4 (-PO 4 ) 6 . The ionic conductivity of this compound is 1.1 Â 10 À7 S cm À1 at 573 K and the activation energy is 0.92 eV, 21 and Mg-ion conduction of this compound at 350 C has been recently conrmed using non-blocking Mg metal electrodes. 22 Interestingly, both the NASICON-type and b-iron sulfate-type structural frameworks have a polymorphic relationship represented by A x B 4 (XO 4 ) 6 frameworks. In case of Na ion migration, the NASICON-type structure has been reported to show higher conductivity than that of the b-iron sulfate-type one for the composition of Na 3 Zr 2 Si 2 PO 12 . 4 The lower activation energy for the Mg ion conduction in the NASICON-type (Mg 0.1 Hf 0.9 ) 4/ 3.8 Nb(PO 4 ) 3 structure than in the b-iron sulfate-type Mg 0.7 (-Zr 0.85 Nb 0.15 ) 4 P 6 O 24 structure has been attributed to the higher structural symmetry of the former structure. 15 NASICON-type MgZr 4 (PO 4 ) 6 has not been synthesised primarily because the stability of the NASICON-type structure is lower for small interstitial cations like Mg 2+ ions. 19 However, the relationship between the structure and ionic conductivity is not yet clear for the case of Mg ions, because there are no experimental ionic conductivity data of both structures with the same composition. In this study, we evaluate the ionic conductivity for both the NASICON-type and b-iron sulfate-type structures with the stoichiometric composition of MgZr 4 (PO 4 ) 6 (MZP) using rst-principles molecular dynamics (FPMD) study based on density functional theory (DFT). So far, MZP has been experimentally prepared only in the b-iron sulfate-type structure, and no report on the synthesis of the NASICON-type structure is known. We also investigated the phase stability of NASICON-type and b-iron sulfate-type structures under various temperature and pressure conditions.

Computational details
Structure relaxation parameters including the lattice parameters and internal coordinates were calculated by rst-principles DFT calculations. In detail, a Vienna ab initio simulation package (VASP) 23-26 based on DFT 27 with a projector augmentedwave (PAW) 28,29 method and a plane-wave basis set was used. Further, we used a generalized gradient approximation (GGA)type exchange-correlation functional developed by Perdew, Burke, and Ernzernhof and later modied for solid materials (PBEsol). 30 The cutoff energy for the plane-wave basis was set to 500 eV and k-point resolution was set to $1000 (i.e., N x , N y , and N z , the number of grids in k x -, k y -, and k z -directions of the reciprocal space were set to satisfy N x Â N y Â N z Â N atom z 1000, where N atom is the number of atoms in each unit cell).
The NASICON-type structure was constructed by replacing Li ions in the structure of LiZr 2 (PO 4 ) 3 taken from the inorganic crystal structure database (ICSD) (ID: 92250, original data reported by Catti et al.) 31 with Mg ions and vacancies, because there are no reports of the experimental synthesis. The crystal structure input for b-iron sulfate-type MZP was extracted from the ICSD (ID: 250452, original data reported by Gobechiya et al.). 32 Table 1 lists DFT-derived relaxed lattice parameters, cell volume, and calculated total electron energies for (a) NASICONtype and (b) b-iron sulfate-type MZP. In the NASICON-type structure, three Mg ions and vacancies reside at the 6b and/or 36f Wyckoff position (36f sites are split around the 6b site). The most stable Mg/vacancy arrangement was determined using a genetic algorithm (GA), as described in our previous paper. [33][34][35] Mg-ion occupation at 6b sites is the most stable aer GA optimisation, and the energy difference among the GAderived structures (2400 samples) is within 60 meV per atom.
In contrast, Mg ions and vacancies occupied the 4e sites of the Wyckoff position at 1 : 1 molar ratio in the b-iron sulfatetype MZP. Three congurations of Mg/vacancy arrangements were calculated for the unit cell of Mg 2 Zr 8 P 12 O 48 , and the most stable one was chosen for following calculations. The volumetric and energetic differences among the three congurations are 0.19Å 3 per atom and 12 meV per atom, respectively, indicating a small dependence on the Mg arrangement.
The obtained Mg/vacancy congurations for NASICON-type and b-iron sulfate-type MZP structures and their lattice parameters were thereaer used as structural inputs for following calculations, unless specied otherwise.
As listed in Table 1, the b-iron sulfate-type MZP is slightly more stable than the NASICON-type one by only 14 meV per atom, and the cell volume per atom difference between b-iron sulfate-type MZP (14.4Å 3 per atom) and NASICON-type (14.8Å 3 per atom) is small. Therefore, both structures are reasonably stable and are expected to be good candidates for comparing the Mg-ion diffusivity.
FPMD calculations were performed to compare the Mg-ion conductivity of the NASICON-type and the b-iron sulfate-type MZP structures. Initially, all the Mg atoms were located in the 6b sites for the NASICON-type structure and 4e sites for the biron sulfate-type structure.
The supercell structure including eight MZP units (i.e., where, the supercell matrix M sc is expressed as follows: FPMD simulations were performed for NASICON-type and biron sulfate-type structure in the temperature range of 1573 to 1973 K with increments of 50 K for 100 ps. The DFT calculation settings used were equivalent to those of the structure relaxation calculations except for the cutoff energy and k-point grids as follows: the cutoff energy for the plane-wave basis was set to 350 eV and 1 Â 1 Â 1 k-point grid (only G point) was employed to reduce the computational cost. For MD simulations, the time step was set to 1 fs, and the NVT canonical ensemble was adopted using a Nosé thermostat. 36 The mean square displacements (MSD) of all elements were calculated, and then the diffusion coefficients of Mg-ion conduction and corresponding migration energies were determined. Above FPMD studies have been performed in our previous studies, showing reasonable accordance between the experimental and the computational results. 10,11,37 We evaluated the thermodynamic stability of the MZP polymorphs by rst principles DFT calculations. Temperature-and/ or pressure-dependent phase stability was also evaluated using DFT-based phonon calculations or cell-volume constraint calculations, respectively. We used the same lattice as that used for FPMD simulations. Structural relaxation calculations were performed for a structure in which the lattice constant was varied from 95 to 103% under the condition of constant lattice volume. The energy-volume relationship is tted with the Murnaghan equation of state: 38,39 EðV where, B 0 is the bulk modulus at zero pressure, B 0 is its rst derivative, E 0 is the minimum energy, and V 0 is the volume at the minimum energy. In addition, the bulk modulus B can be expressed as, where, V and p are the volume and pressure, respectively. The enthalpy-pressure diagram was constructed by combining eqn (1) and (2). Phonon frequencies and their contribution to free energy changes were calculated using the PHONOPY code. 40 Density functional perturbation theory (DFPT), 41,42 also implemented in VASP, was used to estimate the real-space force constants of the supercells. The other calculation conditions were the same as those used for the structure relaxation calculations.
From these calculations, the free energy was calculated by changing the temperature from 0 to 2000 K and the pressure was changed from 0 to 20 GPa. Table 1 DFT-derived relaxed lattice parameters, cell volume, and calculated total electron energies for the most stable structure of (a) NASI-CON-type and (b) b-iron sulfate-type MZP   .1 A 2 ). Therefore, Zr, P, and O ions remain near their original sites during the MD run, and these MSD values represent the magnitude of thermal vibrations. Similar MSD proles were obtained for the NASICON-type structure, in which only Mg ions are found to diffuse in the lattice. Diffusion path of Mg ions at 1573 K and 1973 K, i.e. both the lowest and the highest simulation temperatures, are visualised through isosurfaces of the Mg-ion probability density distribution from the FPMD simulation for (Fig. 3(a and b) and ESI Fig. S1(a) †) the NASICON-type and (Fig. 3(c and d) and ESI Fig. S1(b) †) b-iron sulfate-type structures. We conrmed the formation of a 3-dimensional pathway for both structures. Moreover, metastable Mg-ion sites are indicated in both (Fig. 3b) NASICON-type and (Fig. 3d) b-iron sulfate-type structures by varying the isosurface level (Li ion probability density). In Fig. 3b and d, the stable and metastable sites are marked with orange and blue circles, respectively. On one hand, the metastable Mg sites in b-iron sulfate-type structure are linked to three neighbouring 4e stable Mg/vacancy sites with a straight diffusion path. On the other hand, the metastable Mg sites (18e sites) are connected to two neighbouring 6b stable sites, which is consistent with the reported Li migration pathway in LiZr 2 (-PO 4 ) 3 and its derivative. 43 Diffusion coefficients were extracted through a simple linear tting with the following expression: where, n, D, and t represent the dimension of Mg-ion diffusion, ionic diffusion coefficient, and time, respectively. The ionic conductivity of Mg is estimated by the Nernst-Einstein relationship: where, R, T, z, F, and c represent the gas constant, temperature, ionic valence, Faraday constant, and ionic concentration, respectively. Since number of Mg 2+ hopping events are limited within the simulation time (100 ps) at ambient temperature according to Boltzmann statistics and attempt frequency, high temperature simulation (>1573 K) is necessary to evaluate migration energy in sufficient accuracy. Therefore, low temperature conductivities are extrapolated from high temperature data as reported in our previous FPMD calculations. 10,11,37 The Arrhenius plots for the NASICON-type and the b-iron sulfate-type structures is shown in Fig. 4, and the migration energies of Mg were evaluated by a simple linear tting scheme. Higher Mg-ion conductivity and lower migration energy are derived for the NASICON-type structure compared to those of the b-iron sulfate-type structure.
The migration energy of Mg in the b-iron sulfate-type structure obtained in this study (0.71 eV) is different from the experimental value (1.6 eV), 21 which may owe to inclusion of grain boundary conduction for experimental results. The b-iron sulfate-type structure has a higher migration energy than the NASICON-type structure (0.63 eV) in present computations. This tendency agrees that the migration energy of b-iron sulfate-type Mg 0.7 (Zr 0.85 Nb 0.15 ) 4 (PO 4 ) 6 (0.92 eV) 21 is higher than that of NASICON-type (Mg 0.1 Hf 0.9 ) 4/3.8 Nb(PO 4 ) 3 (0.68 eV). 18 This result suggests that the NASICON-type structure has better ion conductivity than the b-iron sulfate-type one.
Since NASICON-type MZP has higher Mg-ion conductivity, the pressure-and temperature-dependence of the phase stability was investigated for NASICON-type and b-iron sulfatetype structures. Fig. 5 displays the difference between the phase stability of b-iron sulfate-type structure and NASICONtype one, DE, as functions of temperature and pressure. DE is described as follows:  where, E(b-iron sulfate-type) and E(NASICON-type) represent the free energies of the corresponding structures at a given pressure and temperature. Since DE is negative for all pressure and temperature conditions, the formation of the NASICONtype structure is thermodynamically disfavoured. However, the energy difference is only 14 meV per atom between two structures (at 0 K and at 0 atm), so that formation of NASICONtype structure may be feasible such as by doping techniques. Also, since the NASICON-type structure is stable under negative pressure condition, thin lm formation technique using selection of appropriate substrate may stabilizes NASICON-type structure. In addition, we computed phase stability of MgTi 4 (-PO 4 ) 6 composition by DFT. NASICON-type structure is more stable than b-iron sulfate-type structure by 7 meV per atom, showing accordance with experimental formation of NASICONtype structure. 44 These results validate the present computational results.

Conclusions
First-principles DFT calculations were performed to clarify the difference between the Mg-ion conductivity in b-iron sulfatetype and NASICON-type MZP. The results of the FPMD calculations suggest that the NASICON-type structure has higher Mgion conductivity than the b-iron sulfate-type one. However, evaluation of the pressure-and temperature-dependent free energies by DFT calculations indicated that the formation of NASICON-type MZP is disfavoured. Therefore, the composition should be appropriately optimised to stabilize the NASICONtype structure.

Conflicts of interest
There are no conicts to declare.