Graphics processing units accelerated semiclassical initial value representation molecular dynamics

A Perspective on the Investigation of Spectroscopy and Kinetics of Complex Molecular Systems with Semiclassical Approaches

In this Perspective we show that semiclassical methods provide a rigorous hierarchical way to study the vibrational spectroscopy and kinetics of complex molecular systems. The time averaged approach to spectroscopy and the semiclassical transition state theory for kinetics, which have been first adopted and then further developed in our group, provide accurate quantum results on rigorous physical grounds and can be applied even when dealing with a large number of degrees of freedom. In spectroscopy, the multiple coherent, divide-and-conquer, and adiabatically switched semiclassical approaches have practically permitted overcoming issues related to the convergence of results. In this Perspective we demonstrate the possibility of studying the semiclassical vibrational spectroscopy of a molecule adsorbed on an anatase (101) surface, a system made of 51 atoms. In kinetics, the semiclassical transition state theory is able to account for anharmonicity and the coupling between the reactive and bound modes. Our group has developed this technique for practical applications involving the study of phenomena like kinetic isotope effect, heavy atom tunneling, and elusive conformer lifetimes. Here, we show that our multidimensional anharmonic quantum approach is able to tackle on-the-fly the thermal kinetic rate constant of a 135 degree-of-freedom system. Overall, semiclassical methods open up the possibility to describe at the quantum mechanical level systems characterized by hundreds of degrees of freedom leading to the accurate spectroscopic and kinetic description of biomolecules and complex molecular systems.

Unsupervised Machine Learning Neural Gas Algorithm for Accurate Evaluations of the Hessian Matrix in Molecular Dynamics

The Hessian matrix of the potential energy of molecular systems is employed not only in geometry optimizations or high-order molecular dynamics integrators but also in many other molecular procedures, such as instantaneous normal mode analysis, force field construction, instanton calculations, and semiclassical initial value representation molecular dynamics, to name a few. Here, we present an algorithm for the calculation of the approximated Hessian in molecular dynamics. The algorithm belongs to the family of unsupervised machine learning methods, and it is based on the neural gas idea, where neurons are molecular configurations whose Hessians are adopted for groups of molecular dynamics configurations with similar geometries. The method is tested on several molecular systems of different dimensionalities both in terms of accuracy and computational time versus calculating the Hessian matrix at each time-step, that is, without any approximation, and other Hessian approximation schemes. Finally, the method is applied to the on-the-fly, full-dimensional simulation of a small synthetic peptide (the 46 atom N-acetyl-l-phenylalaninyl-l-methionine amide) at the level of DFT-B3LYP-D/6-31G* theory, from which the semiclassical vibrational power spectrum is calculated.

Representing molecular ground and excited vibrational eigenstates with nuclear densities obtained from semiclassical initial value representation molecular dynamics

We present in detail and validate an effective Monte Carlo approach for the calculation of the nuclear vibrational densities via integration of molecular eigenfunctions that we have preliminary employed to calculate the densities of the ground and the excited OH stretch vibrational states in the protonated glycine molecule [Aieta et al., Nat Commun 11, 4348 (2020)]. Here, we first validate and discuss in detail the features of the method on a benchmark water molecule. Then, we apply it to calculate on-the-fly the ab initio anharmonic nuclear densities in the correspondence of the fundamental transitions of NH and CH stretches in protonated glycine. We show how we can gain both qualitative and quantitative physical insight by inspection of different one-nucleus densities and assign a character to spectroscopic absorption peaks using the expansion of vibrational states in terms of harmonic basis functions. The visualization of the nuclear vibrations in a purely quantum picture allows us to observe and quantify the effects of anharmonicity on the molecular structure, also to exploit the effect of IR excitations on specific bonds or functional groups, beyond the harmonic approximation. We also calculate the quantum probability distribution of bond lengths, angles, and dihedrals of the molecule. Notably, we observe how in the case of one type of fundamental NH stretching, the typical harmonic nodal pattern is absent in the anharmonic distribution.

Anharmonic calculations of vibrational spectra for molecular adsorbates: A divide-and-conquer semiclassical molecular dynamics approach

The vibrational spectroscopy of adsorbates is becoming an important investigation tool for catalysis and material science. This paper presents a semiclassical molecular dynamics method able to reproduce the vibrational energy levels of systems composed by molecules adsorbed on solid surfaces. Specifically, we extend our divide-and-conquer semiclassical method for power spectra calculations to gas-surface systems and interface it with plane-wave electronic structure codes. The Born-Oppenheimer classical dynamics underlying the semiclassical calculation is full dimensional, and our method includes not only the motion of the adsorbate but also those of the surface and the bulk. The vibrational spectroscopic peaks related to the adsorbate are accounted together with the most coupled phonon modes to obtain spectra amenable to physical interpretations. We apply the method to the adsorption of CO, NO, and H₂O on the anatase-TiO₂ (101) surface. We compare our semiclassical results with the single-point harmonic estimates and the classical power spectra obtained from the same trajectory employed in the semiclassical calculation. We find that CO and NO anharmonic effects of fundamental vibrations are similarly reproduced by the classical and semiclassical dynamics and that H₂O adsorption is fully and properly described in its overtone and combination band relevant components only by the semiclassical approach.

