skip to main content


Title: Enhanced Sampling Path Integral Methods Using Neural Network Potential Energy Surfaces with Application to Diffusion in Hydrogen Hydrates
Abstract

The high cost of ab initio molecular dynamics (AIMD) simulations to model complex physical and chemical systems limits its ability to address many key questions. However, new machine learning‐based representations of complex potential energy surfaces have been introduced in recent years to circumvent computationally demanding AIMD simulations while retaining the same level of accuracy. As these machine learning methods gain in popularity over the next decade, it is important to address the appropriate way to develop and integrate them with well‐established simulation methods. This paper details the parameterization and training, using accurate electronic structure calculations, of artificial neural network potentials (NNPs) to model the intermolecular interactions in a hydrogen clathrate hydrate, wherein hydrogen molecules are confined inside the cavities of the 3D crystalline framework of hydrogen‐bonded water molecules. This new NNP is used in conjunction with new path‐integral based enhanced sampling methods, for inclusion of nuclear quantum effects and promotion of free energy barrier crossing, in order to determine properties that are important for understanding the diffusion of hydrogen gas in the clathrate system. These simulations demonstrate the influence of cage occupancy on the free energy barriers that determine the diffusivity of hydrogen gas through the network of large cages.

 
more » « less
NSF-PAR ID:
10361832
Author(s) / Creator(s):
 ;  ;  ;  
