A planar pentacoordinate oxygen in the experimentally observed [Be₅O₆]²⁻ dianion

Abstract

Small multiply charged anions (SMCAs) are exceptionally challenging to generate in gas-phase experiments due to the spontaneous detachment of excess electrons. The [Be₅O₆]²⁻ dianion, first produced in 2006 via electrospray ionization and initially proposed by a concurrent computational study to adopt a linear O–Be-alternating structure, stands as a rare experimentally observed SMCA. In this study, by applying our recently developed electron-compensation strategy, we designed a starlike D_5h [O©Be₅O₅]²⁻ cluster featuring a planar pentacoordinate oxygen (ppO), which intriguingly shares the molecular formula [Be₅O₆]²⁻. Remarkably, this ppO isomer is not only 55.8 kcal mol⁻¹ more stable than the previously reported linear isomer but also represents the global energy minimum on the [Be₅O₆]²⁻ potential energy surface. By adhering to the principles of the electron-compensation strategy, all Be atoms in the ppO isomer are electronically compensated and geometrically shielded by peripheral O atoms, resulting in a well-defined electronic structure. This is evidenced by a positive first vertical detachment energy of 2.44 eV, which effectively prevents the spontaneous loss of excess electrons. Thus, our work serendipitously uncovered and elaborately rationalized an experimentally unprecedented ppO within the previously generated SMCA [Be₅O₆]²⁻, marking a significant milestone in the field.

---

Introduction

Clusters featuring non-classical bonding are typically generated and observed in the gas phase. In general, only the thermodynamically most stable isomer, known as the global energy minimum (GEM), survives the gas-phase annealing process and is subsequently detected through combined spectroscopy. In contrast, although significantly smaller amounts of local energy minima (LEMs) may persist, their spectroscopic signals are often overshadowed by those of the GEMs. For instance, nearly all clusters featuring non-classical planar hyper-coordination have been confirmed to be GEMs. Notable examples include [CB₄]⁺,^1,2 [CAl₄]⁻,³ [CAl₄Na]⁻,⁴ [CAl₃X]^−/0 (X = Si, Ge),^5,6 [CAl₄H]⁻,⁷ C₂Al₄,⁸ and [C₅Al₅]⁻ (ref. 9) with planar tetracoordinate carbon, [NAl₄]^−/0 (ref. 10 and 11) with a planar tetracoordinate nitrogen, [Al₄X]^−/0,¹² [Al₃X₂]⁻ (X = Si, Ge),¹³ and [Cu₃Si₃]⁻ (ref. 14) with a planar tetracoordinate silicon or germanium, as well as [CoB₈]^−/0, [RuB₉]⁻, [TaB₁₀]⁻, and [NbB₁₀]⁻ with planar octa-, nona-, and decacoordinate transition metals.^15,16

Consequently, GEMs are highly favored in computational predictions of clusters with non-classical bonding due to their significantly greater compatibility with gas-phase experiments compared to LEMs. This preference is particularly evident in the design of clusters featuring non-classical planar hyper-coordination. Over the past two decades, hundreds of GEMs have been computationally predicted in this field,^17–23 most of which exhibit planar pentacoordinate configurations involving H and typical second-row non-metals (ppX, where X = H, C, N, F, etc.).^23–35 However, no cluster with a ppX has been experimentally observed to date. Furthermore, stable planar pentacoordinate oxygen (ppO) remains entirely unexplored, even among computationally predicted GEMs.

From the examples above, it is also evident that most experimentally observed clusters are monoanions, likely because they are characterized using photoelectron detachment spectroscopy, which favors monoanions. To the best of our knowledge, no small multiply charged anion (SMCA) with non-classical planar hyper-coordination has been reported, as SMCAs are exceptionally challenging to generate in gas-phase experiments due to the spontaneous detachment of excess electrons. Combined with the experimental gap concerning ppX, a significant breakthrough would be the observation of ppX in SMCAs.