Improved semiclassical dynamics through adiabatic switching trajectory sampling

We introduce an improved semiclassical dynamics approach to quantum vibrational spectroscopy. In this method, a harmonic-based phase space sampling is preliminarily driven toward non-harmonic quantization by slowly switching on the actual potential. The new coordinates and momenta serve as initial conditions for the semiclassical dynamics calculation, leading to a substantial decrease in the number of chaotic trajectories to deal with. Applications are presented for model and molecular systems of increasing dimensionality characterized by moderate or high chaoticity. They include a bidimensional Henon-Heiles potential, water, formaldehyde, and methane. The method improves accuracy and precision of semiclassical results and it can be easily interfaced with all pre-existing semiclassical theories.

Reduced rovibrational coupling Cartesian dynamics for semiclassical calculations: Application to the spectrum of the Zundel cation

We study the vibrational spectrum of the protonated water dimer, by means of a divide-and-conquer semiclassical initial value representation of the quantum propagator, as a first step in the study of larger protonated water clusters. We use the potential energy surface from the work of Huang et al. [J. Chem. Phys. 122, 044308 (2005)]. To tackle such an anharmonic and floppy molecule, we employ fully Cartesian dynamics and carefully reduce the coupling to global rotations in the definition of normal modes. We apply the time-averaging filter and obtain clean power spectra relative to suitable reference states that highlight the spectral peaks corresponding to the fundamental excitations of the system. Our trajectory-based approach allows for the physical interpretation of the very challenging proton transfer modes. We find that it is important, for such a floppy molecule, to selectively avoid initially exciting lower energy modes, in order to obtain cleaner spectra. The estimated vibrational energies display a mean absolute error (MAE) of ∼29 cm^-1 with respect to available multiconfiguration time-dependent Hartree calculations and MAE ∼ 14 cm^-1 when compared to the optically active experimental excitations of the Ne-tagged Zundel cation. The reasonable scaling in the number of trajectories for Monte Carlo convergence is promising for applications to higher dimensional protonated cluster systems.

Vibrational investigation of nucleobases by means of divide and conquer semiclassical dynamics

In this work, we report a computational study of the vibrational features of four different nucleobases employing the divide-and-conquer semiclassical initial value representation molecular dynamics method. Calculations are performed on uracil, cytosine, thymine, and adenine. Results show that the overall accuracy with respect to experiments is within 20 wavenumbers, regardless of the dimensionality of the nucleobase. Vibrational estimates are accurate even in the complex case of cytosine, where two relevant conformers are taken into account. These results are promising in the perspective of future studies on more complex systems, such as nucleotides or nucleobase pairs.

Semiclassical vibrational spectroscopy with Hessian databases

We report on a new approach to ease the computational overhead of ab initio "on-the-fly" semiclassical dynamics simulations for vibrational spectroscopy. The well known bottleneck of such computations lies in the necessity to estimate the Hessian matrix for propagating the semiclassical pre-exponential factor at each step along the dynamics. The procedure proposed here is based on the creation of a dynamical database of Hessians and associated molecular geometries able to speed up calculations while preserving the accuracy of results at a satisfactory level. This new approach can be interfaced to both analytical potential energy surfaces and on-the-fly dynamics, allowing one to study even large systems previously not achievable. We present results obtained for semiclassical vibrational power spectra of methane, glycine, and N-acetyl-L-phenylalaninyl-L-methionine-amide, a molecule of biological interest made of 46 atoms.

"Divide and conquer" semiclassical molecular dynamics: A practical method for spectroscopic calculations of high dimensional molecular systems

We extensively describe our recently established "divide-and-conquer" semiclassical method [M. Ceotto, G. Di Liberto, and R. Conte, Phys. Rev. Lett. 119, 010401 (2017)] and propose a new implementation of it to increase the accuracy of results. The technique permits us to perform spectroscopic calculations of high-dimensional systems by dividing the full-dimensional problem into a set of smaller dimensional ones. The partition procedure, originally based on a dynamical analysis of the Hessian matrix, is here more rigorously achieved through a hierarchical subspace-separation criterion based on Liouville's theorem. Comparisons of calculated vibrational frequencies to exact quantum ones for a set of molecules including benzene show that the new implementation performs better than the original one and that, on average, the loss in accuracy with respect to full-dimensional semiclassical calculations is reduced to only 10 wavenumbers. Furthermore, by investigating the challenging Zundel cation, we also demonstrate that the "divide-and-conquer" approach allows us to deal with complex strongly anharmonic molecular systems. Overall the method very much helps the assignment and physical interpretation of experimental IR spectra by providing accurate vibrational fundamentals and overtones decomposed into reduced dimensionality spectra.