Publisher / Repository:
Wiley Blackwell (John Wiley & Sons)
Date Published:
Journal Name:
Advanced Theory and Simulations
Volume:
4
Issue:
4
ISSN:
2513-0390
Format(s):
Medium: X
Sponsoring Org:
National Science Foundation
More Like this
  1. Largely due to superior properties compared to traditional materials, the use of polymer matrix composites (PMC) has been expanding in several industries such as aerospace, transportation, defense, and marine. However, the anisotropy and nonhomogeneity of these structures contribute to the difficulty in evaluating structural integrity; damage sites can occur at multiple locations and length scales and are hard to track over time. This can lead to unpredictable and expensive failure of a safety-critical structure, thus creating a need for non-destructive evaluation (NDE) techniques which can detect and quantify small-scale damage sites and track their progression. Our research group has improved upon classical microwave techniques to address these needs; utilizing a custom device to move a sample within a resonant cavity and create a spatial map of relative permittivity. We capitalize on the inevitable presence of moisture within the polymer network to detect damage. The differing migration inclinations of absorbed water molecules in a pristine versus a damaged composite alters the respective concentrations of the two chemical states of moisture. The greater concentration of free water molecules residing in the damage sites exhibit highly different relative permittivity when compared to the higher ratio of polymer-bound water molecules in the undamaged areas. Currently, the technique has shown the ability to detect impact damage across a range of damage levels and gravimetric moisture contents but is not able to specifically quantify damage extent with regards to impact energy level. The applicability of machine learning (ML) to composite materials is substantial, with uses in areas like manufacturing and design, prediction of structural properties, and damage detection. Using traditional NDE techniques in conjunction with supervised or unsupervised ML has been shown to improve the accuracy, reliability, or efficiency of the existing methods. In this work, we explore the use of a combined unsupervised/supervised ML approach to determine a damage boundary and quantification of single-impact specimens. Dry composite specimens were damaged via drop tower to induce one central impact site of 0, 2, or 3 Joules. After moisture exposure, Entrepreneur Dr, Raleigh, North Carolina 27695, U.S.A. 553 each specimen underwent dielectric mapping, and spatial permittivity maps were created at a variety of gravimetric moisture contents. An unsupervised K-means clustering algorithm was applied to the dielectric data to segment the levels of damage and define a damage boundary. Subsequently, supervised learning was used to quantify damage using features including but not limited to thickness, moisture content, permittivity values of each cluster, and average distance between points in each cluster. A regression model was trained on several samples with impact energy as the predicted variable. Evaluation was then performed based on prediction accuracy for samples in which the impact energies are not known to the model.

     
    more » « less
  2. Proton transport in aqueous systems occurs by making and breaking covalent bonds, a process that classical force fields cannot reproduce. Various attempts have been made to remedy this deficiency, by valence bond theory or instantaneous proton transfers, but the ability of such methods to provide a realistic picture of this fundamental process has not been fully evaluated. Here we compare an ab initio molecular dynamics (AIMD) simulation of an excess proton in water to a simulation of a classical H3O+ in TIP3P water. The energy gap upon instantaneous proton transfer from H3O+ to an acceptor water molecule is much higher in the classical simulation than in the AIMD configurations evaluated with the same classical potential. The origins of this discrepancy are identified by comparing the solvent structures around the excess proton in the two systems. One major structural difference is in the tilt angle of the water molecules that accept an hydrogen bond from H3O+. The lack of lone pairs in TIP3P produces a tilt angle that is too large and generates an unfavorable geometry after instantaneous proton transfer. This problem can be alleviated by the use of TIP5P, which gives a tilt angle much closer to the AIMD result. Another important factor that raises the energy gap is the different optimal distance in water-water vs H3O+-water H-bonds. In AIMD the acceptor is gradually polarized and takes a hydronium-like configuration even before proton transfer actually happens. Ways to remedy some of these problems in classical simulations are discussed.

     
    more » « less
  3. Background

    Metamodels can address some of the limitations of complex simulation models by formulating a mathematical relationship between input parameters and simulation model outcomes. Our objective was to develop and compare the performance of a machine learning (ML)–based metamodel against a conventional metamodeling approach in replicating the findings of a complex simulation model.

    Methods

    We constructed 3 ML-based metamodels using random forest, support vector regression, and artificial neural networks and a linear regression-based metamodel from a previously validated microsimulation model of the natural history hepatitis C virus (HCV) consisting of 40 input parameters. Outcomes of interest included societal costs and quality-adjusted life-years (QALYs), the incremental cost-effectiveness (ICER) of HCV treatment versus no treatment, cost-effectiveness analysis curve (CEAC), and expected value of perfect information (EVPI). We evaluated metamodel performance using root mean squared error (RMSE) and Pearson’s R2on the normalized data.

    Results

    The R2values for the linear regression metamodel for QALYs without treatment, QALYs with treatment, societal cost without treatment, societal cost with treatment, and ICER were 0.92, 0.98, 0.85, 0.92, and 0.60, respectively. The corresponding R2values for our ML-based metamodels were 0.96, 0.97, 0.90, 0.95, and 0.49 for support vector regression; 0.99, 0.83, 0.99, 0.99, and 0.82 for artificial neural network; and 0.99, 0.99, 0.99, 0.99, and 0.98 for random forest. Similar trends were observed for RMSE. The CEAC and EVPI curves produced by the random forest metamodel matched the results of the simulation output more closely than the linear regression metamodel.

    Conclusions

    ML-based metamodels generally outperformed traditional linear regression metamodels at replicating results from complex simulation models, with random forest metamodels performing best.

    Highlights

    Decision-analytic models are frequently used by policy makers and other stakeholders to assess the impact of new medical technologies and interventions. However, complex models can impose limitations on conducting probabilistic sensitivity analysis and value-of-information analysis, and may not be suitable for developing online decision-support tools. Metamodels, which accurately formulate a mathematical relationship between input parameters and model outcomes, can replicate complex simulation models and address the above limitation. The machine learning–based random forest model can outperform linear regression in replicating the findings of a complex simulation model. Such a metamodel can be used for conducting cost-effectiveness and value-of-information analyses or developing online decision support tools.

     
    more » « less
  4. ABSTRACT

    In the interstellar medium (ISM), the formation of complex organic molecules (COMs) is largely facilitated by surface reactions. However, in cold dark clouds, thermal desorption of COMs is inefficient because of the lack of thermal energy to overcome binding energies to the grain surface. Non-thermal desorption methods are therefore important explanations for the gas-phase detection of many COMs that are primarily formed on grains. Here, we present a new non-thermal desorption process: cosmic ray sputtering of grain ice surfaces based on water, carbon dioxide, and a simple mixed ice. Our model applies estimated rates of sputtering to the three-phase rate equation model nautilus-1.1, where this inclusion results in enhanced gas-phase abundances for molecules produced by grain reactions such as methanol (CH3OH) and methyl formate (HCOOCH3). Notably, species with efficient gas-phase destruction pathways exhibit less of an increase in models with sputtering compared to other molecules. These model results suggest that sputtering is an efficient, non-specific method of non-thermal desorption that should be considered as an important factor in future chemical models.

     
    more » « less
  5. This data set for the manuscript entitled "Design of Peptides that Fold and Self-Assemble on Graphite" includes all files needed to run and analyze the simulations described in the this manuscript in the molecular dynamics software NAMD, as well as the output of the simulations. The files are organized into directories corresponding to the figures of the main text and supporting information. They include molecular model structure files (NAMD psf or Amber prmtop format), force field parameter files (in CHARMM format), initial atomic coordinates (pdb format), NAMD configuration files, Colvars configuration files, NAMD log files, and NAMD output including restart files (in binary NAMD format) and trajectories in dcd format (downsampled to 10 ns per frame). Analysis is controlled by shell scripts (Bash-compatible) that call VMD Tcl scripts or python scripts. These scripts and their output are also included.

    Version: 2.0

    Changes versus version 1.0 are the addition of the free energy of folding, adsorption, and pairing calculations (Sim_Figure-7) and shifting of the figure numbers to accommodate this addition.


    Conventions Used in These Files
    ===============================

    Structure Files
    ----------------
    - graph_*.psf or sol_*.psf (original NAMD (XPLOR?) format psf file including atom details (type, charge, mass), as well as definitions of bonds, angles, dihedrals, and impropers for each dipeptide.)

    - graph_*.pdb or sol_*.pdb (initial coordinates before equilibration)
    - repart_*.psf (same as the above psf files, but the masses of non-water hydrogen atoms have been repartitioned by VMD script repartitionMass.tcl)
    - freeTop_*.pdb (same as the above pdb files, but the carbons of the lower graphene layer have been placed at a single z value and marked for restraints in NAMD)
    - amber_*.prmtop (combined topology and parameter files for Amber force field simulations)
    - repart_amber_*.prmtop (same as the above prmtop files, but the masses of non-water hydrogen atoms have been repartitioned by ParmEd)

    Force Field Parameters
    ----------------------
    CHARMM format parameter files:
    - par_all36m_prot.prm (CHARMM36m FF for proteins)
    - par_all36_cgenff_no_nbfix.prm (CGenFF v4.4 for graphene) The NBFIX parameters are commented out since they are only needed for aromatic halogens and we use only the CG2R61 type for graphene.
    - toppar_water_ions_prot_cgenff.str (CHARMM water and ions with NBFIX parameters needed for protein and CGenFF included and others commented out)

    Template NAMD Configuration Files
    ---------------------------------
    These contain the most commonly used simulation parameters. They are called by the other NAMD configuration files (which are in the namd/ subdirectory):
    - template_min.namd (minimization)
    - template_eq.namd (NPT equilibration with lower graphene fixed)
    - template_abf.namd (for adaptive biasing force)

    Minimization
    -------------
    - namd/min_*.0.namd

    Equilibration
    -------------
    - namd/eq_*.0.namd

    Adaptive biasing force calculations
    -----------------------------------
    - namd/eabfZRest7_graph_chp1404.0.namd
    - namd/eabfZRest7_graph_chp1404.1.namd (continuation of eabfZRest7_graph_chp1404.0.namd)

    Log Files
    ---------
    For each NAMD configuration file given in the last two sections, there is a log file with the same prefix, which gives the text output of NAMD. For instance, the output of namd/eabfZRest7_graph_chp1404.0.namd is eabfZRest7_graph_chp1404.0.log.

    Simulation Output
    -----------------
    The simulation output files (which match the names of the NAMD configuration files) are in the output/ directory. Files with the extensions .coor, .vel, and .xsc are coordinates in NAMD binary format, velocities in NAMD binary format, and extended system information (including cell size) in text format. Files with the extension .dcd give the trajectory of the atomic coorinates over time (and also include system cell information). Due to storage limitations, large DCD files have been omitted or replaced with new DCD files having the prefix stride50_ including only every 50 frames. The time between frames in these files is 50 * 50000 steps/frame * 4 fs/step = 10 ns. The system cell trajectory is also included for the NPT runs are output/eq_*.xst.

    Scripts
    -------
    Files with the .sh extension can be found throughout. These usually provide the highest level control for submission of simulations and analysis. Look to these as a guide to what is happening. If there are scripts with step1_*.sh and step2_*.sh, they are intended to be run in order, with step1_*.sh first.


    CONTENTS
    ========

    The directory contents are as follows. The directories Sim_Figure-1 and Sim_Figure-8 include README.txt files that describe the files and naming conventions used throughout this data set.

    Sim_Figure-1: Simulations of N-acetylated C-amidated amino acids (Ac-X-NHMe) at the graphite–water interface.

    Sim_Figure-2: Simulations of different peptide designs (including acyclic, disulfide cyclized, and N-to-C cyclized) at the graphite–water interface.

    Sim_Figure-3: MM-GBSA calculations of different peptide sequences for a folded conformation and 5 misfolded/unfolded conformations.

    Sim_Figure-4: Simulation of four peptide molecules with the sequence cyc(GTGSGTG-GPGG-GCGTGTG-SGPG) at the graphite–water interface at 370 K.

    Sim_Figure-5: Simulation of four peptide molecules with the sequence cyc(GTGSGTG-GPGG-GCGTGTG-SGPG) at the graphite–water interface at 295 K.

    Sim_Figure-5_replica: Temperature replica exchange molecular dynamics simulations for the peptide cyc(GTGSGTG-GPGG-GCGTGTG-SGPG) with 20 replicas for temperatures from 295 to 454 K.

    Sim_Figure-6: Simulation of the peptide molecule cyc(GTGSGTG-GPGG-GCGTGTG-SGPG) in free solution (no graphite).

    Sim_Figure-7: Free energy calculations for folding, adsorption, and pairing for the peptide CHP1404 (sequence: cyc(GTGSGTG-GPGG-GCGTGTG-SGPG)). For folding, we calculate the PMF as function of RMSD by replica-exchange umbrella sampling (in the subdirectory Folding_CHP1404_Graphene/). We make the same calculation in solution, which required 3 seperate replica-exchange umbrella sampling calculations (in the subdirectory Folding_CHP1404_Solution/). Both PMF of RMSD calculations for the scrambled peptide are in Folding_scram1404/. For adsorption, calculation of the PMF for the orientational restraints and the calculation of the PMF along z (the distance between the graphene sheet and the center of mass of the peptide) are in Adsorption_CHP1404/ and Adsorption_scram1404/. The actual calculation of the free energy is done by a shell script ("doRestraintEnergyError.sh") in the 1_free_energy/ subsubdirectory. Processing of the PMFs must be done first in the 0_pmf/ subsubdirectory. Finally, files for free energy calculations of pair formation for CHP1404 are found in the Pair/ subdirectory.

    Sim_Figure-8: Simulation of four peptide molecules with the sequence cyc(GTGSGTG-GPGG-GCGTGTG-SGPG) where the peptides are far above the graphene–water interface in the initial configuration.

    Sim_Figure-9: Two replicates of a simulation of nine peptide molecules with the sequence cyc(GTGSGTG-GPGG-GCGTGTG-SGPG) at the graphite–water interface at 370 K.

    Sim_Figure-9_scrambled: Two replicates of a simulation of nine peptide molecules with the control sequence cyc(GGTPTTGGGGGGSGGPSGTGGC) at the graphite–water interface at 370 K.

    Sim_Figure-10: Adaptive biasing for calculation of the free energy of the folded peptide as a function of the angle between its long axis and the zigzag directions of the underlying graphene sheet.

     

    This material is based upon work supported by the US National Science Foundation under grant no. DMR-1945589. A majority of the computing for this project was performed on the Beocat Research Cluster at Kansas State University, which is funded in part by NSF grants CHE-1726332, CNS-1006860, EPS-1006860, and EPS-0919443. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1548562, through allocation BIO200030. 
    more » « less