# Enhanced structural description of sodium vanadium phosphate glasses: A combined experimental and molecular dynamics study

S.D. Wansi Wendji, R. Piotrowski, C. Massobrio, M. Boero, C. Tugène, F. Shuaib, D. Hamani, P.-M. Geffroy, P. Thomas, A. Bouzid, et al.

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

## To cite this version:

S.D. Wansi Wendji, R. Piotrowski, C. Massobrio, M. Boero, C. Tugène, et al., Enhanced structural description of sodium vanadium phosphate glasses: A combined experimental and molecular dynamics study. Journal of Non-Crystalline Solids, 2025, 655, pp.123420. 10.1016/j.jnoncrysol.2025.123420. hal-05180861

HAL Id: hal-05180861

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

Submitted on 23 Jul 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Enhanced structural description of sodium vanadium phosphate glasses: A
combined experimental and molecular dynamics study

S.D. Wansi Wendji<sup>a,b,1</sup>, R. Piotrowski<sup>c,1</sup>, C. Massobrio<sup>b,d</sup>, M. Boero<sup>b,d</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. Bouzid<sup>c</sup>, O. Masson<sup>c</sup>, G. Delaizir<sup>c</sup>, G. Ori<sup>a,b,2</sup>

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

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

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

<sup>d</sup>Université de Strasbourg, CNRS, Laboratoire ICube, UMR 7357, Strasbourg, France

**Abstract**

Structural and bonding insights into sodium vanadium phosphate (NVP) glasses are crucial for optimizing their performance as cathode materials in sodium-ion batteries. This study quantitatively assesses the structural features and bonding characteristics of two NVP glass compositions: $(Na_2O)_\alpha-(V_xOy)_{(1-2\alpha)}-(P_2O_5)_\alpha$ ($\alpha$ = 0.375, 0.285). We combine experimental characterization (differential scanning calorimetry, X-ray diffraction and X-ray photoelectron spectroscopy,) and atomistic modeling (classical molecular dynamics (CMD), and Born-Oppenheimer molecular dynamics (BOMD)). This work provides a quantitative analysis of the different $VO_n$ units in the two NVP glass models, superseding previous knowledge based largely on CMD simulations. Our results show that the account of the electronic structure, inherent in BOMD simulations is essential for capturing the V-O bonding environment. This includes the splitting of the V-O peak in the pair distribution function due to both short V=O and longer V-O bonds, a higher degree of polymerization in the phosphate network and a more significant role for $V^{5+}$ as a network former.

**Keywords:** Sodium vanadium phosphate glasses, X-ray diffraction, X-ray photoelectron spectroscopy, molecular dynamics, Born-Oppenheimer

**1. Introduction**

The incorporation of alkali ions, such as lithium and sodium, into transition metal oxide (TMO) containing glasses, like vanadium phosphate (VP) and vanadium tellurite (VT) glasses, significantly enhances their functional properties, making these materials highly attractive for solid-state batteries and electrochemical storage devices [1–6]. This is due to their unique combination of properties, including high glass-forming ability, broad vitrification range, low glass transition and melting temperatures, and high thermal expansion coefficients [7–9]. TMO glasses, such as those involving vanadium and iron oxides, are characterized by multi-electron reactions, enhancing energy density performance [5, 7, 8, 10]. Alkali ions in these glasses act as modifiers by disrupting covalent bonds thereby allowing

for mobile charge carriers and increasing ionic conductivity, crucial for energy storage efficiency [11]. While initial studies focused on lithium oxide due to the dominance of lithium-ion batteries [6, 12, 13], there is growing interest in sodium-based systems, such as sodium vanadium phosphate (NVP) glasses, driven by sodium's abundance and lower cost [14–18].

Previous experimental studies on NVP systems have predominantly concentrated on optical, magnetic, electrical, and spectroscopic properties, with limited attention to detailed structural characterization [19–24]. This gap has been partially addressed by classical molecular dynamics (CMD) simulations using empirical potentials, which have been employed to model materials relevant to energy storage applications [25–27], including glasses and glass-ceramics [10, 28–30]. However, CMD struggles to accurately capture local transition metal bonding environments, particularly for vanadium oxides. Advances in high-performance computing, code efficiency, and parallelization now enable the

<footnote n="1"><sup>1</sup>These authors contributed equally to this work.</footnote>

<footnote n="2"><sup>2</sup>Corresponding author e-mail: guido.ori@cnrs.fr</footnote>
study of complex glasses and amorphous materials us-
ing first-principles molecular dynamics (FPMD) [31-
34], whether through Car-Parrinello (CPMD) [35–37]
or Born-Oppenheimer (BOMD) [38] approaches. This
has been recently exemplified in our study of a proto-
typical VP glass, performed by combining CMD with
BOMD [38].

Our results for 50VₓOᵧ–50P₂O₅ glass highlighted the
local coordination environment of vanadium sites, ad-
dressing the controversial presence of four-fold coor-
dinated V⁴⁺ ions. Experimental techniques such as
XANES [39–41], as well as NMR and EPR [42], have
suggested tetrahedral coordination for V⁴⁺ in certain
glass systems. However, studies on crystalline vana-
dium oxides and amorphous Na₂O–VₓOᵧ–P₂O₅ sys-
tems indicate that only V⁵⁺ ions adopt tetrahedral co-
ordination [43–47], leaving the coordination state of
V⁴⁺ unresolved. In our previous BOMD study on
50VₓOᵧ–50P₂O₅ glass, we clearly demonstrated the ab-
sence of V⁴⁺ ions in four-fold coordination, showing the
predominance of V⁵⁺ in this configuration, this result
being largely in contrast with the CMD models.

In the present work, we extend this combined
CMD/BOMD approach to NVP glasses. In particu-
lar, we study two (Na₂O)α–(VₓOᵧ)(1-2α)–(P₂O₅)α, α =
0.375, 0.285 glasses both experimentally and by MD
simulations. These glasses are synthesized using melt-
quenching technique and characterized experimentally
through density measurements, X-ray photoelectron
spectroscopy (XPS), differential scanning calorimetry
(DSC) and X-ray diffraction (XRD). X-ray diffraction,
in particular, provided total structure factors and pair
correlation functions, which are essential for evaluat-
ing the quantitative accuracy of our atomistic models.
On the atomic-scale simulation side, we first generated
NVP glass models via the melt-quenching approach by
CMD, followed by a BOMD refinement to obtain more
quantitative insights into the glass structure and proper-
ties. Notably, we employed the recently updated BMP-
shrm potential in our CMD simulations [25, 48, 49],
which outperforms other empirical potentials for VP
glasses [38]. The further refinement of CMD-derived
VP models through BOMD equilibration at 300K al-
lowed capturing the vanadium’s local environment by
accounting for electronic structure effects not within the
reach of CMD calculations.

The paper is organized as follows: Section 2 details the
experimental and computational methods. Results are
presented in Section 3 for thermal behaviour and Vana-
dium sites speciation, Section 4 for structure factors,
Section 5 for total and partial pair correlation functions,
coordination numbers, and structural units, and Sec-

tion 6 for bond angle distributions, Q" distributions and
preference ratios. Section 7 presents a detailed analysis
of the results, highlighting general remarks and revisit-
ing the structure of NVP glasses with respect to short-
and intermediate-range orders. Finally, concluding re-
marks are presented in Section 8.

## 2. Experimental and computational methodologies

2.1. NVP glass samples synthesis and characterization

Glasses were prepared using sodium carbonate Na₂CO₃ (Strem Chemicals, 99.5%), vanadium pentoxide V₂O₅ (Strem Chemicals, 99.5%), and ammonium phosphate dibasic (NH₄)₂HPO₄ (Sigma-Aldrich, 98%) as raw materials. The precursors were weighed in appropriate amounts, mixed in an agate mortar, and heated in a platinum crucible to produce two glass compositions, (Na₂O)α–(VₓOᵧ)(1-2α)–(P₂O₅)α, with α = 0.375 and 0.285. These are referred to as NVP25 and NVP43, respectively, based on the total vanadium content, as shown in Table 1. A preliminary heat treatment at 300°C for 2h was required in order to decompose Na₂CO₃ and (NH₄)₂HPO₄. The temperature was then slowly raised at the heating rate of 5°C/min from room temperature up to 900°C. The melt was kept at this temperature for 1h with a homogenization each 15 min. Finally, the melt was quenched at room temperature between two stainless steel plates to get glass sample. Density was measured on the glass powder using helium pycnometer (Micromeritics AccuPyc II 1340) in a 1 cm³ cell. The thermal properties of the two glasses were determined using differential scanning calorimetry (DSC, model TA Instrument AQ20) from 40° to 550°C at 10 °C/min, allowing to determine the glass transition and crystallization onset temperatures (Tg and Tc, respectively), together with the thermal stability that is defined according to the Dietzel criterion as ΔT=Tc-Tg. XPS measurements to determine the speciation of V sites were performed using a Kratos Axis Ultra DLD spectrometer with a monochromatic Al Kα source (12 mA, 15 kV). The working function of the instrument was calibrated to give a binding energy of the metallic gold Au 4f₅/₂ of 84.0 eV, and the dispersion was adjusted to give a binding energy of metallic copper Cu 2p₃/₂ at 932.6 eV. The Kratos charge neutralization system was used for each analysis, and still completed so as to set the V 2p₃/₂ (V⁵⁺) signal at 517.2 eV. The pressure in the analysis chamber was around 5×10⁻⁹ Torr. Survey spectra were obtained over an analysis area of approximately 700x300 µm and a resolution of 160 eV (Pass Energy).
<table><caption>Table 1: Glass composition and V⁴⁺ content of the two NVP synthesized and simulated models by CMD and BOMD.</caption><thead><tr><th>Systems</th><th></th><th>Na₂O<br>(%)</th><th>VₓOᵧ [V₂O₄ + V₂O₅]<br>(%)</th><th>P₂O₅<br>(%)</th><th>V⁴⁺<br>(%)</th></tr></thead><tbody><tr><td rowspan="3">NVP25</td><td>Exp.</td><td>37.50</td><td>25.00 [1.25 + 23.75]</td><td>37.50</td><td>5.0<sup>b</sup></td></tr><tr><td>CMD</td><td>37.50</td><td>25.00 [2.78 + 22.22]</td><td>37.50</td><td>11.1</td></tr><tr><td>BOMD</td><td>37.5</td><td>25.00 [2.78 + 22.22]<sup>c</sup></td><td>37.5</td><td>11.1<sup>c</sup></td></tr><tr><td rowspan="3">NVP43</td><td>Exp.</td><td>28.50</td><td>43.00 [5.59 + 37.41]</td><td>28.50</td><td>13.0<sup>b</sup></td></tr><tr><td>CMD</td><td>28.57</td><td>42.86 [6.43 + 36.43]</td><td>28.57</td><td>16.7</td></tr><tr><td>BOMD</td><td>28.57</td><td>42.86 [7.86 + 35.00]<sup>c</sup></td><td>28.57</td><td>18.3<sup>c</sup></td></tr></tbody></table>

