MD Model Potential

We are going to run all the simulations in the HPC facility of IIT Kanpur. For this you need to login into your HPC account by typing

ssh -X username@hpc2010.hpc.iitk.ac.in

Then enter your password, which is already given to you.

The files required for running the simulation are already present in the home folder of your account. 

For today's tutorial go to the folder T1_MD by following command

cd T1_MD/

Introduction

The aim of this tutorial is to run a molecular dynamics (NVT) simulation on a model system and to generate the free energy surface from the probability distribution.

The model system:

The model system that we are going to use for our simulation is a one-dimensional double well potential with barrier equal to 0.002 a.u. and having two minima at -1 and 1 (shown below). The potential has form

This function can be plotted  by  the gnuplot software.

For opening gnuplot, type gnuplot in the terminal. Once gnuplot terminal appears, then you need to type the following:

load "01-plot.gpi"

The following plot will appear on your screen.

To close gnuplot, type q in the gnuplot terminal and press ENTER.

Now, we will perform MD simulation (NVT) on this model potential energy surface at a temperature of 65 K (kBT<<barrier height). 

Enter the directory MD_65K by typing

cd MD_65K/

In this folder, we have an input file  01-input_md that contains following input parameters:

Now, we are going to submit the run for MD simulation by the following command:

qsub sge_submit.sh

You can check the status of your run by typing: 

qstat -u username

After the run is over, the files TRAJ, AVGTEMP and ENERGIES will be generated as output containing the following information:

• TRAJ: step number, position, velocity, potential energy

• AVGTEMP: step number, instantaneous temperature

• ENERGIES: step number, potential energy, kinetic energy, total energy

• free_energy: coordinate, probability, free energy

Before analysing the results, we are going to perform the same MD simulation at a higher temperature (650 K;  kBT>barrier height).

Now go to the directory MD_650K 

cd ../MD_650K

Now, we are going to submit the run for MD simulation by the following command:

qsub sge_submit.sh

Once the runs are over, we will transfer the data to your system. 

For this we need to open a new terminal and make a new directory Workshop_Tutorial

mkdir Workshop_Tutorial

cd Workshop_Tutorial

Now copy the data from HPC to our system by the following command

scp  -rp username@hpc2010.hpc.iitk.ac.in:/home/username/T1_MD .

Analysis of Results:

Now, we can see the trajectory of the particle moving on the potential energy surface.  This can be  visualized from the TRAJ file, by plotting the coordinates versus potential energy together with the function f(x).  

In case of simulation at 65 K, 

• we need to go to the directory MD_65K.

• Open gnuplot

• Load the script  02-plot.gpi

Here we observe the particle is stuck in one of the well and the barrier crossing did not take place.

In order to plot the probability as a function of x, load the script 01-prob.gpi in gnuplot. 

To plot free energy as a function of x, load the script 02-prob.gpi in gnuplot

In the case of low temperature, it can be seen that the free energy surface is mimicking the bottom of one of the wells, but not the full potential energy landscape.

In case of simulation at 650 K, 

• we need to go to the directory MD_650K.

• Open gnuplot

• Load the script  02-plot.gpi

In order to plot the probability as a function of x, load the script 01-prob.gpi in gnuplot.

The free energy as a function of x can be reconstructed in the similar way.

To plot free energy as a function of x, load the script 02-prob.gpi in gnuplot