Hi,
Sorry for taking so long to get back.
On Wed, Mar 7, 2012 at 11:23 AM, Jan-Philip Gehrcke <jgehrcke.googlemail.com
> wrote:
> Hello!
>
> Within an MMPBSA analysis for a complex comprised of a charged receptor
> and a charged ligand, I would like to use a non-linear PB solver.
>
> I am using MMPBSA.py (of AmberTools 1.5), which itself is using PBSA in
> order to derive the electrostatic energy of solvation. As MMPBSA.py has
> no input variable that allows for selection of linear or non-linear
> solver, I did the following:
>
> 1) call to MMPBSA with -make-mdins
> 2) sed -ie 's/npbopt = 0/npbopt = 1/g' _MMPBSA_pb.mdin
> 3) call to MMPBSA with -use-mdins
>
This is the best (only) way of doing this short of changing the MMPBSA
source code itself. I've been tempted to add this option to the MMPBSA.py
input file for awhile, but having never performed non-linear PB
calculations myself, I've resisted.
> The third step fails with the following output:
>
> calculating complex contribution...
> PB Info in pb_read(): npopt has been overwritten with inp
> PB Info in pb_read(): eneopt has been reset to be 1 for nonlinear PB
> equation
> PB Bomb in pb_read(): cutnb=0 cannot be used with eneopt=1
> 0.000000000000000E+000
>
>
>
> MMPBSA.py, I believe, evaluates the energy "EEL" of the PBSA output. For
> eneopt=2 (default), this -- according to the AT manual -- is "the
> Coulombic energy". Now, the output above indicates that eneopt=1 is
> enforced in case of the non-linear solver. In this case, EEL "includes
> both the reaction field energy and the Coulombic energy".
>
I did not know this (I haven't studied the PB options in depth). I imagine
that it's difficult (impossible?) to extract the "gas phase" coulomb energy
in the context of the non-linear PB equation, so that option is not
supported.
> Does it make sense to evaluate this kind of EEL in terms of MMPBSA?
The only effect this has is that DELTA G solv and DELTA G gas as reported
are no longer properly decomposed. DELTA G gas is computed as the sum of
the gas phase interactions (bond/angle/dihedral, EEL, and VDW), whereas
DELTA G solv is composed of EPB/EGB and ESURF/ECAVITY. Since EPB and EEL
are combined into what MMPBSA.py thinks is a "gas phase" term, any attempt
by MMPBSA.py to separate the gas and solvation terms will fail.
Outside of this, though, the total DELTA G should be perfectly fine, since
everything is mixed together in the end in any case, anyway.
If
> it does, besides changing npbopt in _MMPBSA_pb.mdin I also have to set
> cutnb to a positive value. Which values would make sense?
>
When eneopt=1, as stated in the manual, cutnb is a cutoff applied _only_ to
VDW interactions (1/r**6), not to coulomb interactions. If you want to set
it as infinite, that's fine, I'm sure, but almost certainly unnecessary.
In typical Amber simulations (explicit solvent), the cutoff is set to 8 Å
(which is the direct space cutoff for electrostatics and the hard cutoff
for VDW interactions). If this cutoff is good enough for dynamics, I would
imagine it would be good enough for non-linear PB.
You could always test the effect an 8 Å vs. a 99 Å cutoff for 1 or 2
snapshots and see if it makes any kind of difference (in either timing or
numerical results).
HTH,
Jason
--
Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
352-392-4032
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Sat Mar 10 2012 - 14:00:02 PST