<sup>b</sup>Estimated by XPS measurements.

<sup>c</sup>Estimated by spin topology analysis at PBE0 level from configurations sampled from BOMD trajectory.

<table><caption>Table 2: Density (d), glass transition (T<sub>g</sub>) and crystallization onset temperature (T<sub>c</sub>) of synthesized NVP25 and NVP43 glasses.</caption><thead><tr><th>Systems</th><th>d<sup>a</sup><br>(g/cm³)</th><th>T<sub>g</sub><sup>b</sup><br>(°C)</th><th>T<sub>c</sub><sup>b</sup><br>(°C)</th><th>ΔT<sub>T<sub>c</sub>-T<sub>g</sub></sub><br>(°C)</th></tr></thead><tbody><tr><td>NVP25</td><td>2.807</td><td>327</td><td>436</td><td>109</td></tr><tr><td>NVP43</td><td>2.939</td><td>297</td><td>377</td><td>80</td></tr></tbody></table>

<sup>a</sup>Estimated by He pycnometer measurements. <sup>b</sup>Estimated by DSC measurements.

High Resolution spectra were obtained over the same analysis area at a resolution of 20 eV (Pass Energy). The spectra were analyzed using CasaXPS software [50] (v. 2.3.25) with curves generally mixing 70% Gaussian and 30% Lorentzian, and Shirley mode for inelastic electron background subtraction.

The experimental total X-ray structure factor S(k) and reduced total pair distribution functions G(r) of our NVPs systems were determined through X-ray total scattering following procedures similar to those employed in [51]. X-ray scattering measurements were conducted at room temperature using a specialized laboratory setup based on a Bruker D8 Advance diffractometer. This instrument was equipped with a silver sealed tube (λ = 0.559422 Å) and a rapid LynxEye XE-T detector. Modifications were made to this setup to maximize the collected intensities, minimize spurious signals from the empty environment, and achieve good counting statistics up to a large scattering vector length of 21.8 Å⁻¹. Approximately 20 mg of each sample's powder were placed in a thin-walled (0.01 mm) borosilicate glass capillary with a diameter of about 0.7 mm to limit absorption effects. The μR values (where R is the capillary radius and μ is the sample's linear attenuation coefficient) were estimated based on precise measurements of the mass and dimensions of the samples, yielding values of 0.17 and 0.22 for NVP25 and NVP43 glasses, respectively.

After sealing, the capillary was mounted on a goniometric head and aligned so that its axis coincided with the diffractometer's goniometer axis. Data acquisition involved multiple scans over the 0–152°, 50–152°, and 100–152° 2θ ranges with a step size of 0.01°. These scans were then merged, resulting in a total equivalent acquisition time of approximately 50 hours per sample. The raw data were corrected, normalized, and Fourier transformed using custom software [52] to obtain the reduced atomic pair distribution functions G(r). Corrections accounted for capillary contributions, empty environment, Compton and multiple scatterings, absorption, and polarization effects. The necessary X-ray mass attenuation coefficients, atomic scattering factors, and Compton scattering functions for data correction and normalization were calculated from tabulated data provided by the DABAX database [53]. Absorption corrections were evaluated using a numerical midpoint integration method, where the sample cross-section was divided into small subdomains, following a method similar to that proposed by A. K. Soper and P. A. Egelstaf [54]. The normalization constant was determined using the high-angle method, ensuring that the coherently scattered intensity oscillated weakly around the sample's average atomic scattering factor at high scattering vector lengths (from k<sub>max</sub>/2 to k<sub>max</sub>).

## 2.2. CMD and BOMD simulations details

The computational procedure employed in the previous work [38] was adopted, consisting of a full CMD thermal cycle followed by a final equilibration by means of BOMD at 300 K. We employed BMP-shrm (detail in the SI) developed by Bertani *et al.* [48] based on the original parameters developed for V⁵⁺ and V⁴⁺ by Ori *et al.* [25] and Pedone *et al.* [55] together with the V–O–V and V–O–P three-body parameters recently reported by Malavasi *et al.* [49]. This choice of potential parameters was intended to exploit the best force field among the
three tested previously [38]. Three replicas of NVP25 and NVP43 with chemical compositions detailed in Table 1 were generated randomly and placed in a cubic cell with dimensions adjusted to match the reference densities values (see Table 2). In our CMD simulations, the short-range interactions were truncated at a distance of 7.0 Å and the long-range interactions were calculated using the Ewald summation method with a precision of $10^{-5}$ eV, up to a cutoff distance of 8.5 Å. The velocity Verlet algorithm was employed with a time step of 1 fs to integrate the equation of motion. The NVP glasses models were obtained by melt-quenching thermal procedure [38] with the Nosé-Hoover thermostat [56–58] used to control the temperature. First, the initial structures were relaxed at 300 K in the NVT ensemble. Each system was heated to 4000 K and kept for 150 ps at this temperature followed by the cooling from 4000 K to 300 K with a constant rate of 2.5 K/ps. The relevant properties are averaged over the structures obtained from the final 150 ps of a 300 ps trajectory at 300 K across different replicas. BOMD calculations for these two glasses were performed at 300 K after the CMD dynamics. BOMD calculations for these two glasses were performed at 300 K after the CMD dynamics. BOMD simulations were performed using the Born-Oppenheimer scheme, with atom-centered Gaussian-type basis functions to describe the orbitals, and an auxiliary plane-wave basis set to re-expand the electron density. We employed double-ζ polarized MOLOPT basis sets [59] for V, P, and O, and used a cutoff energy of 800 Ry for the plane waves at the Γ point only. Using the CMD-derived structures as a starting point, we performed a BOMD simulation at 300 K for approximately 5 ps using Goedecker–Teter–Hutter pseudopotentials to describe core-valence interactions [60] together with the semilocal Perdew-Burke-Ernzerhof (PBE) [61, 62] exchange-correlation functional within a spin unrestricted DFT formalism. A final optimization of the electronic structure was obtained employing PBE0 hybrid functional [63] containing 25 % nonlocal Hartree-Fock exchange, with a cutoff energy of 800 Ry. For the BOMD analysis, the last 3 ps of the trajectories were used and the results averaged over the three replicas. Final optimization using the PBE0 hybrid functional yields a more accurate representation of the electronic structure, energy bandgap, and spin topology [64], consistent with prior electronic analyses [38]. The CMD simulations were performed using the DL-POLY code (version 4.10.0) [65], while the BOMD simulations utilized the CP2K package [66] with the DZVP-MOLOPT-SRGTH basis set [67].

### 3. Characteristic temperatures, thermal stability and Vanadium sites speciation

The experimental and simulated chemical compositions of the two NVP glasses are reported in (Table 1). The densities of the synthesized glasses, NVP25 and NVP43, are 2.807 and 2.939 g/cm³, respectively (Table 2). These values indicate that NVP43 has a greater mass per unit volume, likely due to its higher vanadium oxide content at the expense of sodium oxide. The DSC thermograms of the two synthesized glasses are presented in Figure 1. DSC thermograms illustrate the thermal behavior of the two glass compositions, while Table 2 lists the characteristic temperature values obtained. NVP25 exhibits a higher glass transition tem-

<figure><img src="image_3.png" /><figcaption>Figure 1: DSC plots of NVP25 (top) and NVP43 (bottom) glass powder showing glass transition temperature ($T_g$) and the crystallization onset ($T_c$).</figcaption></figure>

perature ($T_g$) of ~327°C compared to NVP43 (~297°C), indicating a more rigid structure, likely due to its higher phosphate content. Phosphorus, being a stronger network former than vanadium, contributes to a more interconnected and rigid glass network. In contrast, NVP43, with a higher vanadium content, forms a less rigid structure. Although NVP25 likely contains more sodium, which typically increases flexibility, its higher $T_g$ em-
phasizes the dominant role of phosphorus in enhanc-
ing network stiffness. This is consistent with the even
higher T<sub>g</sub> of ~397°C in VP50 glass, further support-
ing the role of phosphorus [68]. NVP25 also shows a
higher crystallization temperature (T<sub>c</sub>) of ~436°C com-
pared to NVP43 (~377°C), indicating greater thermal
stability, with a larger ΔT=T<sub>c</sub>-T<sub>g</sub> = (109°C vs 80°C).
The lower crystallization intensity in NVP25 suggests
reduced or less extensive crystallization, aligning with
its increased stability. In contrast, NVP43 exhibits a
more pronounced crystallization due to its higher vana-
dium content. Overall, the higher T<sub>g</sub> and T<sub>c</sub> of NVP25
reflect a more thermally stable and rigid structure, while
NVP43 shows increased flexibility and crystallization.
The presence of a mixed valence for the vanadium (spe-
ciation of vanadium sites) in the synthesized NVP25
and NVP43 glasses was assessed in the synthesized
samples using XPS measurements. Figure 2 (Left)
presents the XPS spectra of the O 1s and V 2p core lev-
els for each sample. These spectra exhibit three distinct
features: (1) O 1s core level at about 531 eV; (2) V 2p<sub>3/2</sub>
at about 517.2 eV; and (3) V 2p<sub>1/2</sub> core level states at
about 524 eV [69]. Deconvolution of the V 2p<sub>3/2</sub> and V
2p<sub>1/2</sub> signals enables quantification of the relative con-
tent of V in different oxidation states, as reported in
Table 1. XPS analysis indicates that vanadium exists
primarily in the +5 oxidation state in both investigated
glasses, with a distinct presence of vanadium in the +4
oxidation state. Specifically, NVP25 and NVP43 show
V<sup>5+</sup>:V<sup>4+</sup> ratios of 95:5% and 87:13%, respectively. The
presence of vanadium in other oxidation states, while
not entirely ruled out, is minimal and negligible.

