Browsing by Author "Yang, Weitao"

Molecular dynamics (MD) simulations is an invaluable tool to investigate chemical and biological processes in atomic details. The accuracy of MD simulations strongly depends on underlying force fields. In conventional molecular mechanics (MM) force fields, the total energy is divided into bond energy, angle energy, dihedral energy, electrostatic interactions and van der Waals interactions. Each of these energy terms is parameterized by fitting to either experimental data or quantum mechanical (QM) calculations. In this dissertation, our aim is to develop accurate force fields for small molecules, water and proteins fully from QM calculations of small fragments. In the framework of conventional MM force fields, we calculated both transferable and molecule-specific atomic polarizabilities of small molecules by electrostatic potential fitting. Atomic polarizabilities are the key physical quantities in induced dipole polarization model. Molecular polarizabilities recovered from our atomic polarizabilities show good agreement with those obtained from QM calculations. We believe the main limitation of conventional MM force fields is the limited form of its Hamiltonian. Going beyond conventional MM force fields, we adopt the many-body expansion method and residue-based systematic molecular fragmentation (rSMF) method to start afresh building force fields for water and proteins, respectively. We used electrostatically embedded two-body expansion as the Hamiltonian of bulk water. QM reference of electrostatically embedded water monomer and dimer at the level of CCSD/aug-cc-pVDZ are parameterized by neural network (NN). Compared with experimental results, our water force fields show good structural and dynamical properties of bulk water. We developed rSMF to partition general proteins into twenty amino acid dipeptides and one peptide bond. The total energy of proteins is the combination of the energy of these small fragments. The QM reference energy of each fragment is parameterized by NN. Our protein force fields compare favorably with full QM calculations for both homogeneous and heterogeneous polypeptides in terms of energy and force errors.

Developing accurate and efficient theoretical approaches to describe the electronic structure has been a long-standing task in quantum chemistry. The main workhouse in quantum chemistry, density functional theory (DFT), has been widely used because of the good accuracy and the affordable computational cost. However, the applicability of commonly used density functional approximations (DFAs) is limited by intrinsic problems such as the delocalization error. Green's function theory that recently has gained increasing attention is shown to outperform the Kohn-Sham DFT approach on many aspects but is also computationally demanding. In this work, DFT and Green's function theory are connected to develop accurate and robust approaches for describing both ground state and excited state properties. For ground state calculations, the renormalized singles (RS) Green's function that captures all singles contributions from the KS Green's function is applied in the GW and the T-matrix approximation to predict accurate quasiparticle (QP) energies. GRSWRS and GRSTRS are shown to outperform over commonly used G0W0 and G0T0 for predicting ionization potentials (IPs) and core-level binding energies (CLBEs). The RS with correlation (RSc) Green's function that also includes higher order contributions in GW is shown to provide further improvements over GRSWRS. The concept of RS has also been used in the multireference DFT approach, which describes strongly correlated systems. We also provide an analytical approach to calculate QP energies of DFAs that can be expressed as a functional of the non-interacting Green's function. For excited state calculations, we combine localized orbital scaling correction (LOSC) with Bethe-Salpeter equation (BSE) to calculate excitation energies of molecular systems. QP energies from LOSC that systematically eliminates the delocalization error are used in BSE, which bypasses the expensive GW calculations. BSE/LOSC is shown to predict accurate excitation energies of valence, charge transfer and Rydberg excitations. We also combine the RS Green's function with BSE. BSE/GRSWRS is shown to provide a comparable accuracy to the computationally expensive BSE/evGW. We show that combining the merit of DFT and Green's function theory leads to accurate and efficient theoretical approaches for describing both the ground state and the excited state.

