Name: multiscale sampling of a heterogeneous water/metal catalyst interface using density functional theory and force-field molecular dynamics
Uploaded: 2019-04-12T14:45:00.000+00:00
Duration: 10 min 52 s
Description: Watch this Scientific Journal Video about Multiscale Sampling of a Heterogeneous Water/Metal Catalyst Interface using Density Functional Theory and Force-Field Molecular Dynamics at JoVE.com

This research was funded by the National Science Foundation through award number CBET-1438325. Fellowship support for CJB through NASA Training Grant NX14AN43H is gratefully acknowledged. Simulations were performed on the Palmetto Supercomputer Cluster, which is maintained by the Cyberinfrastructure Technology Group at Clemson University. We thank Dr. Paul J. Meza-Morales for testing the protocol.

1. Generate the adsorbate structure
<ol>
	<li>Create a VASP POSCAR file comprising a supercell with periodic boundary conditions, as you would if you were performing simulations of adsorbates on metal surfaces under vacuum. The supercell should include an initial guess of the adsorbate structure and the metal surface as well as vacuum space above the adsorbate for adding H2O molecules. Details on how to do this are provided in previous work35,44. 
	NOTE: It is important that the height of the vacuum space be at least 12 &#197; above the top of the adsorbate.</li>
	<li>Relax the structure and minimize its energy using the VASP code. Details on how to do this are provided in previous work35,44. This will produce a file called CONTCAR, which will be used in the next section.</li>
</ol>
2. Add explicit H2O molecules
<ol>
	<li>Add N water molecules to the vacuum space in the CONTCAR created in step 1.2 using the MCPliQ code, where N = &#961;V, &#961; is the density of water, and V is the volume of the vacuum space above the adsorbate. 
	NOTE: &#961; should be taken as the density of water as determined by the TIP3P/CHARMM water model at the simulation temperature. V will be refined in the next step.
	<ol>
		<li>Specify the following information in the MCPliQ master_input.txt file: The number of H2O molecules to add (N) by changing the first argument in line 28, the path to the water.txt file by changing the second argument in line 28, and the minimum and maximum height of the supercell that can be occupied by water molecules by changing the Minimum z-coordinate in line 11 and the Maximum z-coordinate in line 12.</li>
		<li>Execute the MCPliQ code by typing mcpliq from the command-line interface to insert water molecules into the CONTCAR file. The code will output one or more files with the file extension of .POSCAR. 
		NOTE: If more than one POSCAR files are produced, they will be named POSCAR_n.POSCAR. Select the file where n is largest.</li>
	</ol>
	</li>
	<li>Generate the LAMMPS input files for an NPT simulation and equilibrate the cell volume using FFMD in the NPT ensemble in LAMMPS.
	<ol>
		<li>Execute the script lmps_bond_angle.py on the .POSCAR file generated in step 2.1.2 by typing lmps_bond_angle.py $filename.POSCAR at the command-line interface, where $filename is the name of the .POSCAR file generated in step 2.1.2. This script creates a file called $filename.POSCAR.bond_angle_info.txt which lists bonds and angles that will be used in the LAMMPS data file.</li>
		<li>Open VMD and select File &#62; New Molecule in the main window to open the Molecule File Browser window. Select VASP_POSCAR from the Determine File Type dropdown menu. Click Browse and navigate to the $filename.POSCAR file. Click Load to open the $filename.POSCAR file.</li>
		<li>Open the Tk Console within VMD by selecting Extensions &#62; Tk Console from the VMD Main window.</li>
		<li>Execute the following command in the Tk console: topo writelammpsdata $WDPATH/data.myadsorbate full, where $WDPATH is the directory on the computer where VMD will write the LAMMPS data file and data.myadsorbate is the name of the LAMMPS data file.</li>
		<li>Delete the Bonds and Angles section at the bottom of the data.myadsorbate file. Then, append the bond and angle lists in the file $filename.POSCAR.bond_angle_info.text into data.myadsorbate. 
		NOTE: The indexes for the O-H bond type and H-O-H angle type for the water molecules in the $filename.POSCAR.bond_angle_info.txt file are both set to 1. Thus, bond and angle types for adsorbates should start counting at 2.</li>
		<li>Edit the data.myadsorbate file by adding the Lennard-Jones parameters to the Pair Coeffs section and the Coulomb parameters to the Atoms section. Lennard-Jones and Coulomb parameters for the H2O molecules, adsorbate atoms, and metal surface atoms need to be added. 
		NOTE: The Lennard-Jones parameters for Pt atoms, water molecules, and adsorbate atoms in this protocol are obtained from the UFF38, TIP3P/CHARMM37, and OPLS-AA39 force fields, respectively. Coulomb parameters for water molecules and adsorbates atoms are obtained from the TIP3P/CHARMM37 and OPLS-AA39 force fields, respectively. Coulomb parameters for Pt atoms are set to 0 in this protocol. Alternatively, calculated partial charges could be used for the Coulomb parameters for adsorbate atoms and Pt atoms.</li>
		<li>Copy the LAMMPS input file input.equil into the directory $WDPATH. Edit the group variable on line 34 to indicate the atom type indexes for the water oxygen and water hydrogen atoms and the group variable on line 35 to indicate the atom type indexes for the Pt and adsorbate atoms.</li>
		<li>Execute the LAMMPS software by typing mpiexec -np XX lmp_mpi &#60; input.equil at the command-line interface, where XX is the number of CPU cores to use, and lmp_mpi is the name of the LAMMPS executable. Doing this will run an energy minimization to refine the H2O configuration, followed by a FFMD simulation performed at constant number of H2O molecules (N), volume (V), and temperature (T) to bring the water to the simulation temperature, followed by a FFMD simulation run at constant N, pressure (P), and temperature (T) to determine the physically correct height of the simulation box. Output files that will be used in section 3 are called data.myadsorbate_npt and log.myadsorbate. 
		NOTE: The duration of the NPT simulation should be long enough to comprise an &#8220;equilibration&#8221; run, where the volume of the supercell comes to steady state, and a &#8220;production&#8221; run, which is used to sample the ensemble averages (here, the height of the supercell). During the equilibration run, the volume of the supercell when plotted against time should level off to a steady state value. Once this occurs, the NPT simulation can be said to be in its production run. Verify equilibration of the NPT simulation by ensuring that the fluctuations in the height of the supercell (lz) are minimal or have converged to a steady value. If large fluctuations occur, then re-generate a H2O configuration by decreasing the timestep on line 92 in the input.equil file and repeating step 2.2.8 or starting again from step 2.1.1 .</li>
	</ol>
	</li>
