[gmx-users] Calculating electric field vector at snapshots

[David van der Spoel]

Aaron Fafarman wrote:
> Thank you for considering my proposed scheme for calculating the
> electric field vector at a single site at every snapshot of an MD
> trajectory. As discussed previously on this list, I have attempted to
> use a virtual site with no mass and with an infinitesimal charge,
> constructed from two atoms in the region of interest, and printing the
> output of the force on the virtual site at every step.  The problem I
> have run into is that the force on the virtual site is zero at every
> step (see gmxdump below ) . I think this is because mdrun has already
> redistributed the forces from the virtual site to the constructing
> particles. Is there a way (a small code modification perhaps) that I
> could capture the information of the force on the virtual sites before
> it's been redistributed?
> 
> 
> excerpt from gmxdump of trajectory from mdrun with integrator = md;
> atom 1796 is the virtual site
> .....
>       f[ 1795]={ 8.67412e+02, -1.41267e+02,  2.31425e+02}
>       f[ 1796]={ 0.00000e+00,  0.00000e+00,  0.00000e+00}
>       f[ 1797]={-3.19446e+02, -1.38371e+03, -8.82301e+02}
> .....
> 
> 
> Just in case there's some other explanation for the zero forces, I've
> also included excerpts of the output of the gmxdump of the .tpr file
> below:
> ....
> atom[  1796]={type=  9, typeB=  9, ptype=   VSite, m= 0.00000e+00, q=
> 1.00000e-04, mB= 0.00000e+00, qB= 1.00000e-04, resnr=  119} grpnrs=[ 0
> 0 0 0 0 0 0 0 0 0 ]}
> ....
> Virtual site 2:
>             nr: 4
>             multinr[division over processors]: 4
>             iatoms:
>                0 type=364 (VSITE2) 1796 210 211
> ....
> 
> Thanks in advance for your help.
> 
> -Aaron
> 
> 
>> Message: 5
>> Date: Thu, 05 Jun 2008 08:53:32 +0200
>> From: David van der Spoel <spoel at xray.bmc.uu.se>
>> Subject: Re: [gmx-users] Calculating electric field vector at
>>        snapshots
>> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>> Message-ID: <48478D6C.3050904 at xray.bmc.uu.se>
>> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>>
>> Aaron Fafarman wrote:
>>> Dear GMX community,
>>>
>>> I would like to calculate the electric field vector at a point in
>>> space in a protein at many snapshots during an MD trajectory. I
>>> haven't found a built-in function for doing this in Gromacs--is there
>>> one? If not, would people please comment on this proposal for how to
>>> do this:
>>>
>>> 1. Include a virtual site referenced to the two atoms (constructing
>>> particles) that define the point in space I would like to interrogate.
>>> Give it only a very small charge, so small that the force from the
>>> virtual site to the constructing particles is negligible.
>>>
>>> 2. Have the calculation output the force on the virtual site at a
>>> determined interval by setting nstfout in the.mdp file
>>>
>>> 3. Divide the force by the small charge, and hence get the field
>>>
>>> It seems like this should work as long as the force on virtual sites
>>> is reported in the output of forces from nstfout, and as long as using
>>> such a small charge on the virtual site will not run into rounding
>>> errors. Also I'm presuming that the coulomb interaction is not
>>> calculated between the constructing atoms and the virtual site--is
>>> that correct? Any thoughts?
>>>
>> I think it will work. Just try it.

It could be that the force is set to zero after distributing it. Check 
the src/mdlib/vsite.c code. You could modify this code, such that it 
does not reset the forces. Be careful that this does not generate other 
problems (you can check by comparing runs with and without this feature 
using gmxcheck -e -e2 -tol 0).

The proper way to do this would be to recompute the non-bonded forces 
(short range only) with LJ parameters set to zero. The contribution to 
potential and field from PME can be extracted without too much trouble.



>>
>>> Thanks,
>>>
>>> -Aaron
>>> _______________________________________________
>>> gmx-users mailing list    gmx-users at gromacs.org
>>> http://www.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
>>
>> --
>> David.
>> ________________________________________________________________________
>> David van der Spoel, PhD, Assoc. 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://folding.bmc.uu.se
>> ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
>>
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.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


-- 
David van der Spoel, Ph.D., Professor of Biology
Molec. Biophys. group, Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205. Fax: +4618511755.
spoel at xray.bmc.uu.se	spoel at gromacs.org   http://folding.bmc.uu.se