Quantum chemistry methods nowadays reach its maturity with various robust ground state correlation methods. However, many problems related to response do not have satisfactory solutions. Chemical reactivity indexes are some static response to external fields and number of particle change. These chemical reactivity indexes have important chemical significance, while not all of them had analytical expressions for direct evaluations. By solving coupled perturbed self-consistent field equations, analytical expressions were obtained and verified numerically. In the particle-particle (pp) channel, the response to the pairing field can describe N±2 excitations, i.e. double ionization potentials and double electron affinities. The linear response time-dependent density-functional theory (DFT) with pairing fields is the response theory in the density-functional theory (DFT) framework to describe $N\pm 2$ excitations. Both adiabatic and dynamic kernels can be included in this response theory. The correlation energy based on this response, the correlation energy of the particle-particle random phase approximation (pp-RPA), can also be proved equivalent to the ladder approximation of the well-established coupled-cluster doubles. These connections between the response theory, ab initio methods, and Green's function theory would be beneficial for further development. Based on RPA and pp-RPA, the theory of second RPA and the second pp-RPA with restrictions can be used to capture single and double excitations efficiently. We also present a novel methods, variational fractional spin DFT, to calculate singlet-triplet energy gaps for diradicals, which are usually calculated through spin-flip response theories.

Density functional theory (DFT) has become the main working horse for performing electronic structure calculations for chemical and physical systems nowadays. The theory is exact in principle, however, a density functional approximation (DFA) to the unknown exchange-correction energy $E_{\rm{xc}}$ in DFT has to be used in practice. Conventional DFAs have gained much success, while they usually possess intrinsic error and fail to describe some critical physical properties. In this dissertation, we focus on the delocalization error, a key concept to understand the systematic error existing in conventional DFAs, and we present the scaling correction methods which are designed to systematically and effectively reduce the delocalization error. The applications and improvements of two recently developed scaling correction methods, namely the the global scaling correction (GSC) method and the localized orbital scaling correction (LOSC) method, are mainly discussed in this dissertation. First, we demonstrate that the scaling corrections is capable of accurately predicting the quasiparticle energies and photoemission spectra from the orbital energies of the (generalized) Kohn-Sham DFT calculations. Second, we present a new method called QE-DFT, which is developed based on the connection between the orbital energies and the quasiparticle energies, to describe the difficult excited-state problems, including the low-lying, Rydberg and charge transfer excitations and conical intersections. We further derive the analytical gradients of the QE-DFT method, and demonstrate the application of QE-DFT for describing the potential energy surface and geometry optimization of excited states. Third, we show the application of LOSC to describe the polymer polarizability, which is a challenging problem for conventional DFAs. Fourth, we present the recent development for both GSC and LOSC methods. Specifically, we developed the analytic and exact second-order correction under the framework of GSC and achieved much improved accuracy compared to the original work of GSC. We also developed a new and robust self-consistent approach for LOSC method to avoid the convergence difficulties in the original LOSC work, which comes from using the approximate LOSC effective Hamiltonian. Finally, we developed the implementation of the scaling correction methods as a library with the supports to multiple programming languages. In summary, we demonstrated with extensive results that the GSC and LOSC are powerful and effective scaling correction methods to conventional DFAs to largely reduce the delocalization error. With the further developments to GSC and LOSC, they should be of great potential for broad application to describing challenging electronic structure problems of complex systems with high accuracy in the future.

The accurate theoretical description of electronic excitations is of crucial importance for understanding many important processes in various areas of science and technology, including specific problems such as photocatalysis, spectroscopy as well as various biochemical reactions. While the current state of the art for an accurate quantum chemistry description of the electronic ground state is relatively mature, the accurate and efficient treatment of electronic excited states is still teeming with both challenges and opportunities. In this treatise, the following different yet related methodologies have been employed to tackle this problem. Firstly, the orbitals and orbital energies obtained with ground state Kohn-Sham density functional theory calculations are adopted to describe the electron addition excitation, offering valuable qualitative insight into the photocatalytic mechanism and product selectivity of CO2 reduction with rhodium nanoparticles. Secondly, the development of the T-matrix method based on the time-dependent pairing density linear response and the particle-particle random phase approximation (pp-RPA) takes the leap from qualitative correctness to quantitative reliability for describing electron addition and removal excitations. Thirdly, the pp-RPA has been used to study neutral excitations, particularly the excited state potential energy surfaces and geometry optimization. Also, an active orbital selection approach has been developed which enables ~100-fold computational time savings and the application of the pp-RPA to large systems including polyacetylenes and polydiacetylenes. Lastly, a generalized theoretical framework from the superconductive Gorkov Green’s function perspective has been proposed which unifies the pp-RPA with the conventional particle-hole (ph) RPA, allowing for a systematic improvement of electron correlation treatment by incorporating higher-order self energy contributions in the respective ph and pp Bethe-Salpeter equations (BSEs). Initial numerical results in the ph channel demonstrate a significant accuracy improvement compared with conventional methods such as time-dependent density functional theory (TDDFT) within the adiabatic approximation.

