Quantum Monte Carlo Simulations of the Vibrational Wavefunction of the Aromatic Cyclo[10]carbon Using a Full Dimensional Permutationally Invariant Potential Energy Surface.

New experimental measurements [Sun et al., Nature 2023, 623, 972] of the cyclic C10 reveal a cumulenic pentagon-like D5h structure at ∼5 K. However, the long-standing presumption that a large zero-point vibrational energy combined with an extremely flat D5h ↔ D10h ↔ D5h isomerization pathway washes out the pentagonal D5h structure and yields a symmetric D10h decagon remains at odds with the experiment. We resolve this issue with our fitting approach based on a bond-order charge-density matrix expressed in permutationally invariant polynomials. We train the model on τHCTH/cc-pVQZ data morphed to reproduce a relativistic all-electron CCSDT(Q)/CBS D5h-D10h potential energy barrier (benchmarked previously by others). Large scale diffusion Monte Carlo simulations in full dimensionality show that the vibrational ground state of C10 has compositional character of more than 96% D5h, fully reflecting the experimental imaging data. Quantum mechanical variational calculations in 1-D further suggest persistence of the D5h symmetry structure at higher temperatures.

Fitting Potential Energy Surfaces by Learning the Charge Density Matrix with Permutationally Invariant Polynomials.

The electronic energy in the Hartree-Fock (HF) theory is the trace of the product of the charge density matrix (CDM) with the one-electron and two-electron matrices represented in an atomic orbital basis, where the two-electron matrix is also a function of the same CDM. In this work, we examine a formalism of analytic representation of a generic molecular potential energy surface (PES) as a sum of a linearly parameterized HF and a correction term, the latter formally representing the electron correlation energy, also linearly parameterized, by expressing the elements of CDM using permutationally invariant polynomials (PIPs). We show on a variety of numerical examples, ranging from exemplary two-electron systems HeH+ and H3+ to the more challenging cases of methanium (CH5+) fragmentation and high-energy tautomerization of formamide to formimidic acid that such a formulation requires significantly fewer, 10-20% of PIPs, to accomplish the same accuracy of the fit as the conventional representation at practically the same computational cost.

A Fermi resonance and a parallel-proton-transfer overtone in the Raman spectrum of linear centrosymmetric N4H+: A polarizability-driven first principles molecular dynamics study.

We present molecular dynamics (MD), polarizability driven MD (α-DMD), and pump-probe simulations of Raman spectra of the protonated nitrogen dimer N4H+, and some of its isotopologues, using the explicitly correlated coupled-cluster singles and doubles with perturbative triples [CCSD(T)]-F12b/aug-cc-pVTZ based potential energy surface in permutationally invariant polynomials (PIPs) of Yu et al. [J. Phys. Chem. A 119, 11623 (2015)] and a corresponding PIP-derived CCSD(T)/aug-cc-pVTZ-tr (N:spd, H:sp) polarizability tensor surface (PTS), the latter reported here for the first time. To represent the PTS in terms of a PIP basis, we utilize a recently described formulation for computing the polarizability using a many-body expansion in the orders of dipole-dipole interactions while generating a training set using a novel approach based on linear regression for potential energy distributions. The MD/α-DMD simulations reveal (i) a strong Raman activity at 260 and 2400 cm-1, corresponding to the symmetric N-N⋯H bend and symmetric N-N stretch modes, respectively; (ii) a very broad spectral region in the 500-2000 cm-1 range, assignable to the parallel N⋯H+⋯N proton transfer overtone; and (iii) the presence of a Fermi-like resonance in the Raman spectrum near 2400 cm-1 between the Σg + N-N stretch fundamental and the Πu overtone corresponding to perpendicular N⋯H+⋯N proton transfer.

Permutationally invariant polynomial representation of polarizability tensor surfaces for linear regression analysis.

A linearly parameterized functional form for a Cartesian representation of molecular dipole polarizability tensor surfaces (PTS) is described. The proposed expression for the PTS is a linearization of the recently reported power series ansatz of the original Applequist model, which by construction is non-linear in parameter space. This new approach possesses (i) a unique solution to the least-squares fitting problem; (ii) a low level of the computational complexity of the resulting linear regression procedure, comparable to those of the potential energy and dipole moment surfaces; and (iii) a competitive level of accuracy compared to the non-linear PTS model. Calculations of CH4 PTS, with polarizabilities fitted to 9000 training set points with the energies up to 14,000 cm-1 show an impressive level of accuracy of the linear PTS model obtained with ~1600 parameters: ~1% versus 0.3% RMSE for the non-linear vs. linear model on a test set of 1000 configurations.

