I'm finishing up a set of improvements to the mdgx code, and I thought I'd
share a bit of what I've done to see what people might think or want to
try. I will confess to not knowing much about AMBER's capabilities for
treating extra points: my understanding is that AMBER is hard-coded to
implement some "famous" cases of extra points with information about the
frames being specified through topology data that would normally apply to
regular atoms. Therefore, when sander or pmemd parse a TIP4P topology,
they get the O-EPW distance from a bond equilibrium value (%FLAG
BOND_EQUIL_VALUE); similarly with TIP5P. I've programmed mdgx to
recognize TIP4P and TIP5P in this manner as well; there are "hard-coded"
extra point rules for dealing with these sorts of things just like in
sander and pmemd. However, I wanted to make mdgx capable of more general
implementations. It's possible to place massless sites in virtually any
geometry, whether the frame is rigid or not, for frames of as little as
two or as many as four sites. I figured, with the efforts that are spent
parameterizing and simulating "all-atom" models, it may be helpful to have
more flexibility in just where charges get placed. From the standpoint of
coding, I also like having a demarcation between what are "extra points"
and what are normal atoms, particularly because this frees us from having
to pass information about the frames through topology data that is really
meant for other things.
So, here's an example of what I'm now able to do with mdgx. This input
file will morph a system of TIP3P waters into TIP4P-Ew:
>>>>
&rule
ResidueName WAT,
ExtraPoint EPW,
EPtype HW,
FrameStyle 2,
FrameAtom1 O,
FrameAtom2 H1,
FrameAtom3 H2,
Vector12 0.106676721,
Vector13 0.106676721,
Charge -1.04844,
&end
&rule
ResidueName WAT,
ExtraPoint O,
FrameStyle 0,
Charge 0.0,
Sigma 3.16435,
Epsilon 0.16275,
&end
&rule
ResidueName WAT,
ExtraPoint H1,
FrameStyle 0,
Charge 0.52422,
Sigma 0.0,
Epsilon 0.0,
&end
&rule
ResidueName WAT,
ExtraPoint H2,
FrameStyle 0,
Charge 0.52422,
Sigma 0.0,
Epsilon 0.0,
&end
<<<<
The first &rule namelist says "add an extra point to residues named WAT,
using the frame style 2 with the frame O, H1, H2, and vector coefficients
as given." The second three rules have frame style 0, which means "no
extra point, just find all residues WAT and change the nonbonded
parameters of the O/H1/H2 atoms as given." mdgx then takes these rules and
modifies the topology at run time, inserting extra points, altering
charges, adding new van-der Waals types, and meticulously updating all the
connectivity. I seem to be getting the right results under NPT dynamics,
so I think it really has made TIP4P-Ew out of TIP3P.
Of course, as I said, mdgx can also just read a LEaP-generated prmtop
specifying TIP4P or TIP5P. The point is that, using the same &rule
syntax, I can now add new sites to proteins or other small molecules
starting from a plain ff99 topology.
If anyone wants to offer some feedback, I'd like to see if people have
interest in this sort of capability. My interest is to find out whether
breaking out of the one-atom, one-site paradigm might help us develop a
force field that is overall more compatible with TIP4P-Ew, one that can
improve hydration free energy estimates even without going to
polarization. I also have some fluctuating-charge polarizable models in
mind, but they're much further off. One might have to take some care with
the accuracy of PME calculations, but I suspect this sort of custom extra
points approach might also serve to test a quick-and-dirty static dipole
or static quadrupole model (the multipoles being approximated as tight
clusters of extra points for proof-of-concept). Lots of possibilities
when you can put charges wherever you want...
Dave Cerutti
_______________________________________________
AMBER-Developers mailing list
AMBER-Developers.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber-developers
Received on Thu May 26 2011 - 16:00:02 PDT