The accurate description of ground and electronic excited states is an important and challenging topic in quantum chemistry. The pairing matrix fluctuation, as a counterpart of the density fluctuation, is applied to this topic. From the pairing matrix fluctuation, the exact electron correlation energy as well as two electron addition/removal energies can be extracted. Therefore, both ground state and excited states energies can be obtained and they are in principle exact with a complete knowledge of the pairing matrix fluctuation. In practice, considering the exact pairing matrix fluctuation is unknown, we adopt its simple approximation --- the particle-particle random phase approximation (pp-RPA) --- for ground and excited states calculations. The algorithms for accelerating the pp-RPA calculation, including spin separation, spin adaptation, as well as an iterative Davidson method, are developed. For ground states correlation descriptions, the results obtained from pp-RPA are usually comparable to and can be more accurate than those from traditional particle-hole random phase approximation (ph-RPA). For excited states, the pp-RPA is able to describe double, Rydberg, and charge transfer excitations, which are challenging for conventional time-dependent density functional theory (TDDFT). Although the pp-RPA intrinsically cannot describe those excitations excited from the orbitals below the highest occupied molecular orbital (HOMO), its performances on those single excitations that can be captured are comparable to TDDFT. The pp-RPA for excitation calculation is further applied to challenging diradical problems and is used to unveil the nature of the ground and electronic excited states of higher acenes. The pp-RPA and the corresponding Tamm-Dancoff approximation (pp-TDA) are also applied to conical intersections, an important concept in nonadiabatic dynamics. Their good description of the double-cone feature of conical intersections is in sharp contrast to the failure of TDDFT. All in all, the pairing matrix fluctuation opens up new channel of thinking for quantum chemistry, and the pp-RPA is a promising method in describing ground and electronic excited states.

Molecular dynamics (MD) simulation is an extremely powerful, highly effective, and widely used approach to understand the nature of chemical processes in atomic details for small molecules, biomolecules, and materials. The accuracy of MD simulation results is highly dependent on force fields. Quantum mechanical (QM) calculation has excellent accuracy, but the computational cost is not affordable for long MD simulations. Therefore, traditional molecular mechanical (MM) force fields, which divide energy into classical bond, angle, dihedral, electrostatic and van der Waals (vdW) terms, or hybrid QM/MM methods, which consider tradeoff between QM accuracy and MM efficiency, are generally utilized in MD simulations. Machine learning (ML) provides capability for generating an accurate potential at QM level without increasing much computational effort, and ML-based potentials had rapid development and widespread applications during the past decade. In this dissertation, we apply MD and ML techniques to develop new methods for simulating on small molecules and proteins. First, we train ML models to increase the accuracy of QM/MM from the semiempirical to the ab initio level. Active learning is performed to efficiently update ML models on the fly with gradient boosting technique, and new data from MD simulations are sampled according to the boundary of reference energies, distance-based clustering, and density-based clustering. Solvation free energies of small molecules obtained from QM/MM ML models show good agreement with experiment. Next, force fields based on neural network (NN) are constructed for QM/MM vdW interaction, which is normally described with Lennard-Jones (LJ) potential. We develop a new QM/MM NN architecture, dubbed QM-NN/MM-NN, and new input features based on center of mass for NN, which better describes non-bonded interactions than other descriptors. NN force fields greatly outperform LJ potentials and show good transferability to different small molecules. In addition, general and transferable NN force fields based on CHARMM force fields, named CHARMM-NN, are constructed for proteins, according to residue-based systematic molecular fragmentation method. NN is based on atom types and new input features that are similar to MM inputs are proposed, which enhances the compatibility of CHARMM-NN with MM MD. The validations on geometric data, relative potential energies and reorganization energies demonstrate that the potential energy minima of CHARMM-NN are very similar to QM, but the simulations of peptides and proteins indicate that the solvent effects and non-bonded interactions should be modeled in future development of NN force fields. Finally, we develop a piecewise approach to run all-atom steered MD (SMD) simulations within small water box, avoiding the huge amounts of computational resources required to run all-atom SMD simulations using a large water box. The robustness of this approach is validated with a small protein NI3C. Compared to coarse-grained SMD, the all-atom SMD simulations on luciferase reveal more atomic resolution details on force-extension plots and the key secondary structures related to mechanical stability in unfolding pathway.

