Reference Energies for Cyclobutadiene: Automerization and Excited States (2204.05098v2)

Abstract: Cyclobutadiene is a well-known playground for theoretical chemists and is particularly suitable to test ground- and excited-state methods. Indeed, due to its high spatial symmetry, especially at the $D_{4h}$ square geometry but also in the $D_{2h}$ rectangular arrangement, the ground and excited states of cyclobutadiene exhibit multi-configurational characters and single-reference methods, such as adiabatic time-dependent density-functional theory (TD-DFT) or equation-of-motion coupled cluster (EOM-CC), are notoriously known to struggle in such situations. In this work, using a large panel of methods and basis sets, we provide an extensive computational study of the automerization barrier (defined as the difference between the square and rectangular ground-state energies) and the vertical excitation energies at $D_{2h}$ and $D_{4h}$ equilibrium structures. In particular, selected configuration interaction (SCI), multi-reference perturbation theory (CASSCF, CASPT2, and NEVPT2), and coupled-cluster (CCSD, CC3, CCSDT, CC4, and CCSDTQ) calculations are performed. The spin-flip formalism, which is known to provide a qualitatively correct description of these diradical states, is also tested within TD-DFT (combined with numerous exchange-correlation functionals) and the algebraic diagrammatic construction [ADC(2)-s, ADC(2)-x, and ADC(3)] schemes. A theoretical best estimate is defined for the automerization barrier and for each vertical transition energy.

• The paper establishes theoretical best estimates for cyclobutadiene, reporting an automerization barrier of 8.93 kcal/mol and key vertical excitation energies.
• The paper benchmarks multiple methods including coupled cluster, configuration interaction, and multi-reference approaches to tackle complex electronic structures.
• The paper provides practical insights on selecting computational methods and basis sets for accurately modeling antiaromatic systems and multi-configurational effects.

This paper presents a comprehensive benchmark paper evaluating various quantum chemistry methods for calculating the automerization barrier (AB) and vertical excitation energies of cyclobutadiene (CBD), a molecule known for its challenging electronic structure involving antiaromaticity, multi-configurational character, and states with significant double excitation contributions. Highly accurate theoretical best estimates (TBEs) are established using extrapolated coupled-cluster singles, doubles, triples, and quadruples (CCSDTQ) energies with the aug-cc-pVTZ basis set, providing valuable reference data.

Key Quantities Calculated:

1. Automerization Barrier (AB): The energy difference between the square ($D_{4h}$, transition state, ground state $^1B_{1g}$) and rectangular ($D_{2h}$, minimum, ground state $1^1A_g$) geometries. The TBE is determined to be 8.93 kcal/mol using geometries optimized at the CASPT2(12,12)/aug-cc-pVTZ level.
2. Vertical Excitation Energies at $D_{2h}$ Geometry: Calculated at the CC3/aug-cc-pVTZ optimized geometry of the $1^1A_g$ state.
   • $^3B_{1g}$ (TBE: 1.433 eV)
   • $^1B_{1g}$ (TBE: 3.125 eV)
   • $2^1A_g$ (double excitation, TBE: 4.038 eV)
3. Vertical Excitation Energies at $D_{4h}$ Geometry: Calculated at the RO-CCSD(T)/aug-cc-pVTZ optimized geometry of the $^3A_{2g}$ state (a minimum on the triplet surface). Energies are relative to the $^1B_{1g}$ ground state at this geometry.
   • $^3A_{2g}$ (TBE: 0.144 eV)
   • $^1A_{1g}$ (double excitation, TBE: 1.500 eV)
   • $^1B_{2g}$ (TBE: 1.849 eV)

Methods Benchmarked and Practical Implementation Insights:

The paper systematically evaluates several families of methods, offering practical guidance:

1. Selected Configuration Interaction (SCI - CIPSI):
   • Implementation: Used via the QUANTUM PACKAGE software. Employs iterative selection of determinants based on perturbative energy contribution, followed by extrapolation to the full CI (FCI) limit ($E_{PT2} \to 0$). State-averaged calculations were used for excited states.
   • Application: Provides near-FCI quality benchmarks for small molecules. Used here to validate the accuracy of the CCSDTQ-based TBEs.
   • Considerations: Computationally very expensive, limiting applicability to smaller systems or active spaces. The extrapolation introduces an error bar.
2. Coupled Cluster (CC - CCSD, CC3, CCSDT, CC4, CCSDTQ):
   • Implementation: Primarily used CFOUR, DALTON, and MRCC codes. EOM formalism used for excited states. For the $D_{4h}$ geometry, a closed-shell $^1A_g$ state was used as the reference, making the true $^1B_{1g}$ ground state appear as a de-excitation with double-excitation character relative to the reference.
   • Application: Offers a systematic hierarchy for improving accuracy.
   • Findings & Takeaways:
     • CCSD struggles with the multi-reference character and double excitations.
     • Inclusion of triples (CC3, CCSDT) is crucial for reasonable accuracy (errors often < 0.1-0.2 eV for single excitations, but larger for doubles or multi-reference cases). CC3 ($\order*{N^7}$) is a cheaper approximation to CCSDT ($\order*{N^8}$).
     • Inclusion of quadruples (CC4 approx. $\order*{N^9}$, CCSDTQ $\order*{N^{10}}$) is necessary for high accuracy (< 0.05 eV), especially for states with double excitation character ($2^1A_g$, $^1A_{1g}$) and the $D_{4h}$ ground state relative to the chosen reference.
     • CCSDTQ results align well with CIPSI, confirming its benchmark quality for these systems.
