This is an old revision of the document!

Molecular Dynamics simulation of a small molecule

TO USE THE FUNCTION LIBRARY (VERSION UP TO DATE) IN THE INTERACTIVE SHELL: you@eulerX ~$ module load courses mmm vmd ; mmm-init

REMEMBER: this is the command to load the module for the cp2k program:

you@eulerX ~$ module load new cp2k

and to submit the job:

you@eulerX ~$ bsub < jobname

Concerning temperature control, in these exercises we will use the NOSE-HOOVER chains method. This has been briefly described in the lecture, and is presented in this paper by Glenn Martyna (1992).

Download the 4.1 exercise into your $HOME folder and unzip it:

you@eulerX ~$ wget http://www.cp2k.org/_media/exercises:2015_ethz_mmm:exercise_4.1.zip
you@eulerX ~$ unzip exercises:2015_ethz_mmm:exercise_4.1.zip
you@eulerX ~$ cd exercise_4.1

All files of this exercise (input and scripts are all commented) can be also downloaded from the wiki: exercise_4.1.zip

You will start from a configuration already computed in the second lecture (inp.a.pdb) which is included in the repository of this exercise as well. Use the file inp.nve for the first simulation, which is a constant energy simulation.

md_part.inp.nve

&MD                                           ! This section defines the whole set of parameters needed perform an MD run.
  ENSEMBLE NVE                                ! The ensemble/integrator that you want to use for MD propagation
  STEPS 100000                                ! The number of MD steps to perform
  TIMESTEP [fs] 1.0                           ! The length of an integration step
  TEMPERATURE 100.0                           ! The temperature in K used to initialize the velocities with init and pos restart velocities
&END MD

Perform a constant energy simulation, 100000 time steps, with a time step of 1 fs.

you@eulerX exercise_4.1$ bsub cp2k.popt -i inp.nve -o out.nve

Using a different input file, modify the time step and the name of the project. Do it for 0.1, 2, 3, 4 fs.
Access the corresponding *.ener files. How is the energy conservation? How is the behavior of potential and kinetic energy, and how the temperature?

- Plot with gnuplot the different energy conservations and discuss them.

Perform now a constant Temperature simulations. The system is in contact with a thermostat, and the conserved quantity includes the thermostat degrees of freedom.

The first simulation is done at 100 K:

you@eulerX exercise_4.1$ bsub cp2k.popt -i inp.100 -o out.100

Then, perform a simulation at 300 K, using the restart file from the previous simulation: inp.300.

you@eulerX exercise_4.1$ bsub cp2k.popt -i inp.300 -o out.300

Now you have the following outputs to study with vmd:

nve_md-pos-1.pdb
md.100-pos-1.pdb
md.300-pos-1.pdb

Open (for example) nve_md-pos-1.pdb with VMD:

vmd nve_md-pos-1.pdb

Open Tk Console (Extensions menu > Tk console). And to define the two dihedrals (PHI and PSI) from there, you can enter:

source "dihedrals.vmd"

You can also pick from the extensions the “RMSD trajectory tool” and use it to align the molecule along the trajectory. Remember to replace “protein” with “all” in the selection, and then use “align”. You will see that now the molecule is well aligned along the path.

Using “Labels” menu, plot now the two dihedral angles graph.

Which differences do you notice between the nve, the 100 K and the 300 K case? Can you explain them?