Updated source for Van der Waals radii.

Implemented new Van der Waals radii in atomprop. This may
have an effect on results from gmx sas or editconf and therefore
a message is printed to notify the users of this.
The genbox program is adapted to have an extra option to scale
the Van der Waals radii, in case it is difficult to add extra
atoms. The default value for this option is 0.57 which gives approximately
the same density as old genbox.

Change-Id: Ib5d86a920ae7732ee8a29505449e6d34e3d3980c

diff --git a/share/top/vdwradii.dat b/share/top/vdwradii.dat
index d8c1ddbd327dfc37458f3e16cb9fa5b88fba01ec..1aef7329ca1cf088c7bd874726b26aaa40a7e2bc 100644 (file)
--- a/share/top/vdwradii.dat
+++ b/share/top/vdwradii.dat
@@ -3,12 +3,17 @@
 ; longest matches are used
 ; '???' or '*' matches any residue name
 ; 'AAA' matches any protein residue name
-???  C     0.15
-???  F     0.12
-???  H     0.04
-???  N     0.110
-???  O     0.105
-???  S     0.16
+; Source: http://en.wikipedia.org/wiki/Van_der_Waals_radius
+; These come from A. Bondi, "van der Waals Volumes and Radii",
+; J. Phys. Chem. 68 (1964) 441-451
+???  H     0.12
+???  C     0.17
+???  N     0.155
+???  O     0.152
+???  F     0.147
+???  P     0.18
+???  S     0.18
+???  Cl    0.175
 ; Water charge sites
 SOL  MW    0
 SOL  LP    0

diff --git a/src/gromacs/gmxana/gmx_genbox.cpp b/src/gromacs/gmxana/gmx_genbox.cpp
index 2dc75e4f49fee386a9b6a5567d7bb6df66419a9e..2f7e68518b75c8cc2c4505522af15da9b7f7f9ab 100644 (file)
--- a/src/gromacs/gmxana/gmx_genbox.cpp
+++ b/src/gromacs/gmxana/gmx_genbox.cpp
@@ -103,7 +103,7 @@ static gmx_bool in_box(t_pbc *pbc, rvec x)
 }
 
 static void mk_vdw(t_atoms *a, real rvdw[], gmx_atomprop_t aps,
-                   real r_distance)
+                   real r_distance,real r_scale)
 {
     int i;
 
@@ -117,6 +117,10 @@ static void mk_vdw(t_atoms *a, real rvdw[], gmx_atomprop_t aps,
         {
             rvdw[i] = r_distance;
         }
+        else 
+        {
+            rvdw[i] *= r_scale;
+        }
     }
 }
 