The V speciation is approached differently in CMD and BOMD simulations. In CMD, V speciation is predetermined by assigning potential parameters corresponding to V<sup>5+</sup> and V<sup>4+</sup> oxidation states[25, 48]. In BOMD the initial configuration specifies only the total number of V sites, while the final distribution is determined by analyzing the local V site spin topology of the converged ground-state energy configuration using PBE0 hybrid functional calculations, following a procedure validated for VP50 glass [38]. Converged spin localization at V sites enabled a classification based on the oxidation states in BOMD simulations. As previously reported for VP50 glass, PBE0-level calculations were necessary to obtain a well-defined spin distribution, facilitating clear distinction between V<sup>5+</sup> (σ<sub>spin</sub> ~ 0) and V<sup>4+</sup> (σ<sub>spin</sub> ~ 1) sites, as shown in Figure 2 (Right). Based on these σ<sub>spin</sub> values, BOMD calculations gives chemical compositions (in mol%) corresponding to 37.5Na<sub>2</sub>O - 25.0V<sub>x</sub>O<sub>y</sub>[22.2V<sub>2</sub>O<sub>5</sub> + 2.8V<sub>2</sub>O<sub>4</sub>] - 37.5P<sub>2</sub>O<sub>5</sub> and 28.5Na<sub>2</sub>O - 43.0V<sub>x</sub>O<sub>y</sub>[35.8V<sub>2</sub>O<sub>5</sub> +

7.2V₂O₄] - 28.5P₂O₅ for NVP25 and NVP43. Despite a slight increase in V⁴⁺ content compared to the experimental data, the BOMD V distribution aligns well with the experimental compositions. Notably, NVP43 exhibits a higher content of V⁴⁺ than NVP25 (see Table 1). Electronic density of states calculations using the PBE0 hybrid functional resulted in open bandgap values of 1.53±0.39 eV for the NVP25 model and 1.94±0.33 eV for the NVP43 model. These values are consistent with previous findings for VP50 [38], identifying the two NVP models as semiconducting glasses with large bandgaps, although slightly smaller than the 2.4–2.8 eV range observed for VP50 glass.

4. Reciprocal space properties: X-ray Total Struc-
ture Factors

Figure 3 shows the reciprocal space X-ray total struc-
ture factors for the NVP25 and NVP43 glasses. The pat-
terns discernible in the plots reflects not only the distinct
peaks but also the presence of a degree of intermediate-
range order, as indicated by a small shoulder around ~1
Å⁻¹ that was not observed in the VP50 glass [38]. This
feature implies additional structural features deserv-
ing further exploration through computational methods.
The X-ray scattering technique is found to be more sen-
sitive to the local arrangement of vanadium sites com-
pared to neutron scattering [38]. The BOMD simula-
tions align closely with the experimental data in terms
of peak positions and intensities, demonstrating their ef-
fectiveness in accurately capturing the glass structures
of NVP. In contrast, CMD simulations exhibit a sig-
nificant overestimation of the first peak intensity and a
pronounced peak at approximately ~2 Å⁻¹. BOMD, on
the other hand, shows less sharp peaks that are more
in line with the experimental data, although both meth-
ods underdescribe the shoulder at ~1 Å⁻¹, which is a
signature of intermediate-range order, with BOMD per-
forming slightly better in this regard. These disparities
emphasize the inherent limitations of the CMD method-
ology, particularly its struggle to quantitatively model
the short- and intermediate-range order of vanadium-
containing glasses, which is critical for understanding
the local environment of vanadium sites. The differ-
ences between CMD and BOMD results diminish at
higher *k* values (beyond approximately 4 Å⁻¹). This
supports the idea that CMD can serve as a useful tool
for preliminary structural modeling, while the consider-
ation of BOMD approaches is essential to obtain more
accurate and reliable structural descriptions of complex
oxide glasses like NVP25 and NVP43. This statement is
<figure><img src="image_4.png" /><figcaption>Figure 2: (Color online) Left: XPS spectra of NVP25 (top) and NVP43 (bottom) glasses in the O 1s and V 2p regions, revealing the presence of V$^{5+}$ and V$^{4+}$ oxidation states. Right: Spin densities distributions of V sites calculated at the PBE0 level ($\rho_{spin}$, in e units) for both NVP25 (top) and NVP43 (bottom) glass models by BOMD. Inset: Averaged values of the energy bandgap (in eV) calculated at the PBE0 level for both for both NVP25 and NVP43 glass models. Bandgap values for both spin-up and spin-down states are averaged over three configurations simulated for each system.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 3: (Color online) X-ray total structure factors for NVP25 and NVP43 glasses, comparing CMD, BOMD, and experimental data.</figcaption></figure>

further reinforced by the analysis of the $R_{\chi}^{S(k)}$ values [70] as defined by Eq. (1) and presented in Table 3.

$$R_{\chi} = 100 \times \left( \frac{\sum_{i} [\text{Exp}^{(i)} - \text{Sim}^{(i)}]^{2}}{\sum_{i} [\text{Exp}^{(i)}]^{2}} \right)^{1/2}, \quad (1)$$

$R_{\chi}^{S(k)}$ quantifies the level of agreement of the theoretical prediction with the experimental data, reminding that a

lower $R_{\chi}^{S(k)}$ values indicate a better agreement. BOMD yields $R_{\chi}^{S(k)}$ values of (10.1±0.9)% for NVP25 and (9.8±1.2)% for NVP43, whereas CMD results in significantly higher values of (14.4±0.2)% and (17.1±1.0)%, respectively (see Table 3, first row). The improved agreement with experimental data achieved by refining the CMD simulation with an equilibration at 300 K through BOMD is consistent with findings for VP50 glass reported in previous studies [38]. When comparing NVP25 and NVP43, notable differences in their structural features emerge. Although NVP25 contains a higher sodium oxide content, which typically suggests a more disordered structure, the peaks for NVP43 are more broadened. This indicates that NVP43 has a more complex structural environment, likely resulting from its higher vanadium content, which contributes to greater disorder within the glass structure.

# 5. Real space properties

## 5.1. Total pair distribution function

Figure 4 presents a comparative analysis of the reduced total X-ray pair distribution function $G(r)$ for NVP25 and NVP43 glasses. The experimental data, which exhibits greater complexity than those of the VP50 glass [38, 71] due in part to additional Na-O interactions, are compared with results from CMD and
BOMD simulations. The experimental $G(r)$ was derived via Fourier transform of the total structure factors, as outlined in Eq. (2).

$$G(r) = \frac{2}{\pi} \int_{0}^{k_{\text{max}}} k [S(k) - 1] \sin(kr) dk \quad (2)$$

This allows for direct comparison with our calculated results for the NVP25 and NVP43 glasses (Eq. (3)).

$$G(r) = 4\pi r \rho [g_{\text{tot}}(r) - 1], \quad (3)$$

where $\rho$ is the number density and $g_{\text{tot}}(r)$ is the unitless total pair correlation function ($g_{\text{tot}}(r) \to 1$ for $r \to \infty$). The positions of the first four peaks' maxima in the 1–4 Å interval are summarized in Table 4, comparing CMD and BOMD values with experimental data. The analysis reveals significant discrepancies between CMD predictions and experimental measurements, while BOMD results demonstrate markedly improved agreement across the examined $r$ range. In the short-range order region (1–3 Å), BOMD simulations accurately capture the experimental features for both glass compositions. The first peak, centered at approximately 1.5–1.6 Å, attributed to P–O and V=O bonds, is well-reproduced by BOMD in terms of position, intensity, and shape. In contrast, CMD simulations consistently overestimate this peak's intensity by about 20–30%, indicating a greater degree of short-range overstructuring, particularly pronounced for NVP43. Both CMD and BOMD shift the experimentally observed peak at about 1.96 Å to lower values of around 1.8 Å. Beyond the short-range order (3–5 Å), CMD shows a consistent shift of approximately 0.5 Å toward larger $r$ values compared to experimental data for both NVP25 and NVP43, indicating systematic errors in predicting second-neighbor distances. Conversely, BOMD results closely align with experimental peak positions, with deviations typically under 0.1 Å (Table 4). BOMD also better reproduces peak shapes and relative intensities more accurately in the 2–4 Å range, while CMD predicts shifted sharper, more intense peaks, particularly overestimating and right-shifting the fourth peak by up to 20–30% and 0.3 Å for both NVP models. The superior performance of BOMD in reproducing the experimental $G(r)$ underscores the critical role of explicitly accounting for electronic structure effects in modeling sodium vanadium phosphate glasses. Furthermore, the limitations of the CMD method are emphasized by the ($R_{\chi}^{G(r)}$) parameter values calculated within the 1–4 Å range (Table 3, second row). CMD consistently yields significantly higher $R_{\chi}^{G(r)}$ values compared to BOMD, indicating poorer agreement with experimental data.

For instance, the NVP43 system exhibits an $R_{\chi}^{G(r)}$ of 96.2 ± 3.1% for CMD versus 65 ± 4.8% for BOMD, partially due to the observed shift between experimental and calculated second and 4th peak's positions. This comparative analysis demonstrates that BOMD simulations provide a more accurate structural description of NVP25 and NVP43 glasses compared to CMD, particularly in the short and medium-range order regions significantly affecting the glasses' physical and chemical properties.

<figure><img src="image_6.png" /><figcaption>Figure 4: (Color online) Reduced total pair distribution function $G(r)$ for NVP25 and NVP43 glasses, illustrating the comparison between the calculated data (CMD and BOMD) and the experimental data in the range 1–4 Å range.</figcaption></figure>

<table><caption>Table 3: Comparison of the agreement between CMD, and BOMD simulations and experimental data using goodness-of-fit $R_{\chi}$ parameters for X-ray total structure factor $S(k)$ and reduced pair distribution function $G(r)$ of NVP25 and NVP43 glasses at 300 K.</caption><thead><tr><th rowspan="2">X-rays</th><th colspan="2">NVP25</th><th colspan="2">NVP43</th></tr><tr><th>CMD</th><th>BOMD</th><th>CMD</th><th>BOMD</th></tr></thead><tbody><tr><td>R<sup>S(k)</sup><sub>χ</sub></td><td>14.4 ± 0.2</td><td>10.1 ± 0.9</td><td>17.1 ± 1.0</td><td>9.8 ± 1.2</td></tr><tr><td>R<sup>G(r)</sup><sub>χ</sub></td><td>87.3 ± 2.3</td><td>62.9 ± 2.0</td><td>96.2 ± 3.1</td><td>65.0 ± 4.8</td></tr></tbody></table>

