Resolving the trifecta of discrepancies regarding the phosphorus suboxide P₄O

Abstract

Unsatisfied with vibrational and mechanistic discrepancies between three previous investigations, we characterized the P₄O potential energy surface with high-level coupled cluster methods for resolution, up to CCSD(T)/cc-pV(T+d)Z. Our computations substantiate the previous vibrational band assignments of Andrews and coworkers to cyclic and terminal P₄O. The antisymmetric P–O stretch of bridge-bonded P₄O was previously assigned to a band at 856 cm⁻¹, and the symmetric stretch was assigned to a lower frequency band at 553 cm⁻¹. In contrast, we find the symmetric P–O stretch of bridge-bonded P₄O to lie higher than the antisymmetric, and formally assign the symmetric stretch to a previously unassigned 770 cm⁻¹ band. Tentative assignments for the now-elusive antisymmetric P–O stretch of bridge-bonded P₄O are discussed herein. To resolve mechanistic discrepancies, we demonstrate that the isomerization pathway between terminal and bridge-bonded P₄O is a two-step process, requiring passage through a local-minimum intermediate on the potential energy surface. It was previously suggested that the antisymmetric P–O stretch band originally assigned to bridge-bonded P₄O (856 cm⁻¹) could instead belong to a D_2d-symmetry P₄O₂ species, but our computed fundamental vibrational frequency (654 cm⁻¹) does not support this suggestion. New insight has revealed our previous substantiation of the vibrational band assigned to cyclic P₄O₂ to be incorrect, and we now assign the positive combination antisymmetric P–O stretch to an unassigned 919 cm⁻¹ band. Focal point analysis with basis sets up to cc-pV(Q+d)Z and correlation treatment to CCSDT(Q) is performed with the addition of zero-point vibrational energies for cyclic (−95.5 kcal mol⁻¹), bridge-bonded (−92.5), and terminal P₄O (−81.7) relative to tetrahedral P₄ plus ³P oxygen as well as D_2d-symmetry P₄O₂ (−62.0) relative to tetrahedral P₄ plus O₂. We provide the energy of oxygen (+58.6 kcal mol⁻¹) using our methods for comparison to previously reported P₄O₂ species with optimal error cancellation.

---

1. Introduction

The recent discovery of phosphorus allotropes as exfoliable van der Waals (vdW) materials has resurged interest in phosphorus oxidation chemistry due to these materials being unstable in ambient conditions.^1–10 Light, heat, water, and oxygen all influence the degradation of phosphorus materials. Visible wavelengths of light contain sufficient energy to split white phosphorus (P₄ → 2P₂),^11,12 as well as instigate reactions of black and violet phosphorus films with oxygen by exciting them above their bandgaps.^1,13 Excess heat will inevitably reorder any phosphorus allotrope into amorphous, red phosphorus in the absence of contaminants. Reactions of water with phosphorus oxides are generally established, and produce various phosphorus acids. Oxidation of phosphorus in excess oxygen produces phosphorus pentoxide, which can exist as P₄O₁₀ adamantane-like molecules in the gas-phase or repeating P₂O₅ subunits in the condensed-phase. Intermediate phosphorus oxides involved in the formation of P₄O₁₀ remain elusive despite extensive investigation.^12,14–30

Andrews was the first to investigate the lower oxides of phosphorus in white phosphorus (P₄) oxidation reactions in a series of articles,^{12,20–24,26} the first published in 1988.¹² Therein, Andrews and Withnall reacted P₄ molecules with oxygen atoms from photolysis of ozone and microwave discharge of molecular oxygen in low-temperature argon matrices. White phosphorus is known to be quite volatile, but surprisingly, the P₄ and O₃ molecules were not reactive upon codeposition. Andrews and Withnall observed broadening and small shifts of the P₄ and O₃ fundamentals, indicating formation of P₄–O₃ complexes but no further reaction. They instigated reactions via 520–1000 nm photolysis, which produced PO, PO₂, and P₄O species. Reactions were spontaneous between P₄ and oxygen atoms from microwave discharge of O₂, and the aforementioned species were also formed along with number of other phosphorus oxides.¹²

Following their identification in 1988,¹² P₄O species were investigated theoretically by Lohr in 1990³¹ at the suggestion of Andrews. Lohr characterized three minima, cyclic (C-), bridge-bonded (B-), and terminal (T-)P₄O (see Section 3.1 Geometries). Of these, T- and B-P₄O were identified by Andrews and Withnall in 1988,¹² and C-P₄O was identified by Mielke, McCluskey, and Andrews in 1990.²² Hatree–Fock (HF) with a 6-31G* basis set was used by Lohr to optimize the geometries and compute vibrational frequencies of the P₄O isomers. Good agreement was found between theoretical vibrational frequencies and those measured by Andrews and coworkers for C- and T-P₄O; however, Lohr was unable to validate the identity of the bands assigned to B-P₄O. Andrews and Withnall assigned vibrational bands at 856 and 553 cm⁻¹ to the antisymmetric and symmetric P–O stretches of B-P₄O.¹² Lohr found the symmetric stretch to have a higher frequency (836 cm⁻¹) than the antisymmetric (734 cm⁻¹), which resulted from a more constrained P–O–P angle than initially thought by Andrews and Withnall. Lohr did not compute isotopic frequencies, but the disagreement between the ¹⁶O frequencies cast doubt on the previous experimental assignment. Lohr suggested the possibility that a D_2d-symmetry P₄O₂ molecule could be responsible for the bands assigned to B-₄O.³¹

In 2023, Yao and coworkers provided theoretical predictions for the reactions P₄ + O₂ and P₄O + O₂ using the B3LYP/6-311++G** method.³² They did not report symmetries for the vibrational modes, but their two highest frequency vibrations of B-P₄O (presumably the symmetric and antisymmetric P–O stretches) were predicted to be 705 and 629 cm⁻¹. These results appear to support previous theoretical predictions,³¹ that the experimentally observed B-P₄O bands likely belong to a different species¹²; however, there remains a difference of over 100 cm⁻¹ between the two theoretical predictions for the location of the B-P₄O vibrational bands.^31,32 Additionally, the first theoretical study found a singular transition state along the potential energy surface (PES) between T- and B-P₄O with a barrier of 13 kcal mol⁻¹,³¹ whereas the second predicted two transition states (2 and 3 kcal mol⁻¹ barriers) with an intermediate along the PES from T- to B-P₄O.³² Simply, the mechanistic discrepancy in the literature is whether the reaction from T- to B-P₄O proceeds through a P₃–PO intermediate or a direct insertion of the oxygen atom into the P–P bond.^31,32

All considered, there currently exists a trifecta^12,31,32 of discrepancies regarding P₄O species. The experimental assignments for the antisymmetric and symmetric stretches of B-P₄O do not agree with the two theoretically predicted band locations. The two theoretical predictions for where these vibrational bands should lie do not agree with one another, nor whether this molecule was present in the previous experiments. Lastly, the two theoretical mechanisms for the isomerization reaction between T- and B-P₄O do not agree with each other. These discrepancies are resolved herein.

2. Methods

