# Structure, bonding and ionic mobility in Na-V-P-O glasses for energy storage applications

Steve Dave Wansi Wendji, Rémi Piotrowski, Antonio Familiari, Carlo Massobrio, Mauro Boero, Christine Tugène, Firas Shuaib, David Hamani, Pierre-Marie Geffroy, Philippe Thomas, et al.

► To cite this version:

Steve Dave Wansi Wendji, Rémi Piotrowski, Antonio Familiari, Carlo Massobrio, Mauro Boero, et al. Structure, bonding and ionic mobility in Na-V-P-O glasses for energy storage applications. Chemical Communications, 2025, 61 (59), pp.10993-10996. 10.1039/d5cc00443h. hal-05178747

HAL Id: hal-05178747

https://hal.science/hal-05178747v1

Submitted on 23 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Structure, bonding and ionic mobility in Na-V-P-O glasses for energy storage applications through experiments, first-principles and machine learning-accelerated molecular dynamics

S.D. Wansi Wendji,<sup>a,b,†</sup> R. Piotrowski,<sup>c,†</sup> A. Familiari,<sup>a,b,d</sup> C. Massobrio,<sup>b,e</sup> M. Boero,<sup>b,e</sup> C. Tugène,<sup>a</sup> F. Shuaib,<sup>c</sup> D. Hamani,<sup>c</sup> P.-M. Geffroy,<sup>c</sup> P. Thomas,<sup>c</sup> A. Pedone,<sup>d</sup> A. Bouzid,<sup>c</sup> O. Masson,<sup>c</sup> G. Delaizir,<sup>c</sup> and Guido Ori*<sup>a,b</sup>

Na-V-P-O glasses are promising cathode materials for sodium ion batteries, and yet a thorough understanding of their atomic scale behavior has so far been elusive. In this work we integrate structural and electrochemical experiments with first-principles and large-scale machine learning-accelerated molecular dynamics to elucidate quantitatively the interplay among structure, bonding, and ion mobility on space and time scales of unprecedented extensions. We unravel the existence of a broad V coordination distribution together with heterogeneous Na-ion mobility featuring percolation channels. Our results are instrumental in the search of NVP glasses optimization for electrochemical applications.

Sodium vanadium phosphate (NVP) oxides have long been a focus of research in energy storage, driven by the electrochemical properties of crystalline NVP phases.<sup>1,2</sup> Nowadays glassy and glass-ceramic counterparts are emerging as valuable alternatives, targeting improvements upon the crystalline forms in terms of synthesis, stability, scalability, and cost.<sup>3-6</sup> In particular, transition metal oxide-containing glasses enable multielectron reactions and impact thermal behavior, enhancing their energy storage potential.<sup>7,8</sup> We present here a novel characterization of NVP glasses that integrates experiments with first-principles and machine learning-accelerated molecular dynamics simulations. Our focus is on the link between structural, bonding, and ion mobility properties in such non-crystalline systems. Notably, this integrated approach is applied to NVP glasses for the first time. While detailed structural data are available for binary systems such as V<sub>x</sub>O<sub>y</sub>-P<sub>2</sub>O<sub>5</sub> and Na<sub>2</sub>O-V<sub>x</sub>O<sub>y</sub> (e.g., total structure factors and correlation functions), similar insights for ternary NVP glasses remain scarce.<sup>9,10</sup>

On the theoretical side, classical MD (CMD) has provided valuable insights, but current force fields struggle to capture the complex V local environment, particularly VO<sub>n</sub> polyhedra, as shown in recent studies of binary VP and ternary NVP glasses.<sup>11,12</sup> In these systems, we find that the account of chemical bonding via first-principles molecular dynamics (FPMD) significantly improves the description of the local V environment. Relevant findings were the presence of distinct single

<footnote><sup>a</sup> Université de Strasbourg, CNRS, Institut de Physique et Chimie des Matériaux de Strasbourg, UMR 7504, F-67034 Strasbourg, France.</footnote>

<footnote><sup>b</sup> ADYNMAT CNRS consortium, F-67034 Strasbourg, France.</footnote>