</ol>
3. Extract the proper height of the supercell
<ol>
	<li>Execute the script get_npt_lz.py on the log.myadsorbate file by typing get_npt_lz.py log.myadsorbate at the command-line interface. This script outputs the average supercell height from the &#8220;production&#8221; run portion of the NPT simulation in the avg_lz.txt file. 
	NOTE: The get_npt_lz.py script assumes that LAMMPS writes the length of the cell z-dimension (lz) to the log.myadsorbate file every 1000 fs (customizable at line 20 of the get_npt_lz.py script), which is the default in the provided input.equil LAMMPS input file. The get_npt_lz.py script detects and discards the first 2 ns (customizable at line 19 of the get_npt_lz.py script) worth of lz values in the log.myadsorbate file, as they comprise the equilibration portion of the simulation, while the remaining 3 ns comprise the &#8220;production&#8221; portion and are thus used by the get_npt_lz.py script to compute the average z-dimension length. In addition to the avg_lz.txt file, the get_npt_lz.py script outputs a file called npt_data.txt, which provides values of lz as a function of the timestep, as well as a file called npt_plot.png, which plots the same data. The plot can be used to verify equilibration of the NPT simulation.</li>
	<li>Reconstruct the supercell using the average height determined in NPT.
	<ol>
		<li>Copy the data.myadsorbate_npt file into a new directory, referred to here as $WD2PATH, and rename it data.myadsorbate.</li>
		<li>Edit the new data.myadsorbate file so that the lz height is equal to the average value output from the get_npt_lz.py script by changing the zlo and zhi arguments in the data.myadsorbate file such that zlo is 0.0 and zhi is the lz value from the avg_lz.txt file produced in step 3.1.</li>
	</ol>
	</li>
</ol>
4. Generate configurations of H2O molecules
<ol>
	<li>Copy the LAMMPS input file input.prod into $WD2PATH. Edit the group variable on line 32 to indicate the atom type indexes for the water oxygen and water hydrogen atoms and the group variable on line 33 to indicate the atom type indexes for the Pt and adsorbate atoms.</li>
	<li>Execute the LAMMPS software by typing mpiexec -np XX lmp_mpi &#60; input.prod into the command-line interface, where XX is the number of CPU cores to use, and lmp_mpi is the name of the LAMMPS executable. Doing this will run a constant NVT simulation on the H2O molecules. The key output file from this simulation is the dump.myadsorbate.lammpstrj file. 
	NOTE: The duration of the NVT simulation should be long enough to comprise an equilibration run, where the energy of the system comes to steady state, and a production run, from which the ensemble averages (here, the spatial positions of the water molecules) are sampled. During the equilibration run, the energy of the system when plotted against time should level off to a steady state value. Once this occurs, the NVT simulation can be said to be in its production run.</li>
</ol>
5. Determine the hydrogen bond lifetime for proper time sampling
<ol>
	<li>Edit the hb_lifetime_dist.py script to specify: the timestep of the first frame of the dump.myadsorbate.lammpstrj file by changing the actualStart variable on line 22, how often frames are written to the LAMMPS trajectory file by changing the timestep variable on line 23, the first and last timesteps the script should consider (i.e., the production portion of the trajectory) by changing the N_first and N_last variables on lines 24 and 25, whether consecutive frames are considered or frames are skipped by changing the nevery variable on line 26, and the number of lines per frame section of the trajectory file by changing the frameLine variable on line 27. Additionally, edit lines 31 through 35 to specify which atom types within the data.myadsorbate file belong to the adsorbate and which atom types belong to the H2O molecules. 
	NOTE: The hb_lifetime_dist.py script analyzes the H2O configurations in the production run and determines if any H2O molecules are hydrogen bonded to the adsorbate. It then counts the simulation time that each hydrogen bond remains intact and reports this information as a distribution of hydrogen bond lifetimes in units of ps. The specific version of the script that is provided with this protocol assumes that LAMMPS writes the configuration of H2O molecules to the dump.myadsorbate.lammpstrj file every 1000 fs, which is the default in the provided input.prod LAMMPS input file. It detects and discards the first 2 ns worth of configurations in the the dump.myadsorbate.lammpstrj file, as they comprise the equilibration portion of the simulation, and uses the remaining 3 ns to calculate hydrogen bond lifetimes.</li>
	<li>Execute the script hb_lifetime_dist.py on the dump.myadsorbate.lammpstrj file by typing hb_lifetime_dist.py at the command-line interface. Doing this will produce a file called distribution_HB_lifetime.dat.</li>
	<li>Plot the data in the distribution_HB_lifetime.dat file to view the distribution of hydrogen bond lifetimes that occurred during the NVT simulation.</li>
	<li>Determine the time increment to use for the time sampling interval based on the calculated hydrogen bond lifetimes. The best choice is the maximum hydrogen bond lifetime; alternatively, a value that would capture the 95% confidence interval can be used.</li>