Geometries and harmonic vibrational frequencies were computed using the coupled cluster singles, doubles, and perturbative triples wavefunctions [CCSD(T)]³³ with the correlation consistent polarized valence triple-ζ basis set including tight d-type functions on phosphorus atoms [cc-pV(T+d)Z].^34–39 Second order vibrational perturbation theory (VPT2)^40,41 was used to compute anharmonicity using the B3LYP-D3BJ/cc-pV(T+d)Z method.^42–45 Fundamental frequencies reported are thus CCSD(T)/cc-pV(T+d)Z harmonic plus δB3LYP-D3BJ/cc-pV(T+d)Z anharmonic frequencies.^46–51 Isotopic (¹⁸O) harmonic vibrational frequencies were also computed, and the fundamental isotopic vibrational frequencies reported are a product of the ¹⁶O fundamental frequency with a ratio of the harmonic frequencies . For Hartree–Fock (HF) and correlation up to CCSD(T), single point energies were computed with the cc-pV(X+d)Z [X = D, T, Q] basis sets. Energies were extrapolated to near the complete basis set limit using focal point analysis^52–55 with additional CCSDT(Q)/cc-pV(D+d)Z single point and CCSD(T)/cc-pV(T+d)Z zero-point vibrational energy corrections.^33,56,57 DLPNO-CCSDT and DLPNO-CCSDT(Q) was utilized for P₄O₂,^58–60 while canonical full T and (Q) was used for the P₄O species.^56,57 MOLPRO 2022.1 was used to optimize geometries and compute harmonic vibrational frequencies.^61–63 ORCA 6.1 was used to obtain effects of anharmonicity, locate transition states via the nudged elastic band code, and for intrinsic reaction coordinate computations.^64–68 Psi4 1.10 was used for DLPNO-CCSDT and -(Q) computations.^58–60,69 MRCC 2022 was used for canonical CCSDT and -(Q) computations.^70,71

3. Results and discussion

3.1. Geometries

In a recent article regarding P₄O₂ species,⁷² we reported a systematic underestimation of antisymmetric P–O vibrational frequencies relative to experiment using the CCSD(T)/cc-pV(T+d)Z method. This behavior was unexpected, as both HF and post-HF methods typically overestimate vibrational frequencies. Additionally, experimental vibrational frequencies are typically redshifted ca. 10 cm⁻¹ or more when measured in an argon matrix,^73,74 which was utilized by Andrews and coworkers.^12,22 We initially suggested that a methodological overestimation of P–O bond lengths may be responsible, since increased bond lengths correspond to lower vibrational frequencies. We are now more confident that this is the case. The computational gold standard, CCSD(T) with a triple-ζ polarization quality basis set, does not precisely reproduce the experimental geometry of the phosphorus monoxide radical diatomic (PO) as desired for accurate computation of vibrational frequencies. Experimentally, the equilibrium P–O bond length is 1.476 Å,⁷⁵ whereas the CCSD(T)/cc-pVTZ value is 1.495 Å. Inclusion of tight d-functions on phosphorus [cc-pV(T+d)Z] improves agreement and yields 1.487 Å, but remains 0.011 Å short. Correlation of all electrons except the phosphorus 1s core with the cc-pwCVQZ basis set achieves 0.001 Å accuracy relative to experiment for the phosphorus monoxide bond length. It appears the overestimation is both due to the incompleteness of the triple-ζ basis set and not correlating the core electrons (Table 1), each affecting the bond length by a similar amount. We also tested this hypothesis on P₂O species, where again an overestimation was observed. For linear P₂O, CCSD(T)-AE/cc-pwCVQZ gave a P–O bond length of 1.466 Å and CCSD(T)/cc-pV(T+d)Z yielded 1.476 Å. For bridged P₂O, CCSD(T)-AE/cc-pwCVQZ gave a P–O bond length of 1.736 Å and CCSD(T)/cc-pV(T+d)Z yielded 1.747 Å. The CCSD(T)-AE/cc-pwCVQZ method is computationally infeasible for the system at hand; however, the systematic 0.01 Å overestimation of P–O bond lengths and accompanying underestimation of P–O vibrational frequencies should be expected for the computationally-feasible CCSD(T)/cc-pV(T+d)Z method.

Table 1 Effect of basis set and correlation of core electrons on the P–O bond length of phosphorus monoxide. FC denotes frozen core and AE denotes correlation of all electrons (exluding the 1s core of phosphorus atoms). shows the difference in bond length between quadruple-ζ and triple-ζ computations, and shows the difference in bond length between AE and FC computations. Bond lengths are in angstrom

Experimental phosphorus monoxide bond length [r(P–O)] = 1.476
Method	r(P–O)	Method	r(P–O)
CCSD(T)-FC/cc-pwCVTZ	1.483	CCSD(T)-FC/cc-pwCVQZ	1.480	−0.03
CCSD(T)-AE/cc-pwCVTZ	1.480	CCSD(T)-AE/cc-pwCVQZ	1.476	−0.04
	−0.003		−0.004

---

Fig. 1 contains the CCSD(T)/cc-pV(T+d)Z optimized geometries of C-P₄O, B-P₄O, T-P₄O, and D_2d-symmetry P₄O₂. C-P₄O is a planar, C_2v ring in which the P–O–P motif at the apex is connected via single bonds. Double P=P bonds lie on the edges of the C-P₄O ring, whereas the phosphorus atoms at the base are connected via a single P–P bond. C-P₄O shares little resemblance with the P₄ tetrahedron, and has been predicted to be formed from the reaction of P₂O + P₂.²² B- and T-P₄O closely resemble the P₄ tetrahedron, and have been predicted to be formed from the reaction P₄ + O.¹² B-P₄O (C_2v) also contains a P–O–P motif akin to C-P₄O, but each of the phosphorus atoms in the P–O–P motif of B-P₄O are connected to both of the phosphorus atoms at the base via single P–P bonds. T-P₄O features a double P=O bond at the apex of a trigonal pyramid of singly-bonded phosphorus atoms. The D_2d-symmetry P₄O₂ structure resembles the adamantane-backbone analogue, P₄O₆, where four of the bridge-bonded oxygen atoms are replaced with P–P bonds.

Fig. 1 CCSD(T)/cc-pV(T+d)Z geometries of cyclic P₄O (C_2v), bridge-bonded P₄O (C_2v), terminal P₄O (C_3v), and D_2d-symmetry P₄O₂ are depicted.

As mentioned, we previously observed an underestimation of P–O–P vibrational frequencies in P₄O₂ species, which stemmed from the overestimation of P–O bond lengths in P–O–P motifs. Considering the difference between the P–O bond lengths in C-P₄O (1.649 Å), B-P₄O (1.692), and D_2d-symmetry P₄O₂ (1.696), it is simply not possible to precisely anticipate the P–O bond lengths in phosphorus suboxides. Nevertheless, we could expect our CCSD(T)/cc-pV(T+d)Z method to overestimate P–O bond lengths ca. 0.01 Å, as shown for the smaller phosphorus oxides.

The following discussion is of the computed CCSD(T)/cc-pV(T+d)Z optimized geometries, and the previously mentioned inconsistencies with experiment should be kept in mind. Recent gas-phase electron diffraction data by Wu and coworkers yielded 2.199 Å P–P bond lengths for the white phosphorus tetrahedron.⁷⁶ Double P=P bonds in C-P₄O have a length of 2.064 Å, and the single P–P bond has a length of 2.166 Å. Longer than the single bonds in both C-P₄O and P₄, the single P–P bonds of B-P₄O have lengths of 2.239 Å (base–base) and 2.248 Å (base–apex). P–P bonds in T-P₄O dramatically vary in length, as the edges of the trigonal pyramid have lengths of 2.134 Å and the base has lengths of 2.300 Å. We attribute this feature to the electronegative oxygen atom pulling electron density from the phosphorus atom at the apex, which in turn draws the base closer to the apex and causes a large lengthening of the basal P–P bonds. D_2d-symmetry P₄O₂ has similarly elongated P–P bonds, with lengths of 2.315 Å.

