===============================================================================
Poisson-Bolztmann Equation Module
The PBEQ module allows the setting up and the numerical solution of
the Poisson-Boltzmann equation on a discretized grid for a solute molecule.
WARNING: The module has not been used in numerous situations and caution
should be the rule. All problems should be reported to
Benoit Roux at rouxb@PLGCN.Umontreal.CA, phone (514) 343-7105 or
Dmitrii Beglov at beglovd@PLGCN.Umontreal.CA,
phone (514) 343-6111x3953.
Wonpil Im at imwo@PLGCN.Umontreal.CA, phone (514) 343-6111x3953.
* Menu:
* Syntax:: Syntax of the PBEQ commands
* Function:: Purpose of each of the commands
* Examples:: Usage examples of the PBEQ module
Syntax
[SYNTAX PBEQ functions]
Syntax:
PBEQ enter the PBEQ module
END exit the PBEQ module
Subcommands:
SOLVe convergence-specifications grid-specifications
membrane-specifications physical variable-specifications
boundary-specifications force-specifications
atoms-selection
ITERate convergence-specifications
ENPB
CAPAcitance
WRITE property [CARD write-range] [UNIT integer]
COOR coordinate-manipulation-command
SCALar scalar-manipulation-command
HELP
RESET
convergence-specifications::=[MAXI integer] [DEPS real] [DOMEga real]
MAXI [2000] : number of iterations
DEPS [0.000002] : parameter of convergence
DOME [1.0] : initial mixing factor
grid-specifications::= [NCEL integer] [DCEL real]
NCEL [65] : number of grid point in 1D
DCEL [0.1] : size of grid unit cell
physical variable-specifications::= [WATR real] [EPSW real] [EPSP real]
[CONC real] [TEMP real]
WATR [0.0] : solvent probe radius
EPSW [80.0] : bulk solvent dielectric constant
EPSW [1.0] : protein interior dielectric constant
CONC [0.0] : salt concentration
TEMP [300.0] : Temperature
membrane-specifications:: [TMEMb real] [HTMEmb real] [ZMEMb real] [EPSM real]
[EPSH real] [VMEMB real]
TMEMB [0.0] : thickness of membrane (along Z)
HTMEMB [0.0] : thickness of headgroup region
ZMEMB [0.0] : membrane position (along Z)
EPSM [1.0] : membrane dielectric constant
EPSH [EPSM] : membrane headgroup dielectric constant (optional)
VMEMB [0.0] : potential difference across membrane (entered in [volts])
boundary-specifications::= [SMOOTH] [SWIN real]
SWIN [0.5] : solute-solvent boundary Smoothing WINdow
force-specifications::= [FORCE] [STEN real] [NPBEQ integer]
STEN [0.0] : surface tension coefficient (in kcal/mol/A^2)
NPBEQ [1] : the frequency for calculating solvation forces
during minimizations and MD simulations
atoms-selection::= a selection of a group of atoms
range::= XFIRST real XLAST real YFIRST real YLAST real ZFIRST real ZLAST real
property::= [PHI] { KCAL/MOL }
{ VOLTS }
[FKAPPA2]
[CHRG]
[EPSX]
[EPSY]
[EPSZ]
[TITLE]
PHI : Electrostatic Potential [ KCAL/MOL ]
[ VOLTS ]
FKAPPA2 : Debye screening factor
CHRG : charges on the lattice
EPSX : X sets of dielectric function
EPSY : Y sets of dielectric function
EPSZ : Z sets of dielectric function
TITLE : formatted title line
General discussion regarding the PBEQ module
1. SOLVE
Prepare grids and solve PB equation for the selected atoms and return the
electrostatic free energy in ?enpb = (1/2)*Sum Q_i PHI_i over the lattice.
The factor of 1/2 is there for the linear response free energy of charging.
The atomic contributions are returned in WMAIN (destroying the radii).
The maximum number of iterations can be specified. The convergence parameters
DEPS and DOMEga should not be modified. The numerical algorithm is pretty
much the same as DelPhi as described in: Klapper et al. Proteins 1, 47 (1986).
2. ITERATE
Pursue the iteration on the grid. SOLVE must have been called first.
3. ENPB
Compute the electrostatic PB energy Sum Q_i PHI_i over the lattice.
Notice that the electrostatic energy is twice as much as the electrostatic
free energy (see above). The value of the electrostatic energy is passed
through the substitution parameter enpb.
4. CAPACITANCE
Compute the capacitance based on the net induced charge in the double
layer. The induced charge beyond the limits of the box are estimated based on
the analytical solution to a planar membrane.
5. WRITE
The WRITE command is used to write out the grid properties. By default,
a binary file of the property will be written for the whole grid. The keyword
CARD implies that a formatted output will be produced. In that case, the
spatial range can be specified for the output. By default, the electrostatic
potential PHI is given in [UNIT CHARGE]/[ANGS]. If specified, the PHI can be
given in [VOLTS] or in [KCAL/MOL].
7. RESET
Resets all assignements of the PBEQ module and free the HEAP array.
Destroys all lists and grids. By default, the grids and arrays remain assigned
when exiting and re-entering the module. This is to allow multiple call to
PBEQ without having to free the HEAP and other arrays if they are going to be
used again.
8. Miscellaneous command manipulations
*note misc: (miscom.doc) are supported within the PBEQ module,
allowing opening and closing of files, streaming of files, label assignments
(e.g., LABEL), and repeated loops (e.g., GOTO), parameter substitutions
(e.g., @1,@2, etc...) control (e.g., IF 1 eq 10.0 GOTO LOOP) and CALC
(e.g., CALC energy = ?enpb).
9. COORMAN and SCALAR commands
*note misc: (corman.doc) and (scalar.doc) are supported within
the PBEQ module, allowing the easy manipulation of charges, radii, rotation
and translations of molecules, etc...
10. SMOOTH
SMOOTH changes the attribute of the solute-solvent boundary. By default
(NO SMOOTH), the boundary is defined by the van der Waals surface or the
molecular surface (with WATR). SMOOTH keyword change the boundary as a region
having +/- SWIN (Smoothing WINdow) from the surface of the solute. Within the
solute-solvent boundary, the dielectric constant and the Debye screening factor
will be changed continuously from EPSP and 0 to EPSW and the screening factor
at bulk solvent.
11. FORCE
This keyword invokes the calculation of the solvation free energy and
forces and must be followed by SMOOTH keyword. The solvation energy is taken
as a sum of electrostatic and nonpolar solvatin energy. The former is
calculated from the PB equation and the latter by using the surface tension
coefficient (STEN) that relates free energy with surface area. Note that the
calculated surface is approximately the van der Waals surface. If membrane is
considered, the surface of the membrane is also approximately included.
The corresponding forces are also calculated and will be used in minimizations
and MD simulations where NPBEQ can be used to specify the frequency for
calculating the solvation forces. Note that SWIN must be equal or greater to
DCEL to get correct solvation free energy and forces.
12. A set of "ATOMIC BORN RADII"
Atomic radii derived from solvant electrostatic charge distribution may be
used (test/data/radius.str). These radii were tested with free energy
perturbation with explicit solvant. (Reference: M. Nina, D. Beglov and B. Roux.
Atomic Radii for Continuum Electrostatics Calculations based on Molecular
Dynamics Free Energy Simulations. J. Phys. Chem., 1997.)
There are no definitive method to modify the optimized radii when using SWIN
at this point. A preliminary attemp is to do as follows (with DCEL of 0.2):
STREAM RADIUS.STR
SCALAR WMAIN ADD {SWIN}
SCALAR WMAIN MULT {FACTOR}
SCALAR WMAIN SET 0.0 SELE TYPE H* END
with SWIN 0.5 and FACTOR 0.927
Examples
This examples are meant to be a partial guide in setting up
an input file for PBEQ. There is one test file, tstpbeq.inp.
Example (1)
-----------
This example shows how to perform two PB calculations, one for a surrounding
dielectric of 80 (water) and one for a surrounding of 1.0 (vaccum). The
difference between the two energies then corresponds to the electrostatic
contribution to the solvation free energy. The salt concentration was zero
in this calculation.
PBEQ
scalar wmain = radius
SOLVE epsw 80.0 conc 0.0 ncel 30 dcel 0.4 maxi 200
set ener80 = ?ENPB
scalar wmain = radius
SOLVE epsw 1.0
set ener1 = ?ENPB
CALC total = @ener80 - @ener1
RESET
END
Example(2)
----------
This example shows how to use a set of atomic Born radii with a smoothing
window.
set sw 0.4
set factor 0.939
PBEQ
stream radius.str
scalar wmain add @sw
scalar wmain mult @factor
scalar wmain set 0.0 sele type H* end
scalar wmain show
SOLVE epsw 80.0 ncel 100 dcel 0.3 -
smooth swin @sw force sten 0.03 npbeq 1
RESET
END
Example(3)
----------
This example shows how to set up a membrane potential.
Note that a non-zero concentration is required for a sensible
system with a membrane potential.
PBEQ
scalar wmain = radius
scalar wmain show
SOLVE epsw 80.0 ncel 150 dcel 0.5 conc 0.150 -
Tmemb 25.0 Zmemb 0.0 epsm 2.0 vmemb 0.100
RESET
END
CHARMM Documentation / rvenable@deimos.cber.nih.gov