# Molecular Dynamics

Submitted
Shpiro, B. ; Fabian, M. D. ; Rabani, E. ; Baer, R. Forces from stochastic density functional theory under nonorthogonal atom-centered basis sets. arXiv:2108.06770 [physics] Submitted. Publisher's VersionAbstract

We develop a formalism for calculating forces on the nuclei within the linear-scaling stochastic density functional theory (sDFT) in a nonorthogonal atom-centered basis-set representation (Fabian et al. WIREs Comput Mol Sci. 2019;e1412. https://doi.org/10.1002/wcms.1412) and apply it to Tryptophan Zipper 2 (Trp-zip2) peptide solvated in water. We use an embedded-fragment approach to reduce the statistical errors (fluctuation and systematic bias), where the entire peptide is the main fragment and the remaining 425 water molecules are grouped into small fragments. We analyze the magnitude of the statistical errors in the forces and find that the systematic bias is of the order of \$0.065\textbackslash,eV/\textbackslashr\A\\$ (\$\textbackslashsim1.2\textbackslashtimes10ˆ\-3\E\_\h\/a\_\0\\$) when 120 stochastic orbitals are used, independently of systems size. This magnitude of bias is sufficiently small to ensure that the bond lengths estimated by stochastic DFT (within a Langevin molecular dynamics simulation) will deviate by less than 1% from those predicted by a deterministic calculation.

2021
Gope, K. ; Livshits, E. ; Bittner, D. M. ; Baer, R. ; Strasser, D. Two pathways and an isotope effect in H3+ formation following double ionization of methanol. Natural Sciences 2021, ntls.10022. Publisher's VersionAbstract

The trihydrogen ion has a central role in creating complex molecules in the interstellar medium. Therefore, its formation and destruction mechanisms in high photon energy environments involving organic molecules are drawing significant experimental and theoretical attention. Here, we employ a combination of time-resolved ultrafast extreme-ultraviolet pump and near-infrared probe spectroscopy applied to the deuterated CH3OD methanol molecule. Similar to other double-ionization studies, the isotopic labeling reveals two competing pathways for forming trihydrogen: A) H+3 + COH+ and B) H+3 + HCO+. We validate our high-level ab initio nonadiabatic molecular dynamic simulations by showing that it closely reproduces the essential features of the measured kinetic energy release distribution and branching ratios of the two pathways of the deuterated system. The success of ab initio simulation in describing single photon double-ionization allows for an unprecedented peek into the formation pathways for the undeuterated species, beyond present experimental reach. For this case, we find that the kinetic energy release of pathway B shifts to lower energies by more than 0.6 eV due to a dynamical isotope effect. We also determine the mechanism for trihydrogen formation from excited states of the dication and elucidate the isotope effect’s role in the observed dynamics.

2020
Livshits, E. ; Luzon, I. ; Gope, K. ; Baer, R. ; Strasser, D. Time-resolving the ultrafast H2 roaming chemistry and H3+ formation using extreme-ultraviolet pulses. Communications Chemistry 2020, 3 49. Publisher's Version
Gope, K. ; Livshits, E. ; Bittner, D. M. ; Baer, R. ; Strasser, D. Absence of Triplets in Single-Photon Double Ionization of Methanol. The Journal of Physical Chemistry Letters 2020, 11, 8108–8113. Publisher's VersionAbstract

Despite the abundance of data concerning single-photon double ionization of methanol, the spin state of the emitted electron pair has never been determined. Here we present the ﬁrst evidence that identiﬁes the emitted electron pair spin as overwhelmingly singlet when the dication forms in low-energy conﬁgurations. The experimental data show that while the yield of the CH2O+ + H3+ Coulomb explosion channel is abundant, the metastable methanol dication is largely absent. According to high-level ab initio simulations, these facts indicate that photoionization promptly forms singlet dication states, where they quickly decompose through various channels, with signiﬁcant H3+ yields on the low-lying states. In contrast, if we assume that the initial dication is formed in one of the low-lying triplet states, the ab initio simulations exhibit a metastable dication, contradicting the experimental ﬁndings. Comparing the average simulated branching ratios with the experimental data suggests a \textgreater3 order of magnitude enhancement of the singlet:triplet ratio compared with their respective 1:3 multiplicities.

Arnon, E. ; Rabani, E. ; Neuhauser, D. ; Baer, R. Efficient Langevin dynamics for "noisy" forces. J. Chem. Phys. 2020, 152, 161103. Publisher's VersionAbstract