3.2. Vibrational frequencies

One vibrational band was previously assigned to C-P₄O by Andrews and Mielke in their matrix isolation experiment reacting P₂ with O₃ in argon.²² This C-P₄O band was located at 826 cm⁻¹ (redshifted 36 cm⁻¹ with ¹⁸O) and was assigned to the b₂ P–O stretching mode. We compute a fundamental frequency of 787 cm⁻¹ with a reasonable intensity (33 km mol⁻¹) and an isotopic shift of 35 cm⁻¹ (Table 2). As discussed previously, our computed value lies untypically lower than the experimental value; however, we are now able to attribute this feature to the overestimation of the P–O bond lengths in P–O–P motifs by the CCSD(T)/cc-pV(T+d)Z method. In order to demonstrate, we shortened the P–O bond lengths of the C-P₄O structure by 0.01 Å, consistent our observed degree of overestimation, and computed harmonic CCSD(T)/cc-pV(T+d)Z vibrational frequencies from the modified, non-stationary C-P₄O geometry. Others have demonstrated the effectiveness of similar approaches.^77,78 The computed b₂ P–O fundamental vibrational frequency from our adjusted geometry was 823 cm⁻¹, which compares much better to the 826 cm⁻¹ value obtained experimentally by Andrews and Mielke. To further demonstrate the rationality of this approach, the non-stationary harmonic frequencies were also computed in the normal mode basis of the equilibrium geometry using the CMA1-B method.^79–81 All off diagonal elements were included, and we found maximum total energy distribution (TED)^82–84 coupling to be 0.2% (1 cm⁻¹).

Table 2 Composite fundamental vibrational frequencies and isotopic shifts in cm⁻¹ and harmonic IR intensities in km mol⁻¹ for cyclic P₄O

Cyclic P4O (C2v)
^16O Freq.	Int.	^18O Freq.	Sym.
a Ref. 22, experiment of Mielke, McCluskey, and Andrews^22
201	0	201	a2
298	3	287	b1
299	1	299	b2
432	8	423	a1
459	9	458	a1
537	3	537	a1
561	9	561	b2
651	34	634	a1
787 (826)^a	33	752 (790)^a	b2

---

An unassigned 676 cm⁻¹ band observed in the same experiment could appear to be the a₁ P–O stretch mode of C-P₄O since it has a reasonable intensity of 34 km mol⁻¹ and our computed frequency for this vibration is 651 cm⁻¹. However, the 676 cm⁻¹ band increases upon annealing and is absent after photolysis, whereas the 826 cm⁻¹ band decreases upon annealing and increases after photolysis. Additionally, our computed isotopic shift is 17 cm⁻¹ compared to the 29 cm⁻¹ shift observed by Andrews and Mielke. These characteristics indicate the 826 and 676 cm⁻¹ bands do not belong to the same molecular species.

There is a concerning discrepancy between our predicted a₁ and b₂ P–O stretching frequencies of B-P₄O and those assigned experimentally (Table 3). Andrews and coworkers assigned an 856 cm⁻¹ band to the b₂ P–O stretching mode of B-P₄O and a 553 cm⁻¹ band to its symmetric counterpart.¹² Understandably, Andrews and coworkers concluded that these bands belonged to the P–O stretching modes of B-P₄O because they lie ca. 100 cm⁻¹ below the P–O stretching modes of P₄O₆, and B-P₄O would have a more constrained P–O–P angle. Lohr was the first to report theoretical vibrational frequencies for B-P₄O,³¹ where he found the energetic ordering of the symmetric and antisymmetric P–O stretching modes (respectively 836 and 734 cm⁻¹) to be opposite of the assignments by Andrews and coworkers. Similar to Lohr, we predict the a₁ P–O stretching mode to lie higher than the b₂ (respectively 740 and 642 cm⁻¹); however, our predicted vibrational frequencies are much lower than those predicted by Lohr due to the severe underestimation of P–O bond lengths by HF theory. As we did for C-P₄O, we reduced the P–O bond lengths of B-P₄O by 0.01 Å and computed harmonic CCSD(T)/cc-pV(T+d)Z vibrational frequencies from the adjusted non-stationary geometry. After the addition of the stationary B3LYP-D3BJ/cc-pV(T+d)Z anharmonicity, the fundamental frequencies for the symmetric and antisymmetric P–O stretches of the adjusted geometry were 758 and 664 cm⁻¹, respectively. We again verified this approach with CMA1-B,^79–81 obtaining a TED^82–84 coupling of 0.1% (1 cm⁻¹). Even after the adjustment, our computed frequencies for B-P₄O are still not in agreement with the assignments of Andrews and coworkers.¹²

Table 3 Composite fundamental vibrational frequencies and isotopic shifts in cm⁻¹ and harmonic IR intensities in km mol⁻¹ for bridge-bonded P₄O

Bridge-bonded P4O (C2v)
^16O Freq.	Int.	^18O Freq.	Sym.
a Ref. 21, experiment of Andrews and Mielke; originally assigned 856 and 822 cm^−1 (ref. 21)
316	4	304	b1
327	<1	325	a1
332	0	332	a2
414	3	414	b1
426	2	423	a1
443	9	443	b2
551	5	550	a1
642	30	617	b2
740 (770)^a	58	712 (742)^a	a1

---

In an attempt to remedy our discrepancy, we have systematically performed a series of geometry optimizations followed by vibrational frequency computations to reveal their dependence on the level of theory (Table 4). As mentioned, HF severely underestimates the P–O bond length (1.653 Å) compared to all other levels of theory, which in turn heftily raises its predicted a₁ P–O vibrational frequency (845 cm⁻¹). CCSD(T)/cc-pV(T+d)Z predicts P–O bond lengths of 1.692 Å and a vibrational frequency of 747 cm⁻¹. B3LYP compares exceptionally well to CCSD(T) with the same basis set in this case, with P–O bond lengths of 1.692 Å and a vibrational frequency of 739 cm⁻¹. Increasing the cardinality of the basis set by one [cc-pV(Q+d)Z] with the CCSD(T) wavefunction predictably shortens P–O bond lengths (1.688 Å), in turn slightly raising the vibrational frequency (750 cm⁻¹); while augmentation [aug-cc-pV(T+d)Z] predictably lengthens P–O bond lengths (1.696 Å), and slightly lowers the vibrational frequency (736 cm⁻¹). Inclusion of relativistic effects [exact two component (X2C)]^85,86 at the CCSD(T)/cc-pV(T+d)Z level of theory does not have a reportable impact on the geometry, and lowers the vibrational frequency a mere 2 cm⁻¹. Correlation of all electrons, excluding the 1s of phosphorus atoms, in the CCSD(T) wavefunction and use of the appropriate polarized weighted core basis set [cc-pwCVTZ] shortens the P–O bond lengths (1.685 Å) and raises the vibrational frequency (751 cm⁻¹). These observations agree with the previous findings that increasing the cardinality of basis set (triple-ζ → quadruple-ζ) and correlation of core electrons decreases P–O bond lengths ca. 0.01 Å. All considered, we believe the bands originally assigned to B-P₄O by Andrews and coworkers belong to a different molecular species; however, we do not agree with the claim of Lohr³¹ that B-P₄O is not evidenced in the experiments of Andrews and coworkers.^12,21