</ol>
6. Sample configurations of liquid H2O molecules
<ol>
	<li>Determine the number of configurations from the production run of the NVT FFMD trajectory for further calculations. The number of configurations should be selected so that the minimum time between configurations is equal to or greater than the time sampling interval identified in section 5.
	<ol>
		<li>Edit the default value for the num_frames variable on line 21 of the lammps_frames.py script to specify the number of configurations to extract.</li>
		<li>Execute the script lammps_frames.py on the file dump.myadsorbate.lammpstrj by typing lammps_frames.py at the command-line interface. Doing this will output a list of simulation times corresponding to the configurations that should be extracted from the dump.myadsorbate.lammpstrj file. These configurations can be used as starting structures in AIMD or QM simulations. 
		NOTE: 1) The lammps_frames.py script automatically detects the LAMMPS log and dump files as well as the production portion of the trajectory within the dump file and divides the number of configurations within the dump file into 10 groups. Alternatively, the user can specify the log file, the dump file, and the number of configurations from the command-line interface using the -l, -d, and -n options, respectively. To do so, the user should type lammps_frames.py -n XX -l $logfilename -d $dumpfilename at the command-line interface, where XX is the desired number of configurations, $logfilename is the name of the LAMMPS logfile, and $dumpfilename is the name of the LAMMPS trajectory (dump) file. The simulation times that are output refer to the median times in each group. 2) If the configurations will be calculated in VASP with the LDIPOLE flag turned on, a small layer of vacuum space should be added to the top of the supercell above the water layer. This will facilitate convergence of the electronic structure in the VASP calculation. Adding an additional 3 &#197; of vacuum space above the H2O molecules has been successful in the simulations discussed below.</li>
	</ol>
	</li>
</ol>

The authors disclose no conflicts of interest.

The method as presented was selected for its ease of implementation, but multiple customizations could be made. For one, the force fields used in the FFMD simulations can be modified. Changing the force field parameters and/or potentials can be done by editing the LAMMPS input and data files. Similarly, solvents other than H2O could be employed. To make this modification, the desired solvent molecule would need to be inserted starting from step 2.1.1, and the LAMMPS input files would need to be edited to incorporate the appropriate potentials and parameters. Inserting the new solvent molecule would also require supplying the internal coordinates of the solvent molecule in a .txt file analogous to the water.txt file.
Another modification that could be made is to modify the area of the surface slab. The results discussed in this manuscript employed 3 Pt x 3 Pt or 4 Pt x 4 Pt surface slabs, which have surface areas less than 120 &#197;2. As the slab surface area increases, the computational expense also increases. Computational expense has the largest impact on section 5 of this protocol. If the data processing steps in section 5 become computationally prohibitive, big data post processing strategies such as those discussed in Li et al. 201845 can be employed.
Possible sources of uncertainty for this procedure include the force field employed, the sampling method, and the sampling frequency. The water structure is determined by the force field that is used, meaning that the choice of force field could influence the specific configurations of H2O molecules. Our group has assessed how the choices of force field for H2O molecules and Pt atoms influence the interaction energies calculated in FFMD and found that the choice of force field contributes less than 0.1 eV to this interaction energy. Another source of uncertainty is the sampling method, which influences the specific configurations that are used to calculate a quantity of interest. Our group has compared the performance of the &#8220;time sampling&#8221; method presented in this protocol with an &#8220;energy sampling&#8221; method, which is biased to lower energy configurations of H2O molecules, on the interaction energies calculated in DFT and found both of these sampling methods give statistically equal values35,46. The sampling frequency can also influence the results. We have assessed how increasing the number of configurations from 10 to 30,000 influences the average interaction energies calculated in FFMD for 40 different C3HxO3 adsorbates and found that the sampling frequency contributes less than 0.1 eV to the average interaction energy44.
The main limitation to this method is that the adsorbates are approximated by the structures under vacuum during the FFMD simulations. In reality, the adsorbates would exhibit conformational changes (bond stretches, angle bends, torsional motions, etc.) due to normal thermal movements, including interactions with solvent molecules. Attempts to include conformational changes of adsorbates into the FFMD simulations would require detailed development of force fields for catalytic surface adsorbates, i.e., which comprise terms that describe bond stretches, angle bends, and torsional terms, amongst others. As a future direction of this protocol, we are developing such force fields for adsorbates at solid surfaces, which we will use to determine the extent to which using rigid adsorbates influences the results.

