1. Introduction

The CHARMM/APBS interface makes most of the APBS functionality available to CHARMM users. The interface was designed as an extension to the CHARMM’s PBEQ module. The CHARMM/APBS module is invoked by APBS keyword inside of the PBEQ section. The module’s keywords are identical or similar to APBS keywords to allow for easy transition between the codes.

The current version of the CHARMM/APBS module supports serial execution only. Parallel capability will be added in later versions.

2. Installing the iAPBS Interface

2.1. Requirements

To compile iAPBS/CHARMM interface two packages are required:

  • APBS (version 1.4 or later)

  • charmm (version c41b2 or later)

2.2. Building iAPBS interface

Building of the iAPBS interface requires compilation and installation of MALOC and APBS libraries. The following describes all steps (assumes bash as login script, modify appropriately for t/csh). The instructions also assume working gcc/gfortran compilers. The use of Intel compilers is recomended for better performance of the compiled executables.

# create a build directory and cd to it, then:
P=`pwd`
export APBS_PREFIX=${P}
git clone https://github.com/Electrostatics/apbs-pdb2pqr
cd apbs-pdb2pqr
git submodule init
git submodule update

cd $APBS_PREFIX/apbs-pdb2pqr/apbs/
mkdir build
cd build
cmake -DCMAKE_INSTALL_PREFIX:PATH=${APBS_PREFIX} \
 -DCMAKE_BUILD_TYPE=Release -DBUILD_DOC=OFF \
 -DBUILD_SHARED_LIBS=OFF -DENABLE_QUIET=ON \
 -DENABLE_iAPBS=ON ..
make install

These steps should build the apbs executable in ${APBS_PREFIX}/bin and all necessary apbs and iapbs libraries in ${APBS_PREFIX}/lib.

2.3. Building CHARMM/APBS module

Currently, building of CHARMM/APBS module is supported on x86_64 Linux platform only (gnu target in the CHARMM installation process). Both GNU and Intel compilers have been tested. The build process is a bit complex at this time, future releases will include more automated compilation.

cd $APBS_PREFIX
# untar charmm source here
cd $APBS_PREFIX/charmm

export APBS_LIBDIR=$APBS_PREFIX/lib
export APBS_LIB="-L${APBS_LIBDIR} -liapbs -lapbs_routines -lapbs_mg -lapbs_generic -lapbs_pmgc -lmaloc -lz"
export pdir=$APBS_PREFIX/apbs-pdb2pqr/apbs/contrib/iapbs/modules/CHARMM/patches
patch -p0 < $pdir/c41b2/apbs.patch
patch -p0 < $pdir/c41b2/Makefile_gnu.patch
patch -p0 < $pdir/c41b2/install.patch

./install.com gnu xxlarge X86_64 APBS keepf

# run tests
export CHARMM_EXE=$APBS_PREFIX/charmm/exec/gnu/charmm
cd $APBS_PREFIX/apbs-pdb2pqr/apbs/contrib/iapbs/modules/CHARMM/examples
$CHARMM_EXE < apbs_elstat.inp
./run.sh

3. Using CHARMM/APBS module

Keyword APBS inside of the PBEQ section in the CHARMM input file enters the CHARMM/APBS module. The module’s keywords are as similar to the original APBS keywords as possible. However, there are several exceptions. Please read the following section carefully and also take a look at the Examples section of this User’s Guide.

The following table lists all CHARMM/APBS keywords with their description. The left side of the table also lists corresponding APBS keywords which CHARMM/APBS keywords mimic very closely. For detailed discussion of APBS keywords please see APBS documentation.

3.1. CHARMM/APBS Module Keywords

Table 1. APBS and CHARMM/APBS input parameters
APBS keyword CHARMM/APBS keyword Description

APBS

Enters the APBS module

mg-auto

MGAUTO

Automatically-configured sequential focusing multigrid calculation

mg-para

MGPARA

Automatically-configured parallel focusing multigrid calculation

mg-manual

MGMANUAL

Manually-configured multigrid calculation

lpbe/npbe

LPBE/NPBE

Linear/full Poisson-Boltzmann equation

bcfl

BCFL [1]

Boundary condition method: 0: zero; 1: sdh; 2: mdh; 4: focus

srfm

SRFM [1]

