Time evolution of angle made by the helix with the normal to the plane using gromacs

I am new to the md simulation community. I have a peptide which is helical in a bilayer. I want to calculate the angle of the peptide relative to the bilayer with gromacs. For that, I created an index with the peptide residue at the N terminal and the C terminal, so that I can then create a vector with the gmx gangle command. The problem is that for all my trials the computer says me : Inconsistency in user input: Number of positions in selection 2 in the first group not divisible by 3 and I don't understand what It means, I tried to find the explanation on internet but I didn't understand the explanation either. Can anyone of you kindly suggest how can I use this gmx gangle command to compute the angle made by the helical's axis with the normal to the plane? Thank you very much in advance.

• Can you please provide additional information on how is the axis of the helix defined? If yes, I can write an answer with the script that you can use. Apr 19 at 16:19
• I defined the helical axis from the initial two carbon alpha atoms, and two final carbon alpha atoms. Thank you very much for your help. Apr 20 at 8:30
• Why two carbon atoms? Apr 20 at 11:29
• I don't have a specific reason for it, but is there any other way around? Shall I use more Carbon alpha, or all carbon alpha? Apr 20 at 11:45
• Can you please add an image for the system that you are studying? Apr 20 at 12:48

1 Answer

set traj   name of trajectory file
set gro    name of corresponding gro file
set output [open filename.dat w]

mol load type trj $$traj waitfor all mol addfile gro$$gro

set numframes [molinfo top get numframes]
for {set frameC 0} {$$frameC <$$numframes} {incr frameC} {
molinfo top set frame $$frameC set sel1 [selection string for helix first atom] set sel2 [selection string for helix second atom] set vec [vecsub [$$sel1 get pos] [$$sel2 get pos]] set nvec [vecnorm$$vec]
set theta [lindex $$nvec 2] puts$$output "$$frameC$$theta"
}


I wrote this script for analyzing a NAMD DCD file and have not tested this for a gromacs trajectory. However, this should give you something to start with.