<footnote><sup>c</sup> Institut de Recherche sur les Céramiques, UMR 7315 CNRS-Université de Limoges, Centre Européen de la Céramique, 12 rue Atlantis 87068 Limoges Cedex, France</footnote>

<footnote><sup>d</sup> Department of Chemical and Geological Sciences, University of Modena and Reggio Emilia, Modena, Italy</footnote>

<footnote><sup>e</sup> Université de Strasbourg, CNRS, Laboratoire ICube, UMR 7357, Strasbourg, France.</footnote>

<footnote>* E-mail: guido.ori@cnrs.fr</footnote>

<footnote>† Supplementary Information available.</footnote>

<footnote>‡ These authors contributed equally to this work.</footnote>
and double/vanadyl bonding fingerprints, the stronger network-forming role of V⁵⁺, and the higher degree of phosphate network polymerization, aligning well with experimental results.

Overall, the present study advances the understanding of NVP glasses by providing a multifaceted characterization of a series of Na₂O-VₓOᵧ-P₂O₅ glasses (Tables 1 and S1 in ESI). Within this purpose, experimental methods include X-ray diffraction, X-ray photoelectron spectroscopy (XPS), and electrochemical impedance spectroscopy (EIS). On the atomic-scale simulation side, we use a full thermal annealing within Born-Oppenheimer MD (f-BOMD) instead of a short equilibration at 300K (s-BOMD), as previously reported,¹¹,¹². To significantly extend the reach of space and time scale we take full advantage of machine learning interatomic potential (MLIP) trajectories based on Gaussian Approximation Potential (GAP), trained on f-BOMD data.¹³⁻¹⁵ f-BOMD enables full relaxation of the glass constituents during the melting and cooling processes, resulting in a more accurate medium-range structure compared to s-BOMD, which only allows local relaxation while retaining most of the CMD medium-range order.

<table><caption>Table 1 Chemical compositions (in mol (%)) and the glass transition temperature ($T_g$), along with density (d) of Na₂O-VₓOᵧ-P₂O₅ systems NVP0, NVP10, and NVP28 glasses. See ESI for simulated compositions.</caption><thead><tr><th>System</th><th>Na₂O<br>(%)</th><th>VₓOᵧ [V₂O₄ + V₂O₅]<sup>a</sup><br>(%)</th><th>P₂O₅<br>(%)</th><th>d<sup>b</sup><br>(g/cm³)</th><th>T<sub>g</sub><sup>c</sup><br>(°C)</th></tr></thead><tbody><tr><td>NVP0<sup>b</sup></td><td>—</td><td>50.0 [17.6 + 32.4]</td><td>50.0</td><td>2.800</td><td>395<sup>d</sup></td></tr><tr><td>NVP10</td><td>10.0</td><td>40.0 [13.2 + 26.8]</td><td>50.0</td><td>2.777</td><td>439</td></tr><tr><td>NVP28</td><td>28.5</td><td>43.00 [5.59 + 37.41]</td><td>28.5</td><td>2.939</td><td>297</td></tr></tbody></table>

a,b,cEstimated by XPS (a), He pycn. (b), and DSC (c) (see ESI). dFrom ref. ⁹.

<figure><img src="image_2.png" /><figcaption>Fig. 1 X-ray S(k) (top) and G(r) (bottom) of NVP10 glass, showing the comparison between the calculated data (CMD, s-BOMD, f-BOMD, and MLIP) and the experimental data.</figcaption></figure>

Indeed, by comparing the performances of the different models, Figure 1 underpins the significant improvements of f-BOMD over CMD and s-BOMD for NVP10, in terms of X-rays total structure factors S(k) and reduced total pair correlation functions G(r). These changes are exemplified, for both S(k) and G(r) (Table 2 and S5 in ESI), by the substantial
<figure><img src="image_3.png" /><figcaption>Fig. 2 Left: Partial pair correlation function $g_{VO}(r)$ for NVP10 glass obtained by CMD, s-BOMD, f-BOMD, and MLIP. Right: Atomistic models by f-BOMD (Na (yellow), V (blue), P (orange), and O (red)) with (top) local spin density isosurfaces (0.05 a.u., in purple) and (bottom) Wannier centers (green) within a V local environment in a VO$_n$ polyhedron.</figcaption></figure>

