Electrostatic potential and forces in OOPSE

[Dan Gezelter]

I've been reading that paper (which is very
interesting, by the way), and am trying to figure
out the best way to include the Wolf potential
in OOPSE.    If I'm reading it correctly, in order
to evaluate the total energy, it would need to work
in a manner very similar to how we do the reaction
field.  The energy on each particle appears to depend
on the total charge contained within the cutoff radius,
right?  So you'd need to make one pass to calculate
the total charge within the radius, and a second pass to
calculate the forces.  I need to read it more
carefully to be sure, but that is a start.

You can effectively turn off the current switching
behavior by not specifying any cutoffGroups
and by setting the switchingRadius and cutoffRadius
to the same value.  Then you could modify the
charge-charge interaction to include the Wolf-screened
potential.

As for the updates in OOPSE-2.0, there has been
a complete rewrite of the C++ side, and the build process,
and the changes we are making in electrostatics are to
allow point-multipole interactions up to
quadrupole-quadrupole.   The new version needs some
testing, but should be closer to seeing daylight
in a week.

Thanks for pointing out the error in the docs.  I'll
fix it.

Best regards,

  --Dan Gezelter

On Mar 2, 2005, at 3:43 PM, Francesco Mercuri wrote:

>
> Dear Professor Gezelter,
> let me first thank you for your help concerning my last message.
> I'm actually trying to introduce a different potential for
> electrostatic interactions in OOPSE, since they are currently
> implemented as a simple charge-charge scheme, which could give
> rise to some inaccuracies. Due to the problems related to the
> implementation of the Ewald summation, we are thinking about
> a modified Wolf-damped Coulomb potential, as described in
> Zahn et al., J. Phys. Chem. B 2002 (106) 10725. As described
> in the paper, this scheme provides Ewald-quality results,
> without the need for considering K points.
> However, I noted that in OOPSE the evaluation of forces related to
> switched potentials (as LJ or Coulomb) proceeds according to two
> subsequent steps: a first term is evaluated as the product
> between the the derivative of the potential, for an atoms pair, and the
> switching function. Then, all these terms are summed up, within
> a single group, and a second term is added to the forces, as the
> product between the total potential (within the same group) and
> the derivative of the switching function (mass-weighted).
> I guess that evaluating forces as the sum of these two terms
> leads to gread advantages in term of performances. But I wonder
> whether an eventual implementation of a shifted-forces potential,
> like the one described in the aforementioned paper, could still
> be possible. Indeed, in this case the potential is no longer
> expressed as a product of a pure Coulomb term and a switching
> function, and its derivatives (forces) have to be evaluated
> every time a charge-charge interaction is computed. Do you
> see any workaround for that? I tried to replace the
> Coulomb*Switch electrostatics with the potential given in eq. (13)
> of the paper cited above, and forces (eq. 12) accordingly,
> in the do_charge_pair subroutine, including the contribution
> for the switching function derivatives (in the do_force_loop)
> only for LJ terms. But it seems to me that substituting a 
> "group-based" correction to forces (as done when adding the switching 
> function derivatives) with a "pair-based" term (as in eq. 12 of Zahn 
> et al) could
> give rise to some inconsistencies. Anyway, do you think that
> this approach could be an efficient way for improving electrostatics?
> May be I'm missing some subtleties in this simple analysis.
> Moreover, in your last message you mentioned the release of the 2.0 
> version of the code. On which time scale do you expect it would be
> released? There will be any improvement concerning electrostatics?
> And a final remark: in eq. (5) of the OOPSE manual (or eq. 3.3 of the
> online version) the denominator of the switching function should
> be (rcut-rsw)^3 instead of (rcut-rsw)^2, right?
> Thank you in advance and best regards,
>
> Francesco Mercuri
>
***********************************************
   J. Daniel Gezelter
   Assistant Professor
   Department of Chemistry and Biochemistry
   251 Nieuwland Science Hall
   University of Notre Dame
   Notre Dame, IN 46556-5670

   phone:  +1 (574) 631-7595
   fax:    +1 (574) 631-6652
   e-mail: gezelter at nd.edu
   web:    http://www.nd.edu/~gezelter
************************************************