Trajectory shows long lines during VMD visualization

GROMACS version:2023
GROMACS modification: Yes/No
Here post your question

I have simulated a small RNA molecule for 1 us. This is a snapshot of the VMD visualization of the trajcectory of the molecule :

I corrected periodicity effect using the following:
gmx trjconv -s md_0_1.tpr -f md_0_1.xtc -o md_0_1_noPBC_whole.xtc -pbc whole
gmx trjconv -s md_0_1.tpr -f md_0_1_noPBC_whole.xtc -o md_0_1_noPBC_nojump.xtc -pbc nojump
gmx trjconv -s md_0_1.tpr -f md_0_1_noPBC_nojump.xtc -o md_0_1_noPBC_mol.xtc -pbc mol -center

first frame:
image

at around 200ns

at around 1us

what could be the reason behind this? I am using a triclinic box with -d 1.0 nm

On checking the minimum distance between periodic images I found the following:

From 200ns onwards, there seems to be some sort of interactions with the periodic images. Does the above visualization show exactly that? I don’t understand I have taken a sufficiently large box with distance from edge as 1.0 nm, why should there be any interactions among the periodic images?

It would be great if someone helped me out here- significance of these long lines in the trajetcory even after correcting pbc effects:

Thanks.

You likely have an issue due to the box shape. You’ve set up a rectangular box that mimics the dimensions of the nucleic acid, but as soon as it rotates, you violate the minimum image convention and get interactions across images (which you have lots of, making the simulation basically useless). Rotation then makes rewrapping the cell problematic and you’re seeing the outcome of that. When building the box, you have to consider what is appropriate if the molecule tumbles through space, which can happen rather quickly.

Thank you for the reply. I checked with dodecahedron box, these lines disappear from the trajectory. A few frames still show less than 1 nm minimum distance with periodic image. Ideally no frame should show any kind of interactions with periodic images but what happens when only a small number of frames show smaller minimum distance. non-bonded cut-off for my system is 1.2nm.

You have potentially a lot of frames with a minimum distance below the nonbonded cutoff of 1.2 nm, not a small number. And in each of these, you will have some amount of atoms that have duplicate forces computed on them.

Your box size should have (at minimum) a 1.2-nm buffer to the box edge, ideally slightly longer. There are long-range ordering effects from water that have to be accounted for, as well, though those effects may be small. But I would consider gmx editconf -d 1.2 the minimum you should consider for constructing your box.

I see. Thanks a lot for the suggestions.

@jalemkul @subha24 Thanks for this helpful information! I too have a very similar issue regarding long lines in the simulation.

Here’s big-picture what I am doing:

  1. I am modeling the self-assembly of a lipid bilayer in water (and 150 mM NaCl) in a 10x10x10 nm^3 box.

  2. I ran the gmx simulation for 300ns. When I visualized the outputted .xtc trajectory, I noticed long lines at all time points except the first. Below is an image at 300ns.


    (All subsequent images have to be GDrive links as this forum only allows me to embed one image per post. Moreover, I can only post two hyperlinks so please copy+paste the text of the URLs provided.)

  3. I then read this thread and noticed how the OP used trjconv to correct for periodicity effects.
    After running…
    gmx trjconv -s md_0_1.tpr -f md_0_1.xtc -o md_0_1_noPBC_whole.xtc -pbc whole
    … I was pleased to no longer find any long lines at any time point, but the membrane seemed to be at opposite ends of the box. Below is an image at 300ns.
    URL-Copy-Pastehttps://drive.google.com/file/d/1WCKfzSB9j49odhDRxnLQcJ88zbSJEuV5/view?usp=sharing

  4. To correct for this opposite-ended membrane, I read @jalemkul your tip from 2009 about using trjconv -center. Thus, I proceeded with running the two additional trjconv functions provided by the OP in this thread (which use -center).
    After running…
    gmx trjconv -s md_0_1.tpr -f md_0_1_noPBC_whole.xtc -o md_0_1_noPBC_nojump.xtc -pbc nojump
    gmx trjconv -s md_0_1.tpr -f md_0_1_noPBC_nojump.xtc -o md_0_1_noPBC_mol.xtc -pbc mol -center
    … I found that the long lines appear again beginning at 70ns.
    Below is an image at 40ns (which does appear to show the membrane nearly centered without any long lines).
    URL-Copy-Pastehttps://drive.google.com/file/d/11AZsB7M9A7hiGvoPF-4YLubW1ioEiDe8/view?usp=sharing

  • Below is an image at 300ns (with very long lines and the membrane once again at opposite ends).
    URL-Copy-Pastehttps://drive.google.com/file/d/1O5_UopSk0ZWhPnlJbs23sJoQWaCirXD5/view?usp=sharing

Considering this simulation takes significant time to run, I wanted to clarify the best path forward.

a) From your suggestion to the OP’s concerns, would you also recommend doing gmx editconf -d 1.2? Also, to clarify, this would have to proceed running the actual simulation, correct, so I would need to run the simulation again?

b) Do you have a sense why the output for -pbc whole doesn’t have these long lines but the raw one and nojump + mol + center one do?

c) Any other relevant tips about bilayer simulation? Basically would like to have the membrane centered in the box with water flanking both sides.

Thanks for the help!

The order of the trjconv commands are important. I’d recommend following this suggested workflow: https://manual.gromacs.org/current/user-guide/terminology.html#periodic-boundary-conditions.

Thanks, Magnus @MagnusL.

From the workflow you suggested, it states:

  1. whole
  2. cluster
    [3) nojump]
  3. center
    [5) put in box]
    [6) fit to other reference structure]

I tried following this protocol (and performing some deviations), but I am still getting some combination of structures that have the bilayer split at opposite ends and/or have long lines.

I understand each circumstance is unique, but perhaps there is a common workflow more specifically for bilayer self-assembly simulations?

I have thus far run the following commands:
gmx trjconv -s box_w_il_peg_chol_hl_water_NaCl_min_md.tpr -f box_w_il_peg_chol_hl_water_NaCl_min_md.xtc -o box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_whole.xtc -pbc whole

gmx trjconv -s box_w_il_peg_chol_hl_water_NaCl_min_md.tpr -f box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_whole.xtc -o box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_cluster.xtc -pbc cluster
Group for centering: 1 Other
Group for output: 0 System

gmx trjconv -s box_w_il_peg_chol_hl_water_NaCl_min_md.tpr -f box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_cluster.xtc -o box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_center.xtc -center

gmx trjconv -s box_w_il_peg_chol_hl_water_NaCl_min_md.tpr -f box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_center.xtc -o box_w_il_peg_chol_hl_water_NaCl_min_md_noPBC_boxmol.xtc -box 10 10 10 -pbc mol

This is the final structure at 300ns.

Any help would be greatly appreciated!

I don’t think you need the -cluster options and I don’t think you should rescale your box.

After some closer thought I think you can manage with:

  1. Making a selection group containing one atom in the center of the bilayer (pick an atom at the end of the tail of one lipid). Use gmx make_ndx, e.g., gmx make_ndx -f confout.gro -o 1_lipid_atom.ndx
  2. gmx trjconv -s topol.tpr -f traj_comp.xtc -o traj_comp_center1.xtc -pbc mol -center -n 1_lipid_atom.ndx. Choose the 1-atom selection group for centering and the whole system for output.
  3. gmx trjconv -s topol.tpr -f traj_comp_center1.xtc -o traj_comp_center1_centerlipids.xtc -pbc mol -center -n index.ndx (assuming you have your lipids as one selection group in index.ndx). Choose the lipids for centering and the whole system for output.

I hope that works.