Intramolecular Proton Transfer in the Hydrogen Oxalate Anion and the Cooperativity Effects of the Low-Frequency Vibrations: A Driven Molecular Dynamics Study.

We report first-principles molecular dynamics (MD) and dipole-driven molecular dynamics (µ-DMD) simulations of the hydrogen oxalate anion at the MP2/aug-cc-pVDZ level of theory. We examine the role of vibrational coupling between the OH stretching bands, that is, the fundamental and a few combination bands spanning the 2900-3100 cm-1 range, and several of the low-frequency bending and stretching fundamental modes. The low-frequency modes between 300 and 825 cm-1 play a crucial role in the proton-transfer motion. Strong involvement of CO2 and CCO bending and the CC stretching vibrations indicate that these large amplitude motions cause the shortening of the O···O distance and thus promote H+ transfer to the other oxygen by bringing it over the 3.4 kcal/mol barrier. Analysis of resonant µ-DMD trajectories shows that the complex spectral feature near 825 cm-1, closely corresponding to both an overtone of two quanta of 425 cm-1 and a combination band of low-frequency CO2 rocking (300 cm-1) and CCO bending (575 cm-1) modes, is involved in the proton transfer. µ-DMD shows that exciting the system at these mode combinations leads to faster barrier activation than exciting at the OH fundamental mode.

On the Cartesian Representation of the Molecular Polarizability Tensor Surface by Polynomial Fitting to Ab Initio Data.

We describe an approach to constructing an analytic Cartesian representation of the molecular dipole polarizability tensor surface in terms of polynomials in interatomic distances with a training set of ab initio data points obtained from a molecular dynamics (MD) simulation or by any other available means. The proposed formulation is based on a perturbation treatment of the unmodified point dipole polarizability model of Applequist [ J. Am. Chem. Soc. 1972, 94, 2952] and is shown here to be, by construction (i) free of short-range or other singularities or discontinuities, (ii) symmetric and translationally invariant, and (iii) nonreliant on a body-fixed coordinate system. Permutational invariance of like nuclei is demonstrated to be readily applicable, making this approach useful for highly fluxional and reactive systems. Derivation of the method is described in detail, adding brief didactic numerical examples of H2 and H2O and concluding with an MD simulation of the Raman spectrum of H5O2+ at 300 K with the polarizability tensor fitted to CCSD(T)/aug-cc-pVTZ data obtained using the HBB-4B potential [ J. Chem. Phys. 2005, 122, 044308].

Analysis of the Proton Transfer Bands in the Infrared Spectra of Linear N2H+···OC and N2D+···OC Complexes Using Electric Field-Driven Classical Trajectories.

In this work, we describe ab initio calculations and assignment of infrared (IR) spectra of hydrogen-bonded ion-molecular complexes that involve a fluxional proton: the linear N2H+···OC and N2D+···OC complexes. Given the challenges of describing fluxional proton dynamics and especially its IR activity, we use electric field-driven classical trajectories, i.e., the driven molecular dynamics (DMD) method that was developed by us in recent years and for similar applications, in conjunction with high-level electronic structure theory. Namely, we present a modified and a numerically efficient implementation of DMD specifically for direct (or "on the fly") calculations, which we carry out at the MP2-F12/AVDZ level of theory for the potential energy surface (PES) and MP2/AVDZ for the dipole moment surfaces (DMSs). Detailed analysis of the PES, DMS, and the time-dependence of the first derivative of the DMS, referred to as the driving force, for the highly fluxional vibrations involving H+/D+ revealed that the strongly non-harmonic PES and non-linear DMS yield remarkably complex vibrational spectra. Interestingly, the classical trajectories reveal a doublet in the proton transfer part of the spectrum with the two peaks at 1800 and 1980 cm-1. We find that their shared intensity is due to a Fermi-like resonance interaction, within the classical limit, of the H+ parallel stretch fundamental and an H+ perpendicular bending overtone. This doublet is also observed in the deuterated species at 1360 and 1460 cm-1.

Assignment of Infrared-Active Combination Bands in the Vibrational Spectra of Protonated Molecular Clusters Using Driven Classical Trajectories: Application to N4H+ and N4D.

