Re: [AMBER] Calculating electrostatic potential in cpptraj

From: Covington, Cody Lance <cody.l.covington.Vanderbilt.Edu>
Date: Sat, 16 Aug 2014 20:11:24 +0000


> Hello all
>
> Is there a way to calculate the electrostatic potential at a given
> atom position or point in space in cpptraj?


​cpptraj cannot (yet) do PME electrostatics. There are a couple actions that do electrostatic calculations, but they do direct space pair-pair interactions only using minimum image convention.


> ​​
> Or any other analysis program?
>

​None in Amber that I know about.


> ​​
> I can calculate it using my own programs, but I want to use PME to get
> long range electrostatics that would go ​​ beyond where periodic
> boundary conditions end. This might seem silly, but I have a
> surfactant system that has a ​​ high ionic strength.
> ​​
>

​​
>
> ​​
> I have also tried to use
> ​​
> image xoffset 0 yoffset 0 zoffset 0
> ​​
> image xoffset 0 yoffset 0 zoffset 1
> ​​
> ... etc
> ​​
> To make the system more like and 'infinite' one, but this is quite
> cumbersome.
>

​It will be _very_ hard to make the system large enough that you can compute an electrostatic potential using just a direct sum-over-pairs (and it will get _very_ slow very quickly). Your best bet really is to use an Ewald-based method. There are a couple options:

1) Write your own PME or Ewald code. The electrostatic potential at any given point is given by the sum of the direct space ESP and the reciprocal space ESP. But I believe that the convolution of the charge grid ​with the electrostatic potential is done in reciprocal space and the result is back-transformed to real space. So I think you would want to back-transform the reciprocal space pair potential before convoluting it with the charge density (if what you're after is the ESP, that is).

2) Use an existing PME or Ewald code to extract what you need. The electric field is easy enough to get by dividing the electrostatic force by the charge of that particle. The ESP seems to me to be a little bit more challenging. Since I don't think you can pull out the reciprocal space pair potential from any of the existing PME/Ewald codes, you would effectively need to find the electrostatic energy of putting a unitary point charge at a particular location. So compute the full electrostatic energy (and only the electrostatic energy), put a unit charge at the location you want the ESP at, and compute another electrostatic energy. I believe the difference is the electrostatic energy of that point charge in the field generated by all of the particles in all periodic boxes. Divide by that charge and you get the ESP.

For #2, the easiest way I can think of doing this (and the way I would use if I had to do it) would be to use the OpenMM capabilities in ParmEd and write a Python script to compute the ESPs. The general steps here would be:

1) Use createSystem to create an OpenMM system with all of the relevant forces

2) Create two new, separate System objects with only a single NonbondedForce object. Populate the new NonbondedForce objects from the NonbondedForce object in the System you created in step 1, making sure to set all of the Lennard-Jones parameters to 0 so the only energies computed come from electrostatics. Add a new particle to the second system (and second NonbondedForce) that will be the unitary charge at the point in space you want to compute the ESP.

3) Create Context objects for each of these systems and compute their total potential energies.

4) Take the difference between the energies computed in (3) and divide by the charge you assigned to the particle you added.

The advantage of (1) is that you'll almost certainly learn a lot more about PME and Ewald. The advantage of (2) is that ParmEd provides Python classes for most of Amber's input and output files (restarts, trajectories, NetCDF restarts and trajectories, prmtop files, etc.) which should make your job easier, and you'll add a powerful tool to your repertoire (and you won't have to program and validate your own PME implementation, of course).

This was quite lengthy, but hopefully it was helpful, Jason

--
Jason M. Swails
BioMaPS,
Rutgers University
Postdoctoral Researcher
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Sat Aug 16 2014 - 13:30:02 PDT
Custom Search