My postdoctoral research consists of utilizing Monte Carlo algorithms to gain a better understanding of varying systems. More specifically, I explore the interfaces of water/oil mixtures (while varying the curvature of the interface), surfactant mesophases, a pressure equation of state for chain molecules, and hydrogen/water mixtures at extreme temperatures and pressures. Our research group’s highly specialized Monte Carlo algorithms and simulation protocols allow for the precise study and prediction of the physical properties of these systems.
We have developed and tested a chain-revised Groot-Warren equation of state (crGW-EOS) to describe systems of homo-oligomeric chains in the framework of dissipative particle dynamics (DPD). To this extent, thermodynamic perturbation theory (TPT) was applied to introduce correction terms that account for the reduction in pressure with an increasing number of bonds at constant bead number density. We have modified this equation of state by introducing a set of switching functions that yields an accurate second virial coefficient in the low-density limit. For chain molecules, the crGW-EOS offers several improvements over the revised Groot-Warren equation of state (rGW-EOS) and the Groot-Warren equation of state (GW-EOS). We have tested the crGW-EOS by using it to predict the pressure of oligomeric systems and the B2 virial coefficient of chain DPD particles for a range of bond lengths. Additionally, we have developed a method for determining the strength of cross-interaction parameters between chains of different composition and sizes and for thermal and athermal mixtures. We explore how system size affects partitioning near the UCST.
Supercritical fluid-fluid phase separations are important to understand for modeling Jovian planets, explaining the phase behavior of gas lasers, and understanding reduction-oxidation chemistry in Earth’s crust and mantle. To explore the supercritical fluid-fluid immiscibility of the binary H2O/H2 mixture, we use Monte Carlo simulations in the NpT Gibbs ensemble. The simulations use molecular mechanics force fields that treat both molecules as rigid and nonpolarizable but quantitatively reproduce the critical points of the neat compounds. We analyze the atom-atom radial distribution functions for H2O in the H2-rich phase, examine the cluster size distribution of these H2O aggregates and the distribution of H-bond energies and geometries.
Knowledge of the interfacial properties of water/oil mixtures is extremely important for the petrochemical industry and for understanding detergency and hydrophobic effects. In 2010, Javadi et al. observed a significant decrease in the surface tension of water when alkanes, including n-hexane, adsorbed from the vapor phase onto the surface. Additionally, Mucic et al. found that the surface tension of water decreases with increasing (partial) pressure of n-hexane. We probe the liquid/vapor interface of water/n-hexane and water/H2S mixtures using configurational-bias Monte Carlo simulations in the osmotic Gibbs ensemble. We study the effects of n-hexane and hydrogen sulfide, at several different partial pressures, on the liquid/vapor surface tension, and analyze the simulation trajectories to determine the n-hexane and hydrogen sulfide loading and to provide molecular-level insights on the structure of the interfaces.
Understanding the effects of surfactants on the structural, transport, and thermodynamic properties of vapor bubbles in aqueous solutions is important for naval research. We use Monte Carlo simulations in the NpT osmotic Gibbs ensemble to study surfactants adsorbing onto liquid/vapor interfaces with different curvature: (i) a vapor bubble in water (with bubble diameters of about 1.7, 2.1, and 2.4 nm), (ii) a free-standing film of water, and (iii) a bulk water phase. The surfactants considered are ethanol and 1-butanol. The water and surfactant molecules are modeled with the TIP4P/2005 and TraPPE-UA force fields, respectively.
Understanding solute uptake into soft microstructured materials, such as bilayers and worm-like and spherical micelles, is of interest in the pharmaceutical, agricultural, and personal care industries. To obtain molecular-level insight on the effects of solutes loading into a lamellar phase, we utilize the Shinoda–Devane–Klein (SDK) coarse-grained force field in conjunction with configurational-bias Monte Carlo simulations in the osmotic Gibbs ensemble. The lamellar phase is comprised of a bilayer formed by triethylene glycol mono-n-decyl ether (C10E3) surfactants surrounded by water with a 50:50 surfactant/water weight ratio. We study the unary adsorption isotherm and the effects on bilayer structure and stability caused by n-nonane, 1-hexanol, and ethyl butyrate at several different reduced reservoir pressures. Additionally, we study the binary adsorption isotherm of n-nonane and 1-hexanol mixtures.
A Langmuir monolayer consists of amphiphilic surfactant molecules that adsorb at a water/air interface to form a unimolecular layer. As well as having multiple applications in templating, chemical sensing, and coatings, these monolayers can serve as a model for lipid bilayers. The precise fluid phase behavior of Langmuir monolayers and other two-dimensional systems is difficult to accurately obtain through experiments. Molecular simulations are needed to probe the vapor-liquid equilibria of these complex systems. In 1994, Siepmann et al. used configurational-bias Monte Carlo simulations in the Gibbs ensemble (CBMC-GE) to study the vapor-liquid coexistence curve (VLCC) of a pentadecanoic acid monolayer, but these simulations were limited to a small system comprised of only 60 surfactant chains and relatively short trajectories due to the computational resources available. We carry out CBMC-GE simulations for much larger systems and using longer trajectories utilizing the TraPPE united-atom force field for the alkyl tails to probe the VLCC of a pentadecanoic monolayer. Furthermore, we investigate modification of the headgroup-headgroup potential to improve the accuracy of the predictions and to allow for extension of the simulations to coexistence of liquid-expanded and liquid-condensed phases.
Helicobactor pylori (HP) are gastrointestinal bacteria that infect the majority of the world’s population with no effective cure. These bacteria survive stomach acidity by producing HP urease. My Ph.D. projects entailed the provision of atomic-level understanding of the structure and function of the enzyme in order to pave the way for developing inhibitors for medical treatment.
The flaps that cover the active site cavity of the protein were found to open wider than previously observed, a finding which has significant implications for the design of HP urease inhibitors.
A possible shuttling process of urea into the active site of the protein was observed when the behavior of HP urease in 10.5M aqueous urea was studied. These findings allow for the targeting of other regions of the protein that may negatively affect the enzymatic hydrolysis of urea, leading to a greater number of potential inhibitor binding sites.
Interactions between aqueous ammonium ions and HP Urease were monitored. A possible path in which ammonium ions are expelled from the active sites was observed. It was found that blockage of this pathway results in another possible binding site for inhibition.
We obtained a structure of Klebsiella aerogenes urease (PDB ID: 1EJX) with the flap region added via homology modeling based on 2UBP (Bacillus pasteurii urease). Our group previously identified a hitherto unobserved wide-open flap state, and this was selected as the target for ligand docking studies. Crystallographic water molecules were stripped from this structure and it was subsequently loaded into Maestro. Next, we prepared a grid centered on the active site, and docked (a) urea, (b) a 36,000 ligand library, (c) several libraries from the ZINC database including biologically relevant compounds, lead-like structures and drug molecules and finally, (d) selected individual compounds from the literature that had been previously identified as promising inhibitors of urease. From the output we generated a list of compounds in their successfully docked poses, sorted by docking score. We are currently analyzing the results with particular focus on the interactions between the binding pocket and the docked ligands. We hope to identify promising candidates to be submitted for experimental assay in order to determine their efficacy as urease inhibitors.
While ab initio calculations and density functional theory have been well validated for a wide variety of organic compounds, there is comparably little supporting work that has been done in the area of transition metal-containing complexes. We have used fifteen density functionals and the MP2, CCSD, and CCSD(T) methods with Dunning double and triple-z basis sets, both with and without diffuse functions, in order to determine their utility in predicting the heats of formation and bond lengths of third row transition metal complexes. Preliminary results reveal that the TPSSTPSS and TPSSKCIS density functionals are quite effective in accurately predicting the heat of formation and equilibrium bond lengths for these metal complexes, and that the B3LYP, BLYP, and correlated methods perform much more poorly for the same compounds. These results suggest levels of computational theory that are suitable for applications such as geometry optimizations of metal-containing enzyme active sites, and for serving as the quantum mechanical segment of QM/MM calculations on larger protein structures.
The exact pathway of amide hydrolysis is currently unknown. We use computational methods in order to probe this reaction mechanism. The current experimental and computational data are in some instances inconsistent with one another. We are furthering our theoretical investigations in order to elucidate the hydrolytic pathway. We are using computational methods in order to further elucidate the various mechanistic pathways of urea hydrolysis by predicting the electronic structures of the ground states, intermediates, and transition states associated with a variety of possible mechanisms. We are using an array of density functionals to characterize these reaction pathways, with the aim of identifying a suitable QM level of theory for studying urease.
We conducted a systematic analysis of the electrostatic component of binding of three HIV-1 Reverse Transcriptase (RT) inhibitors—nevirapine (NVP), efavirenz (EFZ), and the recently-approved rilpivirine (RPV)—to wild-type (WT) and mutant variants of RT.
Electrostatic charge optimization was applied to determine how suited each molecule’s charge distribution is for binding WT and individual mutants of HIV-1 RT. We also determined the contributions of chemical moieties on each molecule toward the electrostatic component of binding, and show that different regions of a drug molecule may be used for recognition by different RT variants. The electrostatic contribution of certain RT residues toward drug binding was also computed to highlight critical residues for each interaction. Finally, the charge distribution of RPV is optimized to promiscuously bind to three RT variants rather than to each one in turn, with the resulting charge distribution being a compromise between the optimal charge distributions to each individual variant. This work demonstrates that even in a binding site considered quite hydrophobic, electrostatics play a subtle yet varying role that must be considered in designing next-generation molecules that recognize rapidly-mutating targets.
This study examines the promiscuity and specificity of the trypsin/bovine pancreatic trypsin inhibitor (BPTI) system through the analysis of the electrostatic component of ΔG. The charge contribution for each residue on the ligand is computationally eliminated and effects on ΔG are observed. Studies are also initiated for residues on the receptor (trypsin) near the binding site. A very positive change in ΔG indicates that that particular residue is important to binding. A parallel study is being conducted to evaluate the efficiency and accuracy of various numerical methods to determine ΔG electrostatic. The most accurate method, the finite difference method, is being compared with boundary element method, generalized Born, boundary-integral-based electrostatic estimation-Coulomb field approximation, and boundary-integral-based electrostatic estimation-preconditioned methods. The goal is to determine the method that optimizes accuracy and efficiency. Another undergraduate at Wellesley is currently completing this project.
My first research project utilized semi-empirical PM3 and density functional theory (DFT) B3LYP computations in order to determine the highest occupied molecular orbital (HOMO), lowest unoccupied molecular orbital (LUMO), and band gap energy patterns of increasing length single-walled carbon nanotubes (SWCNTs). Semi-empirical PM3 was first used to determine the band gap energy pattern. Changing the length in a carbon nanotube makes an oscillating logarithmic band gap energy pattern with a period of three points. Each point differed by one ring of carbon hexagons around the circumference. Both the HOMO and LUMO energies also create an oscillating pattern. In these patterns, at a length where HOMO energy is at a local maximum, LUMO and band gap energy is usually at a local minimum. To verify this result, the change in energy of the HOMO, LUMO, and band gap energies were graphed. It was found that in most cases these energies intercept the x-axis at the same point. Although the DFT and semi-empirical computations produced different values for the energies, a similar oscillating pattern was observed.