The thesis contains two directions in the simulations of biomolecular systems. The first part (Chapter 2 - Chapter 4) mainly focuses on the simulations of electron transfer processes in condensed phase; the second part (Chapter 5 - Chapter 6) investigates the conformational sampling of polysaccharides and proteins. Electron transfer (ET) reaction is one of the most fundamental processes in chemistry and biology. Because of the quantum nature of the processes and the complicated roles of the solvent, calculating the accurate kinetic and dynamic properties of ET reactions is challenging but extremely useful. Based on the Marcus theory for thermal ET in weak coupling limit, we combined the rigorous ab initio quantum mechanical (QM) method and well-established molecular mechanical (MM) force field and developed an approach to directly calculate a key factor that affects the ET kinetics: the redox free energy. A novel reaction order parameter fractional number of electrons (FNE) was used to characterize the ET progress and to drive the QM/MMMD sampling of the nonadiabatic free energy surface. This method was used for two aqueous metal cations, iron and ruthenium in solution, and generated satisfactory results compared to experiments. In order to further reduce the computational cost, a QM/MM-minimum free energy path (MFEP) method is implemented and combined with the FNE in the calculation of redox free energies. The calculation results using QM/MM-MFEP+FNE generated identical results as the direct QM/MM-MD method for the two metal cations, demonstrating the consistency of the two different sampling strategy. Furthermore, this new method was applied to the calculation of organic molecules and enhanced the computational efficiency 15-30 times than the direct QM/MM-MD method, while maintaining high accuracy. Finally, I successfully extended the QM/MM-MFEP+FNE method to a series of redox proteins, azurin and its mutants, and obtained very accurate redox free energy differences with relative error less than 0.1 eV. The new method demonstrated its excellent transferability, reliability and accuracy among various conditions from aqueous solutions to complex protein systems. Therefore, it shows great promises for applications of the studies on redox reactions in biochemistry. In the studies of force-induced conformational transitions of biomolecules, the large time-scale difference from experiments presents the challenge of obtaining convergent sampling for molecular dynamics simulations. To circumvent this fundamental problem, an approach combining the replica-exchange method and umbrella sampling (REM-US) is developed to simulate mechanical stretching of biomolecules under equilibrium conditions. Equilibrium properties of conformational transitions can be obtained directly from simulations without further assumptions. To test the performance, we carried out REM-US simulations of atomic force microscope (AFM) stretching and relaxing measurements on the polysaccharide pustulan, a (1→6)-β-D-glucan, which undergoes well-characterized rotameric transitions in the backbone bonds. With significantly enhanced sampling convergence and efficiency, the REMUS approach closely reproduced the equilibrium force-extension curves measured in AFM experiments. Consistent with the reversibility in the AFM measurements, the new approach generated identical force-extension curves in both stretching and relaxing simulations, an outcome not reported in previous studies, proving that equilibrium conditions were achieved in the simulations. In addition, simulations of nine different polysaccharides were performed and the conformational transitions were reexamined using the REM-US approach. The new approach demonstrated consistent and reliable performance among various systems. With fully converged samplings and minimized statistical errors, both the agreement and the deviations between the simulation results and the AFM data were clearly presented. REM-US may provide a robust approach to modeling of mechanical stretching on polysaccharides and even nucleic acids. However, the performance of the REM-US in protein systems, especially with explicit solvent model, is limited by the large system size and the complex interactions. Therefore, a Go-like model is employed to simulate the protein folding/unfolding processes controlled by AFM. The simulations exquisitely reproduced the experimental unfolding and refolding force extension relationships and led to the full reconstruction of the vectorial folding pathway of a large polypeptide, the 253-residue consensus ankyrin repeat protein, NI6C. The trajectories obtained in the simulation captured the critical conformational transitions and the rate-limiting nucleation event. Together with the AFM experiments, the coarse-grained simulations revealed the protein folding and unfolding pathways under the mechanical tension.