reduction in the $R_\chi$ parameter measuring the quantitative agreement with respect to experimental data. In reciprocal space, f-BOMD enhances $S(k)$ by providing a broadened first peak and a minor contribution around $1 \AA^{-1}$, more accurately capturing intermediate-range order. In direct space and referring to $G(r)$, f-BOMD performs better when compared to experimental data for the first and second peak, by revealing effectively the presence of single V–O and double V=O bonds, respectively. This holds true also for the fourth peak, corresponding to polyhedral connections involving V–V and V–P distances. Figure 2 provides in a comparative fashion the partial pair correlation function $g_{VO}(r)$ for V–O pairs obtained from the different computational schemes. The s-BOMD model features well separated single and double V–O bonds, with peaks at ~1.6 Å (V=O) and ~1.83 Å (V–O), consistent with previous results$^{11,12}$, improving upon the single broadened peak obtained with CMD. This separation is better defined within f-BOMD showing a well-defined minimum between peaks, and a V–O bond distance (~1.88 Å) closer to experimental data (~1.92 Å).

<table><caption>Table 2 Comparison of the agreement between CMD, s-BOMD, f-BOMD, and MLIP models (small (~400 atoms) and large (~3200 atoms)) and experimental data using goodness-of-fit $R_\chi$ parameters for X-ray total structure factor $S(k)$ of NVP0, NVP10, and NVP28 glasses at 300 K, as well as for neutron total structure factor for NVP0.</caption><thead><tr><th rowspan="2">Method</th><th rowspan="2">Size</th><th colspan="2">NVP0</th><th>NVP10</th><th>NVP28</th></tr><tr><th>R<sub>χ</sub><sup>Sχ(k)</sup></th><th>R<sub>χ</sub><sup>SV(k)</sup></th><th>R<sub>χ</sub><sup>Sχ(k)</sup></th><th>R<sub>χ</sub><sup>SV(k)</sup></th></tr></thead><tbody><tr><td>CMD</td><td>small</td><td>11.3±0.1</td><td>9.8±0.1</td><td>17.2±0.1</td><td>17.1±1.0</td></tr><tr><td>s-BOMD</td><td>small</td><td>8.4±0.2</td><td>8.4±0.2</td><td>14.4±0.2</td><td>9.8±1.2</td></tr><tr><td>f-BOMD</td><td>small</td><td>-</td><td>-</td><td>12.8±0.1</td><td>-</td></tr><tr><td rowspan="2">MLIP</td><td>small</td><td>8.7±0.1</td><td>9.9±0.1</td><td>11.5±0.5</td><td>10.4±0.5</td></tr><tr><td>large</td><td>7.1±0.1</td><td>8.1±0.1</td><td>11.9±0.1</td><td>10.3±0.1</td></tr></tbody></table>
<table><caption>Table 3 Average V coordination number ($n_V$), bond lengths ($r_{ij}$ (in Å), taken as the position of the first maximum of the $g_{\alpha\beta}(r)$), and distribution of individual $n_V(l)$ structural units of V of l-fold coordinated computed for NVP10 glass, comparing CMD, s-BOMD, f-BOMD, and MLIP models.</caption><thead><tr><th></th><th>Exp.</th><th>CMD</th><th>s-BOMD</th><th>f-BOMD</th><th>MLIP</th></tr></thead><tbody><tr><td>$n_V$</td><td>4.20-5.40</td><td>4.33</td><td>4.73</td><td>4.60</td><td>4.93</td></tr><tr><td>$r_{VO}$</td><td>1.59-1.78</td><td>-</td><td>1.59</td><td>1.57</td><td>1.56</td></tr><tr><td></td><td>1.8-2.3(1.92)</td><td>1.78</td><td>1.82</td><td>1.86</td><td>1.86</td></tr><tr><td colspan="6">Structural units</td></tr><tr><td>$l=3$</td><td>VO<sub>3</sub></td><td>5.1±2</td><td>-</td><td>-</td><td>-</td></tr><tr><td>$l=4$</td><td>VO<sub>4</sub></td><td>60.6±2</td><td>38.4±2</td><td>44.7±0.1</td><td>19.4±3.3</td></tr><tr><td>$l=5$</td><td>VO<sub>5</sub></td><td>30.1±2</td><td>50.2±3</td><td>50.3±0.2</td><td>64.8±1.6</td></tr><tr><td>$l=6$</td><td>VO<sub>6</sub></td><td>4.2±1</td><td>11.4±5</td><td>5.0±0.1</td><td>15.9±2.5</td></tr></tbody></table>

