# [gmx-users] Which variable represents the potential energy ?

David spoel at xray.bmc.uu.se
Mon Dec 8 08:22:00 CET 2003

On Sun, 2003-12-07 at 23:07, masakatsu-ito at aist.go.jp wrote:
> Dear David,
>
> > > So I am trying to use my code instead of update() for integrating
> > > EOM. That is, my REMD framework would pass the atom positions to
> > > do_force() function of GROMACS, and would receive the the potential
> > > energy and force vector.  Then my framework integrates equations of
> > > motion with the Gaussian thermostatting technique.
>   ;
> > > So I have a favor to ask of you and any of GROMACS people.  Would you
> > > kindly tell me which part of do_md() is responsible for dealing with
> > > the periodic boundary condition? I think I have to shift an atom
> > > position when the atom happens to cross the boundary of a box. But I
> > > cannot find which function shift atom positions in do_md().
> >
> > part of it is in force.c
> > part of it is in update.c
> >
> > I doubt whether your code should be doing any PBC. Why do you think
> > that?
>
> Because I would like to minimize the modification of GROMACS code.
> But if I can find out which variables of GROMACS define a box for
> periodic boundary condition, it would be easier that my REMD framework
> shift an atom position.
I'm afraid not...
>
> So could you please tell me if box_size[] variables in do_md() defines
> the box for PBC? Can I shift i'th atom position (x[3*i], x[3*i+1],
> x[3*i+2]) in a way like the following?
>
>     x[3*i    ] = x[3*i] -
>              floor( x[3*i] / box_size[0] ) * box_size[0];    // for x
>     x[3*i + 1] = x[3*i + 1] -
>              floor( x[3*i +1] / box_size[1] ) * box_size[1]; // for y
>     x[3*i + 2] = x[3*i + 2] -
>              floor( x[3*i +2] / box_size[2] ) * box_size[2]; // for z
>
First, we don't shift atoms, but charge groups. Second, boxes may be
triclinic, Third shifting is done differently for non bonded forces than
for bonded forces.

> > > The other question is about sum_epot(). Would you kindly tell me why
> > > sum_epot() is not called immediately after calling do_force() in
> > > do_md() function? Should I call another function between do_force()
> > > and sum_epot() to get the potential energy?
> > >
> > because there are corrections which depend on virial (long range
> > corrections) which should be summed first, in particular when running in
> > parallel
>
> Then, does that long range correction work if I remove
> the update() procedure in do_md() function?
> And, is update() function unnecessary for calculating the potential
> energy and force vector?
If you use constraints you need to call update. It affects the forces
and the pressure, but not the energies.

It seems simplest to implement your code as another integrator in
update. The code for the real integrators is very simple.

>
>
> Thanks a lot in advances.
> (I feel sorry that I ask too many questions.)
> Masakatsu
>
> Masakatsu Ito , Ph.D
>
> Grid Technology Research Center
> National Institute of Advanced Industrial Science and Technology
> Tsukuba Central 2, Tsukuba, Ibaraki 305-8568 Japan
> Phone : +81-29-861-5730  Fax : +81-29-861-5301
> E-mail masakatsu-ito at aist.go.jp
>
> _______________________________________________
> gmx-users mailing list
> gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
--
David.
________________________________________________________________________
David van der Spoel, PhD, Assist. Prof., Molecular Biophysics group,
Dept. of Cell and Molecular Biology, Uppsala University.
Husargatan 3, Box 596,  	75124 Uppsala, Sweden
phone:	46 18 471 4205		fax: 46 18 511 755
spoel at xray.bmc.uu.se	spoel at gromacs.org   http://xray.bmc.uu.se/~spoel
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++