AMBER: mbondi radii for gaff atom types

From: Sara J <sara_compchem.yahoo.com>
Date: Wed, 28 May 2008 08:03:52 -0700 (PDT)

Hi,
I am trying to run a PB(and GB) calculation on a small molecule with GAFF atom types. In my tleap script, I have added the following lines to generate the prmtop file with GAFF atom types and to use modified Bondi radii:
 source leaprc.gaff
 set default PBradii mbondi2

My output from the prmtop file looks like:

%FLAG AMBER_ATOM_TYPE
%FORMAT(20a4)
c3 c3 c3 c3 c3 c3 c3 oh hc hc hc hc hc hc hc hc hc hc h1 h1
hc hc hc ho
%FLAG RADIUS_SET
%FORMAT(1a80)
H(N)-modified Bondi radii (mbondi2)
%FLAG RADII
%FORMAT(5E16.8)
  1.70000000E+00 1.70000000E+00 1.70000000E+00 1.70000000E+00 1.70000000E+00
  1.70000000E+00 1.70000000E+00 1.50000000E+00 1.20000000E+00 1.20000000E+00
  1.20000000E+00 1.20000000E+00 1.20000000E+00 1.20000000E+00 1.20000000E+00
  1.20000000E+00 1.20000000E+00 1.20000000E+00 1.20000000E+00 1.20000000E+00
  1.20000000E+00 1.20000000E+00 1.20000000E+00 1.20000000E+00


I am bit worried to see the radii for ho and hc atom type. When I look through the code of leap (/src/leap/src/leap/unitio.c ;
line 5212), I see the following:
 if( GDefaults.iGBparm == 1 || GDefaults.iGBparm == 2 ) {
                            if( sAtomType(aAtomA)[0] == 'c' ||
                                    sAtomType(aAtomA)[0] == 'c' ) dGBrad = 1.3;
                            if( sAtomType(aAtomA)[0] == 'o' ||
                                    sAtomType(aAtomA)[0] == 'o' ) dGBrad = 0.8;
                            if( sAtomType(aAtomA)[0] == 's' ||
                                    sAtomType(aAtomA)[0] == 's' ) dGBrad = 0.8;
                            if( (sAtomType(aAtomA)[0] == 'n' ||
                                    sAtomType(aAtomA)[0] == 'n') &&
                                    GDefaults.iGBparm == 2) dGBrad = 1.3;
(note that I have changed the atom types to lower case to be compatible with gaff type)
According to this code, the 'h' bonded to 'c' should have 1.3 and bonded to 'o' should have 0.8 radii, but I dont see that in my prmtop. I would really appreciate if someone can point me in the right direction. Am I generating the right prmtop file? If yes, why aren't the radii correct?

Thanks for the response.

Sara



      
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber.scripps.edu
To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
      to majordomo.scripps.edu
Received on Sun Jun 01 2008 - 06:07:11 PDT
Custom Search