2. Generate a solvated protein system

2.1. Generate a topology

Using the modified PDB file (chain A of 4AKE with crystal waters removed), generate a topology file for the CHARMM27 force field together with the TIP3P water model using the pdb2gmx tool:

cd top
gmx pdb2gmx -f ../coord/4ake_a.pdb -o protein.pdb -p 4ake.top -i protein_posre.itp -water tip3p -ff charmm27

Note

Total charge -4.000e (in the next step we will add ions to neutralize the system; we need a net-neutral system to properly handle electrostatics)

2.2. Solvate the protein

2.2.1. Adding water

Create a simulation box with editconf and add solvent with solvate:

cd ../solvation
gmx editconf -f ../top/protein.pdb -o boxed.pdb -c -d 0.5 -bt dodecahedron
gmx solvate -cp boxed.pdb -cs spc216 -p ../top/4ake.top -o solvated.pdb

Attention

In order to reduce the system size and make the simulations run faster we are choosing a very tight box (minimum protein-edge distance 0.5 nm, -d 0.5); for simulations you want to publish this number should be 1.2...1.5 nm so that the electrostatic interactions between copies of the protein across periodic boundaries are sufficiently screened.

solvate updates the number of solvent molecules (“SOL”) in the topology file (check the [ system ] section in top/system.top) [2].

2.2.2. Adding ions

Ions can be added with the genion program in Gromacs.

First, we need a basic TPR file (an empty file is sufficient, just ignore the warnings that grompp spits out by setting -maxwarn 10), then run genion (which has convenient options to neutralize the system and set the concentration (check the help!); genion also updates the topology’s [ system ] section if the top file is provided [2]; it reduces the “SOL” molecules by the number of removed molecules and adds the ions, e.g. “NA” and “CL”).

touch ions.mdp
gmx grompp -f ions.mdp -p ../top/4ake.top -c solvated.pdb -o ions.tpr
printf "SOL" | gmx genion -s ions.tpr -p ../top/4ake.top -pname NA -nname CL -neutral -conc 0.15 -o ionized.pdb

The final output is solvation/ionized.pdb. Check visually in VMD (but note that the dodecahedral box is not represented properly.)

Footnotes

[1]Often you would actually want to retain crystallographic water molecules as they might have biological relevance. In our example this is likely not the case and by removing all of them we simplify the preparation step somewhat. If you keep them, pdb2gmx in the next step will actually create entries in the topology for them.
[2](1, 2) The automatic modification of the top file by solvate and genion can become a problem if you try to run these commands multiple times and you get error messages later (typically from grompp) that the number of molecules in structure file and the topology file do not agree. In this case you might have to manually delete or adjust the corresponding lines in :file”system.top file.