Efficient Boltzmann-sampling using first-principles methods is challenging for extended systems due to the steep scaling of electronic structure methods with the system size. Stochastic approaches provide a gentler system-size dependency at the cost of introducing "noisy" forces, which serve to limit the efficiency of the sampling. In the first-order Langevin dynamics (FOLD), efficient sampling is achievable by combining a well-chosen preconditioning matrix S with a time-step-bias-mitigating propagator (Mazzola et al., Phys. Rev. Lett., 118, 015703 (2017)). However, when forces are noisy, S is set equal to the force-covariance matrix, a procedure which severely limits the efficiency and the stability of the sampling. Here, we develop a new, general, optimal, and stable sampling approach for FOLD under noisy forces. We apply it for silicon nanocrystals treated with stochastic density functional theory and show efficiency improvements by an order-of-magnitude.

2019
Cytter, Y. ; Rabani, E. ; Neuhauser, D. ; Preising, M. ; Redmer, R. ; Baer, R. Transition to metallization in warm dense helium-hydrogen mixtures using stochastic density functional theory within the Kubo-Greenwood formalism. Physical Review B 2019, 100. Publisher's VersionAbstract

Abstract The Kubo-Greenwood (KG) formula is often used in conjunction with Kohn-Sham (KS) density functional theory (DFT) to compute the optical conductivity, particularly for warm dense mater. For applying the KG formula, all KS eigenstates and eigenvalues up to an energy cutoff are required and thus the approach becomes expensive, especially for high temperatures and large systems, scaling cubically with both system size and temperature. Here, we develop an approach to calculate the KS conductivity within the stochastic DFT (sDFT) framework, which requires knowledge only of the KS Hamiltonian but not its eigenstates and values. We show that the computational effort associated with the method scales linearly with system size and reduces in proportion to the temperature unlike the cubic increase with traditional deterministic approaches. In addition, we find that the method allows an accurate description of the entire spectrum, including the high-frequency range, unlike the deterministic method which is compelled to introduce a high-frequency cut-off due to memory and computational time constraints. We apply the method to helium-hydrogen mixtures in the warm dense matter regime at temperatures of \textbackslashsim60\textbackslashtext\kK\ and find that the system displays two conductivity phases, where a transition from non-metal to metal occurs when hydrogen atoms constitute \textbackslashsim0.3 of the total atoms in the system.

Luzon, I. ; Livshits, E. ; Gope, K. ; Baer, R. ; Strasser, D. Making Sense of Coulomb Explosion Imaging. J. Phys. Chem. Lett. 2019, 10, 1361–1367. Publisher's VersionAbstract

A multifaceted agreement between ab initio theoretical predictions and experimental measurements, including branching ratios, channel-speciﬁc kinetic energy release, and three-body momentum correlation spectra, leads to the identiﬁcation of new mechanisms in Coulomb-explosion (CE) induced two- and three-body breakup processes in methanol. These identiﬁed mechanisms include direct nonadiabatic Coulomb explosion responsible for CO bond-breaking, a long-range “ inverse harpooning” dominating the production of H2+ + HCOH+, a transient proton migration leading to surprising energy partitioning in three-body fragmentation and other complex dynamics forming products such as H2O+ and H3+. These mechanisms provide general concepts that should be useful for analyzing future time-resolved Coulomb explosion imaging of methanol as well as other molecular systems. These advances are enabled by a combination of recently developed experimental and computational techniques, using weak ultrafast EUV pulses to initiate the CE and a high-level quantum chemistry approach to follow the resulting ﬁeld-free nonadiabatic molecular dynamics.

2017
Luzon, I. ; Jagtap, K. ; Livshits, E. ; Lioubashevski, O. ; Baer, R. ; Strasser, D. Single-photon Coulomb explosion of methanol using broad bandwidth ultrafast EUV pulses. Phys. Chem. Chem. Phys. 2017, 19, 13488–13495.Abstract

Single-photon Coulomb explosion of methanol is instigated using the broad bandwidth pulse achieved through high-order harmonics generation. Using 3D coincidence fragment imaging of one molecule at a time, the kinetic energy release (KER) and angular distributions of the products are measured in different Coulomb explosion (CE) channels. Two-body CE channels breaking either the C–O or the C–H bonds are described as well as a proton migration channel forming H2O+, which is shown to exhibit higher KER. The results are compared to intense-field Coulomb explosion measurements in the literature. The interpretation of broad bandwidth single-photon CE data is discussed and supported by ab initio calculations of the predominant C–O bond breaking CE channel. We discuss the importance of these findings for achieving time resolved imaging of ultrafast dynamics.

Arnon, E. ; Rabani, E. ; Neuhauser, D. ; Baer, R. Equilibrium configurations of large nanostructures using the embedded saturated-fragments stochastic density functional theory. J. Chem. Phys. 2017, 146, 224111.Abstract