<table><caption>Table 4: Peak positions in the total pair distribution function: Comparison of calculated values (CMD and BOMD) with experimental measurements for NVP25 and NVP43 systems.</caption><thead><tr><th rowspan="2">Peak positions</th><th colspan="3">NVP25</th><th colspan="3">NVP43</th></tr><tr><th>Exp.</th><th>CMD</th><th>BO</th><th>Exp.</th><th>CMD</th><th>BO</th></tr></thead><tbody><tr><td>1<sup>st</sup> peak</td><td>1.51</td><td>1.50</td><td>1.50/1.62</td><td>1.59</td><td>1.51</td><td>1.59</td></tr><tr><td>2<sup>nd</sup> peak</td><td>1.96</td><td>1.80</td><td>1.84</td><td>1.96</td><td>1.79</td><td>1.78</td></tr><tr><td>3<sup>rd</sup> peak</td><td>2.51</td><td>2.48</td><td>2.56</td><td>2.57</td><td>2.46</td><td>2.58</td></tr><tr><td>4<sup>th</sup> peak</td><td>3.26</td><td>3.52</td><td>3.39</td><td>3.23</td><td>3.53</td><td>3.37</td></tr></tbody></table>
## 5.2. Partial pair correlation function

Among the partial pair correlation functions for VP50 glass, the V–O correlation showed the most significant improvements when moving from the initial CMD model to BOMD calculations [38]. In this work, the $g_{VO}(r)$ partial pair correlation functions (see Figure 5) for NVP25 and NVP43 glasses reveal further discrepancies between CMD and BOMD, particularly in the short-range order. BOMD shows a clear splitting of the first peak into two distinct maxima at 1.6 Å and 1.8–1.9 Å, reflecting a better description of both V=O and V–O bonds. CMD, on the other hand, produces an averaged broader peak at 1.7–1.8 Å, with only a small shoulder at 1.6 Å for NVP25, indicating an oversimplified representation of vanadium-oxygen interactions. Beyond 3 Å, both methods follow similar trends, although CMD consistently predicts slightly larger *r* values. The discrepancies are more pronounced in NVP43, where CMD struggles to accurately capture the increased complexity of Na-V-P-O interactions, showing its limitations with higher vanadium content. These results emphasize the fact that BOMD better captures vanadium’s local coordination environment. Having established the superior performance of BOMD in describing the V–O local environment, the following analysis of all partial pair correlation functions in this work is based on BOMD data rather than CMD ones. Turning to the other partial pair correlation functions, including the breakdown of the V–O correlation into contributions from $V^{5+}$–O and $V^{4+}$–O, Figure 6 presents the BOMD-derived functions for $g_{OO}(r)$, $g_{NaO}(r)$, $g_{V^{4+}O}(r)$, $g_{V^{5+}O}(r)$, and $g_{PO}(r)$ for NVP25 and NVP43 glasses. Table 5 reports the bond lengths $r_{a\beta}$, corresponding to the positions of the maxima in the partial pair correlation functions $g_{a\beta}(r)$. In NVP25, the P–O correlation shows a main peak at 1.50 Å and a shoulder at 1.62 Å. The $V^{5+}$–O correlations exhibit peaks at 1.62 and 1.84 Å, corresponding to V=O and V–O bonds, respectively. $V^{4+}$–O has a main peak at 1.89 Å and a weak shoulder around 1.62 Å, associated with the low concentration of $V^{4+}$. The combined contributions at 1.62 Å account for the splitting of the first peak in the total pair distribution function $G(r)$ (Figure 4). For NVP43, the P–O correlation shows a single peak at 1.53 Å. The $V^{5+}$–O correlation has a less-defined first peak at 1.62 Å (V=O), with a second peak at 1.75 Å, which is shorter than the corresponding peak in NVP25 (1.84 Å). No significant differences are observed in the Na$^{+}$–O and O–O correlations between the two systems. The shorter $V^{5+}$–O distances compared to $V^{4+}$–O, as shown in Table 5, are consistent with a well-known feature in vanadophosphate crystalline phases: the increase in V–

O single bond length with decreasing vanadium oxidation state [46, 47]. Both $V^{5+}$ and $V^{4+}$ also exhibit at least one characteristic vanadyl bond. BOMD calculations successfully capture this nontrivial bonding pattern of $VO_n$ polyhedra, showing good agreement with experimental trends. Specifically, BOMD reproduces short vanadyl bonds alongside longer single V–O bonds, with bond lengths corresponding to oxidation state as follows: $V^{5+}$–O < $V^{4+}$–O. The corresponding distances are 1.75–1.84 Å and 1.89–1.90 Å for NVP25 and NVP43, respectively. These results improve upon previous findings based solely on CMD simulations with empirical force fields, as reported for LiVP systems [6], highlighting the critical need for first-principles approaches to accurately describe such complex V-containing glass systems.

The data for NVP25 and NVP43 closely align, demonstrating that although the vanadium content is higher in NVP43, the local bonding environments around vanadium and sodium remain similar across both systems.

<figure><img src="image_7.png" /><figcaption>Figure 5: (Color online) Partial pair correlation function $g_{VO}(r)$ showing results from CMD and BOMD. On the right side, a zoom in on the first peak between 1 and 3 Å.</figcaption></figure>

## 5.3. Coordination number and structural units

In this section, the average coordination numbers for Na, V, P, and O (BOMD data) are compared with experimental values (see Table 5). Phosphorus has a well-defined average coordination number of 4.0 for both systems, which aligns very well with experimental data. The total coordinations for vanadium predicted by BOMD are 4.76 and 4.75, falling within the experimental range of 4.20–5.40. For sodium, the calculated average coordination number is slightly higher than the experimental reference in both systems. The running average coordination number for cation-oxygen pairs, shown in Figure S1 (Supplementary Information), reveals distinct structural roles within the glass network. The Na$^{+}$–O coordination environment does not exhibit a well-defined plateau, differing from other
<figure><img src="image_8.png" /><figcaption>Figure 6: Partial pair correlation functions for NVP25 and NVP43 glasses, displaying CMD and BOMD results for g<sub>OO</sub>(r), g<sub>NaO</sub>(r), and g<sub>PO</sub>(r), respectively, from left to right.</figcaption></figure>

<table><caption>Table 5: Average coordination numbers $n_i$ bond lengths $r_{a\beta}$ (in Å) (taken as the position of the first maximum of the pair correlation functions $g_{a\beta}(r)$), calculated for the NVP25 and NVP43 models and averaged over the three replicas of BOMD data. We report the distances obtained from the first maximum of $g_{a\beta}(r)$. The average statistical error across the replicas is 0.02 Å.</caption><thead><tr><th></th><th>Exp.<sup>a</sup></th><th>NVP25</th><th>NVP43</th></tr></thead><tbody><tr><td>n<sub>Na</sub></td><td>4.0-6.0</td><td>6.58</td><td>6.81</td></tr><tr><td>n<sub>P</sub></td><td>4.00</td><td>4.00</td><td>4.00</td></tr><tr><td>n<sub>V</sub></td><td>4.20-5.40</td><td>4.76</td><td>4.75</td></tr><tr><td>n<sub>V<sup>5+</sup></sub></td><td>-</td><td>4.69</td><td>4.62</td></tr><tr><td>n<sub>V<sup>4+</sup></sub></td><td>-</td><td>5.37</td><td>5.38</td></tr><tr><td>n<sub>O</sub></td><td>2.00-4.00</td><td>1.71</td><td>1.79</td></tr><tr><td>r<sub>NaO</sub></td><td>2.30-2.60</td><td>2.33</td><td>2.34</td></tr><tr><td>r<sub>PO</sub></td><td>1.52-1.60</td><td>1.50/1.62</td><td>1.53</td></tr><tr><td>r<sub>VO</sub></td><td>1.58</td><td>1.62</td><td>1.62</td></tr><tr><td></td><td>1.90</td><td>1.84</td><td>1.78</td></tr><tr><td>r<sub>V<sup>5+</sup>O</sub></td><td>1.59-1.78<sup>b</sup></td><td>1.62</td><td>1.62</td></tr><tr><td></td><td>1.80-2.35<sup>b</sup></td><td>1.84</td><td>1.75</td></tr><tr><td>r<sub>V<sup>4+</sup>O</sub></td><td>1.59-1.74<sup>b</sup></td><td>1.61</td><td>1.64</td></tr><tr><td></td><td>1.80-2.12<sup>b</sup></td><td>1.89</td><td>1.90</td></tr><tr><td>r<sub>OO</sub></td><td>2.52-2.79</td><td>2.56</td><td>2.58</td></tr><tr><td>r<sub>VV</sub></td><td>3.50</td><td>3.50</td><td>3.47</td></tr><tr><td>r<sub>PV</sub></td><td>3.24</td><td>3.20</td><td>3.29</td></tr></tbody></table>

cations. Na<sup>+</sup> features a high average coordination number in the 6–7 range and a bond length of approximately 2.33 Å, interacting with neighboring oxygens primarily through ionic bonds. In contrast, V<sup>4+</sup>, V<sup>5+</sup>, and P interact with neighboring bridging oxygens pre-

dominantly via covalent bonds [38]. This is in favor of the glass-modifier role of Na<sup>+</sup>, contrasting with the glass-former roles of P and V. Previous studies on amorphous systems containing V, including both pure vanadium oxides and vanadophosphates containing alkali ions [26, 45, 46, 71, 75–78], have reported vanadium coordination numbers ranging from 4 to 6. This indicates that vanadium can act as either a glass network former (i.e. low coordination) or a modifier (i.e. high coordination), depending on the glass composition and in accordance with Zachariasen's model [79]. In the present work for both NVP models, the coordination number of V<sup>4+</sup> is higher than the one of V<sup>5+</sup> and phosphorus, yet the V<sup>4+</sup>-O bond distance remains comparable to V<sup>5+</sup>-O, while being significantly shorter than Na<sup>+</sup>-O bonds. This confirms that Na<sup>+</sup> is a network modifier with high coordination, while V<sup>4+</sup> behaves as an intermediate-to-weak network former, with an average coordination between 5 and 6. In contrast, V<sup>5+</sup> predominantly acts as a strong-to-intermediate glass network former (average coordination of 4 to 5), in part similar to phosphorus (coordination of 4) within the glass structure. To further obtain insights in the local environment of V and Na, we examined the distribution of VO<sub>n</sub> and Na<sup>+</sup>O<sub>n</sub> structural units. Table 6 gives the BOMD-derived values for both systems. The VO<sub>n</sub> units referred to total V are predominantly VO<sub>5</sub> (approximately 50%), followed by VO<sub>4</sub> (about 35%), with a notable presence of VO<sub>6</sub> units (10–13%). These results align with findings from the VP50 system reported in [38]. A breakdown of the vanadium oxidation states in VO<sub>n</sub> units indicates that V<sup>5+</sup> acts as a strong glass network former, predominantly existing as VO<sub>4</sub> (40–42%) and VO<sub>5</sub> (51–54%) units, while V<sup>4+</sup> displays an intermediate behavior, primarily existing as VO<sub>5</sub> (49–57%) and VO<sub>6</sub> (41–44%) units.