@@ -379,7 +383,8 @@ enum {
 
 static char *insert_mols(const char *mol_insrt, int nmol_insrt, int ntry, int seed,
                          t_atoms *atoms, rvec **x, real **r, int ePBC, matrix box,
-                         gmx_atomprop_t aps, real r_distance, real rshell,
+                         gmx_atomprop_t aps, 
+                         real r_distance, real r_scale, real rshell,
                          const output_env_t oenv,
                          const char* posfn, const rvec deltaR, int enum_rot,
                          gmx_bool bCheckAllPairDist)
@@ -424,7 +429,7 @@ static char *insert_mols(const char *mol_insrt, int nmol_insrt, int ntry, int se
     srenew(atoms_insrt.resinfo, atoms_insrt.nres);
 
     /* initialise van der waals arrays of insert molecules */
-    mk_vdw(&atoms_insrt, r_insrt, aps, r_distance);
+    mk_vdw(&atoms_insrt, r_insrt, aps, r_distance, r_scale);
 
     /* With -ip, take nmol_insrt from file posfn */
     if (posfn != NULL)
@@ -532,7 +537,8 @@ static char *insert_mols(const char *mol_insrt, int nmol_insrt, int ntry, int se
 
 static void add_solv(const char *fn, t_atoms *atoms, rvec **x, rvec **v, real **r,
                      int ePBC, matrix box,
-                     gmx_atomprop_t aps, real r_distance, int *atoms_added,
+                     gmx_atomprop_t aps, 
+                     real r_distance, real r_scale, int *atoms_added,
                      int *residues_added, real rshell, int max_sol,
                      const output_env_t oenv)
 {
@@ -580,7 +586,7 @@ static void add_solv(const char *fn, t_atoms *atoms, rvec **x, rvec **v, real **
     rm_res_pbc(atoms_solvt, x_solvt, box_solvt);
 
     /* initialise van der waals arrays of solvent configuration */
-    mk_vdw(atoms_solvt, r_solvt, aps, r_distance);
+    mk_vdw(atoms_solvt, r_solvt, aps, r_distance, r_scale);
 
     /* calculate the box multiplication factors n_box[0...DIM] */
     nmol = 1;
@@ -637,7 +643,7 @@ static void add_solv(const char *fn, t_atoms *atoms, rvec **x, rvec **v, real **
 
 static char *read_prot(const char *confin, t_atoms *atoms, rvec **x, rvec **v,
                        real **r, int *ePBC, matrix box, gmx_atomprop_t aps,
-                       real r_distance)
+                       real r_distance, real r_scale)
 {
     char *title;
     int   natoms;
@@ -661,7 +667,7 @@ static char *read_prot(const char *confin, t_atoms *atoms, rvec **x, rvec **v,
             title, atoms->nr, atoms->nres);
 
     /* initialise van der waals arrays of configuration 1 */
-    mk_vdw(atoms, *r, aps, r_distance);
+    mk_vdw(atoms, *r, aps, r_distance, r_scale);
 
     return title;
 }
@@ -919,7 +925,7 @@ int gmx_genbox(int argc, char *argv[])
 #define NFILE asize(fnm)
 
     static int      nmol_ins   = 0, nmol_try = 10, seed = 1997, enum_rot;
-    static real     r_distance = 0.105, r_shell = 0;
+    static real     r_distance = 0.105, r_shell = 0, r_scale = 0.57;
     static rvec     new_box    = {0.0, 0.0, 0.0}, deltaR = {0.0, 0.0, 0.0};
     static gmx_bool bReadV     = FALSE, bCheckAllPairDist = FALSE;
     static int      max_sol    = 0;
@@ -936,6 +942,8 @@ int gmx_genbox(int argc, char *argv[])
           "Random generator seed"},
         { "-vdwd",   FALSE, etREAL, {&r_distance},
           "Default van der Waals distance"},
+        { "-vdwscale", FALSE, etREAL, {&r_scale},
+          "HIDDENScale factor to multiply Van der Waals radii from the database in share/gromacs/top/vdwradii.dat. The default value of 0.57 yields density close to 1000 g/l for proteins in water." },
         { "-shell",  FALSE, etREAL, {&r_shell},
           "Thickness of optional water layer around solute" },
         { "-maxsol", FALSE, etINT,  {&max_sol},
@@ -978,7 +986,7 @@ int gmx_genbox(int argc, char *argv[])
         /*generate a solute configuration */
         conf_prot = opt2fn("-cp", NFILE, fnm);
         title     = read_prot(conf_prot, &atoms, &x, bReadV ? &v : NULL, &r, &ePBC, box,
-                              aps, r_distance);
+                              aps, r_distance, r_scale);
         if (bReadV && !v)
         {
             fprintf(stderr, "Note: no velocities found\n");
@@ -1019,7 +1027,8 @@ int gmx_genbox(int argc, char *argv[])
     if (bInsert)
     {
         title_ins = insert_mols(opt2fn("-ci", NFILE, fnm), nmol_ins, nmol_try, seed,
-                                &atoms, &x, &r, ePBC, box, aps, r_distance, r_shell,
+                                &atoms, &x, &r, ePBC, box, aps, 
+                                r_distance, r_scale, r_shell,
                                 oenv, opt2fn_null("-ip", NFILE, fnm), deltaR, enum_rot,
                                 bCheckAllPairDist);
     }
@@ -1032,7 +1041,7 @@ int gmx_genbox(int argc, char *argv[])
     if (bSol)
     {
         add_solv(opt2fn("-cs", NFILE, fnm), &atoms, &x, v ? &v : NULL, &r, ePBC, box,
-                 aps, r_distance, &atoms_added, &residues_added, r_shell, max_sol,
+                 aps, r_distance, r_scale, &atoms_added, &residues_added, r_shell, max_sol,
                  oenv);
     }
 

diff --git a/src/gromacs/gmxlib/atomprop.c b/src/gromacs/gmxlib/atomprop.c
index 9d8aa60f85042b00ef68fa09947c6952a112d939..ed2ffbe666e6bc56565a063407535a106347499b 100644 (file)
--- a/src/gromacs/gmxlib/atomprop.c
+++ b/src/gromacs/gmxlib/atomprop.c
@@ -49,6 +49,7 @@
 #include "index.h"
 #include "strdb.h"
 #include "copyrite.h"
+#include "statutil.h"
 
 typedef struct {
     gmx_bool    bSet;
@@ -62,7 +63,7 @@ typedef struct {
 } aprop_t;
 
 typedef struct gmx_atomprop {
-    gmx_bool          bWarned;
+    gmx_bool          bWarned,bWarnVDW;
     aprop_t           prop[epropNR];
     gmx_residuetype_t restype;
 } gmx_atomprop;
@@ -284,7 +285,8 @@ gmx_atomprop_t gmx_atomprop_init(void)
 
     gmx_residuetype_init(&aps->restype);
     aps->bWarned = FALSE;
-
+    aps->bWarnVDW = FALSE;
+    
     return (gmx_atomprop_t)aps;
 }
 
@@ -329,6 +331,18 @@ void gmx_atomprop_destroy(gmx_atomprop_t aps)
     sfree(ap);
 }
 
+static void vdw_warning(FILE *fp)
+{
+    if (NULL != fp)
+    {
+        fprintf(fp,"NOTE: From version 5.0 %s uses the Van der Waals radii\n",
+                ShortProgram());
+        fprintf(fp,"from the source below. This means the results may be different\n");
+        fprintf(fp,"compared to previous GROMACS versions.\n");
+        please_cite(fp,"Bondi1964a");
+    }
+}
+
 gmx_bool gmx_atomprop_query(gmx_atomprop_t aps,
                             int eprop, const char *resnm, const char *atomnm,
                             real *value)
@@ -368,6 +382,11 @@ gmx_bool gmx_atomprop_query(gmx_atomprop_t aps,
     j = get_prop_index(&(ap->prop[eprop]), ap->restype, resname,
                        atomname, &bExact);
 
+    if (!ap->bWarnVDW)
+    {
+        vdw_warning(stdout);
+        ap->bWarnVDW = TRUE;
+    }
     if (j >= 0)
     {
         *value = ap->prop[eprop].value[j];

diff --git a/src/gromacs/gmxlib/copyrite.c b/src/gromacs/gmxlib/copyrite.c
index 1a9a61bb7cd2497971e8704f24229a2c1b164543..4915f238536998063d08d705ba720999aa06617f 100644 (file)
--- a/src/gromacs/gmxlib/copyrite.c
+++ b/src/gromacs/gmxlib/copyrite.c
@@ -485,6 +485,11 @@ void please_cite(FILE *fp, const char *key)
           "Solvation energy in protein folding and binding",
           "Nature",
           319, 1986, "199-203" },
+        { "Bondi1964a",
+          "A. Bondi",
+          "van der Waals Volumes and Radii",
+          "J. Phys. Chem.",
+          68, 1964,"441-451" },
         { "Eisenhaber95",
           "Frank Eisenhaber and Philip Lijnzaad and Patrick Argos and Chris Sander and Michael Scharf",
           "The Double Cube Lattice Method: Efficient Approaches to Numerical Integration of Surface Area and Volume and to Dot Surface Contouring of Molecular Assemblies",