Table 4 Geometric parameters (bond lengths in Å and angles in degrees) and antisymmetric P–O–P harmonic vibrational frequencies (in cm⁻¹) for bridge-bonded P₄O at nine selected levels of theory. Results obtained at the level of theory used throughout this investigation are in bold [CCSD(T)/cc-pV(T+d)Z]

Bridge-bonded P4O (C2v)
Level of theory	rP–O	θP–O–P	ωe
HF/cc-pV(T+d)Z	1.653	97.6	845
B3LYP-D3BJ/cc-pV(T+d)Z	1.692	95.3	739
MP2/cc-pV(T+d)Z	1.694	93.0	742
CCSD/cc-pV(T+d)Z	1.682	95.5	773
CCSD(T)/cc-pV(T+d)Z	1.692	94.8	747
CCSD(T)-X2C/cc-pV(T+d)Z	1.692	94.8	745
CCSD(T)/aug-cc-pV(T+d)Z	1.696	94.9	736
CCSD(T)/cc-pV(Q+d)Z	1.688	94.7	750
CCSD(T)-AE/cc-pwCVTZ	1.685	94.8	751

---

Andrews and Mielke left vibrational bands at 830, 802, and 770 cm⁻¹ unassigned, and assigned bands at 687 (P₄O₃), 668 (P₄O₄), 641 (P₄O₆), 634 (P₄O₇), and 592 cm⁻¹ (P₂-OPOPO₂) in their P₄ + O₃ experiment.²¹ Considering our adjusted fundamental frequency of 758 cm⁻¹ and an identical isotopic shift to Andrews and Mielke (28 cm⁻¹), we believe the unassigned 770 cm⁻¹ band belongs the a₁ P–O vibrational frequency of B-P₄O and not the originally assigned 553 cm⁻¹ band. This 770 cm⁻¹ band can also be seen in the original P₄ + O₃ experiment¹² and the P₄O₆ decomposition experiment.²¹ We are less confident in the correct assignment for the b₂ vibrational mode of B-P₄O, as there are a number of candidate bands in the region. The isotopic shift, or lack thereof, of the 687 cm⁻¹ band would take it out of contention. The 668 cm⁻¹ band assigned to P₄O₄ had an isotopic shift of 39 cm⁻¹, which is quite different from our computed shift of 25 cm⁻¹. Further experiments or more advanced theoretical methods are required to identify the exact location of the b₂ mode of B-P₄O, but we believe it should lie ca. 650–700 cm⁻¹ and not around the original 553 cm⁻¹ assignment.

Our T-P₄O fundamental vibrational frequencies compare well to the experiments of Andrews and Withnall,¹² with the exception of the vibrational modes containing basal P–P stretching character (Table 5). Similar to the overestimation of P–O bond lengths in P–O–P motifs by CCSD(T)/cc-pV(T+d)Z, we believe the basal P–P bonds of T-P₄O are also overestimated. The phosphorus diatomic experimental P=P bond length is 1.893 Å, while the CCSD(T)/cc-pV(T+d)Z length is 1.910 Å. The white phosphorus tetrahedron has experimental P–P bond lengths of 2.199 Å, while the CCSD(T)/cc-PV(T+d)Z length is 2.208 Å. Therefore, in the simplest cases of double P=P and single P–P bonds, CCSD(T)/cc-pV(T+d)Z overestimates the lengths by 0.017 and 0.008 Å, respectively. We do not have enough evidence to perform a systematic adjustment of the base P–P bond lengths of T-P₄O for the most accurate computation of vibrational frequencies as we did for P–O bond lengths in P–O–P motifs of C- and B-P₄O; nevertheless, we believe the vibrational assignments of Andrews and Withnall are correct in their entirety for T-P₄O, as did Lohr.

Table 5 Composite fundamental vibrational frequencies and isotopic shifts in cm⁻¹ and harmonic IR intensities in km mol⁻¹ for terminally-bonded P₄O

Terminally-bonded P4O (C3v)
^16O Freq.	Int.	^18O Freq.	Sym.
a Ref 12, experiment of Andrews and Withnall.^12
247 (243)^a	10	240 (237)^a	e
321 (393)^a	1	321 (393)^a	e
396 (441)^a	<1	391 (441)^a	a1
494 (501)^a	6	494 (501)^a	e
555 (603)^a	<1	550 (603)^a	a1
1259 (1241)^a	143	1218 (1202)^a	a1

---

Lohr³¹ suggested a D_2d-symmetry P₄O₂ species could be responsible for the 856 cm⁻¹ vibrational band originally assigned to B-P₄O by Andrews and Withnall¹² due to the antisymmetric P–O stretches transforming as the degenerate e representation. We compute the e-symmetry antisymmetric P–O stretch to be 654 cm⁻¹ at the CCSD(T)/cc-pV(T+d)Z level of theory (Table 6). After reducing the P–O bond lengths 0.01 Å, we find it to be 676 cm⁻¹. The negative (b₂) and positive (a₁) combinations of symmetric P–O stretches were respectively 730 and 724 cm⁻¹ pre-adjustment and 746 and 741 cm⁻¹ post-adjustment. We find no evidence for the formation of this D_2d-symmetry P₄O₂ molecule in the experiments of Andrews and coworkers;^12,21,22 nevertheless, it is a minimum on the PES and Yao predicted it to be formed from the reaction P₄ + O₂.³²

Table 6 Composite fundamental vibrational frequencies and isotopic shifts in cm⁻¹ and harmonic IR intensities in km mol⁻¹ for D_2d-symmetry P₄O₂

D2d-symmetry P4O2
^16O Freq.	Int.	^18O Freq.	Sym.
233	2	224	e
348	0	344	a1
354	0	354	b1
399	2	398	e
422	1	419	b2
505	0	501	a1
654	23	629	e
724	0	696	a1
730	91	702	b2

---

An erroneous substantiation of a vibrational band assigned to cyclic-P₄O₂ was likely to have been made by us in a recent article.⁷² Andrews and coworkers assigned a 867 cm⁻¹ vibration to the positive combination antisymmetric P–O stretch vibration of cyclic-P₄O₂,²² and we computed a fundamental frequency of 893 cm⁻¹. Their observed isotopic shift for this band was 37 cm⁻¹, and we computed a shift of 38 cm⁻¹. Considering how CCSD(T)/cc-pV(T+d)Z has continued to overestimate bridged P–O bond lengths, we reduced the P–O bond length of cyclic-P₄O₂ 0.01 Å and computed harmonic CCSD(T)/cc-pV(T+d)Z frequencies. After applying the anharmonic correction, we find the adjusted vibrational frequency to be 924 cm⁻¹. Bands at 919 and 937 cm⁻¹ went unassigned and had respective isotopic shifts of 40 and 37 cm⁻¹. We think either of these bands could be responsible for the cyclic-P₄O₂ positive combination antisymmetric P–O stretch vibration instead of the previously assigned 867 cm⁻¹ band. The 919 cm⁻¹ band lies closer to our now adjusted vibrational frequency, but the isotopic shift of the 937 cm⁻¹ band lies closer to that which we compute. Considering the great agreement between the 919 cm⁻¹ band and our adjusted vibrational frequency (924 cm⁻¹), we believe it is a better assignment for the cyclic-P₄O₂ species and not the originally assigned 867 cm⁻¹ band. We also will not discount the possibility of the 937 cm⁻¹ band, but further experimental work is required for absolute verification.