The Kohn-Sham density functional theory has been the most popular method in electronic structure calculations. Although the theory is in principle exact, approximations are needed for the exchange-correlation energy to make practical calculations possible. This dissertation focuses on one aspect to seek for better approximations: optimization. Optimization has many different meanings in quantum chemistry. Specifically, we carry out energy optimization for a few known energy forms to get insight about what self-consistency can achieve. One kind of them is designed from linear response theories (two-electron addition energies from the particle-particle random phase approximation and spin-flip excitation energies from the spin-flip linear-response time-dependent density functional theory). The post-self-consistency versions are already known to produce excellent results for states that involve static correlations, while the energy optimization further improves the accuracy, especially when the total energy is a quantity of interest. The optimization technique used is the generalized optimized effective potential method, which is a non-local counterpart to the well-established optimized effective potential method. The other part is related to the local orbital scaling correction method, which mostly deals with delocalization errors. We designed a robust self-consistent workflow to avoid the difficult exact-Hamiltonian derivation, and observed reliable prediction on electron densities, total energies, and energy-level alignments. Machine-learning can also be considered as an optimization technique, which is used to minimize the errors of a density functional approximation in the provided database. Special treatment was carried out to transfer electron densities in 3D to neural network compatible input variables (translational and rotational invariant; discrete numbers). The outcome is a finite-range non-local density functional correction to low cost functionals, which is able to achieve similar level of accuracy as compared to accurate and expensive functionals, while keeping the computational cost low. In summary, "optimization" over existing functional approximations is a reliable way to get an improvement, and sometimes new physics can be discovered. On the other hand, "optimization" can also be done by "machines", which may become more and more popular in the future.

Hartree-Fock (HF) theory is one of the fundamental theories in electronic structure calculations. It can provide the basic description of the system, however, as the zeroth order approximation, it cannot describe the electronic structure accurately. Nevertheless, it is the starting reference for many beyond-HF wavefunction methods. Another fundamental theory in electronic structure is the density functional theory (DFT). It is usually ignored that DFT has an important connection with the HF theory. This connection is the single excitation (SE) contribution. In HF, SE does not contribute because of the Brillouin's theorem. On the other hand, SE contributes to DFT and plays a significant role in beyond-DFT many-body methods. One example is the correlation energy calculated by random phase approximation (RPA). Here we would like to discuss SE in different scenarios and extend the topic. First, the optimized effective potential (OEP) does not perform accurately in minimizing the RPA total energy functional. To overcome this problem, we developed a generalized optimized effective potential (GOEP) method. This method can describe the dissociation of weakly interacting diatomic systems accurately. From the analysis of the energy structure, the GOEP absorbs the SE contribution in contrast with the OEP method. We also notice that by performing GOEP for RPA, the physical density of the system is no longer the reference density, which is not traditionally recognized in DFT. This conclusion can be generally applied to any exchange-correlation functional that depends explicitly on the external potential. And we have shown that this physical density performs better than both the reference density and the original density of the starting reference calculated from HF or DFT. Second, the $GW$ approximation is widely used in different research areas. Especially, the $G_0W_0$ method can improve the ionization potential for both molecule and solids calculated from the ground state DFT calculation. However, the $G_0W_0$ method has a strong starting-point dependence. We have recognized that this starting-point dependence largely originates from the lack of SE contribution to the single-particle Green's function. We developed an effective and simple solution by using a subspace diagonalization of the HF Hamiltonian with the DFT density matrix to construct the renormalized singles Green's function and replace the reference Green's function $G_0$. Our method works extremely well for molecules and we are still testing it for solid states. Besides, we have developed a new method to extract excitation energies directly from the quasi-particle energies based on the $GW$ approximation. Starting from the $(N-1)$-electron system, we are able to calculate molecular excitation energies with orbital energies at the $GW$ level. We have demonstrated that this method can accurately capture low-lying local excitations as well as charge transfer excitations in many molecular systems. This provides a new perspective in applying the single-particle Green's function.

