Delphi

Delphi is a program to solve the Poisson-Boltzmann equation to obtain potentials and solvation free energies for heterogenous systems.

The home page of Delphi is at

http://wiki.c2b2.columbia.edu/honiglab_public/index.php/Software:DelPhi

A manual of the program (pdf) is found in http:/honiglab.cpmc.columbia.edu/delphi/doc/delphi_manual.pdf.
A similar but better manual in web format is found in http://honiglab.cpmc.columbia.edu/delphi/doc/io.html.

On milleotto:
module add pgi32/7.1
delphi

It is not available on the local system.

Running delphi

The program is run with the command:
delphi inputfile

if inputfile is missing, the name fort.10 is assumed.

This is a sample input file:
exdi=80.0
scale=2.5
indi=1.0
perfil=90.0
prbrad=1.4
bndcon=4
in(crg,file="delphi.crg")
in(siz,file="delphi.siz")
in(pdb)
energy(s,g)

Three files (at least) are read in:
A pdb file (default fort.13)
An atom radii file (siz, default name fort.11)
An atom charge file (crg, derault name fort.12)

For small molecules, scale=2.0 may give very poorly converged solvation energies with errors of 29 kJ/mol for water and 54 kJ/mol for imidazoleH+. Scale=10 gives accurate results. For drug-like molecules, the error is up to 7 kJ/mol.

Moreira et al, J. Mol. Struct 729 (2005) 11 have examined the dependence of Delphi results on some parameters. They recommend the following settlings:
conv = 10^-3
perfil = 90%
scale = 2.5
coloumbic boundary conditions (bndcon=4).

Parameters

exdi dielectric constant outside the protein (default = 80)
indi dielectric constant inside the protein (default = 2.0 !)
scale - lattice spacing in number of grids/Angstrom (default 1.2)
perfil - percentage of the lattice that the largest of the x,y or z linear dimensions of the molecule will fill (default = 80 %).
prbrad radius (in Å) of the water probe molecule (default = 1.4)
bndcon type of boundary condition:
1 - potential is zero
2 - Debue-Huckel potential of the equivalent dipole (default)
3 - focussing from a previous phi-map calculation
4 - Debue-Huckel potential of all charges
salt - salt concentration (default 0)

Functions.

center(0.2,3,2)
Will offset the molecule by 0.2 grids in the x direction, 3 in the y and 2 in the z,
You can also let center open a pdb file and center on the centre of geometry of that file: center(file="file").
acenter(1.0,5.6,7.0)
Centres the molecule at x=1.0 Å, y=5.6 Å, z=7.0 Å
energy(s,g)
Writes out a certain energy, which can be:
G or GRID for the grid energy,
S or SOL or SOLVATION for the corrected reaction field energy
C or COULOMBIC or COU for the coulombic energy
AS or ANASURF or ANALYTICALSURFACE for the analytical surface energy
AG or ANAGRID or ANALYTICALGRID for the analytical grid energy
read or in
Reads an input file. It comes with several specifiers, namely
SIZ: for radius files
CRG: for charge file
PDB: for the pdb file
FRC: for the file use to determine site potentials
PHI: for the phimap used in focussing
write or out
Writes an output file. The specifiers are:
PHI : for phimaps
FRC : for site potentials
EPS : for epsmaps
MODPDB: for modified pdb files
UNPDB: for unformatted pdb file
UNFRC : for unformatted frc files
For example, out(modpdb, file="test.out") writes a modified pdb file called "test.out".

The atom size file

The file can start with any number of comment lines, starting with !
Then a line defining the file format should come. It should be either (exactly for the 16 first characters:)
atom__res_radius
which gives the format
atom(a6), residue(a3), radius (f8.3)
aaaaaarrrffff.fff

or

atom__resnumbc_radius_
which gives the format
atom(a6), residue(a3), residue_number(a4), subunit(a1), radius (f8.3)
aaaaaarrrnnnncffff.fff

Then, the rest of the file is read in this fixed format.
Blank fields are considered as wild-card characters.
No warnings are given for errors.

The atom charge file

The file can start with any number of comment lines, starting with !
Then a line is read but ignored.
The rest of the file is read in the following fixed format (same as second atom size format)

atom(a6), residue(a3), residue_number(a4), subunit(a1), charge(f8.3)
aaaaaarrrnnnncffff.fff

Blank fields are considered as wild-card characters.

Note that only three decimals are read for the charges (but we have changed that to six decimals)!

Warning!