This perspective contrasts with earlier findings, based mainly on CMD simulations, pointing to both V<sup>5+</sup> and V<sup>4+</sup> as intermediate former or modifier based on vanadium concentration. Previously, V<sup>4+</sup> was thought to exhibit low coordination (four-fold) [26, 41, 80]. However, BOMD provides an entirely different perspective: both V<sup>5+</sup> and V<sup>4+</sup> form structural polyhedra that contribute to the glass network. V<sup>5+</sup> exhibits low to intermediate coordination (4–5), reinforcing its role as a strong-to-intermediate glass network former. On the contrary V<sup>4+</sup> primarily shows higher coordination (5–6), indicating an intermediate-to-weak role, with both vanadium oxidation states displaying similar bond distances with oxygen. For Na<sup>+</sup>O<sub>n</sub>, the coordination numbers for Na<sup>+</sup> fluctuate between 4 and 9 across both systems, with Na<sup>+</sup>O<sub>6</sub> and Na<sup>+</sup>O<sub>7</sub> identified as the

<footnote><sup>a</sup>Estimated from X-ray and neutron diffraction measurements on V<sub>x</sub>O-P<sub>2</sub>O<sub>5</sub> [72, 73] and Na<sub>2</sub>O-V<sub>x</sub>O<sub>y</sub> [74] glasses. <sup>b</sup>Characteristic bond lengths intervals analysed from about 520 vanadophosphate crystalline phases [8, 45, 46].</footnote>
dominant structural units.

<table><caption>Table 6: Distribution of the individual $n_α(l)$ structural units where an atom $α$ (V or Na) is *l*-fold coordinated to oxygen atom, computed for the both NVP25 and NVP43 at BOMD level. These values have been calculated including neighbours within the first minimum in $g_{αβ}(r)$. The average is taken over the three replicas with given error bars. The cutoffs of 2.4 and 3.2 Å were used for the V–O and Na–O bonds respectively. Only fraction larger than 0.1 % are reported.</caption><thead><tr><th></th><th></th><th>NVP25</th><th>NVP43</th></tr></thead><tbody><tr><td>V<sup>tot</sup></td><td></td><td></td><td></td></tr><tr><td>l = 4</td><td>O<sub>4</sub></td><td>36.5 ± 7.5</td><td>35.3 ± 1.9</td></tr><tr><td>l = 5</td><td>O<sub>5</sub></td><td>50.6 ± 8.4</td><td>54.1 ± 5.6</td></tr><tr><td>l = 6</td><td>O<sub>6</sub></td><td>12.9 ± 0.9</td><td>10.6 ± 3.7</td></tr><tr><td>V<sup>5+</sup></td><td></td><td></td><td></td></tr><tr><td>l = 4</td><td>O<sub>4</sub></td><td>40.1 ± 7.2</td><td>42.1 ± 2.6</td></tr><tr><td>l = 5</td><td>O<sub>5</sub></td><td>50.8 ± 8.6</td><td>53.6 ± 3.5</td></tr><tr><td>l = 6</td><td>O<sub>6</sub></td><td>9.1 ± 1.4</td><td>4.4 ± 1.0</td></tr><tr><td>V<sup>4+</sup></td><td></td><td></td><td></td></tr><tr><td>l = 4</td><td>O<sub>4</sub></td><td>7.3 ± 10.3</td><td>2.7 ± 1.8</td></tr><tr><td>l = 5</td><td>O<sub>5</sub></td><td>48.6 ± 7.1</td><td>56.7 ± 15.9</td></tr><tr><td>l = 6</td><td>O<sub>6</sub></td><td>44.0 ± 3.8</td><td>40.7 ± 17.7</td></tr><tr><td>Na<sup>+</sup></td><td></td><td></td><td></td></tr><tr><td>l = 5</td><td>O<sub>5</sub></td><td>12.7 ± 3.7</td><td>8.5 ± 0.2</td></tr><tr><td>l = 6</td><td>O<sub>6</sub></td><td>33.2 ± 1.9</td><td>31.7 ± 2.7</td></tr><tr><td>l = 7</td><td>O<sub>7</sub></td><td>36.6 ± 5.3</td><td>34.8 ± 2.2</td></tr><tr><td>l = 8</td><td>O<sub>8</sub></td><td>14.7 ± 0.1</td><td>17.2 ± 4.8</td></tr><tr><td>l = 9</td><td>O<sub>9</sub></td><td>1.9 ± 0.2</td><td>5.6 ± 1.7</td></tr></tbody></table>

## 6. Bond angle distributions, network connectivity, and Na-glass former association preference

Figure 7 shows the bond angle distributions (BAD) for NVP25 and NVP43 glass systems, centered around various oxygen-containing triads. The V–O–V, V–O–P, and P–O–P distributions provide insights into the connectivity between different polyhedra in the glass network. The V–O–P distribution shows a broad main peak centered around 130–140° with a sharp, less intense peak centered at 90°, indicating a variable form of interconnectivity between the vanadium and phosphorus polyhedra. This broad distribution is indicative of a range of V–O–P bond angles representative of the structural flexibility in the glass network. The slight differences between NVP25 and NVP43 in the main peak in this distribution is a sign that the increase in vanadium content affects the V–P connectivity to some extent. The V–O–V BAD displays a broader distribution compared to V–O–P, with a main peak around 140°. This broadness reflects the variety of V–O–V connections in the glass network, including edge-sharing and corner-sharing between vanadium polyhedra. The similarity between NVP25 and NVP43 in this distribution indicates that

the V–V connectivity is not significantly affected by the composition change. The P–O–P distribution shows a sharp peak at about 125°, which is characteristic of P–O–P linkages in phosphate glasses. The sharpness of this peak, compared to the V–O–V and V–O–P distributions stand for a more rigid and well-defined local structure around phosphorus atoms. The O–P–O BAD has a sharp peak centered at approximately 109°, indicative of the tetrahedral geometry of PO₄ units and persisting across NVP25 and NVP43, as shown by the overlapping peaks. The O–V<sup>5+</sup>–O BAD displays a broad distribution ranging from approximately 70° to 120°, with a notable concentration around 108°. This proves the presence of different geometric configurations among the polyhedral units. Notably, the dominance of tetrahedral angles highlights a substantial contribution of V<sup>5+</sup>O₄ tetrahedra to the structural framework of the glass network. In contrast, the O–V<sup>4+</sup>–O distribution shows a primary peak around 90°, a clear sign of the presence of octahedral and pyramidal geometries within the glass network. A secondary, smaller peak is observed near 170°, rather than at 180°, indicating distorted octahedra and pyramids. The absence of a significant peak around 109° implies a low proportion of V<sup>4+</sup>O₄ tetrahedra, which is consistent with the V<sup>4+</sup>O₆ structural unit analysis.

The network connectivity in NVP25 and NVP43 glasses reveals a complex interplay between phosphate and vanadium polyhedra, as demonstrated by $Q^n$ distribution analysis. $Q^n$ refers to a structural unit where n bridging oxygen (BO) atoms are directly linked to a network-forming ion, such as phosphorus in phosphate glasses. This analysis provides valuable insights into the degree of polymerization and the structural organization of these glass systems [11, 26]. Here, we come back to an explicit comparison between CMD and BOMD data, with the intent of highlighting the most significant differences. One of the most striking improvements brought by BOMD is in the description of the phosphate network connectivity. For NVP25, both CMD and BOMD agree on the predominance of $Q^2$ species, but BOMD shows a higher proportion of $Q^3$ units. This suggests a more cross-linked phosphate network and a greater degree of three-dimensional connectivity than was evident from CMD alone. The differences are even more pronounced for NVP43. While CMD suggests an almost equal distribution of $Q^1$, $Q^2$, and $Q^3$ species, BOMD indicates a clear predominance of $Q^2$ and $Q^3$ units. This suggests that the phosphate network in NVP43 maintains a higher degree of polymerization than previously assumed from CMD, and contradicts earlier views that increasing vanadium content would significantly depolymerize the phosphate net-
<figure><img src="image_9.png" /><figcaption>Figure 7: (Color online) Bond angle distributions for NVP25 (green line) and NVP43 (red line) models. The left panel shows angles centered around oxygen: P–O–P, V–O–P, and V–O–V (top to bottom). The right panel displays angles centered around network formers: O–P–O, O–V<sup>5+</sup>–O, and O–V<sup>4+</sup>–O (top to bottom).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 8: (Color online) Comparison of the $Q^n$ distribution, where n represents the number of bridging oxygen atoms directly linked to a network former P, between the CMD and BOMD simulations for NVP25 and NVP43. The data presented is averaged over three replicas, with the standard error bar indicated.</figcaption></figure>

work. A comparison with previously reported data for VP50 glass [38], where BOMD models show a majority of $Q^4$ (55–60%) and a notable fraction of $Q^3$ (37–

<figure><img src="image_11.png" /><figcaption>Figure 9: (Color online) Comparison of the $Q^n$ distribution, where n represents the number of non-bridging oxygen ([NB]) atoms directly linked to a network former $V^{4+}$ and $V^{5+}$, between the CMD and BOMD simulations for NVP25 and NVP43. The data presented is averaged over three replicas, with the standard error bar indicated.</figcaption></figure>

42%) units, further supports the role of Na ions as network modifiers, reducing the extent and interconnection of the phosphate network.

Turning to the vanadium network, BOMD provides a substantially different picture in terms of $Q''_{[NB]}$ (with n indicating non-bridging oxygen ([NB]) in this case)
compared to CMD, particularly for V⁵⁺ species. In both NVP25 and NVP43, CMD predicts a dominance of Q⁰[NB] environments for V⁵⁺, with a largely interconnected nature of such units. However, BOMD reveals a much more diversified distribution of V⁵⁺ environments, with significant proportions of Q¹[NB] and Q²[NB] species. This BOMD result indicates that V⁵⁺ plays a more active role in network formation than previously understood, potentially contributing to the overall complex connectivity of the glass structure. For V⁴⁺ species, while both methods agree on a preference for Q⁰[NB] environments, BOMD is in favor of a more uniform distribution across different Qⁿ[NB] species, especially in NVP43.