Electron transfer (ET) and excitation energy transfer (EET) are widely observed in various research areas, for example, electronic devices, biomolecular systems, light harvesting systems, etc. Here we discuss two topics. First, we study circularly polarized light (CPL) induced ET process. CPL induced coherent electron transfer (ET) process can be affected by the direction of the CPL. In particular, the yield on the acceptor through the CPL-induced ET can be asymmetrical. Previous study suggests that this yield asymmetry on the acceptor is related with the initial angular momentum polarization on the donor, which can be created by different directions of the CPL. Here we further investigate how the CPL affects the yield asymmetry by studying the yield asymmetry dependence on the molecular energetics, CPL field, and the environment perturbation. We have built a simple 4-state Hamiltonian with one ground state, two degenerate excited donor states and one acceptor and provided the optimal choice of parameters to maximize the yield asymmetry. Both analytical and numerical results suggest that the yield asymmetry is mainly created by the phase and population difference between two excited states under L- and R-CPL. Among different parameters, a slow dephasing rate is most important for observing a large yield asymmetry. With a 200 fs dephasing rate, the yield asymmetry can be as large as 5\%. One should perform the experiment in low temperature to slow the dephasing rate. Second, we study the hole length along the DNA base pairs. Knowing the length of the hole is important to understand the mechanism of the ET process through DNA. The length of the hole is determined by the competition between the delocalization (ionization) and localization (solvation) effect of the DNA. We have revisited the previous work in our group (10.1021/jp0132329) with more advanced quantum chemistry methods to calculate the hole length in AT and GC pairs. Localized scaling orbital correction (LOSC) is used to calculate the ionization potential (IP), which is usually underestimated by traditional DFT calculations because of the size-consistent problem. And we will use the Poisson-Boltzmann equation to calculate the solvation energy contribution. From the LOSC calculation, the IP is in agreement with the high-level ab-initio calculations. This result shows that the LOSC calculation overcomes the problem in traditional DFT. And we have also shown that the IP decreases faster for the increasing length of the GC pair as compared to the AT pair.

This thesis contains four projects in Quantum Mechanics/Molecular Mechanics (QM/MM) applications and methodology developments for biological systems. The first part (chapter 2 - chapter 3) is mainly focused on enzymatic reaction mechanism studies; the second part (chapter 4 - chapter 5) is to develop new methods to effectively calculate solvation free energy; the third part (chapter 6) is about methodology development to visualize non-covalent interaction in fluctuating biological environment; the fourth part (chapter 7) is to further apply QM/MM method on understanding proton coupled electron transfer reaction mechanism.

Catalytic mechanism is the key to understand how enzyme facilitates chemical reactions. A comprehensive understanding of the mechanism would shed light on developing new drug candidates or find similar protein targets. In the first two sections, we explored two enzyme systems, 4-oxalocrotonate tautomerase (4-OT) and anhydro-N-acetylmuramic acid kinase (AnmK). 4-OT is an essential enzyme in the degradative metabolism pathway occurring in the Krebs cycle. The proton transfer process catalyzed by 4-OT was studied to elaborate its catalytic mechanism. With QM/MM simulation, we demonstrated that the enzyme works under half-of-the-sites occupation, i.e. only three of its six active sites are occupied by the substrates. Two sequential proton transfers occur: one proton from the C3 position of 2o4hex is initially transferred to the nitrogen atom of the general base, Pro1. Subsequently, the same proton is shuttled back to the position C5 of 2o4hex to complete the proton transfer process in 4-OT. During the catalytic reaction, conformational changes (i.e., 1-carboxyl group rotation) of 2o4hex cannot proceed in the natural hexametric structure. We further found that the docking process of 2o4hex can influence the specific reactant conformations and an alternative substrate (2-hydroxymuconate) may serve as reactant under a different reaction mechanism than 2o4hex. The other enzyme, AnmK, plays an important role in the cell wall recycling process in Escherichia coli. The catalytic process involves four entangled steps: water nucleophile attack, 1, 6-anhydro bond breaking, phosphorylation, and saccharide conformational change. This is a great challenge to present enzymatic simulation methods. We proposed a comprehensive scheme to tackle this problem, and our results indicate a proton shuttle network may exist, which position a water molecule well for the nucleophile attack.