We investigate the utility of the driven molecular dynamics (DMD) approach to complex molecular vibrations by applying it to linear clusters with several degenerate vibrational modes and infrared (IR) intense combination bands. Here, the prominent features in N4H+ and N4D+ IR spectra, reported and described by others previously, have been characterized for the first time by DMD using recently published high-level potential and dipole moment surfaces. Namely, the calculations closely correlate the parallel proton stretch vibration in N4H+, at 750 cm-1, with the one observed experimentally at 743 cm-1. Second, the intense IR-active combination bands found in experimental spectra within 900-1100 cm-1 have been properly recovered by DMD at 950 cm-1 as strongly IR-active and confirmed as consisting of H+ asymmetric stretch and N2···N2 intermolecular symmetric stretch modes. Furthermore, we show that certain combination bands involving overtone transitions may be recovered by DMD using a hard-driving regime, such as the 1409 cm-1 band measured in N4H+, revealed by DMD at 1375 cm-1, and assigned to a progressive combination of the parallel H+ stretch and two quanta of N2···N2 stretch, in agreement with quantum mechanical studies reported previously by others.

Deconstructing Prominent Bands in the Terahertz Spectra of H7O3+ and H9O4+: Intermolecular Modes in Eigen Clusters.

We report experimental vibrational action spectra (210-2200 cm-1) and calculated IR spectra, using recent ab initio potential energy and dipole moment surfaces, of H7O3+ and H9O4+. We focus on prominent far-IR bands, which postharmonic analyses show, arise from two types of intermolecular motions, i.e., hydrogen bond stretching and terminal water wagging modes, that are similar in both clusters. The good agreement between experiment and theory further validates the accuracy of the potential and dipole moment surfaces, which was used in a recent theoretical study that concluded that infrared photodissociation spectra of the cold clusters correspond to the Eigen isomer. The comparison between theory and experiment adds further confirmation of the need of postharmonic analysis for these clusters.

Driven molecular dynamics studies of the shared proton motion in the H5O2+·Ar cluster: the effect of argon tagging and deuteration on vibrational spectra.

We report IR spectra of H5O2(+) and H5O2(+)·Ar and their deuterium isotopologues using ab initio molecular dynamics. The trajectories were propagated as microcanonical (NVE) ensembles at energies corresponding to temperatures 50 and 100 K. The potential energy surface is calculated on-the-fly at the MP2/aug-cc-pVDZ level of theory. The calculations show that adding an argon atom to H5O2(+) introduces symmetry breaking in the Zundel core ion, causes blueshift in the shared proton vibration by about 200 cm(-1), and leads to the splitting of the OH stretch vibrations into four bands. Driven molecular dynamics (DMD) method is used to assign the spectrum by coupling the dipole moment to an external electric field oscillating at frequency ω. The broad feature at 1100 cm(-1) in the H5O2(+)·Ar spectrum is ascribed to the large amplitude shared proton vibration coupled with torsion and wag modes. MD MP2 simulations predict the H/D redshift in the shared proton vibration and water bending vibration to be about 280 and 460 cm(-1), respectively, in good agreement with experimental observations.

Ab Initio Studies of Structural and Vibrational Properties of Protonated Water Cluster H7O3(+) and Its Deuterium Isotopologues: An Application of Driven Molecular Dynamics.

In this work, we present infrared (IR) spectra of H7O3(+) and its deuterium isotopomers calculated by direct molecular dynamics (MD) simulations at the B3LYP/6-31+G** computational level. The calculated spectra obtained at 100, 300, and 500 K were compared to available experimental observations, and spectral features were assigned using normal-mode analysis (NMA) and driven molecular dynamics (DMD). Spectral peaks at 2410 and 2540 cm(-1) were assigned to asymmetric and symmetric stretches of the bridging hydrogen (BH) using NMA. The weak spectral features at 2166 and 2275 cm(-1) were assigned to a combination band of BH asymmetric stretch, H2O in phase wagging, OO stretch, and H3O(+) rocking vibrations by DMD simulations. Our observation of BH stretch vibrations as low as 2166 cm(-1) is in good agreement with the assignment of the low-resolution spectrum obtained by Schwarz at 2200-2300 cm(-1) [Schwarz, H. A. J. Chem. Phys. 1977, 67, 5525-5534] and vibrational predissociation spectrum by Lee et al. ∼2300 cm(-1) [Okumura, M.; Yeh, L. I.; Myers, J. D.; Lee, Y. T. J. Chem. Phys. 1990, 94, 3416-3427].

Calculation of the vibrational spectra of H5O2(+) and its deuterium-substituted isotopologues by molecular dynamics simulations.