To evaluate the preference of Na⁺ ions for different network-forming species (P, V⁵⁺, V⁴⁺), we calculated the normalized preference ratio r⁰Na⁺/β [81] defined as:

$$r_{\alpha/\beta}^{Na^{+}} = \frac{N_{\beta} \frac{CN_{Na^{+}-\alpha}}{CN_{Na^{+}-\beta}}}{N_{\alpha}} \quad (4)$$

where N_α is the number of glass-former species α and CN⁰Na⁺-α is the partial coordination number of Na⁺ around species α.

The calculated preference ratios are given in Table 7. The unit ratio reflects the statistical distribution of Na⁺ ions around the α and β species, which is strictly governed by their relative concentrations. Values less than 1 indicate a preferential coordination of Na⁺ ions by species β, whereas values larger than 1 suggest a stronger affinity of species α for coordinating Na⁺ ions. The results show that Na⁺ has a stronger affinity for P compared to V in both systems ($r⁰Na⁺/V⁵⁺ < 1$). There is a slight preference for V⁴⁺ over V⁵⁺ ($r⁰Na⁺/V⁴⁺/V⁵⁺ > 1$), although both V species are less favored than P. The preference for P over V is more pronounced in NVP25 compared to NVP43, as indicated by the lower $r⁰Na⁺/V⁵⁺$ value for NVP25.

<table><caption>Table 7: Sodium-glass former connection preference ratio for NVP25 and NVP43</caption><thead><tr><th></th><th>R<sup>Na⁺</sup><sub>V⁵⁺/P</sub></th><th>R<sup>Na⁺</sup><sub>V⁴⁺/P</sub></th><th>R<sup>Na⁺</sup><sub>V⁴⁺/P</sub></th><th>R<sup>Na⁺</sup><sub>V⁴⁺/V⁵⁺</sub></th></tr></thead><tbody><tr><td>NVP25</td><td>0.60</td><td>0.60</td><td>0.65</td><td>1.09</td></tr><tr><td>NVP43</td><td>0.66</td><td>0.62</td><td>0.82</td><td>1.31</td></tr></tbody></table>

# 7. Discussion

## 7.1. General remarks

This study focussed on two sodium vanadium phosphate (NVP) glass compositions, (Na₂O)α- (VₓOᵧ)(1-2α)- (P₂O₅)α (α = 0.375, 0.285), using a

combined approach of experimental characterization and atomic-scale modeling with both classical and Born-Oppenheimer molecular dynamics (CMD and BOMD, respectively). This approach provides a comprehensive understanding of the structural and bonding properties of these glasses. We underscored the limitations of classical MD and currently available empirical force fields, which, while valuable tools, fail to accurately capture the local structure, particularly regarding the complex interactions involving vanadium oxides. By explicitly incorporating electronic structure effects, BOMD simulations demonstrate a significantly enhanced ability to improve upon original CMD models. Despite the shortness of the dynamical equilibration at 300 K using BOMD, one can obtain structural features of NVP glasses of unprecedented quality, particularly for short- and medium-range order distances, as demonstrated by their close agreement with experimental X-ray diffraction data and pair distribution function analysis. One of the key findings is the important role of phosphorus in enhancing the rigidity of the NVP glass network. Sodium typically acts as a network modifier, disrupting the glass network, and potentially decreasing rigidity. However, NVP25 (with a higher content of sodium oxide) exhibits a higher glass transition temperature compared to NVP43. This is a key observation to understand the dominating influence of phosphorus, a stronger glass network former than vanadium, in dictating the overall rigidity of the glass structure. BOMD provides valuable insight into the speciation of vanadium sites within the glass structure. The analysis reveals that vanadium exists predominantly in the +5 oxidation state in both glasses, with a smaller but not negligible presence of vanadium in the +4 oxidation state. This is consistent with experimental XPS measurements and highlights the ability of BOMD to accurately capture the oxidation states of transition metals in these glassy systems.

## 7.2. Revisiting NVP glasses structure: short- and intermediate range-orders

A closer examination of the short-range and intermediate-range orders in NVP glasses, using BOMD-derived data, sheds light on the interplay between phosphate and vanadium polyhedra within the glass network. The analysis of partial pair correlation functions, coordination numbers, bond angle distributions, and Qⁿ distributions allows obtaining a more precise understanding of the structural features underlying the properties of these glasses. The splitting of the V-O peak in the pair distribution function, revealed by BOMD simulations but not by CMD, reflects
<figure><img src="image_12.png" /><figcaption>Figure 10: Atomistic models of NVP25 (left) and NVP43 (right) glasses simulated using BOMD. Color legend (online): Na (yellow spheres), V<sup>5+</sup>O<sub>n</sub> polyhedra (blue), V<sup>4+</sup>O<sub>n</sub> polyhedra (light blue), PO<sub>4</sub> polyhedra (orange), O (red spheres).</figcaption></figure>

the presence of both V=O and V–O bonds, indicat-
ing a more complex local environment around vana-
dium atoms than the one accounted for by classical
force fields. Additionally, BOMD simulations shows
a distinct difference in the average coordination num-
bers of V<sup>4+</sup> and V<sup>5+</sup>, with V<sup>4+</sup> exhibiting higher aver-
age coordination and an intermediate-to-weak network
former character, while V<sup>5+</sup>, with lower average coordi-
nation, acts predominantly as a strong-to-intermediate
network former. The analysis of bond angle distribu-
tions provides further evidence for the structural com-
plexity of NVP glasses. The broad distributions ob-
served for V–O–V and V–O–P angles, compared to
the sharper peak for O–P–O, are indicative of a greater
degree of flexibility and a range of possible configu-
rations for vanadium polyhedra within the glass net-
work. This flexibility, coupled with the varying roles
of V<sup>4+</sup> and V<sup>5+</sup>, contributes to the unique properties
of NVP glasses. Finally, the Q<sup>n</sup> distribution analysis,
particularly the differences observed between CMD and
BOMD predictions, highlights the importance of ac-
curate modeling in understanding network connectiv-
ity. Figure 10 shows the degree of complexity of the
glass networks in NVP25 and NVP43, with the inter-
connected phosphate and vanadate polyhedra, as well
as the open structure facilitated by the incorporation of
Na ions. These BOMD insights allow us to revisit our
understanding of NVP glasses in several key points: i)
The phosphate network appears to maintain a higher de-
gree of polymerization than previously thought, even at
higher vanadium contents. This could explain the ob-
served stability and favorable properties of these glasses

in various applications [5, 6, 11, 15, 82]; ii) V<sup>5+</sup> emerges
as a more active network former, potentially creating
additional linkages in the glass structure; iii) BOMD re-
veals a more diverse range of local environments, espe-
cially for vanadium species. This structural heterogene-
ity is crucial for understanding the glass’s behavior un-
der various conditions, such as its response to ion inser-
tion/extraction in battery applications.; iv) The BOMD
results provide a more complete view of how changing
vanadium content affects the glass structure, leading to
the conclusions that the relationship between composi-
tion and network connectivity is more complex than pre-
viously understood.

It remains true that certain aspects still require further
efforts, such as extending the thermal relaxation at the
BOMD level and enhancing the comparison with exper-
imental data on the total structure factor at the level of
the intermediate range order. This feature is still poorly
described by both CMD (to a greater extent) and BOMD
(to a lesser extent). We also remark that the signifi-
cant enhancement in the local V–O bonding environ-
ment achieved with BOMD enables a better quantifica-
tion of oxygen speciation, including bridging oxygens
(BO) and non-bridging oxygens (NBO) within the glass
network. Understanding these factors is crucial as they
significantly impact sodium ion mobility [80, 83–85].
Achieving such improved insight into network connec-
tivity will have important implications for interpreting
the macroscopic properties of NVP glasses, including
thermal stability, conductivity, and ion transport.
## 8. Conclusions

This study provides a comprehensive analysis of two sodium vanadium phosphate (NVP) glasses by combining experimental characterization with atomic-scale modelling. It emphasizes the importance of completing a classical molecular dynamics trajectory with the use of first-principles molecular dynamics (of the Born-Oppenheimer kind in this case) to describe quantitatively the structural and bonding features of these glassy systems. This is crucial even at the relatively low temperature (T= 300 K), proving that the account of electronic structure effects is able to impact the network well beyond a mere structural rearrangement. Although classical molecular dynamics using the Bertani *et al.* force field [25, 48] provides a valuable ground for structural modeling improving upon other available force fields [38], it is unable to describe correctly the local environment around vanadium atoms, particularly the V-O bonding. The present findings reveal that NVP25, which has a higher phosphate content, exhibits a higher glass transition temperature compared to NVP43, underscoring the significant role of phosphorus as a strong network former in enhancing the rigidity of the glass structure. BOMD simulations also provide a clearer picture of vanadium speciation, accurately identifying the presence of both V⁴⁺ and V⁵⁺ oxidation states, along with the splitting of V-O bonds. These details are beyond the reach of CMD. V⁵⁺ primarily works as a strong-to-intermediate glass network former, while V⁴⁺ behaves more like a intermediate-to-weak network former fully in line with a number experimental results [43, 44, 86] but in contrast with previous models based on X-ray absorption near edge structure (XANES) [41, 42, 87] data or CMD calculations [6, 26, 41]. BOMD reveals a complex interplay of short- and intermediate-range order within the glass network, indicating a rich local environment around vanadium atoms. Further investigation into the roles of bridging and non-bridging oxygens could provide deeper insights into their influence on ionic conductivity in NVP glasses.

## Acknowledgement