Note that Delphi does not interpret H-atom names in the format 1HD2 LEU correctly (i.e. as HD21), instead it is interpreted litteraly as "1HD2" and if the same name is not given in the size and charges files, a "default" charge and radii are used for those atoms, without any crash!

Use changepdb format 0.

Warning 2!

Strangely enough the program does not allow radii of 0 (as in Parse for non-polar hydrogens). If they are set to 0.1 (or even 0.001), it works. In our test case, there is no difference with an radius of 0.1 and 0.001, or even with 1.0.

charged atom number 2radius changed from zero to 1.000000

Problems

We have not been able to run focusing calculations for big system.
The problem has been traced to the handling of phimap.

Example of a focusing file (that does not work):
exdi=80.0
scale=2.5
indi=1.0
perfil=90.0
prbrad=1.4
BNDCON=3
in(crg,file="delphi.crg")
in(siz,file="delphi.siz")
in(pdb)
in(phi,unit=14)
energy(s,g)

Changes made to the source files

1. We decided to let the program read in specific radii and charges for all atoms. Therefore, we had to increase five parameters in the file qdiffpar5.h:

        parameter (nclist = 50000)
        parameter (nrlist = 50000)
c
c nrmax= maximum number of entries in radius file
        parameter (nrmax = 50000)
c ncmax= maximum number of entries in charge file
        parameter (ncmax = 50000)

2. Second, we forced the program to read 6 decimals for the charges (needed to get the correct charge on the residues with standard Amber charges):
In file rdhcrg.f, we changed the following line (from F8.3):
202     format(A6,A3,A4,A1,F11.6)

3. To read more than 10000 atoms, I changed (19/1-10) in getatm2.f:
        parameter(natmax=100000)

Sample output file

__________DelPhi V. 4 Release 1.1 ______________
|                                                |
| A program to solve the PB equation             |
| in 3D, using non-linear form, incorporating    |
| many dielectric regions, multisalt ionic       |
| strength, different probe radii, periodic      |
| and focussing boundary conditions, utilizing   |
| stripped optimum successive over-relaxation    |
| and an improved algorithm for mapping the      |
| Mol. Surface to the finite-Difference grid     |
| Recompiled on Linux and PC                     |
|    January 2002 -------- Walter Rocchia       |
|__________________           ___________________|
                    DelPhi V. 4

program started on Mon Jan 18 2010
              at 15:19:06
opening parameter file delphi.in
atom radii read from file
delphi.siz

! Atomic radii for delphi automatically generated by WtDelph
! PARSE radii
reading pK style radius file
# of radius parameter records:           3

atomic charges read from file
delphi.crg

! Atomic charges for delphi automatically generated by WtDel
! Amber charges from prmtop file
# of charge parameter records:           3
assigning charges and radii...

opening formatted file:delphi.pdb
You are not reading from an objectfile! assuming     having only molecules, and one dielectric medium
number of atoms read in =            3

atomic coordinates, charges and radii written to file
delphi.mpb

Direct mapping of epsilon: (0/1)(n/y)           1
time to read in and/or assign rad/chrg= 0.0000000E+00
object number :              1, is a molecule