In Figure 2, present a spin topology analysis of the NVP10 model obtained from f-BOMD, characterized by spin localization on V sites corresponding to paramagnetic $V^{4+}$ and allowing for a precise speciation between $V^{5+}$ and $V^{4+}$ sites. Additionally, we provide an atomistic perspective on the unique bonding characteristics of VO$_n$ polyhedra, based on maximally localized Wannier functions (WFCs).¹¹ Within the spin-unrestricted DFT-BOMD framework, this method concur to identify different bonding types: single V–O bonds, defined by two WFCs shared between connected atoms, and double/vanadyl V=O bonds, characterized by more than two WFCs shared between connected atoms. Accurate assessment of the local V environment enables better partitioning bridging and non-bridging oxygen contents (BO and NBO, respectively) in the glass network, a key factor impacting Na ion dynamics.¹⁶⁻¹⁸ In Table 3, we also report the V coordination number ($n_V$), calculated by integrating the first peak of the $g_{VO}(r)$, as well as the distribution of individual $n_V(l)$ structural units for V with l-fold coordination in the NVP10 glass. The f-BOMD data provide a novel perspective, distinct from CMD results and building upon s-BOMD data, showing a higher $n_V$ (4.6) and a broad coordination distribution (4 to 6), with a dominant presence of VO$_5$ polyhedra (~50%). These findings strongly underscore the role of V as a network former in the glass structure.¹¹,¹²

To enhance computational efficiency without compromising the accuracy achieved through BOMD, we leveraged accelerated- MD by a MLIP of GAP type, fitted over the NVP10 model for the f-BOMD datasets (Fig. S2 in ESI). Our new MLIP shows a remarkable overall performance in terms of energy, forces, and virials amounting to 5.8 meV/atom, 0.4 eV/Å and 16.9 meV/atom respectively when calculated with respect to DFT accuracy for testing sets. Close behaviors are found on training and testing sets. Figure 1 compares the performance of the newly developed MLIP for NVP10 glass, showing that the MLIP model closely matches f-BOMD results for both $S(k)$ and $G(r)$, accurately describing the NVP10 structure and bonding distances. This is further substantiated by the comparable $R_χ$ values for the MLIP and f-BOMD in ~400-atom models, with additional improvements observed when considering a significantly larger system (~3200 atoms, clearly a dimension beyond the capability of FPMD/BOMD when aiming at the production of a significant time trajectory). MLIP also reproduces accurately V–O and V=O bond distances, $n_V$, and $n_V(l)$ when compared to f-BOMD data. The MLIP approach allows reducing drastically the computational cost from about 220 days for f-BOMD to 4 days with MLIP (Table S6 in ESI).
Figure 3 allows a comparison between the total correlation functions $T(r)$ for NVP0 and $G(r)$ for NVP28 glasses, highlighting the robustness and transferability of our MLIP to compositions beyond the training set (see also structure factors $S(k)$ comparison in Fig. S4 and S5 in ESI). For NVP0, key improvements include the precise description of the first peak near

<figure><img src="image_4.png" /><figcaption>Fig. 3 Comparison of experimental and calculated total correlation functions data for NVP0 ($T(r)$) and NVP28 ($G(r)$) glasses.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 4 Left: MSD of Na ions in NVP28 glass at 1200 K, simulated by MLIP, highlighting highly mobile (green) and trapped ions (orange). Center: Snapshot with Na ions color-coded by MSD value and $PO_n/VO_n$ units as transparent polyhedra. Right: Arrhenius plot of ionic conductivity ($\log(\sigma)$ vs. $1000/T$), comparing this work (green circles: sim. and exp.) with exp. data for other Na-ion NVP glasses (grey and orange symbols.¹⁹,²⁰)</figcaption></figure>