Surface calculation method: 0: mol; 1: smol; 2: spl2

pdie

PDIE [2.0]

Solute dielectric

sdie

SDIE [78.54]

Solvent dielectric

sdens

SDENS [10.0]

Vacc sphere density

srad

SRAD [1.4]

Solvent radius

swin

SWIN [0.3]

Cubic spline window

temp

TEMP [298.15]

Temperature (in K)

gamma

GAMMA [0.105]

Surface tension for apolar energies/forces (in kJ/mol/A2)

chgm

CHGM [1]

Charge discretization method: 0: spl0; 1: spl2

calcenergy

CALCE [1]

Energy calculation flag: 0: Do not perform energy calculation; 1: Calculate total energy only; 2: Calculate per-atom energy components

calcforce

CALCF [0]

Atomic forces calculation: 0: Do not perform force calculation; 1: Calculate total force only; 2: Calculate per-atom force components

write pot

WPOT

Writes electrostatic potential data to iapbs-pot.dx in DX format

write charge

WCHG

Writes charge data to iapbs-charge.dx in DX format

write smol

WSMOL

Writes molecular surface data to iapbs-smol.dx in DX format

write kappa

WKAPPA

Writes the ion-accessibility kappa map to iapbs-kappa.dx in DX format

write diel

WDIEL

Writes dielectric maps to iapbs-diel[x,y,z].dx in DX format

read charge

RCHG

Reads charge data from iapbs-charge.dx in DX format

read kappa

RKAPPA

Reads the ion-accessibility kappa map from iapbs-kappa.dx in DX format

read diel

RDIEL

Reads dielectric maps from iapbs-diel[x,y,z].dx in DX format

ionq

IONQ{1|2}

Counterion charges (in e)

ionc

IONC{1|2}

Counterion concentrations (in M)

ionr

IONR{1|2}

Counterion radii (in A)

dime

DIMX [65] DIMY [65] DIMZ [65]

Grid dimensions (in x, y and z)

cmeth

CMET [1]

Centering method: 0: Center on a point 1: Center on a molecule

center

CNTX CNTY CNTZ

Grid center if CMET 0

ccmeth

CCME [1]

Coarse grid centering method: 0: Center on a point 1: Center on a molecule

ccenter

CCNX CCNY CCNZ

Coarse grid center if CCME 0

fcmeth

FCME [1]

Fine grid centering method: 0: Center on a point 1: Center on a molecule

fcenter

FCNX FCNY FCNZ

Fine grid center if FCME 0

grid

GRDX GRDY GRDZ

Grid spacings

glen

GLNX GLNY GLNZ

Grid side lengths

cglen

CGLX CGLY CGLZ

Coarse grid side lengths

fglen

FGLX FGLY FGLZ

Fine grid side lengths

pdime

PDIX PDIY PDIZ

Grid of processors to be used in calculation

ofrac

OFRA [0.1]

Overlap fraction between processors

DEBUG [0]

Debuging flag

UPDATE [1]

How often are solvation forces calculated during minimization or MD. Defaults to 1 which means solvation forces are updated every step.

UMETHOD [1]

Method by which are solvation forces updated, when using UPDATE keyword. 0: no solvation forces are included between updates 1: forces from previous APBS calculation are used between updates

SFORCE

Requests calculation of solvation forces (this keyword must be present when doing MD or minimization with APBS calculated solvation forces).

Note
Notes
  • Values in square brackets [] are defaults.

  • To disable creation of io.mc file before attempting any extensive minimization or MD simulation with the APBS module please set the MCSH_HOME environment variable to /dev/null (export MCSH_HOME=/dev/null).

  • Each occurrence of the APBS keyword triggers APBS calculation. This not may be desirable in all cases, for example when setting up a molecular dynamics simulation with APBS calculated solvation forces. In such a case keyword SKIP in the APBS section will prevent from starting APBS calculation, the section will just initialize all APBS parameters. The actual APBS calculation will be performed when DYNAMICS section starts the simulation.

  • When the SFORCE keyword is specified the APBS calculation is initialized (all parameters are set) with the APBS keyword inside of the PBEQ section but the actual electrostatic calculation starts during minimization or MD step.

  • When performing solvation forces calculation, for example during minimization or molecular dynamics, the keywords SRFM and CHGM must be both set to 2. Also, no counterions must be present.

