[gmx-users] Cannot compute velocity of COM of a group of atoms

Sat Aug 11 01:44:06 CEST 2018

Hi,

It doesn't. If the input trajectory doesn't have any velocities, there's no
way to compute any from e.g. forces, and in general they would not match
the original trajectory anyway.

On the original question, the alternative is to compute the displacement of
the COM over time and thus the velocity.

Mark

On Wed, Aug 8, 2018, 17:35 Justin Lemkul <jalemkul at vt.edu> wrote:

>
>
> On 8/8/18 10:33 AM, ARNAB MUKHERJEE wrote:
> > Hi,
> >
> > Thank you very much for your suggestion.
> >
> > By "nsteps" I meant the # of step interval between which I would record
> the
> > velocity. So my simulation run has 10 million steps (200 ns) and I chose
> > nstvout = 1000.
> >
> > But I was wondering is it possible to do a mdrun -rerun (changing the
> .mdp
> > file), and compute the velocity for the old trajectory (.trr) file,
> rather
> > than running a new simulation again?
>
> The -rerun function is for re-evaluating energies of provided
> configurations. I don't think it provides any reliable estimate of
> velocities.
>
> -Justin
>
> > Thanks and Regards,
> >
> > Arnab
> >
> > On Wed, Aug 8, 2018 at 4:15 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
> >
> >>
> >> On 8/7/18 1:47 PM, ARNAB MUKHERJEE wrote:
> >>
> >>> Dear all,
> >>>
> >>> I have simulated a system of DNA and Protein, and I want to calculate
> the
> >>> velocity of the center of mass of protein as a function of time. So I
> used
> >>> the following command :
> >>>
> >>> gmx traj -f traj_comp.trr -s md_run-E-Field.tpr -n index.ndx -ov
> >>> test-vel.xvg -com
> >>>
> >>> I have pasted below the file that it generates :
> >>>
> >>> # GROMACS:      gmx traj, VERSION 5.0.6
> >>> # Executable:   /applic/applications/GROMACS/5.0.6/bin/gmx
> >>> # Library dir:  /applic/applications/GROMACS/5.0.6/share/gromacs/top
> >>> # Command line:
> >>> #   gmx traj -f traj_comp.trr -s md_run-E-Field.tpr -n index.ndx -ov
> >>> test-vel.xvg -com
> >>> # gmx is part of G R O M A C S:
> >>> #
> >>> # Georgetown Riga Oslo Madrid Amsterdam Chisinau Stockholm
> >>> #
> >>> @    title "Center of mass velocity"
> >>> @    xaxis  label "Time (ps)"
> >>> @    yaxis  label "Velocity (nm/ps)"
> >>> @TYPE xy
> >>> @ view 0.15, 0.15, 0.75, 0.85
> >>> @ legend on
> >>> @ legend box on
> >>> @ legend loctype view
> >>> @ legend 0.78, 0.8
> >>> @ legend length 2
> >>> @ s0 legend "Protein X"
> >>> @ s1 legend "Protein Y"
> >>> @ s2 legend "Protein Z"
> >>> "test-vel.xvg" 24L, 731C                                      24,1
> >>> Bot
> >>>
> >>> So it doesn't have the velocities. Is it due to the fact that in my
> .mdp
> >>> file, I had set nstvout  = 0 ? I have pasted below the initial part of
> my
> >>> .mdp file
> >>>
> >>> title       = NVT equilibration with position restraint on all solute
> >>> (topology modified)
> >>> ; Run parameters
> >>> integrator  = md        ; leap-frog integrator
> >>> nsteps      = 5000000   ; 1 * 500000 = 500 ps
> >>> ;nsteps          = 5000
> >>> dt          = 0.02      ; 1 fs
> >>> ; Output control
> >>> nstxout     = 0 ; save coordinates every 10 ps
> >>> nstvout     = 0 ; save velocities every 10 ps
> >>> nstcalcenergy   = 50
> >>> nstenergy   = 1000  ; save energies every 1 ps
> >>> nstxtcout       = 2500
> >>> ;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0
> ps
> >>>                                ; nstxout-compressed replaces nstxtcout
> >>> ;compressed-x-grps  = System  ; replaces xtc-grps
> >>> nstlog      = 1000  ; update log file every 1 ps
> >>> ; Bond parameters
> >>> continuation    = no   ; first dynamics run
> >>> constraint_algorithm = lincs ; holonomic constraints
> >>> constraints = none          ; all bonds (even heavy atom-H bonds)
> >>> constrained
> >>> ;lincs_iter = 2                 ; accuracy of LINCS
> >>> lincs_order = 4                 ; also related to accuracy
> >>>
> >>> Is there a way I can still compute the velocities, or do I need to run
> the
> >>> simulation again, putting nstvout = nsteps ?
> >>>
> >> Yes, you need a non-zero value for nstvout, but you shouldn't set it to
> >> the same value as nsteps, because then you'll only get one frame. You
> need
> >> some sensible output frequency that allows you to collect relevant
> >> statistics of the quantity of interest.
> >>
> >> -Justin
> >>
> >> --
> >> ==================================================
> >>
> >> Justin A. Lemkul, Ph.D.
> >> Assistant Professor
> >> Virginia Tech Department of Biochemistry
> >>
> >> 303 Engel Hall
> >> 340 West Campus Dr.
> >> Blacksburg, VA 24061
> >>
> >> jalemkul at vt.edu | (540) 231-3129
> >> http://www.thelemkullab.com
> >>
> >> ==================================================
> >>
> >> --
> >> Gromacs Users mailing list
> >>
> >> * Please search the archive at http://www.gromacs.org/Support
> >> /Mailing_Lists/GMX-Users_List before posting!
> >>
> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >>
> >> * For (un)subscribe requests visit
> >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> >> send a mail to gmx-users-request at gromacs.org.
> >>
>
> --
> ==================================================
>
> Justin A. Lemkul, Ph.D.
> Assistant Professor
> Virginia Tech Department of Biochemistry
>
> 303 Engel Hall
> 340 West Campus Dr.
> Blacksburg, VA 24061
>
> jalemkul at vt.edu | (540) 231-3129
> http://www.thelemkullab.com
>
> ==================================================
>
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>