# [gmx-users] manual eq. 4.74-4.75 (dihedral restraints) head scratcher

Berk Hess gmx3 at hotmail.com
Sun Apr 18 17:48:09 CEST 2010

> Date: Sat, 17 Apr 2010 14:58:22 -0400
> From: jalemkul at vt.edu
> To: gmx-users at gromacs.org
> Subject: Re: [gmx-users] manual eq. 4.74-4.75 (dihedral restraints)	head	scratcher
>
>
>
> Daniel L. Ensign wrote:
> > Hello gmx-users, you rock and rollers,
> >
> > Equations 4.74 and 4.75 in my copy of the manual have (please pardon my
> > pseudo-LaTeX):
> >
> > (4.74)
> > \phi' = (\phi - \phi_0) MOD 2\pi
> >
> > (4.75)
> > V(\phi') = \frac{1}{2}k ( \phi' - \phi_0 - \Delta \phi )^2 if \phi' >
> > \Delta \phi
> > or
> > V(\phi') = 0 if \phi' \leq \Delta \phi
> >
> > but should there be absolute values around all of the \phi-\phi_0? Which
> > way is it in the code -- with the absolute distance between \phi and
> > \phi_0 or the directed distance?
> >
>
> It looks like absolute values are considered.  From src/gmxlib/dihres.c:
>
>      /* phi can jump if phi0 is close to Pi/-Pi, which will cause huge
>       * force changes if we just apply a normal harmonic.
>       * Instead, we first calculate phi-phi0 and take it modulo (-Pi,Pi).
>       * This means we will never have the periodicity problem, unless
>       * the dihedral is Pi away from phiO, which is very unlikely due to
>       * the potential.
>       */
>      dp = phi-phi0;
>      if (fabs(dp) > dphi) {
>        /* dp cannot be outside (-2*pi,2*pi) */
>        if (dp >= M_PI)
>      dp -= 2*M_PI;
>        else if(dp < -M_PI)
>      dp += 2*M_PI;
>
>
> > Also, as far as I can tell (and some mornings I definitely don't read
> > too good) neither the manual nor
> > http://wiki.gromacs.org/index.php/Dihedral_Restraints define the fields
> > in [ dihedral_restraints ], although the latter does name them. There, I
> > see
> >
> > [ dihedral_restraints ]
> > ; ai   aj    ak    al  type  label  phi  dphi  kfac  power
> >     5    7     9    15     1      1  180     0     1      2
> >
> > ai, aj, ak, al = atom numbers, obviously
> > type = ?, but I'm guessing there's only one type anyway
>
> Probably so.
>
> > label = what is this one?
>
> Looks to be bookkeeping.  The code doesn't seem to use it other than to print
> debug information, but I could be wrong since I haven't surfed around it very long.
>

power and label have been copied from distance restraints, but not used at all.
I have already removed them for the next release.

Berk

> > phi means phi0 ?
>
> Yes.
>
> > dphi means Delta phi, I guess
>
> Yes.
>
> > kfac is the force constant, probably
>
> Indirectly.  This term is equivalent to the fac value in distance restraints.
> Since the force constant is specified in the .mdp file, different restraints
> would otherwise have to be restrained with equivalent force constants.  The
> value of kfac is multiplied by the value of dihre_fc in the .mdp file, so that
> different restraints could have different force constants.
>
> > power = what is this one? Does 2 give me harmonic constraints?
> >
>
> Not a clue on this one.  Also doesn't seem to be used in the code, but maybe
> it's somewhere outside of dihres.c.
>
> > I'd be happy to translate any answers given to the wiki, assuming I get
> > answered and that I'm allowed to edit the wiki.
> >
>
> The wiki link you posted above doesn't actually exist any more, and actually
> points to an empty page, even though the content is still online:
>
> http://www.gromacs.org/index.php?title=Documentation/How-tos/Dihedral_Restraints
>
> The old wiki doesn't actually exist, but has been merged into the Gromacs site,
> so if you register as a user you can make contributions.  If you have troubles,
> I might get around to updating this page with the above information, unless
> there is more information to be considered, or I'm wrong :)
>
> -Justin
>
> > May the Force be with you (as long as it's calculated in hand-tuned
> > assembly loops),
> > Dan
> >
>
> --
> ========================================
>
> Justin A. Lemkul
> Ph.D. Candidate
> ICTAS Doctoral Scholar
> MILES-IGERT Trainee
> Department of Biochemistry
> Virginia Tech
> Blacksburg, VA
> jalemkul[at]vt.edu | (540) 231-9080
> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
>
> ========================================
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php

_________________________________________________________________
Express yourself instantly with MSN Messenger! Download today it's FREE!
http://messenger.msn.click-url.com/go/onm00200471ave/direct/01/
-------------- next part --------------
An HTML attachment was scrubbed...