An ab initio Langevin dynamics approach is developed based on stochastic density functional theory (sDFT) within a new embedded fragment formalism. The forces on the nuclei generated by sDFT contain a random component natural to Langevin dynamics and its standard deviation is used to estimate the friction term on each atom by satisfying the fluctuation–dissipation relation. The overall approach scales linearly with system size even if the density matrix is not local and is thus applicable to ordered as well as disordered extended systems. We implement the approach for a series of silicon nanocrystals (NCs) of varying size with a diameter of up to 3nm corresponding to Ne = 3000 electrons and generate a set of configurations that are distributed canonically at a fixed temperature, ranging from cryogenic to room temperature. We also analyze the structure properties of the NCs and discuss the reconstruction of the surface geometry.

2010
Paul, A. K. ; Adhikari, S. ; Baer, M. ; Baer, R. H2 + photodissociation by an intense pulsed photonic Fock state. Phys. Rev. A 2010, 81, 013412.Abstract

We study the photodissociation of the H+2 molecule by ultrashort Fock-state electromagnetic pulses (EMPs). We use the Born-Oppenheimer treatment combined with an explicit photon number representation via diabatic electrophoton potential surfaces for simplification of the basic equations. We discuss the issue of the number of photon states required and show that six photon states enable good accuracy for photoproduct kinetic energies of up to 3 eV. We calculate photodissociation probabilities and nuclear kinetic-energy (KE) distributions of the photodissociation products for 800 nm, 50-TW/cm2 pulses. We show that KE distributions depend on three pulse durations of 10, 20, and 45 fs and on various initial vibrational states of the molecule. We compare the Fock-state results to those obtained by “conventional,” i.e., coherent-state, laser pulses of equivalent electric fields and durations. The effects of the quantum state of EMPs on the photodissociation dynamics are especially strong for high initial vibrational states of H+2. While coherent-state pulses suppress photodissociation for the high initial vibrational states of H+2, the Fock-state pulses enhance it.

2009
Paul, A. K. ; Adhikari, S. ; Mukhopadhyay, D. ; Halasz, G. J. ; Vibok, A. ; Baer, R. ; Baer, M. Photodissociation of H2+ upon Exposure to an Intense Pulsed Photonic Fock State. J. Phys. Chem. A 2009, 113, 7331–7337.Abstract

Producing and controlling nonclassical light states are now the subject of intense experimental efforts. In this paper we consider the interaction of such a light state with a small molecule. Specifically, we develop the theory and apply it numerically to calculate in detail how a short pulse of nonclassical light, such as the high intensity Fock state, induces photodissociation in H2+. We compare the kinetic energy distributions and photodissociation yields with the analogous results of quasi-classical light, namely a coherent state. We find that Fock-state light decreases the overall probability of dissociation for low vibrational states of H2+ as well as the location of peaks and line shapes in the kinetic energy distribution of the nuclei.

2008
Granot, R. S. ; Baer, R. A spline for your saddle. J. Chem. Phys. 2008, 128, 184111–9.Abstract

Pinpointing extrema on a multidimensional hypersurface is an important generic problem with a broad scope of application in statistical mechanics, biophysics, chemical reaction dynamics, and quantum chemistry. Local minima of the hypersurface correspond to metastable structures and are usually the most important points to look for. They are relatively easy to find using standard minimizing algorithms. A considerably more difficult task is the location of saddle points. The saddle points most sought for are those which form the lowest barriers between given minima and are usually required for determining rates of rare events. We formulate a path functional minimum principle for the saddle point. We then develop a cubic spline method for applying this principle and locating the saddle point(s) separating two local minima on a potential hypersurface. A quasi-Newton algorithm is used for minimization. The algorithm does not involve second derivatives of the hypersurface and the number of potential gradients evaluated is usually less than 10% of the number of potential evaluations. We demonstrate the performance of the method on several standard examples and on a concerted exchange mechanism for self-diffusion in diamond. Finally, we show that the method may be used for solving large constrained minimization problems which are relevant for self-consistent field iterations in large systems.

Granot, R. S. ; Baer, R. A tight-binding potential for helium in carbon systems. J. Chem. Phys. 2008, 129, 214102.Abstract

The presence of helium in carbon systems, such as diamonds and fullerenes is of interest for planetary sciences, geophysics, astrophysics, and evolution biology. Such systems typically involve a large number of atoms and require a fast method for assessing the interaction potential and forces. We developed a tight-binding approach, based on density functional calculations, which includes a many-body potential term. This latter term is essential for consolidating the density functional results of helium in bulky diamond and Helium passing through a benzene ring which is important for helium-fullerene applications. The method is simple to apply and exhibits good transferability properties.

2006
Livshits, E. ; Baer, R. Time-dependent density-functional studies of the D2 Coulomb explosion. J. Phys. Chem. A 2006, 110, 8443–8450. Publisher's VersionAbstract

