[gmx-developers] Implementation of fluctuating charge into GROMACS: Pressure calculation

[David van der Spoel]

On 2010-10-16 22.16, Lee-Ping Wang wrote:
> Hi there,
>
> I’ve gotten the QTPIE interaction to work with PBC using a
> custom-designed shift function; the matrix elements are taken smoothly
> to zero at the cutoff radius. This ensures that there is no force
> discontinuity, which I’ve found causes some artifacts in the dynamics
> and energy minimization.
>
> Implementing the Ewald summation has been a bit problematic, because as
> you said, there needs to be a pairwise decomposition of the
> electrostatic energy in computing the matrix elements. I think I have
> something that works, but the forces aren’t implemented yet, and it’s
> expensive because it requires a double loop in the k-space summation
> (the normal Ewald summation only needs a single loop). Hopefully, the
> shift function will perform well enough and I won’t worry so much about
> the Ewald summation.
>
> I am trying to build the QTPIE forces into the “fshift” vector to get
> the pressures right, but I don’t have any way of telling whether my
> calculated pressures are correct. I suspect there is a problem, since
> the box keeps on blowing up if I turn pressure coupling on, but I need
> some reference method to compare against; how do you all make sure that
> the pressure is being computed properly? Do you simply repeat energy
> calculations at different box sizes and compute the pressure by finite
> difference?

You can run NVT and monitor the pressure for one thing. But the way to 
go is to put something really simple in a box, e.g. two water molecules, 
such that you can test your implementation analytically.

>
> Thanks,
>
> - Lee-Ping
>
> *From:* gmx-developers-bounces at gromacs.org
> [mailto:gmx-developers-bounces at gromacs.org] *On Behalf Of *Gerrit Groenhof
> *Sent:* Wednesday, October 13, 2010 2:34 AM
> *To:* Discussion list for GROMACS development
> *Subject:* Re: [gmx-developers] rlist larger than L/2 in simple neighbor
> search
>
> On 10/12/2010 06:49 PM, Lee-Ping Wang wrote:
>
> Hi there,
>
> For the past few months, I’ve been adding a charge equilibration model
> to GROMACS (to be specific, QTPIE – a refinement of QEq). For the model
> to work properly, each atom must “see” each other atom in the
> electrostatic interaction (it is a screened Coulomb interaction). This
> is easy enough to do with pbc = no; I just set all the cutoffs to zero,
> and add a new neighbor-list that includes excluded atom pairs.
>
>
> Since we'd be interested in that, I would appreciate if you can keep us
> informed of your progress.
>
> Now I would like to run the simulation with PBC turned on, which means
> either implementing Ewald summation or using cutoffs for the
> electrostatic interaction. When cutoffs are turned on (with or without
> PBC), the model goes haywire even with very long cutoff lengths. I
> haven’t implemented Ewald summation for the screened Coulomb
> interaction, and I’m not sure if it’s worth the effort if an easier fix
> is possible.
>
>
> If I recall correctly, the screened Coulomb here means the overlap
> integral between s-type Slater functions. This interaction behave as a
> normal coulomb at large distances, so that coudl explain the issues with
> increasing the cut-off. But an Ewald expression will be difficult as
> well, since the all pairs-wise interactions are required to compute the
> charge flow parameters. I have so far only seen it applied to molecules.
> For non bonded system like a box of water, a direct charge flow between
> atoms far apart is unlikely though. Could one not simply scale also the
> charge transfer variables to zero with distance as well?
>
> As a quick fix, I changed the simple neighbor-searching code to set rcut
> = 1000000 when QTPIE is turned on.
>
> I am worried about the other possible consequences of my change; will
> this adversely affect the molecular dynamics or create artifacts in some
> unexpected way?
>
> You will get the ordinary artefacts caused by cutoff, like ordering etc.
>


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