How to do MDs

Amber

Prepare the system (Protein(s) only)

Prepare your protein by removing all hydrogens (like -ignh for pdb2gmx) and feed it to tleap.

$ cat build.leap

source leaprc.protein.ff14SB
source leaprc.RNA.OL3
source leaprc.DNA.OL15
source leaprc.water.tip3p
mol = loadpdb mol.pdb
solvatebox mol TIP3PBOX 10
addions mol K+ 154 Cl- 91
saveamberparm mol ionized.parm7 ionized.rst7

Points to be caution and ALWAYS save you some time:

1. Properly distinguish all the components using chain, e.g. chain A for protein subunit A, chain B for DNA strand 1, chain C for DNA strand 2 etc.
2. Consider a fix.tcl (see below)
3. Remember to calculate number of ions using this
4. You may do a trial run of tleap to confirm the net charge of your solute
5. If a PDB file can go through pdb2gmx, it will be probably fine for tleap

$ cat fix.tcl
mol new mol.pdb
set sel [atomselect top "resname ILE and name CD"]
$sel set name CD1
set sel [atomselect top "resname ASN and name OC1"]
$sel set name O
set sel [atomselect top "resname ASN and name OC2"]
$sel set name OXT
set sel [atomselect top all]
$sel writepdb mol.pdb

The ionized.parm7 and ionized.rst7 are now ready for simulations!

MD parameters

Some sample parameters can be found here and the run commands are

export prmtop=ionized.parm7
pmemd.cuda -O -i 01_min.in -o 01_min.out -p $prmtop -c ionized.rst7 -r 01_min.rst -ref ionized.rst7
pmemd.cuda -O -i 02_heat.in -o 02_heat.out -p $prmtop -c 01_min.rst -r 02_heat.rst -x 02_heat.nc -ref 01_min.rst
pmemd.cuda -O -i 03_eq.in -o 03_eq.out -p $prmtop -c 02_heat.rst -r 03_eq.rst -x 03_eq.nc -ref 02_heat.rst
pmemd.cuda -O -i 05_Pull.in -o 05_Pull.out -p $prmtop -c 03_eq.rst -r 05_Pull.rst -x 05_Pull.nc
pmemd.cuda -O -i 06_Prod.in -o 06_Prod.out -p $prmtop -c 05_Pull.rst -r 06_Prod.rst -x 06_Prod.nc

Gromacs

Prepare the system (Protein(s) only)

Justin's tutorials are always the good starting points (here), especially the "Lysozyme in water" one (here).

Simple steps to prepare a system, and it usually works!

# Create a cubic box of 2 nm on each side and place the system to the center and align its principle axes to the reference axes
# -princ usually helps to reduce your system size
gmx editconf -f conf.gro -o editconf.gro -bt triclinic -d 2 -c -princ
# Solvate the box
gmx solvate -cp editconf.gro -o solvate.gro -p topol.top
# Add ions to make it neutral and of 0.15 M NaCl
# If you want KCl, add -pname K
gmx grompp -f ions.mdp -c solvate.gro -o ions.tpr -p topol.top
gmx genion -s ions.tpr -o ionized.pdb -conc 0.15 -neutral -p topol.top

Things to note:

1. You can find ions.mdp in Justin's tutorial or here (it does not even matter).
2. Biological systems might require surprisingly large boxes, if you are interested, read this.
3. You can consider using -d 1.2 if you are confident that your system does not extend its dimension a lot.

Standard MD (Protein(s) only)

A standard MD includes steps as follows:

1. Energy minimization
2. Annealing
3. Pre-equilibrium
4. Production run

Here are some standard run scripts.

gmx mdrun -pin on -ntmpi 1 -ntomp 12 -v -deffnm 

Things to note:

1. In step2-3, position restraints have been added to favor the initial structures (usually crystal)
2. You should consider release the restraints step-wisely (usually people are just too lazy to do so)

Trajectory manipulation

Usually you will have multiple trajectories from production runs, you would like to combine them and remove the PBC effects in order to calculate something or make videos.

gmx trjcat -f *.xtc -o xtc -dt 100 -n ndx

Notes:

1. -settime if you did not continue

gmx trjconv -f pre.trr -o unwrap.xtc -s tpr -n ndx -pbc mol (-dt 100)
gmx trjconv -f unwrap.xtc -o fit.xtc -s tpr -n ndx -fit rot+trans

Notes:

1. Usually protein for least squares fit and centering
2. Try "-pbc nojump" if "mol" does not work
3. "-pbc nojump" is good for dimer; "-pbc mol" is good for multimer
4.  !!!Never use "-center" with "-pbc mol"!!!

Split

echo 0 | gmx trjconv -f ed.trr -s ../ed.tpr -split 20
for ii in $(seq 0 10); do
   jj=$(($ii + 1))
   echo 0 | gmx trjconv -f trajout$ii.xtc -s ../ed.tpr -skip 20 -o $jj.gro
done
rm -f trajout$ii.xtc \#*

Anton

Alt 1:

Load dtr in vmd

pbc unwrap -all
animate write trr

gmx trjcat -o xtc
gmx trjconv -pbc nojump -o unwrap.xtc

Alt 2:

Load dtr in vmd

animate write trr

Load trr in vmd

pbc unwrap -all
animate write trr

gmx trjcat -o unwrap.xtc

Alt 3:

Load dtr in vmd

animate write trr

gmx trjcat -f *.trr -o xtc -cat
gmx trjconv -f xtc -o xtc -timestep 120
gmx trjconv -s tpr -pbc res -center -o unwrap.xtc