Real-time first principle simulations are presented of the D2 Coulomb explosion dynamics detonated by exposure to very intense few-cycle laser pulse. Three approximate functionals within the time-dependent density functional theory (TDDFT) functionals are examined for describing the electron dynamics, including time-dependent Hartree-Fock theory. Nuclei are treated classically with quantum corrections. The calculated results are sensitive to the underlying electronic structure theory, showing too narrow kinetic energy distribution peaked at too high kinetic energy when compared with recent experimental results (Phys. Rev. Lett. 2003, 91, 093002). Experiment also shows a low energy peak which is not seen in the present calculation. We conclude that while Ehrenfest-adiabatic-TDDFT can qualitatively account for the dynamics, it requires further development, probably beyond the adiabatic approximation, to be quantitative.

2005
Shemesh, D. ; Baer, R. ; Seideman, T. ; Gerber, B. R. Photoionization dynamics of glycine adsorbed on a silicon cluster:“On-the-fly” simulations. The Journal of chemical physics 2005, 122, 184704.Abstract

Dynamics of glycine chemisorbed on the surface of a silicon cluster is studied for a process that involves single-photon ionization, followed by recombination with the electron after a selected time delay. The process is studied by “on-the-fly” molecular dynamics simulations, using the semiempirical parametric method number 3 (PM3) potential energy surface. The system is taken to be in the ground state prior to photoionization, and time delays from 5 to 50 fs before the recombination are considered. The time evolution is computed over 10 ps. The main findings are (1) the positive charge after ionization is initially mostly distributed on the silicon cluster. (2) After ionization the major structural changes are on the silicon cluster. These include Si–Si bond breaking and formation and hydrogen transfer between different silicon atoms. (3) The transient ionization event gives rise to dynamical behavior that depends sensitively on the ion state lifetime. Subsequent to 45 fs evolution in the charged state, the glycine molecule starts to rotate on the silicon cluster. Implications of the results to various processes that are induced by transient transition to a charged state are discussed. These include inelastic tunneling in molecular devices, photochemistry on conducting surfaces, and electron-molecule scattering.

2004
Kallush, S. ; Band, Y. B. ; Baer, R. Rotational Aspects of short-pulse population transfer in diatomic molecules. Chem. Phys. Lett. 2004, 392, 23.Abstract

A fully-selective population transfer scheme for diatomic molecules using short-duration (

Baer, R. ; Siam, N. Real-time study of the adiabatic energy loss in an atomic collision with a metal cluster. The Journal of chemical physics 2004, 121, 6341–6345.
2003
Neuhauser, D. ; Baer, R. ; Kosloff, R. Quantum soliton dynamics in vibrational chains: Comparison of fully correlated, mean field, and classical dynamics. J. Chem. Phys. 2003, 118, 5729–5735.Abstract

The dynamics of a chain of vibrational bonds which develop a classical solitary compression wave is simulated. A converged fully correlated quantum mechanical calculation is compared with a time dependent mean field approach (TDSCF) and with a classical simulation. The dynamics were all generated from the same Hamiltonian. The TDSCF and classical calculations show a fully developed solitary wave with the expected dependence of group velocity on amplitude. The full quantum calculations show a solitary-like wave which propagates for a while but then degrades. The robustness of the compression wave depends on the initial preparation. Evidence of partial recurrence of the wave has also been observed. (C) 2003 American Institute of Physics.

Anderson, S. M. ; Neuhauser, D. ; Baer, R. Trajectory-dependent cellularized frozen Gaussians, a new approach for semiclassical dynamics: Theory and application to He–naphtalene eigenvalues. The Journal of chemical physics 2003, 118, 9103–9108.Abstract

A semiclassical cellular method is proposed. Signals generated by semiclassical techniques generally deteriorate over time as trajectories become chaotic. One approach to remedy this problem has been to have each trajectory weighted by an entire cell of nearby trajectories (Filinov transform). But even in this approach the exponential part of the propagator typically becomes large and positive over time. Here the cellularization (Filinov) parameter is subject to constraints which make it time dependent and trajectory dependent. It also depends on dimensionality, so it ends up as a matrix. Physically, the Filinov transform is done differently in different directions associated with the stability matrix for the phase-essentially a more confined integration in directions where the matrix diverges and a wider integration in other directions. This squelches the contribution from any part of a trajectory that becomes excessively chaotic. A trajectory-dependent cellurized frozen Gaussian is applied here within the Herman-Kluk semiclassical approach. It is tested by looking at a single-particle three-dimensional problem, He attached to a rigid immovable naphtalene, where it is shown to be more accurate than the original HK approach, without the divergence of the correlation function common in the usual cellular dynamics (HK) formulation, and is able to separate a low-lying excited state from the ground state. (C) 2003 American Institute of Physics.

2002
Szichman, H. ; Baer, R. A five-dimensional quantum mechanical study of the H+ CH 4 –> H 2+ CH 3 reaction. The Journal of chemical physics 2002, 117, 7614–7623.