The structure and bonding of V O –P O (VP50) glass is studied by classical and Born–Oppenheimer molecular dynamics (CMD and BOMD). By testing the performance of three different empirical CMD force fields sets and through a thorough comparison with experimental data, our study showcases the effectiveness of BOMD in addressing the shortcomings of CMD and achieves a significantly improved and superior quantitative description of VP50 glass’ structures. BOMD allows us to achieve an unprecedented agreement with experimental data in terms of both reciprocal space (neutron and X-ray structure factors) and real space (total pair correlation functions) properties. The key improvement is ascribed to a better description of the local electronic and bonding environment around both P and, especially, V sites that cannot be obtained by empirical force fields. A clear signature of single and double bonds is found with BOMD together with an overall better description of VO coordinating polyhedra distribution that constitute the network of this glass. Our study is enriched by a thorough analysis of bond angle distributions around VO units, order and connectivity parameters and local bonding features based on Wannier functions formalism. As a byproduct of this work, we assessed which CMD scheme aligns more closely with the BOMD data and experimental findings. Our analysis revealed that including three-body potential parameters associated with the local environment of V sites significantly improves the performance of CMD. However, while enhanced, CMD still falls short of achieving the full accuracy demonstrated by BOMD in describing VP50 glass. All together, the results found define the computational grounds for a deep understating of VP amorphous glasses in conjunctions with depicting the necessary requirements for the development of interatomic potentials aiming to a quantitative comprehension and design of VP-based amorphous materials.