3.3. Terminal to bridge-bonded P₄O reaction coordinate

Yao and coworkers recently reported that the T- to B-P₄O isomerization reaction proceeds through an intermediate minimum on the potential energy surface,³² while Lohr appears to have found only the second transition state (TS2).³¹ Substantiating the prediction of Yao and coworkers, we also found a T- to B-P₄O isomerization reaction that proceeds through an intermediate (IM). Fig. 2 contains CCSD/cc-pV(T+d)Z optimized stationary points along the reaction coordinate (included in the SI) with CCSD(T)/cc-pV(T+d)Z single point energies. CCSD/cc-pV(T+d)Z hessian analysis was performed on the stationary points to verify that saddle points contained one imaginary vibrational frequency and that the minima only featured real vibrational modes. Intrinsic reaction coordinate computations using CCSD/cc-pV(T+d)Z were also undertaken to connect the saddle points to the minima via steepest descent along the potential energy surface. We found no evidence for a direct isomerization reaction, containing only one saddle point, between T- and B-P₄O.

Fig. 2 CCSD(T)/cc-pV(T+d)Z//CCSD/cc-pV(T+d)Z stationary points along the T-P₄O → B-P₄O reaction coordinate. TS denotes transition state and IM denotes intermediate.

3.4. Energetics

Focal point analysis^52–55 was used to extrapolate the energies of C-, B-, and T-P₄O relative to tetrahedral P₄ plus ³P oxygen (Tables 7–9). HF is found to be wholly inadequate for accurate prediction of energetics for all three P₄O isomers, as the respective δMP2/cc-pV(Q+d)Z values for C-, B-, and T-P₄O are −45.89, −48.82, and −59.23 kcal mol⁻¹. MP2 performs better for C- and B-P₄O than for T-P₄O, with respective δCCSD extrapolated values of +5.62 and +7.47 kcal mol⁻¹ for C- and B-P₄O and +13.81 kcal mol⁻¹ for T-P₄O. As typical, CCSD(T) predicts lower energies than CCSD, and extrapolated δCCSD(T) values for C-, B-, and T-P₄O are −4.90, −3.61, and −5.55 kcal mol⁻¹. Inclusion of perturbative quadruple excitations [CCSDT(Q)] with the cc-pV(D+d)Z basis set has a larger effect on the energies than what is desirable for our purpose, considering the computational rigor. CCSDT(Q)/cc-pV(D+d)Z corrections for C-, B-, and T-P₄O are −0.78, −0.28, and −0.64 kcal mol⁻¹.

Table 7 Incremented focal point energies in kcal mol⁻¹ for cyclic P₄O. Total energies are relative to tetrahedral P₄ plus ³P oxygen and are reported as Net/Extrapolation + ΔZPVE

Cyclic P4O
Total = −97.66 + 2.18 = −95.5 kcal mol^−1
	ΔHF	+δ MP2	+δ CCSD	+δ CCSD(T)	+δ CCSDT	+δ CCSDT(Q)	Net
cc-pV(D+d)Z	−38.52	−46.34	+6.19	−4.11	−0.21	−0.78	[−83.76]
cc-pV(T+d)Z	−48.05	−44.58	+5.58	−4.83	[−0.21]	[−0.78]	[−92.87]
cc-pV(Q+d)Z	−50.03	−45.89	+5.60	−4.87	[−0.21]	[−0.78]	[−96.18]
Extrapolation	[−50.53]	[−46.85]	[+5.62]	[−4.90]	[−0.21]	[−0.78]	[−97.66]

---

Table 8 Incremented focal point energies in kcal mol⁻¹ for bridge-bonded P₄O. Total energies are relative to tetrahedral P₄ plus ³P oxygen and are reported as Net/Extrapolation + ΔZPVE

Bridge-bonded P4O
Total = −94.65 + 2.14 = −92.5 kcal mol^−1
	ΔHF	+δ MP2	+δ CCSD	+δ CCSD(T)	+δ CCSDT	+δ CCSDT(Q)	Net
cc-pV(D+d)Z	−32.75	−40.99	+5.67	−2.23	−0.21	−0.28	[−70.79]
cc-pV(T+d)Z	−45.61	−46.08	+7.17	−3.27	[−0.21]	[−0.28]	[−88.28]
cc-pV(Q+d)Z	−47.19	−48.82	+7.35	−3.47	[−0.21]	[−0.28]	[−92.62]
Extrapolation	[−47.20]	[−50.82]	[ + 7.47]	[−3.61]	[−0.21]	[−0.28]	[−94.65]

---

Table 9 Incremented focal point energies in kcal mol⁻¹ for terminal P₄O. Total energies are relative to tetrahedral P₄ plus ³P oxygen and are reported as Net/Extrapolation + ΔZPVE

Terminal P4O
Total = −84.04 + 2.35 = −81.7 kcal mol^−1
	ΔHF	+δ MP2	+δ CCSD	+δ CCSD(T)	+δ CCSDT	+δ CCSDT(Q)	Net
cc-pV(D+d)Z	−15.28	−52.11	+12.05	−4.15	−0.10	−0.64	[−60.24]
cc-pV(T+d)Z	−28.96	−56.23	+13.50	−5.18	[−0.10]	[−0.64]	[−77.61]
cc-pV(Q+d)Z	−30.32	−59.23	+13.68	−5.39	[−0.10]	[−0.64]	[−82.00]
Extrapolation	[−30.15]	[−61.41]	[ + 13.81]	[−5.55]	[−0.10]	[−0.64]	[−84.04]

---

Similarly, focal point analysis was also used to extrapolate the energy of the D_2d-symmetry P₄O₂ species (−62.0 kcal mol⁻¹) relative to tetrahedral P₄ plus O₂ (Table 10). Behaving similarly to the previously investigated P₄O₂ species, δMP2 raises its energy ca. 20 kcal mol⁻¹, but δCCSD lowers it by nearly equal magnitude (ca. 18 kcal mol⁻¹). CCSD(T)/cc-pV(T+d)Z overestimates the net extrapolated energy by about 1 kcal mol⁻¹. For the interested investigator, we provide the energy of oxygen (+58.6 kcal mol⁻¹) using single point energy computations and zero-point vibrational energy contributions at the same levels of theory as used in our focal point analysis for comparison to previously reported P₄O₂ species with optimal error cancellation.⁷² The experimental value for oxygen is 59.0 kcal mol⁻¹.⁸⁷

Table 10 Incremented focal point energies in kcal mol⁻¹ for P₄O₂ (D_2d). Total energies are relative to tetrahedral P₄ plus O₂ and are reported as Net/Extrapolation +ΔZPVE

P4O2 (D2d)
Total = −63.97 + 1.98 = −62.0 kcal mol^−1
	ΔHF	+δ MP2	+δ CCSD	+δ CCSD(T)	+δ CCSDT^a	+δ CCSDT(Q)^b	Net
a DLPNO-CCSDT – CCSD(T).b DLPNO-CCSDT(Q) – DLPNO-CCSDT.
cc-pV(D+d)Z	−43.39	+24.54	−16.74	+1.55	−1.87	+0.78	[−35.12]
cc-pV(T+d)Z	−62.34	+20.88	−18.13	+2.17	[−1.87]	[+0.78]	[−58.52]
cc-pV(Q+d)Z	−64.72	+18.66	−18.54	+2.31	[−1.87]	[+0.78]	[−63.39]
Extrapolation	[−64.76]	[+17.04]	[−18.84]	[+2.41]	[−1.87]	[+0.78]	[−65.25]