In this work, we report the design of a starlike D_5h [O©Be₅O₅]²⁻ cluster (1a in Fig. 1) featuring a planar pentacoordinate oxygen (ppO). By applying our recently developed electron-compensation strategy,³⁶ this dianion exhibits a positive first vertical detachment energy of 2.44 eV, indicating the avoidance of spontaneous electron detachment. Notably, a literature survey revealed that the corresponding [Be₅O₆]²⁻ dianion was generated in 2006 via electrospray ionization,³⁷ but a concurrent computational study³⁸ incorrectly proposed a linear O–Be-alternating structure (0 in Fig. 1). Our ppO isomer is not only 55.8 kcal mol⁻¹ more stable than such a linear isomer but also represents the GEM on the [Be₅O₆]²⁻ potential energy surface. Thus, we have identified a ppO within an experimentally observed SMCA, marking the first experimental observation of a cluster with a ppX and the first SMCA with planar hyper-coordination.

Fig. 1 Optimized structures and relative energies (ΔEs) of low-lying isomers of Be₅O₆²⁻. The interatomic distances are given in Å, while the ΔEs are given in kcal mol⁻¹ and computed using the energy of linear isomer 0 as the reference. Green and red balls denote Be and O atoms, respectively.

Results and discussion

Design of [O©Be₅O₅]²⁻

The design of 1a originated from our systematic effort to create a beryllium-based starlike structure featuring a ppO, which coincidentally resulted in the molecular formula [Be₅O₆]²⁻. Specifically, we began by designing [X©Be_nO_n] model structures (see Scheme 1) using our recently proposed “electron-compensation” strategy. By following this approach, all Be atoms in the model structures are electronically compensated and geometrically shielded by peripheral O atoms, resulting in stable Be_nO_n skeletons. Though these skeletons lack valence electrons for centripetal bonding, the Be atoms on such skeletons are ready to accept the donation from the central atom, that is, they are electronically suited to accommodate a central atom (X) with a fully filled valence shell, such as F⁻, O²⁻, N³⁻, Ne, and others. During the evaluation of the compatibility of these candidate atoms with Be_nO_n skeletons for n = 4–6, we discovered that an oxygen atom could perfectly fit the Be₅O₅ skeleton at a molecular charge of −2.00 |e|, yielding a starlike dianion with D_5h symmetry, [O©Be₅O₅]²⁻ (1a). At the B2PLYP³⁹-D3 (ref. 40) (BJ)⁴¹/aug-cc-pVTZ level of theory, 1a was confirmed as an energy minimum, with the lowest vibrational frequency (ν_min) at 114 cm⁻¹. The distances between the central O atom and each peripheral Be atom are uniformly 1.765 Å, only 0.030 Å longer than the sum (1.735 Å) of the average radii of Be (1.061 Å) and O (0.674 Å) atoms in tetrahedrally bonded crystals.⁴² Thus, the central O atom can be considered to be coordinated by five Be atoms, confirming the formation of a non-classical ppO.

Scheme 1 The X©Be_nO_n (n = 4–6) model structures for designing the clusters with planar hypercoordinations.

Stability consideration

