Variational Monte Carlo reactions

On the basis of Monte Carlo simulations [40] and molecular orbital calculations [26a], hydrogen bonding was proposed as the key factor controlling the variation of the acceleration for Diels-Alder reactions in water. Experimental differences of rate acceleration in water-promoted cycloadditions were recently observed [41]. Cycloadditions of cyclopentadiene with acridizinium bromide, acrylonitrile and methyl vinyl ketone were investigated in water and in ethanol for comparison (Scheme 3). Only a modest rate acceleration of 5.3 was found with acridizinium bromide, which was attributed to the absence of hydrogenbonding groups in the reactants. The acceleration factor reaches about 14 with acrylonitrile and 60 with methyl vinyl ketone, which is the best hydrogen-bond acceptor [41]. 

In recent years, there have been many attempts to combine the best of both worlds. Continuum solvent models (reaction field and variations thereof) are very popular now in quantum chemistry but they do not solve all problems, since the environment is treated in a static mean-field approximation. The Car-Parrinello method has found its way into chemistry and it is probably the most rigorous of the methods presently feasible. However, its computational cost allows only the study of systems of a few dozen atoms for periods of a few dozen picoseconds. Semiempirical cluster calculations on chromophores in solvent structures obtained from classical Monte Carlo calculations are discussed in the contribution of Coutinho and Canuto in this volume. In the present article, we describe our attempts with so-called hybrid or quantum-mechanical/molecular-mechanical (QM/MM) methods. These concentrate on the part of the system which is of primary interest (the reactants or the electronically excited solute, say) and treat it by semiempirical quantum chemistry. The rest of the system (solvent, surface, outer part of enzyme) is described by a classical force field. With this, we hope to incorporate the essential influence of the in itself uninteresting environment on the dynamics of the primary system. The approach lacks the rigour of the Car-Parrinello scheme but it allows us to surround a primary system of up to a few dozen atoms by an environment of several ten thousand atoms and run the whole system for several hundred thousand time steps which is equivalent to several hundred picoseconds. 

For the calculation using Variflex, the number of a variational transition q uantum s tates, N ej, w as given b y t he v ariationally d etermined minimum in Nej (R), as a function of the bond length along the reaction coordinate R, which was calculated by the method developed by Wardlaw-Marcus [6, 7] and Klippenstein [8]. The basis of their methods involves a separation of modes into conserved and transitional modes. With this separation, one can evaluate the number of states by Monte Carlo integration for the convolution of the sum of vibrational quantum states for the conserved modes with the classical phase space density of states for the transitional modes. 

This energy expression can be used to build up the respective variational functional to get the molecular orbitals [above]. A crucial step in the general self-consistent reaction field procedure is the estimation of the solvent charge density needed to obtain the response function G(r,r ) and the reaction potential. The use of Monte Carlo or molecular dynamics simulations of the system consisting the solute and surrounding solvent molecules has been proposed to find the respective solvent static and polarization densities. 

Both continuum model and Monte Carlo approach suggest the solvent to be in equilibrium with the reacting system at each point along the RP. The effects of non-equilibrium solvation were tested by variational transition state (VTST see Chap. 1) calculations of the microsolvated Sj 2 reaction (supermolecule approach)... 

---