Glyala is one of the simplest molecules that exhibits some important features common to larger biomolecules. In particular, it has more than one long-lived conformation, which we will identify in this exercise by mapping out its potential energy surface.

The conformations of glyala dipeptide are characterized by the dihedral angles of the backbone. Below, we color carbons in green, hydrogens in white, oxygen in red and nitrogen in blue, i.e. the torsional angle $\phi$ is N-C-C-N , while $\psi$ is C-N-C-C along the backbone.

 tar -xvf glyala-epot.tar.gz 

With this knowledge at hand, we will fix the dihedral angles and perform geometry optimization for all remaining degrees of freedom.

We have prepared a CP2K input file water.inp for running a MD simulation of liquid water using the force field from the first exercise (parametrized by Praprotnik et al.). Download water.tar.gz

Repeat the MD using initial temperatures 200 and 400 K. In order not to overwrite any of your previous files, it is advisable to run the new simulations in different folders.

Next we are going to analyze the trajectories in order to calculate the radial distribution function (rdf, $g(r)$) as a function of temperature.

VMD comes with an extension for exactly this purpose: In the VMD Main window open “Extensions → Analysis” click on “Radial Pair Distribution function $g(r)$”. In the appearing window use “Utilities → Set unit cell dimensions” to let VMD know the simulation box you used. After that use Selection 1 and 2 to define the atomic types that you want to calculate the rdf for, for example “element H”.

Then we will calculate diffusion coefficient. It is a proportionality constant between the molar flux due to molecular diffusion and the gradient in the concentration of the species (or the driving force for diffusion), which is defined by: $6D=\lim_{t\to\infty} \ \frac{\delta <r^2(t)>}{\delta t}$

To evaluate this expression, all that is needed is to evaluate at each point in time in the calculation the average of the square of the distance that each atom has traveled since the start of the production phase of the dynamics, and examining the slope of this function at a long time. By storing the initial coordinates, it is straightforward to evaluate the square of the distance. However, some care is needed due to the use of periodic boundary conditions: the program stores x, the coordinates, but in many programs, during the dynamics, if any atom has its x, y, or z coordinate become larger than box size or smaller than zero, then it is moved back into the other side of the box. This has the effect of making the raw distance traveled meaningless. The value of D is obtained from the slope, at a long time, of the right-hand side of the above equation (you need to divide by six to obtain D, take the slope, and also be careful with units).

VMD comes with an extension for exactly this purpose: In the VMD Main window open “Extensions → Analysis” click on “RMSD Trajectory Tool”. In the appearing window use “all” to let VMD know the molecule you want to track. Tick “Plot”, and press “RMSD”, you will have the RMSD plot for the water system.

Now, we will move to a more realistic system - Glyala in water. We will preformed a MD of glyala in water and save the trajectory.

The initial geometry provided in the PDB file is a glyala molecule solvated by 73 water molecules. The geometry is not equilibrated. You need first to equilibrate the system at 300K. When the system is equilibrated, you need to analysis the result.

 Selection 1, Selection 2 

 element O and not same residue as element C