---

4. Conclusion

Three P₄O species were characterized using the CCSD(T)/cc-pV(T+d)Z method for resolution of existing vibrational and mechanistic discrepancies in the literature. Vibrational bands assigned to C- and T-P₄O by Andrews and coworkers are substantiated, while we believe the symmetric P–O stretch band of B-P₄O assigned to 553 cm⁻¹ is instead a previously unassigned 770 cm⁻¹ band. The molecular species responsible for the vibrational bands previously assigned to B-P₄O remains elusive, but we predict they do not belong to a D_2d-symmetry P₄O₂ molecule. We found no evidence for a direct isomerization pathway between T- and B-P₄O, instead finding a two-step mechanism proceeding through a P₃-PO intermediate. Our previous vibrational band assignment for the positive combination antisymmetric P–O stretch of cyclic-P₄O₂ has been revised, and we now believe it to be a previously unassigned 919 cm⁻¹ band.

Further experimental investigation into these species is encouraged. Current vibrational assignments rely on matrix-isolation spectroscopy, where definitive characterization can be difficult due to the presence of many different species. Mass-selected gas-phase spectroscopy could provide unambiguous identification; however, low production yields may preclude direct absorption measurements. Action-based spectroscopic techniques could circumvent this limitation, necessitating complementary theoretical characterization of charged P₄O species.

Conflicts of interest

There are no conflicts to declare.

Data availability

The data supporting this article has been included as part of the supplementary information (SI). Supplementary information: Tables S1–S4, containing the computed fundamental, harmonic, and isotopic frequencies at the reported levels of theory. Cartesian coordinates for the optimized geometries are also included. See DOI: https://doi.org/10.1039/d6cp00698a.

Acknowledgements

E. J. P and Y. A. acknowledge the Air Force Office of Scientific Research (AFOSR) grant number FA9550-23-1-0375, the Gordon and Betty Moore Foundation, GBMF12246, and partial support from NSF Grant No. 2152159 (NRT-QuaNTRASE). The research of the H. F. S. group is supported by the U.S. Department of Energy, Office of Basic Energy Sciences, Chemical Sciences Division, Fundamental Interactions, Computational and Theoretical Chemistry (CTC) program, contract DE-SC0018412. E. J. P, A. G. H., and M. E. L. thank Dr. Rollin A. King (Bethel University) and Dr. Wesley D. Allen (University of Georgia) for their guidance.

References