Graphics processing units in bioinformatics, computational biology and systems biology

Several studies in Bioinformatics, Computational Biology and Systems Biology rely on the definition of physico-chemical or mathematical models of biological systems at different scales and levels of complexity, ranging from the interaction of atoms in single molecules up to genome-wide interaction networks. Traditional computational methods and software tools developed in these research fields share a common trait: they can be computationally demanding on Central Processing Units (CPUs), therefore limiting their applicability in many circumstances. To overcome this issue, general-purpose Graphics Processing Units (GPUs) are gaining an increasing attention by the scientific community, as they can considerably reduce the running time required by standard CPU-based software, and allow more intensive investigations of biological systems. In this review, we present a collection of GPU tools recently developed to perform computational analyses in life science disciplines, emphasizing the advantages and the drawbacks in the use of these parallel architectures. The complete list of GPU-powered tools here reviewed is available at http://bit.ly/gputools.

Semiclassical "Divide-and-Conquer" Method for Spectroscopic Calculations of High Dimensional Molecular Systems

A new semiclassical "divide-and-conquer" method is presented with the aim of demonstrating that quantum dynamics simulations of high dimensional molecular systems are doable. The method is first tested by calculating the quantum vibrational power spectra of water, methane, and benzene-three molecules of increasing dimensionality for which benchmark quantum results are available-and then applied to C_{60}, a system characterized by 174 vibrational degrees of freedom. Results show that the approach can accurately account for quantum anharmonicities, purely quantum features like overtones, and the removal of degeneracy when the molecular symmetry is broken.

On-the-Fly ab Initio Semiclassical Calculation of Glycine Vibrational Spectrum

We present an on-the-fly ab initio semiclassical study of vibrational energy levels of glycine, calculated by Fourier transform of the wavepacket correlation function. It is based on a multiple coherent states approach integrated with monodromy matrix regularization for chaotic dynamics. All four lowest-energy glycine conformers are investigated by means of single-trajectory semiclassical spectra obtained upon classical evolution of on-the-fly trajectories with harmonic zero-point energy. For the most stable conformer I, direct dynamics trajectories are also run for each vibrational mode with energy equal to the first harmonic excitation. An analysis of trajectories evolved up to 50 000 atomic time units demonstrates that, in this time span, conformers II and III can be considered as isolated species, while conformers I and IV show a pretty facile interconversion. Therefore, previous perturbative studies based on the assumption of isolated conformers are often reliable but might be not completely appropriate in the case of conformer IV and conformer I for which interconversion occurs promptly.

The importance of the pre-exponential factor in semiclassical molecular dynamics

This paper deals with the critical issue of approximating the pre-exponential factor in semiclassical molecular dynamics. The pre-exponential factor is important because it accounts for the quantum contribution to the semiclassical propagator of the classical Feynman path fluctuations. Pre-exponential factor approximations are necessary when chaotic or complex systems are simulated. We introduced pre-exponential factor approximations based either on analytical considerations or numerical regularization. The approximations are tested for power spectrum calculations of more and more chaotic model systems and on several molecules, for which exact quantum mechanical values are available. The results show that the pre-exponential factor approximations introduced are accurate enough to be safely employed for semiclassical simulations of complex systems.

Improved scaling of time-evolving block-decimation algorithm through reduced-rank randomized singular value decomposition

When the amount of entanglement in a quantum system is limited, the relevant dynamics of the system is restricted to a very small part of the state space. When restricted to this subspace the description of the system becomes efficient in the system size. A class of algorithms, exemplified by the time-evolving block-decimation (TEBD) algorithm, make use of this observation by selecting the relevant subspace through a decimation technique relying on the singular value decomposition (SVD). In these algorithms, the complexity of each time-evolution step is dominated by the SVD. Here we show that, by applying a randomized version of the SVD routine (RRSVD), the power law governing the computational complexity of TEBD is lowered by one degree, resulting in a considerable speed-up. We exemplify the potential gains in efficiency at the hand of some real world examples to which TEBD can be successfully applied and demonstrate that for those systems RRSVD delivers results as accurate as state-of-the-art deterministic SVD routines.

Integrated PK-PD and agent-based modeling in oncology

Mathematical modeling has become a valuable tool that strives to complement conventional biomedical research modalities in order to predict experimental outcome, generate new medical hypotheses, and optimize clinical therapies. Two specific approaches, pharmacokinetic-pharmacodynamic (PK-PD) modeling, and agent-based modeling (ABM), have been widely applied in cancer research. While they have made important contributions on their own (e.g., PK-PD in examining chemotherapy drug efficacy and resistance, and ABM in describing and predicting tumor growth and metastasis), only a few groups have started to combine both approaches together in an effort to gain more insights into the details of drug dynamics and the resulting impact on tumor growth. In this review, we focus our discussion on some of the most recent modeling studies building on a combined PK-PD and ABM approach that have generated experimentally testable hypotheses. Some future directions are also discussed.