In this work, we present infrared spectra of H(5)O(2)(+) and its D(5)O(2)(+), D(4)HO(2)(+), and DH(4)O(2)(+) isotopologues calculated by classical molecular dynamics simulations on an accurate potential energy surface generated from CCSD(T) calculations, as well as on the BLYP DFT potential energy surface sampled by means of the Car-Parrinello algorithm. The calculated spectra obtained with internal energies corresponding to a temperature of about 30 K are in overall good agreement with those from experimental measurements and from quantum dynamical simulations.

Ab Initio Molecular Dynamics Simulations of the Infrared Spectra of H3O2(-) and D3O2(.).

We present the infrared spectra of H3O2(-) and D3O2(-) calculated using MP2 direct molecular dynamics approach at temperatures of 100 and 300 K. The spectral peaks were assigned using the normal-mode analysis, instantaneous normal-mode analysis, isotopic substitution, polarized infrared absorptions, and analysis of the position-position correlation function. Our results predict the bridging hydrogen stretch between 600 and 900 cm(-1) and bridging hydrogen bend vibrations between 1250 and 1650 cm(-1). We also examine two DFT methods (B3PW91 and B3LYP) and report on the differences between them and the MP2 spectra.

Vibrational analysis of the H5O2+ infrared spectrum using molecular and driven molecular dynamics.

Standard molecular and driven molecular dynamics are used to analyze prominent spectral features in the H5O2+ infrared spectrum. In the driven method, the molecular Hamiltonian is augmented with a time-dependent term, mu x epsilon(0) sin(omegat), where mu is the dipole moment of H5O2+, epsilon0 is the electric field, and omega is the frequency. The magnitude of the electric field determines whether the driving is mild (the harmonic limit) or strong (anharmonic motion and mode coupling). We analyze the spectrum in the wavenumber range from 600 to 1900 cm(-1), where recent experimental measurements are available for H5O2+. On the basis of the simulations, we have assigned the broad feature around 1000 cm(-1) to the proton transfer coupled with the torsion motion. Intense absorption near 1780 cm(-1) is assigned to the H2O monomer bend coupled with proton transfer.

All-Atom Calculation of the Normal Modes of Bacteriorhodopsin Using a Sliding Block Iterative Diagonalization Method.

Conventional normal-mode analysis of molecular vibrations requires computation and storage of the Hessian matrix. For a typical biological system such storage can reach several gigabytes posing difficulties for straightforward implementation. In this work we discuss an iterative block method to carry out full diagonalization of the Hessian while only storing a few vectors in memory. The iterative approach is based on the conjugate gradient formulation of the Davidson algorithm for simultaneous optimization of L roots, where in our case 10 < L < 300. The procedure is modified further by automatically adding a new vector into the search space for each locked (converged) root and keeping the new vector orthogonal to the eigenvectors previously determined. The higher excited states are then converged with the orthonormality constraint to the locked roots by applying a projector which is carried out using a read-rewind step done once per iteration. This allows for convergence of as many roots as desired without increasing the computer memory. The required Hessian-vector products are calculated on the fly as follows, Kp = dgp/dt, where K is the mass weighted Hessian, and gp is the gradient along p. The method has been implemented into the TINKER suite of molecular design codes. Preliminary results are presented for the normal modes of bacteriorhodopsin (bR) up to 300 cm(-)(1) and for the high frequency range between 2840 and 3680 cm(-)(1). There is evidence of a highly localized, noncollective mode at ∼1.4 cm(-)(1), caused by long-range interactions acting between the cytoplasmic and extracellular domains of bR.

Normal mode analysis using the driven molecular dynamics method. II. An application to biological macromolecules.

The driven molecular-dynamics (DMD) method, recently proposed by Bowman, Zhang, and Brown [J. Chem. Phys. 119, 646 (2003)], has been implemented into the TINKER molecular modeling program package. The DMD method yields frequencies and normal modes without evaluation of the Hessian matrix. It employs an external harmonic driving term that can be used to scan the spectrum and determine resonant absorptions. The molecular motions, induced by driving at resonant frequencies, correspond to the normal-mode vibrations. In the current work we apply the method to a 20-residue protein, Trp-cage, that has been reported by Neidigh, Fesinmeyer, and Andersen [Nature Struct. Biol. 9, 425 (2002)]. The structural and dynamical properties of this molecule, such as B-factors, root-mean square fluctuations, anisotropies, vibrational entropy, and cross-correlations coefficients, are calculated using the DMD method. The results are in very good agreement with ones calculated using standard normal-mode analysis method. Thus, the DMD method provides a viable alternative to the standard Hessian-based method and has considerable potential for the study of large systems, where the Hessian-based method is not feasible.