1.8 Å⁻¹ in $S(k)$ (Fig. S4) and all four characteristic peaks in X-rays $T(r)$. Notably, MLIP captures the second peak in $T(r)$ at ~1.9 Å (V–O bond) and the low-$k$ region of $S(k)$, reflecting its improved ability to model short-range V–O bonds and intermediate-range order. For NVP28, a high-$Na_2O$ composition relevant to energy applications, MLIP outperforms both s-BOMD and CMD, particularly for the first peak in $S(k)$ (Fig. S5), the low-$k$ region, and the 4th peak of $G(r)$ near 3.3 Å, thereby highlighting realistic $VO_n/PO_n$ polyhedral interconnections.

Fig. S6 in ESI is indicative of the performance of our MLIP when describing NVP thermal relaxation, leading to an estimated
glass transition temperature ($T_g$) in the range 540-650K when comparing 400- and 3200-atom models, a large improvement with respect to the typical overestimation of CMD (~1500K), and fairly in line with respect to the experimental value of ~570K. 18,21 The accurate structural and thermal description of NVP28 by our MLIP allows the reliable assessment of the ionic conductivity ($\sigma_{ion}$) via Na diffusion coefficients obtained from mean square displacements (MSD). 22-24 MSD analysis (Fig. 4) reveals a highly heterogeneous dynamical pattern for Na ions, with a temperature- and composition-dependent fraction of them exhibiting negligible mobility and vibrating within local cages, while some atoms do feature significantly higher mobility. This heterogeneity correlates with partial Na ordering, forming Na-rich regions and percolation channels. Such ordering, especially at intermediate range distance, is reflected by a peak in the Faber-Ziman partial structure factor at ~1 Å (Fig. S7). 25 At T= 473 K, the NVP28 glass exhibits a total experimental conductivity of 3.3 × 10⁻⁵ S/cm (Fig. 4), determined by EIS. The nearly ideal semicircle of the Nyquist plot (Fig. S8), though slightly depressed 26, suggests a coexistence of minor ionic (Na⁺) and dominant electronic conductivity, driven by small polaron hopping of V⁴⁺ and V⁵⁺ sites. This is consistent with the relatively low ionic conductivity contribution (~5%) from Na⁺ ions estimated via MLIP (Fig. 4 and S6 in ESI). Notably, Fig. 4 also includes data for two additional NVP glasses with similar Na₂O content but differing by V/P ratios and V speciation, where mixed ionic and electronic conductivity has been reported, with ionic contributions reaching up to ~50%. 19,20 While the comparable total conductivity across these results is consistent, the structural differences underscore the need for deeper investigations. This is where the new MLIP, with its ability to provide an accurate description of the structure and a detailed partitioning of BO/NBO roles in Na-ion dynamics. This study establishes a comprehensive understanding of the structural and dynamic properties of Na-V-P-O glasses, leveraging a novel combination of experiments, first-principles simulations, and ML-accelerated MD. By providing unprecedented insights into the medium-range order, vanadium speciation, and Na-ion transport pathways, this work lays the foundation for the design of high-performance glass-based materials for energy storage applications. The demonstrated accuracy and scalability of the MLIP approach open avenues for extending these methods to other complex amorphous systems, accelerating innovation in glass science and energy technology.

This work was supported by Agence Nationale de la Recherche (AMSES project, ANR-20-CE08-0021) and EUR-QMAT-ANR. The HPC centers of GENCI (A0160807670, A0160912441, A0160910832 and AD010914978) and of the University of Strasbourg (HPC Équipe@Meso and CPER Alsac calcul/Big Data) are acknowledged for computational resources.

## Conflicts of interest

There are no conflicts to declare.

## Data availability

The data supporting this article is included as part of the SI, and simulated data are available at the NOMAD repository²⁷.

## Notes and references

1. F. Kong, X. Liang, L. Yi, X. Fang, Z. Yin, Y. Wang, R. Zhang, L. Liu, Q. Chen, M. Li et al., *Energy*, 2021, **219**, 119513.