Solvation free energy strongly relates with the molecule solubility and is an important indicator to screen potential drug candidates. To improve the computational speed and simulation accuracy on evaluating solvation free energy, we developed a new approach to combine λ dynamics with metadynamics to compute free energy surface with respect to λ, we named it λ-metadynamics. Particularly, the λ-metadynamics method extends metadynamics to a single virtual variable λ, i.e., the coupling parameter between solute and solvent, to compute absolute solvation free energy as an exemplary application. We demonstrated that λ-metadynamics simulations can recover the potential of mean force surface with respect to λ compared to the benchmark results from traditional λ-dynamics with umbrella sampling. The solvation free energy results for five small organic molecules from λ-metadynamics simulations using the same filling scheme show that the statistical errors are within ±0.5 kcal/mol. Further, a detailed exploration over λ-metadynamics found out the intrinsic problems associated with metadynamics on high second derivative free energy surfaces. We introduced new techniques, such as transformation of virtual variable, region segmentation, mirror filling procedure and independent trajectory average, integrating with the original λ-metadynamics to improve its power. The revised λ-metadynamics is applied to investigate the solvation free energy of 20 molecular sets, and its result suggests a minimum QM basis set (3-21G) with point charge model perform better than more complex ones.

Non-covalent interaction is prevalent, such as hydrogen bond, van der Waals, steric clashes, and it plays a central role in many chemical and biological systems. Although the atomic distance may indicate strength of a specific interaction, a quantitative description is very difficult. We extended a recently introduced non-covalent interaction index (NCI) into its ensemble averaged counterpart, averaged non-covalent interaction index (aNCI). NCI is suitable for quantitatively characterize non-covalent interactions for any static system, while aNCI is capable to determine and visualize non-covalent interactions for fluctuating biological systems. We further applied aNCI on various systems including solute-solvent and ligand-protein non-covalent interactions. For water and benzene molecules in aqueous solution, solvation structures and the specific hydrogen bond patterns were visualized clearly. For the Cl-+CH3Cl SN2 reaction in aqueous solution, charge reorganization influences over solvation structure along SN2 reaction were revealed. For ligand-protein systems, aNCI can recover several key fluctuating hydrogen bond patterns that have potential applications for drug design. Therefore, aNCI, as a complementary approach to the original NCI method, can extract and visualize non-covalent interactions from thermal noise in fluctuating environments.

Proton-coupled electron transfer is a reaction mechanism involves simultaneous transfer of both proton(s) and electron(s). It is pervasive in redox reactions and has significant importance in photosynthesis, respiration, and electrochemical process such as hydrogen-ion discharge. We applied the fractional number of electron approach and QM/MM scheme to study the basic principle about how electron transfer (ET) and proton transfer (PT) coupled together. This approach was employed in the oxidation of an intra-molecular hydrogen bonded phenol (PhOH...N<), with two well-defined reactions coordinates (PT: proton position; ET: energy gap). The computed absolute redox potential agrees well with the experimental data. We found that at a low applied overpotential (< 5.0 eV), proton transfer occurs prior to electron, and this gives a stepwise proton-electron transfer pathway (PT-ET). With a high overpotential (> 7.0 eV), the reaction occurs in a reverse sequence with proton motion triggered by electron transfer (ET-PT). When the overpotential is between 5.0 and 7.0 eV, the reaction proceeds by concerted proton-electron transfer (CPET) mechanism. Our work suggests that the mechanistic details of PCET process can be changed by altering electrode potentials.