How to calculate lipid SCD with gmx order?

GROMACS version: 2020.2

Hi, I am trying to use GROMACS to analyze the deuterium order parameter of POPC in a NAMD trajectory. I kept getting the error message “There are inconsistent shifts over periodic boundaries in a molecule type consisting of 59504 atoms. The longest distance involved in such interactions
is 6.236 nm which is above half the box length. Either you have excessively
large distances between atoms in bonded interactions or your system is
I’m trying to figure out what went wrong in my process and would appreciate some help with this.
Here is my process:

  1. created a tpr file with pdb, top files and a dummy mdp file (the top file was created from the pdb file with the vmd topo pluggin):
    gmx grompp -p -c system.pdb -f system.mdp -maxwarn 10 -o system.tpr
  2. created an index file for POPC sn-1:
    selected carbons on POPC sn-1:
    r POPC & a C22

    r POPC & a C218
    deleted irrelevant groups:
    del 0-29
  3. convert dcd to trr with VMD
  4. converted trr to xtc and tried to fix pbc error:
    gmx trjconv -f traj.trr -s system.pdb -pbc whole -o traj.xtc
    gmx trjconv -f traj.xtc-s -pbc nojump -o traj2.xtc
    gmx trjconv -f traj2.xtc -ur compact -o traj3.xtc
    gmx trjconv -f traj3.xtc -center -o traj4.xtc
    (I tried these steps in different orders and also tried to adjust the box size with ‘-box x, x,x’ but nothing seemed to work)
  5. gmx order -f traj4.xtc -s system.tpr -n sn1.ndx -od SCD_POPC_sn1.xvg

Thanks so much!


Unless you’ve got a really good reason to want to do so, you shouldn’t use gmx order for this. It’s designed for united-atom force fields and I’m guessing that you are likely using CHARMM36, based on the atom-names and the fact you’re using NAMD. Gmx order also doesn’t work correctly for unsaturated carbon atoms.

If you are using an all-atom force fields there are lots of other tools you can use straight away. Take a look at for several examples, such as the calc_op.tcl VMD script.