The authors acknowledge High Performance Computing Centers support from GENCI (Grand Equipment National de Calcul Intensif; grants n. A0160807670, A0160912441, A0160910832, and A0160913426), the University of Strasbourg (Pole HPC Équipe@Meso (Program Investissements d'Avenir) and the CPER

Alsac calcul/Big Data), and the région Nouvelle Aquitaine via the CANaMIAS project AAPR2021-2020-11779110. We acknowledge financial support from the Agence Nationale de la Recherche (ANR) within the framework of the project AMSES (ANR-20-CE08-0021) and the EUR-QMAT-ANR under the contract ANR-11-LABX-0058-NIE and ANR-17-EURE-0024 within the Investissement d'Avenir program ANR-10-IDEX-0002-02.

## Data availability

NVP models obtained by BOMD as well as XRD experimental data are available at the European Center of Excellence Novel Materials Discovery (CoE-NOMAD) repository [88].

## References

[1] G. Linsley, A. Owen, F. Hayatee, Electronic conduction in vanadium phosphate glasses, J. Non-Cryst. Solids 4 (1970) 208-219.

[2] G. Khattak, A. Mekki, L. Wenger, X-ray photoelectron spectroscopy (XPS) and magnetic susceptibility studies of vanadium phosphate glasses, J. Non-Cryst. Solids 355 (2009) 2148-2155.

[3] J. Kjeldsen, A. Rodrigues, S. Mossin, Y. Yue, Critical V₂O₅/TeO₂ ratio inducing abrupt property changes in vanadium tellurite glasses, J. Phys. Chem. B 118 (2014) 14942-14948.

[4] L. McDonald, C. Siligardi, M. Vacchi, A. Zieser, M. Affatigato, Tellurium vanadate glasses: V⁴⁺ colorimetric measure and its effect on conductivity, Front. Mater. 7 (2020) 103.

[5] C. Mugoni, R. Rosa, R. Giovanardi, M. Affatigato, M. Gualtieri, C. Siligardi, S. Andronenko, S. Misra, Synthesis and characterization of (68-x) CuO-xV₂O₅-32TeO₂ (x=0-68 mol%) and (35-x) CuO-xV₂O₅-65TeO₂ (x=0-35 mol%) glasses: Conduction mechanism, structure and epr study, Mater. Chem. Phys. 266 (2021) 124488.

[6] G. Broglia, C. Mugoni, J. Du, C. Siligardi, M. Montorsi, Lithium vanado-phosphate glasses: Structure and dynamics properties studied by molecular dynamics simulations, J. Non-Cryst. Solids 403 (2014) 53-61.

[7] F. Kong, X. Liang, L. Yi, X. Fang, Z. Yin, Y. Wang, R. Zhang, L. Liu, Q. Chen, M. Li, et al., Multi-electron reactions for the synthesis of a vanadium-based amorphous material as lithium-ion battery cathode with high specific capacity, Energy 219 (2021) 119513.

[8] E. Boivin, J. Chotard, C. Masquelier, L. Croguennec, Towards reversible high-voltage multi-electron reactions in alkali-ion batteries using vanadium phosphate positive electrode materials, Molecules 26 (2021) 1428.

[9] R. Giovanardi, M. Montorsi, G. Ori, J. Cho, T. Subhani, A. Boccaccini, C. Siligardi, Microstructural characterisation and electrical properties of multiwalled carbon nanotubes/glass-ceramic nanocomposites, J. Mater. Chem. 20 (2010) 308-313.

[10] K. Nagamine, T. Honma, T. Komatsu, Selective synthesis of Lithium Ion-Conductive β-LiVOPO₄ crystals via Glass-Ceramic Processing, J. Am. Ceram. Soc. 91 (2008) 3920-3925.

[11] R. Brow, C. Click, T. Alam, Modifier coordination and phosphate glass networks, J. Non-Cryst. Solids 274 (2000) 9-16.
[12] S. Ibrahim, M. Marzouk, G. El Komy, Structural characteristics and electrical conductivity of vanadium-doped lithium ultraphosphate glasses, Silicon 9 (2017) 403–410.

[13] T. Pietrzak, E. Jerzye, M. Wasiucionek, J. Nowiński, Nanocrystallisation in vanadate phosphate and lithium iron vanadate phosphate glasses, Phys. Chem. Glasses-Eur. J. Glass Sci. Technol. B 57 (2016) 113–124.

[14] C. Delmas, Sodium and sodium-ion batteries: 50 years of research, Adv. Energy Mater. 8 (2018) 1703137.

[15] J.-Y. Hwang, S.-T. Myung, Y.-K. Sun, Sodium-ion batteries: present and future, Chem. Soc. Rev. 46 (2017) 3529–3614.

[16] M. Slater, D. Kim, E. Lee, C. Johnson, Sodium-ion batteries, Adv. Funct. Mater. 23 (2013) 947–958.

[17] H. Pan, Y.-S. Hu, L. Chen, Room-temperature stationary sodium-ion batteries for large-scale electric energy storage, Energy Environ. Sci. 6 (2013) 2338–2360.

[18] T. Liu, Y. Zhang, Z. Jiang, X. Zeng, J. Ji, Z. Li, X. Gao, M. Sun, Z. Lin, M. Ling, et al., Exploring competitive features of stationary sodium ion batteries for electrochemical energy storage, Energy Environ. Sci. 12 (2019) 1512–1533.

[19] M. El-Desoky, N. Wally, E. Sheha, B. Kamal, Impact of sodium oxide, sulfide, and fluoride-doped vanadium phosphate glasses on the thermoelectric power and electrical properties: structure analysis and conduction mechanism, J. Mater. Sci.: Mater. Electron. 32 (2021) 3699–3712.

[20] F. H. ElBatal, M. A. Marzouk, A. M. Abdelghany, UV-visible and infrared absorption spectra of gamma irradiated V2O5-doped in sodium phosphate, lead phosphate, zinc phosphate glasses: A comparative study, J. Non Cryst. Solids 357 (2011) 1027–1036.

[21] M. Saad, W. Stambouli, N. Sdiri, H. Elhouichet, Effect of mixed sodium and vanadium on the electric and dielectric properties of zinc phosphate glass, Mater. Res. Bull. 89 (2017) 224–231.

[22] A. H. Hammad, A. M. Abdelghany, Optical and structural investigations of zinc phosphate glasses containing vanadium ions, J. Non Cryst. Solids 433 (2016) 14–19.

[23] E. E. Assem, I. Elmehasseb, Structure, magnetic, and electrical studies on vanadium phosphate glasses containing different oxides, J. Mater. Sci. 46 (2011) 2071–2076.

[24] N. Sehrawat, P. Sharma, M. Bala, A. Ohlan, S. Dahiya, R. Punia, A. Maan, Deciphering the thermal, physical, structural, and optical characteristics of sodium-doped vanadophosphate glasses, Opt. Mater. 155 (2024) 115811.

[25] G. Ori, M. Montorsi, A. Pedone, C. Siligardi, Insight into the structure of vanadium containing glasses: a molecular dynamics study, J. Non-Cryst. Solids 357 (2011) 2571–2579.

[26] M. Montorsi, G. Broglia, C. Mugoni, Structural insight into transition metal oxide containing glasses by molecular dynamic simulations, Molecular Dynamics Simulations of Disordered Materials: From Network Glasses to Phase-Change Memory Alloys (2015) 181–213.

[27] J. Du, M. Montorsi, S. Barbi, X. Lu, Rare earth and transition metal containing glasses, Atomistic Simulations of Glasses: Fundamentals and Applications (2022) 367–438.

[28] S. Marijan, M. Razum, T. Klaser, P. Mošner, L. Koudelka, Ž. Skoko, J. Pisk, L. Pavić, Tailoring structure for improved sodium mobility and electrical properties in V2O5-Nb2O5-P2O5 Glass(es)-(Ceramics), J. Phys. Chem. Solids 181 (2023) 111461.

[29] K. Nagamine, T. Honma, T. Komatsu, A fast synthesis of Li₃V₂(PO₄)₃ crystals via glass-ceramic processing and their battery performance, Journal of Power Sources 196 (2011) 9618–9624.

[30] G. Delaizir, V. Seznec, P. Rozier, C. Surcin, P. Salles, M. Dollé, Electrochemical performances of vitreous materials in the sys-

tem Li₂O-V₂O₅-P₂O₅ as electrode for lithium batteries, Solid State Ionics 237 (2013) 22–27.

[31] T. Kühne, M. Iannuzzi, M. Del Ben, V. Rybkin, P. Seewald, F. Stein, T. Laino, J. Hutter, et al., CP2K: An electronic structure and molecular dynamics software package - quickstep: Efficient and accurate electronic structure calculations, J. Chem. Phys. 152 (2020) 194103.

[32] S. Bauer, P. Benner, T. Bereau, V. Blum, M. Boley, C. Carbogno, C. R. Catlow, G. Dehm, S. Eibl, R. Ernstorfer, Á. Fekete, M. Scheffler, et al., Roadmap on data-centric materials science, Model. Simul. Mat. Sci. Eng. 32 (2024) 063301.

[33] F. Mouvet, J. Villard, V. Bolnykh, U. Rothlisberger, Recent advances in first-principles based molecular dynamics, Acc. Chem. Res. 55 (2022) 221–230.

[34] V. Gavini, S. Baroni, V. Blum, D. Bowler, A. Buccheri, J. Chełlikowsky, S. Das, et al., Roadmap on electronic structure codes in the exascale era, Model. Simul. Mater. Sci. Eng. 31 (2023) 063301.

[35] R. Car, M. Parrinello, Unified approach for molecular dynamics and density-functional theory, Phys. Rev. Lett. 55 (1985) 2471–2474.

[36] T.-L. Pham, M. Guerboub, A. Bouzid, M. Boero, C. Massobrio, Y.-H. Shin, G. Ori, Unveiling the structure and ion dynamics of amorphous Na₃OHxCl antiperovskite electrolytes by first-principles molecular dynamics, J. Mater. Chem. A Mater. Energy Sustain. 11 (2023) 22922–22940.

[37] T.-L. Pham, M. Guerboub, S. W. Wendj, A. Bouzid, C. Tugène, M. Boero, C. Massobrio, Y.-H. Shin, G. Ori, Structural properties of amorphous Na₃OCI electrolyte by first-principles and machine learning molecular dynamics (2024). arXiv:2404.11442.

[38] S. Wansi Wendji, C. Massobrio, M. Boero, C. Tugène, E. Levchenko, F. Shuaib, R. Piotrowski, D. Hamani, G. Delaizir, P.-M. Geffroy, et al., Quantitative assessment of the structure and bonding properties of 50VₓOᵧ-50P₂O₅ glass by classical and born-oppenheimer molecular dynamics, J. Non-Cryst. Solids 634 (2024) 122967.

[39] J. Kornatowski, B. Wichterlova, M. Rozwadowski, W. Baur, Simultaneous occurrence of differently coordinated framework heteroatoms in one zeolite: MFI type vanadium silicalite, KVS-5., in: Stud. Surf. Sci. Catal., volume 84, Elsevier, 1994, pp. 117–124.

[40] L. Bogomolova, V. Jachkin, N. Krasil'nikova, EPR study of vanadium-containing amorphous silica formed by sol-gel method, J. Non-Cryst. Solids 241 (1998) 13–26.

[41] X. Lu, L. Deng, S. Saslow, H. Liu, C. Benmore, B. Parruzot, J. Reiser, S. Kim, J. Ryan, J. Vienna, et al., Vanadium oxidation states and structural role in aluminoborosilicate glasses: an integrated experimental and molecular dynamics simulation study, J. Phys. Chem. B 125 (2021) 12365–12377.

[42] G. Tricot, H. Vezin, Description of the intermediate length scale structural motifs in sodium vanado-phosphate glasses by magnetic resonance spectroscopies, J. Phys. Chem. C 117 (2013) 1421–1427.

[43] F. Benzi, E. Paris, S. Della Longa, C. Mugoni, C. Siligardi, G. Giuli, V k-edge XANES full multiple scattering study of v-bearing phosphate glasses, in: Proc. in Phys., Synchrotron Radiation Sci. and Appl., volume 220, Springer International Publishing, 2021, pp. 219–231.

[44] F. Benzi, G. Giuli, S. Della Longa, E. Paris, Vanadium k-edge XANES in vanadium-bearing model compounds: a full multiple scattering study, J. Synchrotron Radiat. 23 (2016) 947–952.

[45] P. Zavalij, M. Whittingham, Structural chemistry of vanadium oxides with open frameworks, Acta Crystallogr. B 55 (1999) 627–663.
[46] M. Schindler, F. C. Hawthorne, W. H. Baur, Crystal chemical aspects of vanadium: polyhedral geometries, characteristic bond valences, and polymerization of (VOₙ) polyhedra, Chem. Mater. 12 (2000) 1248–1259.

[47] S. Boudin, A. Guesdon, A. Leclaire, M.-M. Borel, Review on vanadium phosphates with mono and divalent metallic cations: syntheses, structural relationships and classification, properties., Int. J. Inorg. Chem. 2 (2000) 561–579.

[48] M. Bertani, M. Menziani, A. Pedone, Improved empirical force field for multicomponent oxide glasses and crystals, Phys. Rev. Mater. 5 (2021) 045602.

[49] G. Malavasi, A. Pedone, The effect of the incorporation of catalase mimetic activity cations on the structural, thermal and chemical durability properties of the 45S5 bioglass®, Acta Mater. 229 (2022) 117801.

[50] N. Fairley, V. Fernandez, M. Richard-Plouet, C. Guillot-Deudon, J. Walton, E. Smith, D. Flahaut, M. Greiner, M. Biesinger, S. Tougaard, et al., Systematic and collaborative approach to problem solving using X-ray photoelectron spectroscopy, Appl. Surf. Sci. Adv. 5 (2021) 100112.

[51] M. Micoulaut, A. Piarristeguy, O. Masson, L.-M. Poitras, R. Es-calier, A. Kachmar, A. Pradel, Quantitative assessment of net-work depolymerization in archetypal superionic glasses and its relationship with ion conduction: A case study on Na₂S-GeS₂, Phys. Rev. B 108 (2023) 144205.

[52] O. Masson, PYTSREDX, a data reduction program to obtain the atomic pair distribution function (PDF) from X-ray total scattering data (2022).

[53] Dabaxfiles, http://ftp.esrf.fr/pub/scisoft/xop2.3/, 2024.

[54] A. K. Soper, P. Egelstaff, Multiple scattering and attenuation of neutrons in concentric cylinders: I. isotropic first scattering, Nucl. Instrum. Methods 178 (1980) 415–425.

[55] A. Pedone, G. Malavasi, M. Menziani, A. Cormack, U. Segre, A new self-consistent empirical interatomic potential model for oxides, silicates, and silica-based glasses, J. Phys. Chem. B 110 (2006) 11780–11795.

[56] S. Nosé, A unified formulation of the constant temperature molecular dynamics methods, J. Chem. Phys 81 (1984) 511–519.

[57] W. Hoover, Canonical dynamics: Equilibrium phase-space distributions, Phys. Rev. A 31 (1985) 1695.

[58] C. Massobrio, I. Essomba, M. Boero, C. Diarra, M. Guerboub, K. Ishisone, A. Lambrecht, E. Martin, I. Morrot-Woisard, G. Ori, et al., On the actual difference between the nosé and the nosé-hoover thermostats: a critical review of canonical temperature control by molecular dynamics, Phys. Status Solidi B 261 (2023) 2300209.

[59] J. VandeVondele, J. Hutter, Gaussian basis sets for accurate calculations on molecular systems in gas and condensed phases, J. Chem. Phys. 127 (2007) 114105.

[60] S. Goedecker, M. Teter, J. Hutter, Separable dual-space gaussian pseudopotentials, Phys. Rev. B 54 (1996) 1703.

[61] M. Ernzerhof, G. Scuseria, Assessment of the perdew.burke.ernzerhof exchange-correlation functional, J. Chem. Phys. 110 (1999) 5029.

[62] M. Seth, T. Ziegler, Range-separated exchange functionals with slater-type functions, J. Chem. Theory Comput. 8 (2012) 901.

[63] C. Adamo, V. Barone, Toward reliable density functional methods without adjustable parameters: The PBE0 model, J. Chem. Phys. 110 (1999) 6158–6170.

[64] A. Saúl, G. Radtke, Density functional approach for the magnetism of β-TeVO₄, Phys. Rev. B 89 (2014) 104414.

[65] I. Todorov, W. Smith, U. Cheshire, The DL POLY 4 user manual, STFC, STFC Daresbury Laboratory, Daresbury, Warrington, Cheshire, WA4 4AD, United Kingdom, version 4 (2011).

[66] T. Kühne, M. Iannuzzi, M. D. Ben, V. Rybkin, P. Seewald, F. Stein, T. Laino, R. Khaliullin, O. Schütt, F. Schiffmann, Cp2k: An electronic structure and molecular dynamics software package - quickstep: Efficient and accurate electronic structure calculations, J. Chem. Phys. 152 (2020) 1194103.

[67] J. VandeVondele, J. Hutter, Gaussian basis sets for accurate calculations on molecular systems in gas and condensed phases, J. Chem. Phys. 127 (2007).

[68] A. Feltz, B. Unger, Redox reactions in condensed oxide systems II. variation of the structure of vanadium phosphate glasses in dependence on the oxidation state of vanadium, J. Non Cryst. Solids 72 (1985) 335–343.

[69] M. Biesinger, L. Lau, A. Gerson, R. C. Smart, Resolving surface chemical states in xps analysis of first row transition metals, oxides and hydroxides: Sc, Ti, V, Cu and Zn, Appl. Surf. Sci. 257 (2010) 887–898.

[70] A. Wright, The comparison of molecular dynamics simulations with diffraction experiments, J. Non-Cryst. Solids 159 (1993) 264–268.

[71] U. Hoppe, N. Wyckoff, M. Schmitt, R. Brow, A. Schöps, A. Hannon, Structure of V₂O₅–P₂O₅ glasses by X-ray and neutron diffraction, J. Non-Cryst. Solids 358 (2012) 328–336.

[72] T. Honma, T. Togashi, N. Ito, T. Komatsu, Fabrication of Na₂FeP₂₇ glass-ceramics for sodium ion battery, J. Ceram. Soc. Jpn. 120 (2012) 344–346.

[73] U. Hoppe, G. Walter, R. Kranold, D. Stachel, Structural specifics of phosphate glasses probed by diffraction methods: a review, J. Non Cryst. Solids 263-264 (2000) 29–47.

[74] U. Hoppe, A. Ghosh, S. Feller, A. Hannon, D. Keen, J. Neusefeind, Structural units of binary vanadate glasses by X-ray and neutron diffraction, J. Non-Cryst. Solids 572 (2021) 121120.

[75] M. Ungureanu, M. Lévy, J. Souquet, Mixed conductivity of glasses in the P₂O₅–V₂O₅–Na₂O system, Ionics (Kiel) 4 (1998) 200–206.

[76] A. Mosset, P. Lecante, J. Galy, J. Livage, Structural analysis of amorphous V₂O₅ by large-angle X-ray scattering, Philos. Mag. B 46 (1982) 137–149.

[77] M. Nabavi, C. Sanchez, J. Livage, Structure and properties of amorphous V₂O₅, Philos. Mag. B 63 (1991) 941–953.

[78] U. Hoppe, R. Kranold, E. Gattef, An X-ray diffraction study of the structure of vitreous V₂O₅, Solid State Commun. 108 (1998) 71–76.

[79] W. Zachariassen, The atomic arrangement in glass, J. Am. Chem. Soc. 54 (1932) 3841–3851.

[80] J. Du, M. Montorsi, S. Barbi, X. Lu, Rare earth and transition metal containing glasses, Atomistic Simulations of Glasses: Fundamentals and Applications (2022) 367–438.

[81] A. Tilocca, A. Cormack, N. de Leeuw, The structure of bioactive silicate glasses: new insight from molecular dynamics simulations, Chem. Mater. 19 (2007) 95–103.

[82] T. Jenkins, J. Alarco, B. Cowie, I. Mackinnon, Validating the electronic structure of vanadium phosphate cathode materials, ACS Appl. Mater. Interfaces 13 (2021) 45505–45520.

[83] A. Karthikeyan, P. Vinatier, A. Levasseur, K. Rao, The molecular dynamics study of lithium ion conduction in phosphate glasses and the role of non-bridging oxygen, J. Phys. Chem. B 103 (1999) 6185–6192.

[84] I. Mandal, S. Mannan, L. Wondraczek, N. N. Gosvami, A. Allu, N. Krishnan, Machine learning-assisted design of na-ion-conducting glasses, J. Phys. Chem. C Nanomater. Interfaces 127 (2023) 14636–14644.

[85] S. Keshri, S. Ganisetti, R. Kumar, A. Gaddam, K. Illath, T. Ajithkumar, et al., Ionic conductivity of Na₃Al₂P₃O₁₂ glass electrolytes-role of charge compensators, Inorg. Chem. 60
(2021) 12893–12905.

[86] G. Giuli, E. Paris, J. Mungall, C. Romano, D. Dingwell, V oxidation state and coordination number in silicate glasses by XAS, Am. Mineral. 89 (2004) 1640–1646.

[87] A. Bianconi, A. Giovannelli, I. Dovoli, S. Stizza, L. Palladino, O. Gzowski, L. Murawski, XANES (X-ray absorption near edge structure) of V in vanadium-iron phosphate glasses, Solid State Commun. 42 (1982) 547–551.

[88] Data deposited and available at NOMAD data repository, 2024. URL: ...