[gmx-users] Elastic network model

[mohan maruthi sena]

Hi all,
            I want to create an elastic network model for a protein with
372 residues. Generally in elastic network  models we consider C-alpha
atoms of the residues, connected to all other atoms which fall under
certain cutoff distance (let say 0.8nm). I wish to consider only bonded
stretching potential and does not want to consider non-bonded, angle &
dihedral terms.  Through my earlier postings regarding to same query, I got
to know that I have to create topology file which should contain
information about bonding information with all other atoms.
When I am performing MD run the atoms are overlapping each other and
ultimately the system gets contracted.
When I use integrator = nm , I am getting negative eigen values and rmsf
calculated from g_anaeig option has all zero values.
Please look in to the topology, mdp file and suggest me a way to build an
elastic network model. I do not understand where I went wrong.

Please find the topology file at the following link:

http://researchweb.iiit.ac.in/~bipin.singh/atn.txt



; RUN CONTROL PARAMETERS
integrator               = md
dt                       = 0.002 ; time step (ps)
nsteps                   = 1000000 ; number of steps
disre                   = simple

; OUTPUT CONTROL OPTIONS
nstenergy                = 500
nstxout                  = 500
nstvout                  = 500
nstfout                  = 500
energygrps               = System
;

; NEIGHBORSEARCHING PARAMETERS
nstlist                   = 0      #Frequency to update neighbourlist

pbc                       = no
;
comm_mode                 = ANGULAR


continuation              = yes
; Temperature coupling
tcoupl                   = v-rescale           ; Couple temperature to
external heat bath according to Berendsen method
tc-grps                  = System  ;Non-Protein ; Use separate heat baths
for Protein and Non-Protein groups
tau_t                    = 0.2      ;0.2         ; Coupling time constant,
controlling strength of coupling
ref_t                    = 300        ; Temperature of heat bath


; GENERATE VELOCITIES FOR STARTUP RUN
gen_vel                  = no    ; Assign velocities to particles by taking
them randomly from a Maxwell distribution
gen_temp                 = 300.0  ; Temperature to generate corresponding
Maxwell distribution
gen_seed                 = 9999   ; Seed for (semi) random number
generation. Different numbers give different sets of velocities

; OPTIONS FOR BONDS
constraints              = none ; No constraints except for those defined
explicitly in the topology



My mdp file is as follows:

; VARIOUS PREPROCESSING OPTIONS
title                    = Position Restrained Molecular Dynamics
cpp                      = /lib/cpp


; RUN CONTROL PARAMETERS
integrator               = md
dt                       = 0.002 ; time step (ps)
nsteps                   = 1000000 ; number of steps
disre                   = simple

; OUTPUT CONTROL OPTIONS
nstenergy                = 500
nstxout                  = 500
nstvout                  = 500
nstfout                  = 500
energygrps               = System


; NEIGHBORSEARCHING PARAMETERS
nstlist                   = 0      #Frequency to update neighbourlist

pbc                       = no

comm_mode                 = ANGULAR


continuation              = yes
; Temperature coupling
tcoupl                   = v-rescale           ; Couple temperature to
external heat bath according to Berendsen method
tc-grps                  = System  ;Non-Protein ; Use separate heat baths
for Protein and Non-Protein groups
tau_t                    = 0.2      ;0.2         ; Coupling time constant,
controlling strength of coupling
ref_t                    = 300        ; Temperature of heat bath


; GENERATE VELOCITIES FOR STARTUP RUN
gen_vel                  = no    ; Assign velocities to particles by taking
them randomly from a Maxwell distribution
gen_temp                 = 300.0  ; Temperature to generate corresponding
Maxwell distribution
gen_seed                 = 9999   ; Seed for (semi) random number
generation. Different numbers give different sets of velocities

; OPTIONS FOR BONDS
constraints              = none ; No constraints except for those defined
explicitly in the topology

refcoord_scaling         = com


Thank you in advance,

with regards,
Mohan