MBG wiki | RecentChanges | Blog | 2021-09-24 | 2021-09-23

Performing a molecular dynamics simulation in a hexagonal cell

Using an orthogonal cell for a periodic boundary simulation of a cylindrically-shaped molecule is one of the least efficient choices : the diagonal distances between neighboring images of the solute are much larger than those aimed for. The result is that a significant unit cell volume (and a few thousand water molecules) are being simulated needlessly. A hexagonal cell is a much better solution (make a drawing of the lattice to convince yourself):

 dna in hexagonal cell

The difference in the number of waters needed for the simulation can be quite significant. With a border of 15 Angstroem and an orthogonal box NAMD reports (for the same solute) :

Info: ****************************
Info: STRUCTURE SUMMARY:
Info: 25117 ATOMS
Info: 17038 BONDS
Info: 9599 ANGLES
Info: 2228 DIHEDRALS
Info: 68 IMPROPERS
Info: 0 EXCLUSIONS
Info: 22 CONSTRAINTS
Info: 176 FIXED ATOMS
Info: 24569 RIGID BONDS
Info: 0 RIGID BONDS BETWEEN FIXED ATOMS
Info: 50254 DEGREES OF FREEDOM
Info: 8647 HYDROGEN GROUPS
Info: 108 HYDROGEN GROUPS WITH ALL ATOMS FIXED
Info: TOTAL MASS = 154535 amu
Info: TOTAL CHARGE = 8.23289e-07 e
Info: *****************************

If a hexagonal cell is used, the numbers are :

Info: ****************************
Info: STRUCTURE SUMMARY:
Info: 15436 ATOMS
Info: 10588 BONDS
Info: 6374 ANGLES
Info: 2228 DIHEDRALS
Info: 68 IMPROPERS
Info: 0 EXCLUSIONS
Info: 22 CONSTRAINTS
Info: 176 FIXED ATOMS
Info: 14894 RIGID BONDS
Info: 0 RIGID BONDS BETWEEN FIXED ATOMS
Info: 30886 DEGREES OF FREEDOM
Info: 5416 HYDROGEN GROUPS
Info: 108 HYDROGEN GROUPS WITH ALL ATOMS FIXED
Info: TOTAL MASS = 96259.7 amu
Info: TOTAL CHARGE = 8.23289e-07 e
Info: *****************************

Running a simulation with a hexagonal cell is not difficult but you will have to get a couple of sin()s and cos()s correct. The easy bit first [ and assuming that (i) you do a DNA run, and, (ii) you've already done the Dickerson tutorial ]. When using VMD's solvate pluggin, instead of

...
#
# Make water box
#
solvate psfgen_ions.psf psfgen_ions.pdb -o hydrated -b 1.80 -t 15.0

...

give

...
#
# Make water box
#
package require hexsolvate
hexsolvate psfgen_ions.psf psfgen_ions.pdb -o hydrated -b 1.80 -t 15.0

...

Now, the geometry exercise : with a hexagonal lattice the cell basis vectors are no longer orthonormal. You will have to determine them by hand and use your findings in the heating-up NAMD script. To aid your quest, here is an example :

If the output from pdbset for your (hexagonal) box is :

 Orthogonal Coordinate limits in output file:
                 Minimum   Maximum    Centre     Range
       On X :     -37.70     36.28     -0.71     73.98
       On Y :     -29.17     29.17      0.00     58.34
       On Z :     -25.62     25.62      0.00     51.24

then, the NAMD script should be given something like :

...

#
# Particle Mesh Ewald parameters. 
#
Pme                     on
PmeGridsizeX            80                      # <===== CHANGE ME
PmeGridsizeY            56                      # <===== CHANGE ME
PmeGridsizeZ            56                      # <===== CHANGE ME

...

cellBasisVector1        73.98   00.00   00.00   # <===== CHANGE ME
cellBasisVector2        00.00   43.75   25.26   # <===== CHANGE ME
cellBasisVector3        00.00   00.00   51.24   # <===== CHANGE ME
cellOrigin               0.00    0.00    0.00   # <===== CHANGE ME

...

Can you work out why ?