Modeling molecular phenomena involved in heterogeneous catalysis under liquid conditions is necessary for understanding catalytic function; however, this remains challenging because it requires a fine balance between chemical accuracy and computational expense. In general, since catalysis involves the breaking and forming of chemical bonds, quantum mechanics must be used to at least some degree; however, long simulations are challenging in quantum mechanics, as they require significant computer resources. Since molecules in the liquid phase are under constant thermal motion, simulations must also include configurational sampling, i.e., they must incorporate multiple spatial arrangements of the liquid molecules, as each different spatial arrangement (i.e., each configuration) has a different energy. This means that multiple configurations of liquid molecules must be simulated for each catalytic species of interest. These needs &#8211; to use quantum mechanics and to perform multiple calculations per catalytic species &#8211; can render modeling in heterogeneous catalysis under liquid phase computationally intractable. The purpose of the method described herein is to enable computationally tractable simulations of phenomena in heterogeneous catalysis under liquid phase.
We are particularly interested in heterogeneously catalyzed reactions that are carried out under liquid water. Water molecules have significant influence on catalytic phenomena, such as interacting with catalytic species (e.g., via dispersion forces and hydrogen bonding)1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23, participating in catalytic reactions1,7,8,9,15,21,22,24,25,26,27, and influencing reaction pathways and/or catalytic rates1,11,12,15,18,23,25,27,28,29,30,31. Modeling of these phenomena has been performed using QM and/or ab initio molecular dynamics (AIMD)1,2,6,7,14,22,25,27,28,32,33,34, force field molecular dynamics (FFMD)35, and quantum mechanics/molecular mechanics (QM/MM)10. In AIMD and FFMD, the atoms in the system are moved pursuant to Newton&#8217;s equations of motion according to the forces acting upon them. In AIMD, the system energy and forces are calculated with quantum mechanics, whereas in FFMD, the system energy and forces are calculated using force fields, which are algebraic expressions that are parameterized based on experimental or QM data. In QM/MM, the portion of the system where the bond breaking and forming occurs is calculated with QM, and the remainder of the system is calculated with MM, which employs force fields. Because they directly employ QM, AIMD and QM/MM are better suited for capturing the bond breaking and forming that occurs in aqueous phase heterogeneous catalysis; however, FFMD is significantly more computationally tractable and thus better suited for generating the configurations of liquid H2O molecules. The method presented in this protocol balances chemical accuracy and computational expense by employing a combination of QM and FFMD.
Specifically, this method uses FFMD simulations for generating configurations of liquid H2O and QM to calculate system energies. FFMD is carried out using LAMMPS.36 The force fields used in FFMD in this work employ Lennard-Jones + Coulomb (LJ+C) potentials, where the LJ parameters have been taken from the TIP3P/CHARMM model37 for H2O, the universal force field38 (UFF) for Pt, and the OPLS-AA force field39 for catalytic species, and the Coulomb parameters have been taken from the TIP3P/CHARMM37 model for H2O and the OPLS-AA force field39 for catalytic species. The Coulomb parameters for Pt atoms have been set to 0. QM calculations are performed using the VASP code40,41,42, which is a density functional theory (DFT) code. Water molecule insertions are performed with a code developed in-house called Monte Carlo Plug-in for Quantum Methods (MCPliQ). File conversions from VASP to LAMMPS in this protocol are performed with the Visual Molecular Dynamics (VMD) software43.
The protocol is intended to generate configurations of liquid water molecules around catalytic species on flat transition metal surfaces at low coverage. Coverage is denoted &#952;&#160;and defined as the number of adsorbates per surface metal atom (i.e., the number of surface adsorbates normalized by the number of metal atoms in the topmost layer of the metal slab in the catalyst model). In this manuscript, low coverage is defined as &#952; &#8804; 1/9 monolayer (ML), where 1 ML means one catalytic species per surface metal atom. The catalyst models should be placed in periodic simulation boxes. The simulation boxes do not have to be cubes. This manuscript demonstrates the use of the protocol for generating configurations of liquid H2O that can be used to calculate quantities of interest in aqueous phase heterogeneous catalysis.
This protocol requires that the user has access to installed and working versions of the VASP, MCPliQ, LAMMPS, and VMD software. More information about VASP (https://www.vasp.at/), LAMMPS (https://Lammps.sandia.gov/), and VMD (https://www.ks.uiuc.edu/Research/vmd/) are available on their respective websites. The MCPliQ software is documented at https://github.com/getman-research-group/JoVE_article, along with all input files and Python scripts mentioned in this protocol. This protocol assumes that the executables and scripts mentioned within will be run on a high-performance research computer and are installed in a directory that is in the user&#8217;s $PATH variable. If an executable or script is placed in a location that is not in the user&#8217;s $PATH, then the path to the executable must be included to execute it. Executables and scripts are executed in steps 2.1.2, 2.2.1, 2.2.8, 3.1, 4.2, 5.2, and 6.1.2. For example, to execute the MCPliQ code in step 2.1.2 from a directory that is not in the user&#8217;s $PATH, the user would type $PATHTOMCPLIQ/mcpliq at the command-line interface instead of mcpliq, where $PATHTOMCPLIQ is the location where the mcpliq executable has been stored (e.g., $PATHTOMCPLIQ might be ~/bin). Before starting this protocol, all executables and scripts should be given executable permissions (e.g., in Linux, this could be done by typing chmod +x mcpliq at the command-line interface from the directory where the mcpliq executable is stored). Further, any modules required by any of the software or scripts should be loaded (these dependencies will be specific to individual installations of the various software and the computer where the simulations will be run).

<table><tbody><tr><td>VASP software</td><td>Computational Materials Physics, Dept. of Physics, University of Vienna</td><td>vasp.5.4.4</td><td>Standard parallel VASP executable in the newest version.</td></tr><tr><td>LAMMPS software</td><td>Sandia National Laboratory</td><td>31Mar17-dp</td><td>Double-precision, parallel LAMMPS executable from 31 March 2017.</td></tr><tr><td>VMD software</td><td>Theoretical and Computational Biophysics Group, University of Illinois at Urbana-Champaign</td><td>1.9.3</td><td>Standard VMD executable in the newest version.</td></tr><tr><td>MCPliQ software</td><td>Getman Research Group, Dept. of Chemical and Biomolecular Engineering, Clemson University</td><td></td><td>Executable and input files for the MCPliQ software availabe from the Getman Research Group GitHub page.</td></tr><tr><td>JoVE article scripts</td><td>Getman Research Group, Dept. of Chemical and Biomolecular Engineering, Clemson University</td><td></td><td>Python scripts for this JoVE manuscript available from the Getman Research Group GitHub page.</td></tr><tr><td>H2O PDB file</td><td>Getman Research Group, Dept. of Chemical and Biomolecular Engineering, Clemson University or RCSB Protein Data Bank</td><td></td><td>PDB file for a water molecule, available from the Getman Research Group GitHub page or at http://www.rcsb.org/ligand/HOH.</td></tr></tbody></table>

multiscale sampling of a heterogeneous water/metal catalyst interface using density functional theory and force-field molecular dynamics

A significant number of heterogeneously-catalyzed chemical processes occur under liquid conditions, but simulating catalyst function under such conditions is challenging when it is necessary to include the solvent molecules. The bond breaking and forming processes modeled in these systems necessitate the use of quantum chemical methods. Since molecules in the liquid phase are under constant thermal motion, simulations must also include configurational sampling. This means that multiple configurations of liquid molecules must be simulated for each catalytic species of interest. The goal of the protocol presented here is to generate and sample trajectories of configurations of liquid water molecules around catalytic species on flat transition metal surfaces in a way that balances chemical accuracy with computational expense. Specifically, force field molecular dynamics (FFMD) simulations are used to generate configurations of liquid molecules that can subsequently be used in quantum mechanics-based methods such as density functional theory or ab initio molecular dynamics. To illustrate this, in this manuscript, the protocol is used for catalytic intermediates that could be involved in the pathway for the decomposition of glycerol (C3H8O3). The structures that are generated using FFMD are modeled in DFT in order to estimate the enthalpies of solvation of the catalytic species and to identify how H2O molecules participate in catalytic decompositions.

One use of this protocol is to calculate the energies of interaction between the liquid water and catalytic species, i.e., &#916;Eint35:
&#8710;Eint=ECatalytic species+H2O+EClean catalyst surface-ECatalytic species-EClean catalyst surface+H2O
where ECatalytic species+H2O is the energy of a configuration of H2O molecules around a catalytic species on a metal surface, EClean catalyst surface is the energy of the clean catalyst surface in vacuum, ECatalytic species is the energy of the catalytic species on a metal surface in vacuum, and EClean catalyst surface+H2O is the energy of the configuration of H2O over the catalyst surface with the catalytic species removed. The positions of the H2O molecules used to calculate ECatalytic species+H2O and EClean catalyst surface+H2O should be identical. All values of E are calculated using the VASP code. The quantity &#916;Eint includes all of the physical and chemical interactions between all of the molecules in the liquid water structure and the catalytic species and gives a reasonable estimate of the enthalpy of solvation of the catalytic species, which is needed to calculate its free energy of solvation and total free energy. Table 1 provides values for &#916;Eint calculated for species on a Pt(111) catalyst surface with chemical formulas equal to CxHyOz in units of eV (1 eV = 96.485 kJ/mol). The values were calculated at coverages &#8804;1/9 ML.35,46&#160;The values reported are the averages taken over 10 configurations of liquid H2O, and the uncertainties are reported as standard deviations. All the values are negative, indicating favorable interactions with water.
Another application of this protocol is to generate starting structures for AIMD. Movie 1 is a movie of an AIMD trajectory that was started from a configuration generated by this protocol. At the start of this movie, a COH adsorbate is shown on a Pt(111) surface under a structure of liquid H2O. One H2O molecule is emphasized, which formed a hydrogen bond with COH. Over the course of the movie, this H2O molecule abstracts the proton from the COH adsorbate and deposits a second hydrogen atom on the Pt(111) surface. The H2O molecule thus helps to catalyze the reaction COH* + * &#8594;&#160;CO* + H*, where the *s indicate catalytic sites. This simulation highlights the main strength and the main purpose of the multiscale sampling method described herein. Numerous configurations of H2O molecules are generated with FFMD, due to its strength in computational tractability. However, a limitation of FFMD is that it cannot capture bond breaking and forming unless a reactive force field is implemented. AIMD uses quantum mechanics to calculate energies and thus can capture bond breaking and forming. However, AIMD is too computationally demanding to generate all of the configurations of H2O molecules necessary to ensure sufficient sampling has been achieved. Thus, this protocol combines the two methods.
The structures of liquid H2O molecules generated by this procedure are dependent on input settings. Setting these improperly can have unintended influences on the water structures. For example, when the intermolecular distances become too small or when other parameters in the molecular dynamics input files are set improperly or take on unphysical values, the water structure can become unreasonable. Under these circumstances, the structure of water will &#8220;blow up&#8221; unintendedly during the FFMD trajectory. Figure 1 shows an example of this. The snapshot on the left-hand side is the starting structure for a FFMD run, and the snapshot on the right-hand side is a snapshot taken within 1 ps of starting the simulation. As can be seen, the H2O molecules have moved far away from the surface. This is caused by improper settings made in the simulation input files and is not a structure that is likely to occur in reality.
<img alt="Figure 1" class="xfigimg" src="/files/ftp_upload/59284/59284fig1v2.jpg" /> 
Figure 1: Example of a negative result. The force field molecular dynamics simulation &#8220;blows up&#8221; due to an unphysical setting or value. Left hand image: The starting geometry of the Pt(111) surface, adsorbate, and liquid water structure. Right hand image: The geometry of the Pt(111) surface, adsorbate, and liquid water structure less than 1 ps later. In the right-hand image, the H2O molecules have separated from the surface due to unphysically large forces. <a href="https://www-jove-com.remotexs.ntu.edu.sg/files/ftp_upload/59284/59284fig1v2large.jpg" target="_blank">Please click here to view a larger version of this figure.</a>
<img alt="Movie 1" src="/files/ftp_upload/59284/59284mov1.jpg" /> 
Movie 1: Ab initio molecular dynamics (AIMD) simulation initiated from a configuration generated in multiscale sampling. A H2O molecule that is originally hydrogen bonded to a COH adsorbate on a Pt(111) surface abstracts the proton from COH and deposits a second hydrogen on the Pt(111) surface. This bond breaking and forming event can be captured by AIMD but not with force field molecular dynamics (FFMD) unless a reactive force field is used. The initial configuration of H2O molecules used in this AIMD simulation was generated using FFMD as described in this manuscript. <a href="https://www-jove-com.remotexs.ntu.edu.sg/files/ftp_upload/59284/movie1.mp4" target="_blank">Please click here to view this video. (Right-click to download.)</a>
<table border="1" fo:keep-together.within-page="1" fo:keep-with-next.within-page="always">
	<tbody>
		<tr>
			<td>Catalytic Species</td>
			<td>&#8710;Eint (eV)</td>
		</tr>
		<tr>
			<td>COH</td>
			<td>-0.70 &#177; 0.07</td>
		</tr>
		<tr>
			<td>CO&#160;</td>
			<td>-0.03 &#177; 0.03</td>
		</tr>
		<tr>
			<td>CH2OH</td>
			<td>-0.64 &#177; 0.12</td>
		</tr>
		<tr>
			<td>CHO-CHOH-CH2OH</td>
			<td>-0.93 &#177; 0.22</td>
		</tr>
		<tr>
			<td>COH-COH-CH2OH</td>
			<td>-0.87 &#177; 0.23</td>
		</tr>
		<tr>
			<td>COH-CHOH-COH</td>
			<td>-1.72 &#177; 0.26</td>
		</tr>
		<tr>
			<td>CHOH-COH-CO</td>
			<td>-1.57 &#177; 0.25</td>
		</tr>
		<tr>
			<td>CHO-CO-CO</td>
			<td>-0.31 &#177; 0.19</td>
		</tr>
	</tbody>
</table>
Table 1: Water-catalytic species interaction energy results. Interaction energies in eV calculated for eight CxHyOz adsorbates on Pt(111). Values reported are the averages taken over multiple configurations of liquid H2O. The uncertainties are the standard deviations of the averages. 1 eV = 96.485 kJ/mol.

Watch this Scientific Journal Video about Multiscale Sampling of a Heterogeneous Water/Metal Catalyst Interface using Density Functional Theory and Force-Field Molecular Dynamics at JoVE.com

Multiscale Sampling of a Heterogeneous Water/Metal Catalyst Interface using Density Functional Theory and Force-Field Molecular Dynamics

The goal of the protocol presented here is to generate and sample trajectories of configurations of liquid water molecules around catalytic species on a flat transition metal surface. The sampled configurations can be used as starting structures in quantum mechanics-based methods.

A significant number of heterogeneously-catalyzed chemical processes occur under liquid conditions, but simulating catalyst function under such conditions ...

multiscale-sampling-heterogeneous-watermetal-catalyst-interface-using

Nanyang Technological University

Research

JoVE Journal

Chemistry

12.7K Views.  Clemson University. De protocolmodellen, thermische en kwantumverschijnselen in vloeibare fase heterogene katalyse. Het is voor zover wij weten de kwantummechanica met volledige bemonstering van een expliciete vloeibare omgeving. Configuraties van vloeibare moleculen gegenereerd uit dit protocol vertegenwoordigen die verwacht onder de werkelijke reactieomstandigheden en kunnen worden gebruikt om moleculair niveau verschijnselen die afhankelijk zijn van de ruimtelijke regelingen van moleculen te verkennen....

Video: Multiscale Sampling of a Heterogeneous Water/Metal Catalyst Interface using Density Functional Theory and Force-Field Molecular Dynamics

In Situ SIMS and IR Spectroscopy of Well-defined Surfaces Prepared by Soft Landing of Mass-selected Ions

Soft landing of mass-selected ions onto surfaces is a powerful approach for the highly-controlled preparation of materials that are inaccessible using conventional synthesis techniques. Coupling soft landing with in situ characterization using secondary ion mass spectrometry (SIMS) and infrared reflection absorption spectroscopy (IRRAS) enables analysis of well-defined surfaces under clean vacuum conditions. The capabilities of three soft-landing instruments constructed in our laboratory are illustrated for the representative system of surface-bound organometallics prepared by soft landing of mass-selected ruthenium tris(bipyridine) dications, [Ru(bpy)3]2+ (bpy = bipyridine), onto carboxylic acid terminated self-assembled monolayer surfaces on gold (COOH-SAMs). In situ time-of-flight (TOF)-SIMS provides insight into the reactivity of the soft-landed ions. In addition, the kinetics of charge reduction, neutralization and desorption occurring on the COOH-SAM both during and after ion soft landing are studied using in situ Fourier transform ion cyclotron resonance (FT-ICR)-SIMS measurements. In situ IRRAS experiments provide insight into how the structure of organic ligands surrounding metal centers is perturbed through immobilization of organometallic ions on COOH-SAM surfaces by soft landing. Collectively, the three instruments provide complementary information about the chemical composition, reactivity and structure of well-defined species supported on surfaces.

In Situ SIMS and IR Spectroscopy of Well-defined Surfaces Prepared by Soft Landing of Mass-selected Ions

Soft landing of mass-selected ions onto surfaces is a powerful approach for the highly-controlled preparation of novel materials. Coupled with analysis by in situ secondary ion mass spectrometry (SIMS) and infrared reflection absorption spectroscopy (IRRAS), soft landing provides unprecedented insights into the interactions of well-defined species with surfaces.

Soft landing of mass-selected ions onto surfaces is a powerful approach for the highly-controlled preparation of materials that are inaccessible using ...

Preparation and 3D Tracking of Catalytic Swimming Devices

We report a method to prepare catalytically active Janus colloids that "swim" in fluids and describe how to determine their 3D motion using fluorescence microscopy. One commonly deployed method for catalytically active colloids to produce enhanced motion is via an asymmetrical distribution of catalyst. Here this is achieved by spin coating a dispersed layer of fluorescent polymeric colloids onto a flat planar substrate, and then using directional platinum vapor deposition to half coat the exposed colloid surface, making a two faced "Janus" structure. The Janus colloids are then re-suspended from the planar substrate into an aqueous solution containing hydrogen peroxide. Hydrogen peroxide serves as a fuel for the platinum catalyst, which is decomposed into water and oxygen, but only on one side of the colloid. The asymmetry results in gradients that produce enhanced motion, or "swimming". A fluorescence microscope, together with a video camera is used to record the motion of individual colloids. The center of the fluorescent emission is found using image analysis to provide an x and y coordinate for each frame of the video. While keeping the microscope focal position fixed, the fluorescence emission from the colloid produces a characteristic concentric ring pattern which is subject to image analysis to determine the particles relative z position. In this way 3D trajectories for the swimming colloid are obtained, allowing swimming velocity to be accurately measured, and physical phenomena such as gravitaxis, which may bias the colloids motion to be detected.

A method to prepare catalytically active Janus colloids that can "swim" in fluids and determine their 3D trajectories is presented.

We report a method to prepare catalytically active Janus colloids that "swim" in fluids and describe how to determine their 3D motion using fluorescence ...

Engineering

Probing the Structure and Dynamics of Interfacial Water with Scanning Tunneling Microscopy and Spectroscopy

Water/solid interfaces are ubiquitous and play a key role in many environmental, biophysical, and technological processes. Resolving the internal structure and probing the hydrogen-bond (H-bond) dynamics of the water molecules adsorbed on solid surfaces are fundamental issues of water science, which remains a great challenge owing to the light mass and small size of hydrogen. Scanning tunneling microscopy (STM) is a promising tool for attacking these problems, thanks to its capabilities of sub-&#197;ngstr&#246;m spatial resolution, single-bond vibrational sensitivity, and atomic/molecular manipulation. The designed experimental system consists of a Cl-terminated tip and a sample fabricated by dosing water molecules in situ onto the Au(111)-supported NaCl(001) surfaces. The insulating NaCl films electronically decouple the water from the metal substrates, so the intrinsic frontier orbitals of water molecules are preserved. The Cl-tip facilitates the manipulation of the single water molecules, as well as gating the orbitals of water to the proximity of Fermi level (EF) via tip-water coupling. This paper outlines the detailed methods of submolecular resolution imaging, molecular/atomic manipulation, and single-bond vibrational spectroscopy of interfacial water. These studies open up a new route for investigating the H-bonded systems at the atomic scale.

Here, we present a protocol to investigate the structure and dynamics of interfacial water at the atomic scale, in terms of submolecular resolution imaging, molecular manipulation, and single-bond vibrational spectroscopy.

Water/solid interfaces are ubiquitous and play a key role in many environmental, biophysical, and technological processes. Resolving the internal ...

Synthesis and Performance Characterizations of Transition Metal Single Atom Catalyst for Electrochemical CO2 Reduction

This protocol presents both the synthesis method of the Ni single atom catalyst, and the electrochemical testing of its catalytic activity and selectivity in aqueous CO2 reduction. Different from traditional metal nanocrystals, the synthesis of metal single atoms involves a matrix material that can confine those single atoms and prevent them from aggregation. We report an electrospinning and thermal annealing method to prepare Ni single atoms dispersed and coordinated in a graphene shell, as active centers for CO2 reduction to CO. During the synthesis, N dopants play a critical role in generating graphene vacancies to trap Ni atoms. Aberration-corrected scanning transmission electron microscopy and three-dimensional atom probe tomography were employed to identify the single Ni atomic sites in graphene vacancies. Detailed setup of electrochemical CO2 reduction apparatus coupled with an on-line gas chromatography is also demonstrated. Compared to metallic Ni, Ni single atom catalyst exhibit dramatically improved CO2 reduction and suppressed H2 evolution side reaction.

Synthesis and Performance Characterizations of Transition Metal Single Atom Catalyst for Electrochemical CO2 Reduction

Here, we present a protocol for the synthesis and electrochemical testing of transition metal single atoms coordinated in graphene vacancies as active centers for selective carbon dioxide reduction to carbon monoxide in aqueous solutions.

This protocol presents both the synthesis method of the Ni single atom catalyst, and the electrochemical testing of its catalytic activity and selectivity ...

In situ FTIR Spectroscopy as a Tool for Investigation of Gas/Solid Interaction: Water-Enhanced CO2 Adsorption in UiO-66 Metal-Organic Framework

In situ infrared spectroscopy is an inexpensive, highly sensitive, and selective valuable tool to investigate the interaction of polycrystalline solids with adsorbates. Vibrational spectra provide information on the chemical nature of adsorbed species and their structure. Thus, they are very useful for obtaining molecular level understanding of surface species. The IR spectrum of the sample itself gives some direct information about the material. General conclusions can be drawn concerning hydroxyl groups, some stable surface species and impurities. However, the spectrum of the sample is &#34;blind&#34; with respect to the presence of coordinatively unsaturated ions and gives rather poor information about the acidity of surface hydroxyls, species decisive for the adsorption and catalytic properties of the materials. Furthermore, no discrimination between bulk and surface species can be made. These problems are solved by the use of probe molecules, substances that interact specifically with the surface; the alteration of some spectral features of these molecules as a result of adsorption provides valuable information about the nature, properties, location, concentration, etc., of the surface sites.
The experimental protocol for in-situ IR studies of gas/sample interaction includes preparation of a sample pellet, activation of the material, initial spectral characterization through the analysis of the background spectra, characterization by probe molecules, and study of the interaction with a particular set of gas mixtures. In this paper we investigate a zirconium terephthalate metal organic framework, Zr6O4(OH)4(BDC)6 (BDC = benzene-1,4-dicarboxylate), namely UiO-66 (UiO refers to University of Oslo). The acid sites of the UiO-66 sample are determined by using CO and CD3CN as molecular probes. Furthermore, we have demonstrated that CO2 is adsorbed on basic sites exposed on dehydroxylated UiO-66. Introduction of water to the system produces hydroxyl groups acting as additional CO2 adsorption sites. As a result, CO2 adsorption capacity of the sample is strongly enhanced.

In situ FTIR Spectroscopy as a Tool for Investigation of Gas/Solid Interaction: Water-Enhanced CO2 Adsorption in UiO-66 Metal-Organic Framework

The use of FTIR spectroscopy for investigation of the surface properties of polycrystalline solids is described. Preparation of sample pellets, activation procedures, characterization with probe molecules and model studies of CO2 adsorption are discussed.

In situ infrared spectroscopy is an inexpensive, highly sensitive, and selective valuable tool to investigate the interaction of polycrystalline solids ...

Vibrational Spectra of a N719-Chromophore/Titania Interface from Empirical-Potential Molecular-Dynamics Simulation, Solvated by a Room Temperature Ionic Liquid

The accurate molecular simulation prediction of vibrational spectra, and other structural, energetic and spectral characteristics, of photo-active metal-oxide surfaces in contact with light-absorbing dyes is an ongoing thorny and elusive challenge in physical chemistry. With this in mind, a molecular-dynamics (MD) simulation were performed by using optimized empirical potentials for a well-representative and prototypical dye-sensitized solar cell (DSC) solvated by a widely-studied room temperature ionic liquid (RTIL), in the guise of a [bmim]+[NTf2]- RTIL solvating an N719-sensitizing dye adsorbed onto 101 anatase-titania. In doing so, important insights were gleaned into how using a RTIL as the electrolytic hole acceptor modulates the dynamical and vibrational properties of a N719 dye, estimating the spectra for the DSC photo-active interface via Fourier transformation of mass-weighted velocity autocorrelation functions from MD. The acquired vibrational spectra were compared with the experiment spectra and those sampled from ab initio molecular dynamics (AIMD); in particular, various empirical-potential spectra generated from MD provide insight into how partial-charge charge parameterization of the ionic liquid affects vibrational spectra prediction. In any event, careful fitting of empirical force-field models has been shown to be an effective tool in handling DSC vibrational properties, when validated by AIMD and an experiment.

A dye-sensitized solar cell was solvated by RTILs; using optimized empirical potentials, a molecular dynamics simulation was applied to compute vibrational properties. The obtained vibrational spectra were compared with experiment and ab initio molecular dynamics; various empirical potential spectra show how partial-charge charge parameterization of the ionic liquid affects vibrational spectra prediction.

The accurate molecular simulation prediction of vibrational spectra, and other structural, energetic and spectral characteristics, of photo-active ...

Development of Heterogeneous Enantioselective Catalysts using Chiral Metal-Organic Frameworks (MOFs)

Substrate size discrimination by the pore size and homogeneity of the chiral environment at the reaction sites are important issues in the validation of the reaction site in metal-organic framework (MOF)&#8211;based catalysts in an enantioselective catalytic reaction system. Therefore, a method of validating the reaction site of MOF-based catalysts is necessary to investigate this issue. Substrate size discrimination by pore size was accomplished by comparing the substrate size versus the reaction rate in two different types of carbonyl-ene reactions with two kinds of MOFs. The MOF catalysts were used to compare the performance of the two reaction types (Zn-mediated stoichiometric and Ti-catalyzed carbonyl-ene reactions) in two different media. Using the proposed method, it was observed that the entire MOF crystal participated in the reaction, and the interior of the crystal pore played an important role in exerting chiral control when the reaction was stoichiometric. Homogeneity of the chiral environment of MOF catalysts was established by the size control method for a particle used in the Zn-mediated stoichiometric reaction system. The protocol proposed for the catalytic reaction revealed that the reaction mainly occurred on the catalyst surface regardless of the substrate size, which reveals the actual reaction sites in MOF-based heterogeneous catalysts. This method for reaction site validation of MOF catalysts suggests various considerations for developing heterogeneous enantioselective MOF catalysts.

Development of Heterogeneous Enantioselective Catalysts using Chiral Metal-Organic Frameworks MOFs

Here, we present a protocol for active site validation of metal-organic framework catalysts by comparing stoichiometric and catalytic carbonyl-ene reactions to find out whether a reaction takes place on the inner or outer surface of metal-organic frameworks.

Substrate size discrimination by the pore size and homogeneity of the chiral environment at the reaction sites are important issues in the validation of ...

Visualizing Solution Structure at Solid-Liquid Interfaces using Three-Dimensional Fast Force Mapping

Amongst the challenges for a variety of research fields are the visualization of solid-liquid interfaces and understanding how they are affected by the solution conditions such as ion concentrations, pH, ligands, and trace additives, as well as the underlying crystallography and chemistry. In this context, three-dimensional fast force mapping (3D FFM) has emerged as a promising tool for investigating solution structure at interfaces. This capability is based on atomic force microscopy (AFM) and allows the direct visualization of interfacial regions in three spatial dimensions with sub-nanometer resolution. Here we provide a detailed description of the experimental protocol for acquiring 3D FFM data. The main considerations for optimizing the operating parameters depending on the sample and application are discussed. Moreover, the basic methods for data processing and analysis are discussed, including the transformation of the measured instrument observables into tip-sample force maps that can be linked to the local solution structure. Finally, we shed light on some of the outstanding questions related to 3D FFM data interpretation and how this technique can become a central tool in the repertoire of surface science.

Here, we present a protocol for using three-dimensional fast force mapping - an atomic force microscopy technique - for visualizing solution structure at solid-liquid interfaces with the subnanometer resolution by mapping the tip-sample interactions within the interfacial region.

Amongst the challenges for a variety of research fields are the visualization of solid-liquid interfaces and understanding how they are affected by the ...

Precise Electrochemical Sizing of Individual Electro-Inactive Particles

Nanoimpact electrochemistry enables the time-resolved in situ characterization (e.g., size, catalytic activity) of single nanomaterial units, providing a means of elucidating heterogeneities that would be masked in ensemble studies. To implement this technique with redox inactive particles, a solution-phase redox reaction is used to produce a steady-state background current on a disk ultramicroelectrode. When a particle adsorbs onto the electrode, it produces a stepwise reduction in the exposed electrode area, which produces, in turn, a stepwise decrease in the current commensurate with the size of the adsorbing species. Historically, however, nanoimpact electrochemistry has suffered from &#34;edge effects,&#34; in which the radial diffusion layer formed at the circumference of the ultramicroelectrodes renders the step size dependent not only on the size of the particle but also on where it lands on the electrode. The introduction of electrocatalytic current generation, however, mitigates the heterogeneity caused by edge effects, thus improving the measurement precision. In this approach, termed &#34;electrocatalytic interruption,&#34; a substrate that regenerates the redox probe at the diffusion layer is introduced. This shifts the rate-limiting step of the current generation from diffusion to the homogeneous reaction rate constant, thus reducing flux heterogeneity and increasing the precision of particle sizing by an order of magnitude. The protocol described here explains the set-up and data collection employed in nanoimpact experiments implementing this effect for improved precision in the sizing of redox in-active materials.

As an analytical technique, nanoimpact electrochemistry, an increasingly important approach to counting and characterizing nanometer-scale, electro-inactive particles, suffers from poor precision due to the heterogeneous current distributions that arise from its use of ultramicroelectrodes. Outlined here is a generalized approach, termed "electrocatalytic interruption," that enhances precision in such measurements.

Nanoimpact electrochemistry enables the time-resolved in situ characterization (e.g., size, catalytic activity) of single nanomaterial units, providing a ...

Using Nitrogen Porosimetry to Characterize Metal-Organic Frameworks

Determining Surface Areas and Pore Volumes of Metal-Organic Frameworks

The surface area and pore volume of a metal-organic framework (MOF) can provide insight into its structure and potential applications. Both parameters are commonly determined using the data from nitrogen sorption experiments; commercial instruments to perform these measurements are also widely available. These instruments will calculate structural parameters, but it is essential to understand how to select input data and when calculation methods apply to the sample MOF. This article outlines the use of the Brunauer-Emmett-Teller (BET) method and Barrett-Joyner-Halenda (BJH) method for the calculation of surface area and pore volume, respectively. Example calculations are performed on the representative MOF UiO-66. Although widely applicable to MOFs, sample materials and adsorption data must meet certain criteria for the calculated results to be considered accurate, in addition to proper sample preparation. The assumptions and limitations of these methods are also discussed, along with alternative and complementary techniques for the MOF pore space characterization.

Author Spotlight: Characterizing Porous Materials for Aiding the Development of Robust Metal-Organic Frameworks with Adsorption Behavior 

This article describes the use of nitrogen porosimetry to characterize metal-organic frameworks, using UiO-66 as a representative material.

The surface area and pore volume of a metal-organic framework (MOF) can provide insight into its structure and potential applications. Both parameters are ...

Multiscale Sampling of a Heterogeneous Water/Metal Catalyst Interface using Density Functional Theory and Force-Field Molecular Dynamics

Samenvatting

Meer video's verkennen

Hoofdstukken in deze video