grid size                  :           3
percent of box to be filled:   90.00000
scale,in grids/A :   1.000000
xmin,xmax     (A): -9.131000      -5.777000
ymin,ymax     (A): -15.74500      -12.92700
zmin,zmax     (A):   15.81600       19.01300
x,y,z range   (A):    3.354000       2.818000       3.197001
system geometric center in (A): -7.454000      -14.33600       17.41450
grid box is centered in (A) : -7.454000      -14.33600       17.41450
object centred at (gu)     : 0.0000000E+00 0.0000000E+00 0.0000000E+00
outer dielectric           :   80.00000
dielectric in medium number           1 :   1.000000
first kind salt concentration (M)   : 0.0000000E+00
valences salt 1 are                   1 and           1
second kind salt concentration (M) : 0.0000000E+00
valences salt 2 are                   0 and           0
ionic strength (M)         : 0.0000000E+00
debye length (A)           :   1000000.
absolute temperature (K)   :   297.3342
ion exclusion radius (A)   :   2.000000
probe radius facing water(A:   1.400000
probe radius, internal (A) :   1.400000
boundary conditions        : dipolar
x,y,z periodic bc. and volt. drop flags: F F F F F F
# of linear iterations     : automatic
# of non-linear iterations :           0
non-linear energy calculat.: F
manual relaxation parameter : F
ionic direct energy contribution: F
concentration map output   : F
spherical charge distbn.   : F
INSIGHT format output      : F
site potential output      : F
modified atom file output : T
map file label             :
qdiffxas: qdiffxs4 with an improved surfacing routine
convergence graph turned off
potential listings turned off

start vw surface at   0.0000000E+00
Starting creating Vand der Waals Epsilon Map
Ending creating Vand der Waals Epsilon Map
fill in re-entrant regions at   0.0000000E+00
boundary points facing continuum solvent=            0
total number of boundary points before elab.=            0
radpmax,radprb    1.400000       1.400000       1.400000
# of vertices =          384 # of edges =          372
# of pairs =            3
time to find all pairs =   0.0000000E+00
# pairs analyzed (atom-atom and atom-object)=            3
# exposed pairs (atom-atom and atom-object)=            3
no. arc points =          357
no. surface atoms =            3 nbur =            0
mkacc time =   0.0000000E+00
writing accessible surface arcs data to ARCDAT
grid for indexing accessible points =    1.400000
bgp added m=           0 bgp removed mr =           0
time to grow re-entrant surface =   0.0000000E+00
no. cavity mid-points inaccessible to solvent =            0
after surface elaboration ibnum=            0
     and               ibnumsurf=            0
scaling boundary grid points
time taken =   0.0000000E+00
            0 points had to be assigned by global comparison
time to turn everything in is 0.0000000E+00
number of charges coming from molecules            3

number of atom coordinates read :            3
total number of assigned charges    :            3
net assigned charge              :   0.0000000E+00
assigned positive charge         :   0.8340000
centred at (gu) :   2.000000       1.738501       2.294499
assigned negative charge         : -0.8340000
centred at (gu) :   1.836999       2.009001       1.801498

number of dielectric boundary points           0
no. dielectric boundary points in salt =            0
no. grid points charged and not internal=           0
iepsmp to db, and charging done at 0.0000000E+00
number of grid points assigned charge           1

   setting boundary conditions

some initial phi values:
midg,midg,1; midg,midg,igrid
   -2.721608       2.937844
midg,1,midg; midg,igrid,midg
    1.738883      -1.198084
1,midg,midg; igrid,midg,midg
   -1.356813      0.4846959

gauss-seidel spectral radius is        0.0000000000

estimated iterations to convergence          11

setup time was (sec)   0.0000000E+00

now iterating at: 15:19:07

   rms-change     max change       #iterations
   0.0000000E+00 0.0000000E+00 at          10 iterations
finished qdiffx linear iterations
at                       : 15:19:07
total time elapsed so far:   0.0000000E+00
# loops                  :           10
mean,max change (kT/e)   :   0.0000000E+00 0.0000000E+00

total grid energy          :        107.4648     kt
self-reaction field energy :      0.0000000E+00 kt
total s.charge,no epsin carrying :    0.0000
corrected reaction field energy:   0.0000000E+00 kt
total reaction field energy :     0.0000000E+00 kt

All required energy terms but grid and self_react.: 0.0000000E+00kt
energy calculations done at 0.0000000E+00

total cpu time was (sec)   0.0000000E+00

DelPhi exited at 15:19:07

Another ending:

total grid energy          :        23944.70     kt
self-reaction field energy :      0.0000000E+00 kt
total s.charge,no epsin carrying :   16.0076
corrected reaction field energy:   -7466.163     kt
total reaction field energy :     -7466.163     kt

All required energy terms but grid and self_react.: -7466.163    kt

The interesting result (the solvation energy) is given in the line:
corrected reaction field energy:   -7466.163     kt
The conversion from kt (RT) to kJ/mol is 2.49434 kJ/mol at 300 K.

Compilation

On milleotto (Samuel 2/9-09):

module add pgi32/7.1
(There are other PGI 32-bit compilors, but this was the only one that worked)

In the makefile for Delphi:

FC=pgf77
CC=pgcc

Then type make

Old stuff

The program is compiled with the Portland Group pgf77 and pgcc compilers.
It runs on docenten.

On docenten:
Files are in /sw/pkg/bio/DELPHI/DELPHI_2004_LINUX
module add pgi32/5.2
delphi

On toto7/wheinim64
Files are in /sw/bio/DELPHI/DELPHI_2004_LINUX.
module add pgi/5.2
delphi

Older (for docenten; Jacob, 30/1-08):
Change in makefile
FC=pgf77
CC=pgcc

Then
make delphi

----------

On docenten in

/home/kongsted/delphitest/DELPHI/DELPHI
 
You will find a correctly compiled version (intel) and the makefile. 

Remember to 

module add intel32

before typing 

make

I like the intel compiler which is the main reason for chosing this. It can be shiftet to f.ex. pgf by using the old makefile and 

module add pgi

and chosing a 32 bit compiler as well in this case.