Next, we compared the thermodynamic stability of 1a with the previously reported linear isomer 0. We reoptimized the structure of 0 at the B2PLYP-D3(BJ)/aug-cc-pVTZ level and refined the energies of 0 and 1a using CCSD(T)/aug-cc-pVTZ single-point calculations based on the optimized geometries [abbreviated as CCSD(T)//B2PLYP-D3(BJ)]. Remarkably, 1a is 55.8 kcal mol⁻¹ lower in energy than 0, indicating that 1a should be significantly more likely to be generated and characterized in the gas phase experiment. To further assess the experimental viability of 1a, we extensively explored the singlet and triplet potential energy surfaces of [Be₅O₆]²⁻ using the stochastic search algorithm^43,44 implemented in the GXYZ 3.0 program.⁴⁵ The low-lying isomers are shown in Fig. 1. As illustrated, 1a was confirmed to be the exclusive GEM. In addition to 1a, we identified three new isomers (1b, 1c, and 1d) that are lower in energy than 0 by 17.0, 10.2, and 0.2 kcal mol⁻¹, respectively. Structurally, 1b adopts a D_3h geometry with two axial beryllium atoms linked by three O–Be–O bridges, 1c features a C_s structure with a BeO group attached to a pyramidal Be₄O₅ moiety, and 1d consists of a Be₃O₃ ring with a linear Be₂O₃ tail. Note that the T1 diagnostic values of CCSD(T) calculations for these [Be₅O₆]²⁻ isomers range from 0.016 to 0.018, lower than the threshold of 0.020, so the results from such single-reference calculations should be reliable. Given that 1a is at least 38.8 kcal mol⁻¹ lower in energy than its isomers, it should be the only experimentally observable isomer, dominating the [Be₅O₆]²⁻ dianion generated in electrospray ionization experiments. Therefore, we focus on 1a in the subsequent analysis.

To evaluate the dynamic stability of 1a,^46–49 we performed Born–Oppenheimer molecular dynamics (BOMD) simulations. The structural evolution was monitored using the root-mean-square deviation (RMSD) of atomic positions. As shown in Fig. 2, the RMSD plots for simulations at 298 K, 500 K, and 1000 K exhibit no irreversible upward jumps, and the fluctuations remain small, with variation ranges of 0.02–0.19 Å, 0.03–0.24 Å, and 0.04–0.34 Å, respectively. The average RMSD values are 0.07 Å, 0.10 Å, and 0.15 Å, indicating that 1a is dynamically rigid against isomerization and dissociation at these temperatures.

Fig. 2 RMSD (in Å) versus simulation time (in ps) for the BOMD simulations (at 298, 500 and 1000 K) of 1a at the PBE/DZVP level.

Electronic structure analyses

The electronic stability of 1a was assessed by examining its vertical detachment energy (VDE) using the outer valence Green's function (OVGF) method at the OVGF/aug-cc-pVTZ level.⁵⁰ The positive first VDE of 2.44 eV indicates that electron detachment from 1a is endothermic, explaining its stable existence in electrospray ionization experiments. Additionally, we observed a large HOMO–LUMO gap of 4.76 eV at the B2PLYP-D3(BJ)/aug-cc-pVTZ level, suggesting that electron excitation from occupied to Rydberg orbitals is relatively difficult. Together, the positive VDE and large HOMO–LUMO gap confirm the electronic robustness of 1a.

To understand the stability of 1a, we analyzed its electronic structure using adaptive natural density partitioning (AdNDP)^51,52 to identify characteristic n-center two-electron (nc-2e) bonds. As shown in Fig. 3, among the 24 valence electrons (including two negative charges) in 1a, there are five 1c–2e O lone pairs (occupation numbers, ONs = 1.99 |e|, orbital A), ten 2c–2e Be–O σ bonds (ONs = 1.99 |e|, orbital B), five 3c–2e Be–O–Be π bonds (ONs = 2.00 |e|, orbital C), and four lone pairs on the central ppO (ONs = 1.91–1.97 |e|, orbitals D–G). Notably, the formation of five 3c–2e Be–O–Be π bonds not only reinforces the Be₅O₅ skeleton but also compensates for the electron deficiency of the five Be atoms through O → Be π-backdonation, aligning with our electron-compensation strategy and stabilizing the starlike structure. Simultaneously, the central ppO satisfies the octet rule (with eight valence lone pair electrons). We also considered an alternative AdNDP scheme with three 6c–2e σ bonds (orbitals E′–G′). Although the ONs for E′–G′ (1.99 |e|) are slightly higher than those for E–G (1.91–1.92 |e|), the differences of 0.07–0.08 |e| are negligible, as they are distributed across five Be atoms (each contributing only 0.014–0.016 |e|). Thus, the AdNDP analysis confirms that the central ppO in 1a is a dianion.

Fig. 3 AdNDP view of chemical bonding in 1a. ON denotes the occupation number, which should be close to 2.00 |e| for doubly occupied orbitals.

To further confirm this, the natural bond orbital (NBO)⁵³ analysis was performed for 1a at the B2PLYP-D3(BJ)/aug-cc-pVTZ level. Consistent with the AdNDP analysis, the ppO possesses the electron configuration of 2 s^1.912p_x^1.94p_y^1.942p_z^1.99, suggesting a dianion. The natural charge on ppO (−1.80 |e|), while that on each Be atom is +1.71 |e|, suggesting strong electrostatic characteristics for ppO–Be interactions. Meanwhile, the Wiberg bond index for each ppO–Be bond is 0.08, indicating negligible covalent interactions.

To further investigate the nature of interactions between the central ppO and the peripheral Be₅O₅ skeleton, energy decomposition analysis with natural orbitals for chemical valence (EDA-NOCV)^54,55 was performed at the B3LYP-D3(BJ)/TZ2P level. Table S1† presents the results of different fragmentation schemes. The scheme involving a singlet O²⁻ and a singlet neutral Be₅O₅ fragment yielded the lowest orbital interaction energy, consistent with the AdNDP analysis, and was therefore selected. As shown in Table 1, the total attractive energy (ΔE_attr, −913.6 kcal mol⁻¹) and the Pauli repulsion (ΔE_Pauli, 475.3 kcal mol⁻¹) result in an interaction energy (ΔE_int) of −441.3 kcal mol⁻¹. Within ΔE_attr, electrostatic interactions (ΔE_elstat, −656.7 kcal mol⁻¹) account for 71.9%, while orbital interactions (ΔE_orb, −256.9 kcal mol⁻¹) contribute only 28.1%, indicating that the interaction between the central O²⁻ and the peripheral Be₅O₅ skeleton is predominantly electrostatic. Furthermore, 79.6% of the total ΔE_orb arises from terms exhibiting donation or backdonation characteristics (ΔE_orb(1)–ΔE_orb(3) and ΔE_orb(5)), whereas ΔE_orb(4), which reflects electron-sharing bond characteristics, contributes only 5.3% of ΔE_orb. This further underscores the dianionic nature of the central ppO.

Table 1 Results of EDA-NOCV calculations for the Be₅O₆²⁻ cluster at the B3LYP-D3(BJ)/TZ2P level using O²⁻ (singlet, 2s²2p⁶) and the Be₅O₅ (singlet) ligand (L) as interacting fragments. Energies are given in kcal mol⁻¹. The deformation densities of individual orbitals are given in the ESI^a

Energy terms	Assignments	Interaction energies for O^2− + Be5O5
a a and b are the percentage contributions to ΔEattr and ΔEorb, respectively.
ΔEint		−441.3
ΔEPauli		475.3
ΔEattr		−913.6
ΔEelstat		−656.7 (71.9%)^a
ΔEorb		−256.9 (28.1%)^a
ΔEorb(1)	L ← O (s) σ backdonation	−102.7 (40.0%)^b
ΔEorb(2)	L ← O (px) σ backdonation	−43.4 (16.9%)^b
ΔEorb(3)	L ← O (py) σ backdonation	−43.3 (16.9%)^b
ΔEorb(4)	L–O (pz) electron-sharing bonds	−13.6 (5.3%)^b
ΔEorb(5)	L → O (px/py) σ donation	−15.1 (5.8%)^b
ΔErest		−38.7 (15.1%)^b

---

Conclusions

In summary, our targeted design of the starlike [O©Be₅O₅]²⁻ cluster, featuring a planar pentacoordinate oxygen (ppO), coincidentally shares the molecular formula [Be₅O₆]²⁻ with a small multiply charged anion generated in 2006 via electrospray ionization experiments. Through extensive exploration of the [Be₅O₆]²⁻ potential energy surface, we revealed that this ppO isomer is not only 55.8 kcal mol⁻¹ more stable than the previously proposed linear O–Be-alternating isomer but also represents the global energy minimum. This strongly suggests that nearly all the experimentally generated [Be₅O₆]²⁻ clusters adopt the starlike ppO geometry. The experimental persistence of this dianion is attributed to the mitigation of electron deficiency of beryllium atoms through O → Be π-backdonation, as well as its unique electronic structure, which requires the two additional electrons (corresponding to the molecular charges) for the central ppO to satisfy the octet rule. Moreover, the VDE of 2.44 eV corresponds to an endothermic electron detachment for 1a. This energetically feasible transition implies that 1a is a promising candidate for observation via experimental photoelectron spectra. These findings not only uncover but also rationalize the experimental observation of an unprecedented ppO within a small multiply charged anion [Be₅O₆]²⁻, marking a significant breakthrough in the field.

Computational

The singlet and triplet [Be₅O₆]²⁻ potential energy surfaces were explored using the stochastic search algorithm.^43,44 Randomly generated structures were initially optimized at the B3LYP/6-31G(d) level and then the ten lowest-energy isomers were re-calculated at the B2PLYP³⁹-D3 (ref. 40) (BJ)⁴¹/aug-cc-pVTZ level. The energies of the eight lowest isomers were further refined at the CCSD(T)/aug-cc-pVTZ level. The relative energies of the isomers were compared at the CCSD(T)/aug-cc-pVTZ level, incorporating B2PLYP-D3(BJ)/aug-cc-pVTZ zero-point energy corrections [abbreviated as CCSD(T)//B2PLYP-D3(BJ)]. The geometries of low-lying isomers are also re-optimized at the PBE0-D3/aug-cc-pVTZ level, which are essentially not different from those obtained at the B2PLYP-D3(BJ)/aug-cc-pVTZ level. The structures optimized at the B2PLYP-D3(BJ)/aug-cc-pVTZ level are shown in the text, while the optimized geometries (in Cartesian coordinates) are given in the ESI. Born–Oppenheimer molecular dynamic (BOMD) simulations^46–49 were carried out to assess dynamic stability at the PBE/DZVP level. Vertical detachment energies (VDEs) were calculated using the outer valence Green's function (OVGF) procedure at the OVGF/aug-cc-pVTZ level.⁵⁰ To better understand the chemical bonding, natural bond orbital (NBO)⁵³ analysis and adaptive nature density partitioning (AdNDP) analysis⁵¹ were performed at the B2PLYP-D3(BJ)/aug-cc-pVTZ and B3LYP/6-31G* levels, respectively. The AdNDP analysis was done using the AdNDP program,⁵² while the energy decomposition analysis with natural orbitals for chemical valence (EDA-NOCV)^54,55 was further conducted at the B3LYP-D3(BJ)/TZ2P level using the ADF 2022 program package.⁵⁶ The stochastic search algorithm was implemented using the GXYZ 3.0 program,⁴⁵ the CCSD(T) calculations were carried out using the MolPro 2012.1 package,⁵⁷ and all other calculations were performed using the Gaussian 16 package.⁵⁸

Data availability

The data supporting this article have been included as part of the ESI.†

Author contributions

The manuscript was written through contributions from all authors. All authors have given approval to the final version of the manuscript.

Conflicts of interest

The authors declare no conflict of interest.

Acknowledgements

This work was supported by the NSFC (Grant No. 22073058), the Natural Science Foundation of Shanxi Province (Grant No. 202303021211017 and 202403021222073) and the HPC of Shanxi University.

Notes and references

1. Z. H. Cui, M. Contreras, Y. H. Ding and G. Merino, J. Am. Chem. Soc., 2011, 133, 13228–13231.
2. S. Becker and H. J. Dietze, Int. J. Mass Spectrom. Ion Processes, 1988, 82, 287–298.
3. X. Li, L. S. Wang, A. I. Boldyrev and J. Simons, J. Am. Chem. Soc., 1999, 121, 6033–6038.
4. X. Li, H. F. Zhang, L. S. Wang, G. D. Geske and A. I. Boldyrev, Angew. Chem., Int. Ed., 2000, 112, 3776–3778.
5. L. S. Wang, A. I. Boldyrev, X. Li and J. Simons, J. Am. Chem. Soc., 2000, 122, 7681–7687.
6. X. Li, H. J. Zhai and L. S. Wang, Chem. Phys. Lett., 2002, 357, 415–419.
7. J. Xu, X. Zhang, S. Yu, Y. H. Ding and K. H. Bowen, J. Phys. Chem. Lett., 2017, 8, 2263–2267.
8. F. Dong, S. Heinbuch, Y. Xie, J. J. Rocca and E. R. Bernstein, Phys. Chem. Chem. Phys., 2010, 12, 2569–2581.
9. C. J. Zhang, P. Wang, X. L. Xu, H. G. Xu and W. J. Zheng, Phys. Chem. Chem. Phys., 2021, 23, 1967–1975.
10. S. K. Nayak, B. K. Rao, P. Jena, X. Li and L. S. Wang, Chem. Phys. Lett., 1999, 301, 379–384.
11. B. B. Averkiev, A. I. Boldyrev, X. Li and L. S. Wang, J. Chem. Phys., 2006, 125, 124305.
12. A. I. Boldyrev, X. Li and L. S. Wang, Angew. Chem., Int. Ed., 2000, 39, 3307–3310.
13. X. Li, L. S. Wang, N. A. Cannon and A. I. Boldyrev, J. Chem. Phys., 2002, 116, 1330–1338.
14. X. H. Yin, H. L. Zeng, X. B. Liu, X. L. Xu, H. G. Xu, G. Merino, W. J. Zheng and Z. H. Cui, Angew. Chem., Int. Ed., 2025, 64, e202415789.
15. C. Romanescu, T. R. Galeev, W. L. Li, A. I. Boldyrev and L. S. Wang, Angew. Chem., Int. Ed., 2011, 50, 9334–9337.
16. T. R. Galeev, C. Romanescu, W. L. Li, L. S. Wang and A. I. Boldyrev, Angew. Chem., Int. Ed., 2012, 124, 2143–2147.
17. R. Keese, Chem. Rev., 2006, 106, 4787–4808 , and references therein.
18. L. M. Yang, E. Ganz, Z. Chen, Z. X. Wang and P. v. R. Schleyer, Angew. Chem., Int. Ed., 2015, 54, 9468–9501 , and references therein.
19. P. Nag and S. R. Vennapusa, Chemistry, 2022, 4, 1723–1756 , and references therein.
20. L. Leyva-Parra, L. Diego, O. Yañez, D. Inostroza, J. Barroso, A. Vásquez-Espinal, G. Merino and W. Tiznado, Angew. Chem., Int. Ed., 2021, 60, 8700–8704.
21. X. B. Liu, W. Tiznado, L. J. Cui, J. Barroso, L. Leyva-Parra, L. H. Miao, H. Y. Zhang, S. Pan, G. Merino and Z. H. Cui, J. Am. Chem. Soc., 2024, 146, 16689–16697.
22. G. Castillo-Toraya, F. Ortíz-Chi, J. Barroso, M. Orozco-Ic, L. Leyva-Parra and G. Merino, Angew. Chem., Int. Ed., 2025, 64, e202500292.
23. V. Vassilev-Galindo, S. Pan, K. J. Donald and G. Merino, Nat. Rev. Chem., 2018, 2, 0114 , and references therein.
24. K. Sarmah, A. J. Kalita, S. K. Purkayastha and A. K. Guha, Angew. Chem., Int. Ed., 2024, 63, e202318741.
25. A. J. Kalita, S. S. Rohman, C. Kashyap, S. S. Ullah, I. Baruah and A. K. Guha, Inorg. Chem., 2020, 59, 17880–17883.
26. M. H. Wang, X. Dong, Z. H. Cui, M. Orozco-Ic, Y. H. Ding, J. Barroso and G. Merino, Chem. Commun., 2020, 56, 13772–13775.
27. L. J. Cui, L. H. Miao, M. Orozco-Ic, L. Li, S. Pan, G. Merino and Z. H. Cui, Angew. Chem., Int. Ed., 2025, 64, e202416057.
28. L. X. Bai, Y. X. Jin, M. Orozco-Ic, G. Merino and J. C. Guo, Chem. Commun., 2024, 60, 14996–14999.
29. R. Sun, B. Jin, J. Zhao, X. Wu and C. Yuan, ChemPhysChem, 2025, 26, e202400882.
30. J. H. Bian, B. Jin, X. F. Zhao, R. Sun, C. Yuan, C. Y. Zhou and Y. B. Wu, RSC Adv., 2021, 11, 15841–15846.
31. X. F. Zhao, J. H. Bian, F. Huang, C. Yuan, Q. Wang, P. Liu, D. B. Li, X. T. Wang and Y. B. Wu, RSC Adv., 2018, 8, 36521–36526.
32. R. Sun, X. F. Zhao, B. Jin, B. Huo, J. H. Bian, X. L. Guan, C. Yuan and Y. B. Wu, Phys. Chem. Chem. Phys., 2020, 22, 17062–17067.
33. R. Sun, B. Jin, B. Huo, C. Yuan, H. J. Zhai and Y. B. Wu, Chem. Commun., 2022, 58, 2552–2555.
34. M. H. Wang, X. Dong, Y. H. Ding and Z. H Cui, Chem. Commun., 2020, 56, 7285–7288.
35. L. X. Bai, C. Y. Gao, J. C. Guo and S. D. Li, Molecules, 2024, 29, 3831.
36. B. Jin, C. Yuan, G. Lu and Y. B. Wu, Chem. Commun., 2022, 58, 13095–13098.
37. K. Franzreb and P. Williams, Chem. Phys. Lett., 2006, 419, 379–384.
38. A. Dreuw, Chem. Phys. Lett., 2006, 419, 385–389.
39. S. Grimme, J. Chem. Phys., 2006, 124, 034108.
40. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
41. S. Grimme, S. Ehrlich and L. Goerigk, J. Comput. Chem., 2011, 32, 1456–1465.
42. P. Pyykkö, J. Phys. Chem. A, 2015, 119, 2326–2337.
43. M. Saunders, J. Comput. Chem., 2004, 25, 621–626.
44. P. P. Bera, K. W. Sattelmeyer, M. Saunders, H. F. Schaefer and P. v. R. Schleyer, J. Phys. Chem. A, 2006, 110, 4287–4290.
45. H. G. Lu and Y. B. Wu, in GXYZ 3.0, Shanxi University, Taiyuan, 2021.
46. The BOMD simulations were performed at the PBE⁴⁷/DZVP⁴⁸ level using CP2K package⁴⁹.
47. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865–3868.
48. J. VandeVondele and J. Hutter, J. Chem. Phys., 2007, 127, 114105.
49. T. D. Kühne and J. Hutter, J. Chem. Phys., 2020, 152, 194103.
50. J. V. Ortiz, V. G. Zakrzewski and O. Dolgounircheva. Conceptual Perspectives in Quantum Chemistry, Kluwer Academic, 1997.
51. D. Y. Zubarev and A. I. Boldyrev, Phys. Chem. Chem. Phys., 2008, 10, 5207–5217.
52. The AdNDP program was downloaded freely at http://ion.chem.usu.edu/%7Eboldyrev/adndp.php .
53. A. E. Reed, L. A. Curtiss and F. Weinhold, Chem. Rev., 1988, 88, 899–926.
54. A. Michalak, M. Mitoraj and T. Ziegler, J. Phys. Chem. A, 2008, 112, 1933–1939.
55. M. P. Mitoraj, A. Michalak and T. A. Ziegler, J. Chem. Theory Comput., 2009, 5, 962–975.
56. G. Te Velde, et al , J. Comput. Chem., 2001, 22, 931–967.
57. H. J. Werner, et al, in MolPro 2012.1, University College Cardiff Consultants Limited, Cardiff U.K., 2012.
58. M. J. Frisch, et al, in Gaussian 16, Revision A.03, Gaussian, Inc., Wallingford CT, U.S.A., 2016.

---

Footnote

† Electronic supplementary information (ESI) available: The EDA results of the Be₅O₆²⁻ cluster obtained using different charged fragments at the B3LYP-D3(BJ)/TZ2P level, the shapes of deformation densities (Δρ) for EDA-NOCV analysis of the Be₅O₆²⁻ cluster, and Cartesian coordinates for the structures reported in this work. See DOI: https://doi.org/10.1039/d5sc02361k

---