3. Multi-Reference Methods (CASSCF, CASPT2, NEVPT2):
   • Implementation: Used MOLPRO. State-averaged CASSCF for vertical excitations, state-specific for AB. Two active spaces explored: minimal (4e, 4o) $\pi$ space and extended (12e, 12o) including $\sigma_{CC}/\sigma^*_{CC}$ orbitals. Specific settings like level shifts (0.3 Ha) and IPEA shifts (0.25 Ha) were used for CASPT2. Both partially-contracted (PC) and strongly-contracted (SC) NEVPT2 were tested.
   • Application: Explicitly designed for multi-configurational systems.
   • Findings & Takeaways:
     • Active space selection is critical. The minimal (4e, 4o) space performs poorly, especially for ionic states ($^1B_{1g}$ at $D_{2h}$, $^1B_{2g}$ at $D_{4h}$) and leads to significant disagreement between CASPT2 and NEVPT2.
     • The larger (12e, 12o) space yields much better agreement between CASPT2 and NEVPT2 and improved accuracy overall. Practical Tip: Large discrepancies between CASPT2 and NEVPT2 often signal an inadequate active space.
     • Even with the larger space, accurately describing ionic states remains challenging. PC-NEVPT2 is theoretically slightly more rigorous than SC-NEVPT2.
     • CASSCF alone lacks dynamic correlation and provides only qualitative results.
4. Spin-Flip Methods (SF-TD-DFT, SF-ADC, SF-EOM-CCSD):
   • Implementation: Used Q-CHEM. Starts from a high-spin triplet reference state to access singlet ground and excited states (including doubles) via spin-flipping operators. Tested various TD-DFT functionals (B3LYP, PBE0, BH&HLYP, M06-2X, CAM-B3LYP, LC-$\omega$PBE08, $\omega$B97X-V, M11) within the Tamm-Dancoff Approximation (TDA). ADC variants SF-ADC(2)-s, SF-ADC(2)-x, SF-ADC(3), and the averaged SF-ADC(2.5) were used. SF-EOM-CCSD and its non-iterative triples corrections (dT, fT) were also assessed.
   • Application: Alternative approach for handling multi-reference character and double excitations without explicit multi-reference wave functions.
   • Findings & Takeaways:
     • SF-TD-DFT: Accuracy is highly dependent on the functional. Functionals with a high percentage (~50%) of short-range exact exchange (BH&HLYP, M06-2X, M11) perform best for AB and most excitations, but errors can still be significant (\SIrange{1}{4}{kcal/mol} for AB, up to \SI{0.5}{eV} or more for some excitations, especially at $D_{4h}$). Basis set convergence is notably faster than WFT methods. M06-2X was often the best performer among functionals.
     • SF-ADC: SF-ADC(2)-s ($\order*{N^5}$) offers a good cost-performance balance. SF-ADC(2)-x ($\order*{N^6}$) performs worse than SF-ADC(2)-s. SF-ADC(3) ($\order*{N^6}$) shows opposite error trends to SF-ADC(2)-s, making their average, SF-ADC(2.5), often highly accurate. Generally outperforms SF-TD-DFT.
     • SF-EOM-CCSD: Shows similar performance to SF-ADC(2)-s, often slightly overestimating excitation energies (~0.2 eV). The (dT) and (fT) triples corrections improve accuracy, with (dT) slightly better for excitations.
     • Spin contamination was generally found to be low for the states studied.

Computational Considerations:

• Basis Sets: Wave function methods typically require at least aug-cc-pVTZ for converged results (errors < ~0.05 eV or 1 kcal/mol). SF-TD-DFT converges faster, often achieving reasonable accuracy with aug-cc-pVDZ.
• Geometry: The choice of geometry and the level of theory used for optimization significantly impacts AB calculations. For vertical excitations, high-quality equilibrium geometries (e.g., CC3 or CCSD(T)) are preferred.
• Cost: Methods scale differently: SF-TD-DFT ~$\order*{N^3}$-$\order*{N^4}$, SF-ADC(2)-s ~$\order*{N^5}$, SF-ADC(2)-x/SF-ADC(3)/SF-EOM-CCSD ~$\order*{N^6}$, CC3 ~$\order*{N^7}$, CCSDT ~$\order*{N^8}$, CC4 ~$\order*{N^9}$, CCSDTQ ~$\order*{N^{10}}$, CASPT2/NEVPT2 depend heavily on active space size.

Practical Recommendations:

• For systems with potential multi-reference character or double excitations like CBD:
  • Standard single-reference methods (low-level CC, standard TD-DFT) are often inadequate.
  • SF-TD-DFT can be a cost-effective starting point, but requires careful functional selection (high short-range exact exchange, e.g., M06-2X) and validation, as accuracy can be limited.
  • SF-ADC(2)-s provides a better balance of cost ($\order*{N^5}$) and accuracy. SF-ADC(2.5) can yield excellent results if SF-ADC(3) is feasible.
  • Multi-reference methods (CASPT2/NEVPT2) are powerful but require careful active space selection, which is non-trivial. Use agreement between methods as a diagnostic.
  • Coupled Cluster: For high accuracy, CCSDT ($\order*{N^8}$) is often a good target, but challenging cases may necessitate CC4 ($\order*{N^9}$) or CCSDTQ ($\order*{N^{10}}$).

This paper provides valuable benchmarks and practical insights for choosing and implementing appropriate computational methods for electronically challenging molecular systems.