1. A. Favron, E. Gaufres, F. Fossard, A.-L. Phaneuf-L’Heureux, N. Y. Tang, P. L. Levesque, A. Loiseau, R. Leonelli, S. Francoeur and R. Martel, Photooxidation and quantum confinement effects in exfoliated black phosphorus, Nat. Mater., 2015, 14, 826–832.
2. Y. Abate, D. Akinwande, S. Gamage, H. Wang, M. Snure, N. Poudel and S. B. Cronin, Recent progress on stability and passivation of black phosphorus, Adv. Mater., 2018, 30, 1704749.
3. L. Zhang, H. Huang, B. Zhang, M. Gu, D. Zhao, X. Zhao, L. Li, J. Zhou, K. Wu and Y. Cheng, et al., Structure and properties of violet phosphorus and its phosphorene exfoliation, Angew. Chem., 2020, 132, 1090–1096.
4. W. Ahmad, A. Abbas, U. Younis, J. Zhang, S. H. Aleithan and Z. Wang, Advancements in optoelectronics: harnessing the potential of 2D violet phosphorus, Adv. Func. Mater., 2024, 34, 2410723.
5. P. C. Debnath, K. Park and Y.-W. Song, Recent advances in black-phosphorus-based photonics and optoelectronics devices, Small Methods, 2018, 2, 1700315.
6. Y. Huang, J. Qiao, K. He, S. Bliznakov, E. Sutter, X. Chen, D. Luo, F. Meng, D. Su and J. Decker, et al., Interaction of black phosphorus with oxygen and water, Chem. Mater., 2016, 28, 8330–8339.
7. Q. Li, H. Yang, X. Zheng, Y. Chen, C. Wang, Y. Han, Y. Guo, X. Zheng and Y. Wei, Thinning Effect of Few-Layer Black Phosphorus Exposed to Dry Oxidation, Nanomaterials, 2025, 15, 974.
8. Y. Kaddar, Z. Mansouri, A. Benyoussef, A. El Fatimy and O. Mounkachi, Elucidating the Oxidation Process and Enhanced Stability of Black Phosphorus through NTCDA Passivation: A Molecular Dynamics Study, Adv. Theory Simul., 2025, 2401205.
9. K. L. Kuntz, R. A. Wells, J. Hu, T. Yang, B. Dong, H. Guo, A. H. Woomer, D. L. Druffel, A. Alabanza and D. Tomanek, et al., Control of surface and edge oxidation on phosphorene, ACS Appl. Mater. Interfaces, 2017, 9, 9126–9135.
10. K. L. Utt, P. Rivero, M. Mehboudi, E. O. Harriss, M. F. Borunda, A. A. Pacheco SanJuan and S. Barraza-Lopez, Intrinsic defects, fluctuations of the local shape, and the photo-oxidation of black phosphorus, ACS Cent. Sci., 2015, 1, 320–327.
11. M. S. Oakley, J. J. Bao, M. Klobukowski, D. G. Truhlar and L. Gagliardi, Multireference Methods for Calculating the Dissociation Enthalpy of Tetrahedral P4 to Two P2, J. Phys. Chem. A, 2018, 122, 5742–5749.
12. L. Andrews and R. Withnall, Matrix Reactions of Oxygen Atoms with P4. Infrared Spectra of P₄O, P₂O, PO, and PO₂, Cheminform, 1988, 19, 5605–5611.
13. M. Ghafariasl, S. Singh, S. Gamage, T. Prusnick, M. Snure and Y. Abate, Photodegradation and Thermal Effects in Violet Phosphorus, Adv. Mater. Interfaces, 2024, 11, 2300794.
14. L. R. Maxwell, S. B. Hendricks and L. S. Deming, The Molecular Structure of P₄O₆, P₄O₈, P₄O₁₀ and As₄O₆ by Electron Diffraction, J. Chem. Phys., 1937, 5, 626–637.
15. A. Chapman, Spectra of phosphorus compounds—III The vibrational assignment and force constants of P4O6 and P4O10. Spectrochim, Act A Mol. Biomol. Spectrosc., 1968, 24, 1687–1696.
16. B. Beagley, D. W. J. Cruickshank, T. G. Hewitt and K. H. Jost, Molecular structures of P4O6 and P₄O₈, Trans. Faraday Soc., 1969, 65, 1219–1230.
17. Z. Mielke and L. Andrews, Infrared spectra of phosphorus oxides (P₄O₆, P₄O₇, P₄O₈, P₄O₉ and P₄O₁₀) in solid argon, J. Phys. Chem., 1989, 93, 2971–2976.
18. R. G. Egdell, M. H. Palmer and R. H. Findlay, Electronic structure of the Group 5 oxides: photoelectron spectra and ab-initio molecular orbital calculations, Inorg. Chem., 1980, 19, 1314–1319.
19. J. L. Rose, T. C. VanCott, P. N. Schatz, M. E. Boyle and M. H. Palmer, Vacuum ultraviolet spectroscopy of Group V oxides: P₄O₆ and As₄O₆, J. Phys. Chem., 1989, 93, 3504–3511.
20. L. Andrews, M. McCluskey, Z. Mielke and R. Withnall, Matrix reactions of P4 and P2 with O3 molecules, J. Mol. Struct., 1990, 222, 95–108.
21. Z. Mielke and L. Andrews, Matrix infrared spectra of the products from photochemical reactions of tetraphosphorus with ozone and decomposition of tetraphosphorus hexoxide, Inorg. Chem., 1990, 29, 2773–2779.
22. Z. Mielke, M. McCluskey and L. Andrews, Matrix reactions of P2 and O3 molecules, Chem. Phys. Lett., 1990, 165, 146–154.
23. M. McCluskey and L. Andrews, Matrix infrared spectra of the products of the phosphorus (P2) and ozone reaction, J. Phys. Chem, 1991, 95, 2988–2994.
24. M. McCluskey and L. Andrews, Infrared spectra of products formed by the photoinduced reactions of phosphorus-oxygen ((P2) x (O2) y) clusters in solid argon, J. Phys. Chem., 1991, 95, 2679–2687.
25. R. C. Mowrey, B. A. Williams and C. H. Douglass, Vibrational Analysis of P₄O₆ and P₄O₁₀, J. Phys. Chem. A, 1997, 101, 5748–5752.
26. C. Bauschlicher, M. Zhou and L. Andrews, A study of the products of the reaction of phosphorus and dioxygen, J. Phys. Chem. A, 2000, 104, 3566–3571.
27. I. R. Beattie, K. M. S. Livingston, G. A. Ozin and D. J. Reynolds, Single-crystal Raman spectra of arsenolite (As₄O₆) and senarmonite (Sb₄O₆). The gas-phase Raman spectra of P4O6, P₄O₁₀, and As₄O₆, J. Chem. Soc. A, 1970, 449–451.
28. Y. Wang, Y. Xie, P. Wei, H. F. Schaefer III, P. V. R. Schleyer and G. H. Robinson, Splitting Molecular Oxygen en Route to a Stable Molecule Containing Diphosphorus Tetroxide, J. Am. Chem. Soc., 2013, 135, 19139–19142.
29. J. S. Francisco, Coupled cluster study of the energetic and spectroscopic properties of OPO x (x= O,+ 1,- 1), J. Chem. Phys., 2002, 117, 3190–3195.
30. X. Chu, W. Qian, B. Lu, L. Wang, J. Qin, J. Li, G. Rauhut, T. Trabelsi, J. S. Francisco and X. Zeng, The triplet hydroxyl radical complex of phosphorus monoxide, Angew. Chem., 2020, 132, 22133–22137.
31. L. L. Lohr, An ab initio characterization of the tetraphosphorus oxide P₄O, J. Phys. Chem., 1990, 94, 4832–4835.
32. L. Yao, Y. Chen and G. Liu, Theoretical study of the oxidation mechanisms of excited P4Ox (x = 1–3) and the anharmonic effect on the main reactions, Chem. Phys. Lett., 2023, 814, 140328.
33. J. F. Stanton, Why CCSD(T) works: a different perspective, Chem. Phys. Lett., 1997, 281, 130–134.
34. B. P. Pritchard, D. Altarawy, B. Didier, T. D. Gibsom and T. L. Windus, A New Basis Set Exchange: An Open, Up-to-date Resource for the Molecular Sciences Community, J. Chem. Inf. Model., 2019, 59, 4814–4820.
35. D. Feller, The role of databases in support of computational chemistry calculations, J. Comput. Chem., 1996, 17, 1571–1586.
36. K. L. Schuchardt, B. T. Didier, T. Elsethagen, L. Sun, V. Gurumoorthi, J. Chase, J. Li and T. L. Windus, Basis Set Exchange: A Community Database for Computational Sciences, J. Chem. Inf. Model., 2007, 47, 1045–1052.
37. T. H. Dunning, Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen, J. Chem. Phys., 1989, 90, 1007–1023.
38. T. H. Dunning, K. A. Peterson and A. K. Wilson, Gaussian basis sets for use in correlated molecular calculations. X. The atoms aluminum through argon revisited, J. Chem. Phys., 2001, 114, 9244–9253.
39. D. E. Woon and T. H. Dunning, Gaussian basis sets for use in correlated molecular calculations. III. The atoms aluminium through argon, J. Chem. Phys., 1993, 98, 1358–1371.
40. I. Mills, Vibration-rotation structure in asymmetric and symmetric top molecules, 1972.
41. P. R. Franke, J. F. Stanton and G. E. Douberly, How to VPT2: Accurate and Intuitive Simulations of CH Stretching Infrared Spectra Using VPT2 + K with Large Effective Hamiltonian Resonance Treatments, J. Phys. Chem. A, 2021, 125, 1301–1324.
42. A. D. Becke, Density-functional thermochemistry. III. The role of exact exchange, J. Chem. Phys., 1993, 98, 5648–5652.
43. C. Lee, W. Yang and R. G. Parr, Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density, Phys. Rev. B, 1988, 37, 785–789.
44. S. H. Vosko, L. Wilk and M. Nusair, Accurate spin-dependent electron liquid correlation energies for local spin density calculations: a critical analysis, Can. J. Phys., 1980, 58, 1200–1211.
45. P. J. Stephens, F. J. Devlin, C. F. Chabalowski and M. J. Frisch, Ab initio calculation of vibrational absorption and circular dichroism spectra using density functional force fields, J. Phys. Chem., 1994, 98, 11623–11627.
46. M. Hochlaf, C. Puzzarini and M. L. Senet, Towards the computations of accurate spectroscopic parameters and vibrational spectra for organic compounds, Mol. Phys., 2015, 113, 1661–1673.
47. C. Levi, J. M. Martin and I. Bar, Fundamental vibrational frequencies and dominant resonances in methylamine isotopologues by ab initio and density functional theory methods, J. Comput. Chem., 2008, 29, 1268–1276.
48. R. P. Bettens, Comparison of fundamental and harmonic frequencies of first-row closed-shell diatomics calculated using full ab initio methods and composite methods, J. Phys. Chem. A, 2004, 108, 1826–1829.
49. J. C. Howard and G. S. Tschumper, Benchmark structures and harmonic vibrational frequencies near the CCSD (T) complete basis set limit for small water clusters:(H2O) n= 2, 3, 4, 5, 6, J. Chem. Theory Comput., 2015, 11, 2126–2136.
50. Y. Fang, F. Liu, S. J. Klippenstein and M. I. Lester, Direct observation of unimolecular decay of CH3CH2CHOO Criegee intermediates to OH radical products, J. Chem. Phys., 2016, 145.
51. S. J. Klippenstein, C. R. Mulvihill and P. Glarborg, Theoretical kinetics predictions for reactions on the NH2O potential energy surface, J. Phys. Chem. A, 2023, 127, 8650–8662.
52. A. G. Csaszar, W. D. Allen and H. F. Schaefer III, In pursuit of the ab initio limit for conformational energy prototypes, J. Chem. Phys., 1998, 108, 9751–9764.
53. M. S. Schuurman, S. R. Muir, W. D. Allen and H. F. Schaefer III, Toward subchemical accuracy in computational thermochemistry: Focal point analysis of the heat of formation of NCO and [H,N,C,O] isomers, J. Chem. Phys., 2004, 120, 11586–11599.
54. D. Feller, The use of systematic sequences of wave functions for estimating the complete basis set, full configuration interaction limit in water, J. Chem. Phys., 1993, 98, 7059–7071.
55. T. Helgaker, W. Klopper, H. Koch and J. Noga, Basis-set convergence of correlated calculations on water, J. Chem. Phys., 1997, 106, 9639–9646.
56. J. Noga and R. J. Bartlett, The full CCSDT model for molecular electronic structure, J. Chem. Phys., 1987, 86, 7041–7050.
57. Y. J. Bomble, J. F. Stanton, M. Kallay and J. Gauss, Coupled-cluster methods including noniterative corrections for quadruple excitations, J. Chem. Phys., 2007, 123, 054101.
58. A. Jiang, Z. L. Glick, D. Poole, J. M. Turney, C. D. Sherrill and H. F. Schaefer III, Accurate and efficient open-source implementation of domain-based local pair natural orbital (DLPNO) coupled-cluster theory using a t1-transformed Hamiltonian, J. Chem. Phys., 2024, 161.
59. A. Jiang, H. F. Schaefer III and J. M. Turney, Linear-Scaling Local Natural Orbital-Based Full Triples Treatment in Coupled-Cluster Theory, J. Chem. Theory Comput., 2025, 21, 2386–2401.
60. A. Jiang, H. F. Schaefer III and J. M. Turney, Linear-scaling quadruple excitations in local pair natural orbital coupledcluster theory, J. Chem. Phys., 2025, 162.
61. H.-J. Werner, P. J. Knowles, G. Knizia, F. R. Manby and M. Schutz, Molpro: a general-purpose quantum chemistry program package, WIREs Comput. Mol. Sci., 2012, 2, 242–253.
62. H.-J. Werner, P. J. Knowles, F. R. Manby, J. A. Black, K. Doll, A. Hesselman, D. Kats, A. Kohn, T. Korona and D. A. Kreplin, et al., The Molpro quantum chemistry package, J. Chem. Phys., 2020, 152, 144107.
63. H.-J. Werner, P. J. Knowles, P. Celani, W. Gyorffy, A. Hesselman, D. Kats, G. Knizia, A. Kohn, T. Korona, D. Kreplin, et al., MOLPRO, version 2022.1, a package of ab initio programs.
64. F. Neese, The ORCA program system, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2012, 2, 73–78.
65. F. Neese, An improvement of the resolution of the identity approximation for the formation of the Coulomb matrix, J. Comput. Chem., 2003, 24, 1740–1747.
66. F. Neese, F. Wennmohs, A. Hansen and U. Becker, Efficient, approximate and parallel Hartree–Fock and hybrid DFT calculations. A ‘chain-of-spheres’ algorithm for the Hartree–Fock exchange, Chem. Phys., 2009, 356, 98–109.
67. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu, J. Chem. Phys., 2010, 132, 154104.
68. S. Grimme, S. Ehrlich and L. Goerigk, Effect of the damping function in dispersion corrected density functional theory, J. Comput. Chem., 2011, 32, 1456–1465.
69. D. G. A. Smith, et al., Psi4 1.4: Open-Source Software for High-Throughput Quantum Chemistry, J. Chem. Phys., 2020, 152(18), 184108.
70. D. Mester, P. R. Nagy, J. Csoka, L. Gyevi-Nagy, P. B. Szabo, R. A. Horvath, K. Petrov, B. Hegely, B. Ladoczki and G. Samu, et al., Overview of developments in the MRCC program system, J. Phys. Chem. A, 2025, 129, 2086–2107.
71. M. Kallay, P. Nagy, D. Mester, L. Gyevi-Nagy, J. Csoka, P. Szabo, Z. Rolik, G. Samu, J. Csontos and B. Hegely, et al., Mrcc, a quantum chemical program suite written by, J. Chem. Phys, 2020, 152, 17407.
72. E. J. Poncelet, M. E. Lahm, A. G. Poncelet, J. M. Turney, M. A. Duncan, H. F. Schaefer III and Y. Abate, Vibrational Signatures of Unrealized Phosphorus Suboxide Intermediates in White Phosphorus Oxidation Reactions, J. Phys. Chem. A, 2026, 130, 121–128.
73. P. R. Schreiner, H. P. Reisenauer, D. Ley, D. Gerbig, C.-H. Wu and W. D. Allen, Methylhydroxycarbene: Tunneling Control of a Chemical Reaction, Science, 2011, 332, 1300–1303.
74. B. Myer, Low Temperature Spectroscopy: Optical Properties of Molecules in Matrices, Mixed Crystals, and Frozen Solutions, American Elsevier Publishing Company, 1971.
75. S. Bailleux, M. Bogey, C. Demuynck, Y. Liu and A. Walters, Millimeter-wave spectroscopy of PO in excited vibrational states up to v= 7, J. Mol. Spectrosc., 2002, 216, 465–471.
76. B. M. Cossairt, C. C. Cummins, A. R. Head, D. L. Lichtenberger, R. J. F. Berger, S. A. Hayes, N. W. Mitzel and G. Wu, On the Molecular and Electronic Structures of As P₃ and P₄, J. Am. Chem. Soc., 2010, 132, 84598465.
77. W. D. Allen and A. G. Csaszar, On the ab initio determination of higher-order force constants at nonstationary reference geometries, J. Chem. Phys., 1993, 98, 2983–3015.
78. F. Paw Lowski, A. Halkier, P. Jorgensen, K. L. Bak, T. Helgaker and W. Klopper, Accuracy of spectroscopic constants of diatomic molecules from ab initio calculations, J. Chem. Phys., 2003, 118, 2539–2549.
79. M. E. Lahm, N. L. Kitzmiller, H. F. Mull, W. D. Allen and H. F. Schaefer III, Concordant Mode Approach for Molecular Vibrations, J. Am. Chem. Soc., 2022, 144, 23271–23274.
80. N. L. Kitzmiller, M. E. Lahm, L. N. Olive Dornshuld, J. Jin, W. D. Allen and H. F. Schaefer III, Convergent Concordant Mode Approach for Molecular Vibrations: CMA-2, J. Chem. Theory Comput., 2024, 20, 10886–10898.
81. L. N. O. Dornshuld, M. E. Lahm, N. L. Kitzmiller, W. D. Allen and H. F. Schaefer III, Concordant Mode Approach (CMA): Vibrational Analysis of New and Upgraded Intermolecular Benchmarks for Noncovalent Bonding, J. Phys. Chem. A, 2026, 130, 3249–3260.
82. W. D. Allen, A. G. Csaszar and D. A. Horner, The puckering inversion barrier and vibrational spectrum of cyclopentene. A scaled quantum mechanical force field algorithm, J. Am. Chem. Soc., 1992, 114, 6834–6849.
83. P. Pulay and F. Torok, Uber die Parameterdarstellung der Kraftkonstantenmatrix, I. Acta Chim. Hung, 1965, 44, 287.
84. G. Keresztury and G. Jalsovszky, An alternative calculation of the vibrational potential energy distribution, J. Mol. Struct., 1971, 10, 304–305.
85. D. Peng and M. Reiher, Exact decoupling of the relativistic Fock operator, Theor. Chem. Acc., 2012, 131, 1081.
86. W. Liu and D. Peng, Exact two-component Hamiltonians revisited, J. Chem. Phys., 2009, 131.
87. P. Wang, S. Gong and Y. Mo, Bond dissociation energy of O₂ measured by fully state-to-state resolved threshold fragment yield spectra, J. Chem. Phys., 2024, 160, 164302.

---