2. E. Boivin, J.-N. Chotard, C. Masquelier and L. Croguennec, *Molecules*, 2021, **26**, 1428.
121 3 S. Gandi, V. K. Katta, D. P. Dutta and B. R. Ravuri, *New J Chem*, 2020, **44**, 2897-2906.

122 4 Z. Wei, D. Wang, X. Yang, C. Wang, G. Chen and F. Du, *Adv. Mater. Interfaces*, 2018, **5**, 1800639.

123 5 J. Zhang, Y. Li, Z. Chen, Q. Liu, Q. Chen and M. Chen, *Energy Environ. Mater.*, 2023, **6**, e12573.

124 6 S. Yan, K. P. Abhilash, L. Tang, M. Yang, Y. Ma, Q. Xia, Q. Guo and H. Xia, *Small*, 2019, **15**, e1804371.

125 7 M. S. Whittingham, *Chem. Rev.*, 2014, **114**, 11414-11443.

126 8 R. Giovanardi, M. Montorsi, G. Ori, J. Cho, T. Subhani, A. Boccaccini and C. Siligardi, *J. Mater. Chem.*, 2010, **20**, 308-313.

128 9 U. Hoppe, N. Wyckoff, M. Schmitt, R. Brow, A. Schöps and A. Hannon, *J. Non-Cryst. Solids*, 2012, **358**, 328-336.

129 10 U. Hoppe, A. Ghosh, S. Feller, A. Hannon, D. Keen and J. Neuefeind, *J. Non-Cryst. Solids*, 2021, **572**, 121120.

130 11 S. Wansi Wendji, C. Massobrio, M. Boero et al., *J. Non-Cryst. Solids*, 2024, **634**, 122967.

131 12 S. Wansi Wendji, R. Piotrowski, C. Massobrio et al., *SSRN preprint*, 2024, **4997219**.

132 13 A. P. Bartók, M. C. Payne, R. Kondor and G. Csányi, *Phys. Rev. Lett.*, 2010, **104**, 136403.

133 14 F. Shuaib, G. Ori, P. Thomas, O. Masson and A. Bouzid, *J. Am. Ceram. Soc.*, 2024, **108**, e20128.

134 15 T.-L. Pham, M. Guerboub, S. Wansi Wendji et al., *arXiv*, 2024, 240411442.

135 16 Y. Onodera, Y. Takimoto, H. Hijiya, Taniguchi et al., *NPG Asia Mater.*, 2019, **11**, year.

136 17 J. Habasaki, C. Leon and K. L. Ngai, *Dynamics of glassy, crystalline and liquid ionic conductors*, Springer, Cham, 2016.

137 18 A. Pedone, M. Bertani, L. Brugnoli and A. Pallini, *J. Non-Cryst Solids X*, 2022, 100115.

138 19 A. Sharma, I. Suzuki, K. Toyooka et al., *J. Phys. Chem. C Nanomater. Interfaces*, 2024, **128**, 9793-9801.

139 20 M. Wasiucionek, J. Garbarczyk, P. Kurek and W. Jakubowski, *Solid State Ion.*, 1994, **70-71**, 346-349.

140 21 F. Lodesani, M. C. Menziani, H. Hijiya, Y. Takato, S. Urata and A. Pedone, *Sci. Rep.*, 2020, **10**, 2906.

141 22 T.-L. Pham, M. Guerboub, A. Bouzid et al., *J. Mater. Chem. A*, 2023, **11**, 22922.

142 23 *Theory and simulation in physics for materials applications*, ed. E. Levchenko, Y. J. Dappe and G. Ori, Springer, Cham, 2020.

144 24 C. Massobrio, *The structure of amorphous materials using molecular dynamics*, Inst. of Phys. Pub., London, 2022.

145 25 S. S. Sørensen, M. M. Smedskjaer and M. Micoulaut, *J. Phys. Chem. B*, 2023, **127**, 10179-10188.

146 26 M. Saad, W. Stambouli et al., *Mater. Res. Bull.*, 2017, **89**, 224-231.

147 27 *Data deposited at the NOMAD data repository*, 2024.