4. Examples of using the CHARMM/APBS module

The following examples demonstrate how APBS functionality can be used from within CHARMM. Using CHARMM/APBS module is very similar to using the PBEQ module so examples provided for the PBEQ module can be easily adapted for CHARMM/APBS. The example files are located in CHARMM/examples in iAPBS distribution directory.

4.1. Solvation energy calculation

This example executes a solvation energy calculation and prints out total solvation energy for the given molecule.

* sp_elstat.inp
* external files: top_all22_prot.inp, par_all22_prot.inp and radius.str
*

!if ?apbs .ne. 1 then stop

stream datadir.def

open read card unit 11 name @0top_all22_prot.inp
read rtf card unit 11
close unit 11

open read card unit 11 name @0par_all22_prot.inp
read para card unit 11
close unit 11

 ... read in or generate structure

coor orient

scalar charge show

PBEQ
   stream @0radius.str
   set factor 0.939
   set sw     0.4
   scalar wmain add  @sw
   scalar wmain mult @factor
   scalar wmain set 0.0 sele type H* end
   scalar wmain show

  APBS mgauto lpbe dimx 65 dimy 65 dimz 65 -
  cglx 30 cgly 30 cglz 30 fglx 15 fgly 15 fglz 15 -
  srfm 2 debug 1 -
  calcene 1 calcforce 1 sforce -
  sele all END

END

skip all excl pbelec pbnp
ener
set elstatenergy = ?ENPB
set apolarEN = ?ENNP
stop

4.2. Molecular dynamics in implicit solvent using APBS

This examples shows how to use the CHARMM/APBS module for performing molecular dynamics simulation in implicit solvent with APBS.

* md.inp
* external files: top_all22_prot.inp, par_all22_prot.inp and radius.str
*

if ?apbs .ne. 1 then stop

stream datadir.def

open read card unit 11 name @0top_all22_prot.inp
read rtf card unit 11
close unit 11

open read card unit 11 name @0par_all22_prot.inp
read para card unit 11
close unit 11

 ... read in or generate structure

coor orient

PBEQ
   set factor 0.939
   set sw     0.4
   stream @0radius.str
   scalar wmain add  @sw
   scalar wmain mult @factor
   scalar wmain set 0.0 sele type H* end

  APBS mgauto lpbe -
  grdx 0.5 grdy 0.5 grdz 0.5 -
  swin @sw srfm 2 -
  calcene 2 calcfor 2 debug 0 -
  sforce -
  sele all END
END

skip none

dynamics leap verlet strt nstep 150 timestep 0.001 -
firstt 100.0 finalt 300.0 teminc 100.0 -
twindh 10.0

stop

4.3. Calculation and visualization of electrostatic potential

When this input file is run three DX files will be created. These can be visualized using several applications (for details please see APBS visualization guide). The DX files will be iapbs-pot.dx, iapbs-smol.dx and iapbs-charge.dx which will contain electrostatic potential, solvent accessible surface and charge information, respectively.

* visualization.inp
* external files: top_all22_prot.inp, par_all22_prot.inp and radius.str
*

!if ?apbs .ne. 1 then stop

stream datadir.def

open read card unit 11 name @0top_all22_prot.inp
read rtf card unit 11
close unit 11

open read card unit 11 name @0par_all22_prot.inp
read para card unit 11
close unit 11

... read in or generate structure

coor orient

PBEQ
   stream @0radius.str
   set factor 0.939
   set sw     0.4
   scalar wmain add  @sw
   scalar wmain mult @factor
   scalar wmain set 0.0 sele type H* end
   scalar wmain show

  APBS mgauto lpbe dimx 65 dimy 65 dimz 65 -
  cglx 30 cgly 30 cglz 30 fglx 15 fgly 15 fglz 15 -
  calcene 1 calcforce 0 -
  ionq1 1.0 ionc1 0.15 ionr1 2.0 ionq2 -1.0 ionc2 0.15 ionr2 2.0 -
  wpot wsmol wchg -
  debug 1 -
  sele all END
END

set elstatenergy = ?ENPB
Note
Notes

To determine the correct size of the numerical grid enveloping the studied molecule (parameters cglen and fglen) you can use psize.py tool from the APBS distribution.