********>Bugfix 24:
Author: Dave Case
Date: 10/01/02

Programs: LEaP

Severity: major for certain combinations of parameters

Description: The assignment of GB radii to atoms does not recognize that
             GAFF atom types are lower-case; hence it gives incorrect results
             when the gaff force field is used.
             (This patch also adds an "mbondi2" option which is intended
             for use in conjunction with igb=2.)

Fix:  apply the following patch to amber7/src/leap/src/leap/unitio.c:

------------------------------------------------------------------------------
*** unitio.c	2002/01/30 19:12:27	7.28
--- unitio.c	2002/10/01 23:56:57
***************
*** 5154,5160 ****
                      &dMass, &dPolar, &dEpsilon, &dRStar, &dEpsilon14,
                      &dRStar14, &iElement, &iHybridization, sDesc);
  
! 		if( GDefaults.iGBparm < 3 ) {   /* Bondi or modified Bondi radii */
  			switch( iElement ){
  				case 1: dGBrad = 1.2; 
  
--- 5154,5161 ----
                      &dMass, &dPolar, &dEpsilon, &dRStar, &dEpsilon14,
                      &dRStar14, &iElement, &iHybridization, sDesc);
  
! 		if( GDefaults.iGBparm < 3  || GDefaults.iGBparm == 6 ) {
! 			/* Bondi or modified Bondi radii */
  			switch( iElement ){
  				case 1: dGBrad = 1.2; 
  
***************
*** 5164,5177 ****
  					   iGBparm=2 adds the modifcation from Biopolymers
  					   56: 275 (2001)    */
  
! 					if( GDefaults.iGBparm > 0 ){
  						aAtomA = aAtomBondedNeighbor(saPAtom->aAtom, 0);
! 						if( sAtomType(aAtomA)[0] == 'C' ) dGBrad = 1.3;
! 						if( sAtomType(aAtomA)[0] == 'O' ) dGBrad = 0.8;
! 						if( sAtomType(aAtomA)[0] == 'S' ) dGBrad = 0.8;
! 						if( sAtomType(aAtomA)[0] == 'N' &&
  							GDefaults.iGBparm == 2) dGBrad = 1.3;
  					}
  					break;
  				case 6: dGBrad = 1.7; break;
  				case 7: dGBrad = 1.55; break;
--- 5165,5189 ----
  					   iGBparm=2 adds the modifcation from Biopolymers
  					   56: 275 (2001)    */
  
! 					if( GDefaults.iGBparm == 1 || GDefaults.iGBparm == 2 ){
  						aAtomA = aAtomBondedNeighbor(saPAtom->aAtom, 0);
! 						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;
  					}
+ 
+ 					if( GDefaults.iGBparm == 6 ){  /* try Alexey's scheme */
+ 						aAtomA = aAtomBondedNeighbor(saPAtom->aAtom, 0);
+ 						if( sAtomType(aAtomA)[0] == 'N' ||
+ 						     sAtomType(aAtomA)[0] == 'n' ) dGBrad = 1.3;
+ 					}
+ 
  					break;
  				case 6: dGBrad = 1.7; break;
  				case 7: dGBrad = 1.55; break;
***************
*** 5276,5282 ****
          ParmSetAtom(uUnit->psParameters, iIndex, sType,
                      &dMass, &dPolar, &dEpsilon, &dRStar, &dEpsilon14,
                      &dRStar14, &iElement, &iHybridization, sDesc);
!                 if( GDefaults.iGBparm < 4){
    		    /* for now, hardwire the Bondi radii  */
  		    switch( iElement ){
  			    case 1:  dGBscreen = 0.85; break;
--- 5288,5294 ----
          ParmSetAtom(uUnit->psParameters, iIndex, sType,
                      &dMass, &dPolar, &dEpsilon, &dRStar, &dEpsilon14,
                      &dRStar14, &iElement, &iHybridization, sDesc);
!                 if( GDefaults.iGBparm < 4 || GDefaults.iGBparm == 6 ){
    		    /* for now, hardwire the Bondi radii  */
  		    switch( iElement ){
  			    case 1:  dGBscreen = 0.85; break;

------------------------------------------------------------------------------
Next, apply the following patch to amber7/src/leap/src/leap/commands.c:

*** commands.c	2002/01/14 22:52:38	7.41
--- commands.c	2002/10/02 00:00:58
***************
*** 2540,2547 ****
  		} else if ( !strcmp( sString, "mgbjsb" ) ) {
  			GDefaults.iGBparm = 5;
  			VP0(("Using radii optimized for Amber charges by Jayaram et al. for mod. GB\n"));
  		} else {
! 			VP0(("gbparm value must be 'bondi', 'amber6', 'mbondi', 'pbamber', 'gbjsb', 'mgbjsb'\n"));
  			setUsage();
  		}
  		return(NULL);
--- 2540,2550 ----
  		} else if ( !strcmp( sString, "mgbjsb" ) ) {
  			GDefaults.iGBparm = 5;
  			VP0(("Using radii optimized for Amber charges by Jayaram et al. for mod. GB\n"));
+ 		} else if ( !strcmp( sString, "mbondi2" ) ) {
+ 			GDefaults.iGBparm = 6;
+ 			VP0(("Using H(N)-modified Bondi radii\n"));
  		} else {
! 			VP0(("gbparm value must be 'bondi', 'amber6', 'mbondi', 'mbondi2', 'pbamber', 'gbjsb', 'mgbjsb'\n"));
  			setUsage();
  		}
  		return(NULL);

------------------------------------------------------------------------------

Temporary workarounds: don't use GB with the gaff force field

Files affected: amber7/src/leap/src/leap/{commands.c,unitio.c}