Project

General

Profile

Feature #994

Na/NA confusion in g_sas (vdwradii.dat)

Added by Sven Steinhauer almost 7 years ago. Updated over 6 years ago.

Status:
Closed
Priority:
Normal
Category:
analysis tools
Target version:
Difficulty:
uncategorized
Close

Description

Dear developers,

I am implanting an implicit force field into a Gromacs 4.5.5 git version from 06/2012. During the development of solvent accessible surface handling, I encountered the problem that g_sas does not distinguish between sodium and a special type of nitrogen atoms, defined in my local vdwradii.dat file and which basically come both from the Amber03.ff. My vdwradii.dat reads like:
? Na 0.1
...
? NA 0.185

During execution of g_sas, I get the following message:
"Warning double different entries ??? NA 0.185 and 0.1 on line 39 in file vdwradii.dat"

I expect
'Na' to be a sodium atomtype (as defined in amber03.ff, element number 11, sodium)
'NA' to be a nitrogen atomtype (as defined in amber03.ff, element number 7, nitrogen)

In atomprop.c, an upper case conversion is happening, which is the source of this atomtype confusion. (see below)

For the moment, I can just comment out one of both atomtypes or implement the radii in a different way as a workaround, but this does not seem to be the proper way of fixing this problem to me.

I guess that removing the conversion line(s) without adapting other code will be fatal to the integrity of the atomtype library. Unfortunately, I do not know well enough the details of the database.

Thank you for your help in advance.

Yours sincerely

Sven Steinhauer

/* file src/gromacs/gmxlib/atomprop.c */
atomprop_query(..) {
..
/**/ upstring(atomname);
strncpy(resname,resnm,MAXQ-1);
/**/ upstring(resname);

j = get_prop_index(&(ap->prop[eprop]),ap->restype, resname, atomname,&bExact);
..
}

/* file string2.c*/
void upstring (char *str) {
int i;

for (i=0; (i < (int)strlen(str)); i++)
str[i] = toupper(str[i]);
}

/* file src/gromacs/gmxlib/atomprop.c */
static int get_prop_index(aprop_t *ap,gmx_residuetype_t restype,
char *resnm,char *atomnm,
gmx_bool *bExact) {
int i,j=NOTFOUND; // j is returned
long int alen,rlen;
long int malen,mrlen;
gmx_bool bProtein,bProtWild;

bProtein  = gmx_residuetype_is_protein(restype,resnm);
bProtWild = (strcmp(resnm,"AAA")==0);
malen = NOTFOUND;
mrlen = NOTFOUND;
// for each atom property
for(i=0; (i&lt;ap-&gt;nprop); i++) {
/* return number of matching characters of resnm,
or NOTFOUND if not at least all characters in char ap->resnm[i] match */
rlen = dbcmp_len(resnm, ap->resnm[i]);
// if not already determined: check, if wildcard or wildprot
if (rlen NOTFOUND) {
if ( (strcmp(ap->resnm[i],"
")0) ||
(strcmp(ap->resnm[i],"???")==0) )
rlen=WILDCARD;
else if (strcmp(ap->resnm[i],"AAA")==0)
rlen=WILDPROT;
}
/* return number of matching characters atomnm,
or NOTFOUND if not at least all characters in char *ap->atomnm[i] match */
alen = dbcmp_len(atomnm, ap->atomnm[i]);
if ( (alen > NOTFOUND) && (rlen > NOTFOUND)) {
if ( ( (alen > malen) && (rlen >= mrlen)) ||
( (rlen > mrlen) && (alen >= malen) ) ) {
malen = alen;
mrlen = rlen;
j = i; //property index
}
}
}
*bExact = ((malen  (long int)strlen(atomnm)) &&
((mrlen (long int)strlen(resnm)) ||
((mrlen WILDPROT) && bProtWild) ||
((mrlen WILDCARD) && !bProtein && !bProtWild)));
if (debug) {
fprintf(debug,"searching residue: %4s atom: %4s\n",resnm,atomnm);
if (j == NOTFOUND)
fprintf(debug," not successful\n");
else
fprintf(debug," match: %4s %4s\n",ap->resnm[j],ap->atomnm[j]);
}
return j;
}

Sven Steinhauer
Université de Liège- Gembloux Agro-BioTech
Passage des Déportés 2
5030 Gembloux
phone +32 81 62 25 32

Associated revisions

Revision cc4cd1e6 (diff)
Added by David van der Spoel over 6 years ago

Fixes #994 confusion in atomprop.

The atomprop code modified the contents of the atom
database files by turning everything into uppercase.
This caused confusion between elements and protein
or nucleic acid atom names.

Change-Id: Ia0ac4636471b7051bf38a714ed16e299bf3d64ab

History

#2 Updated by Mark Abraham over 6 years ago

  • Status changed from New to Feedback wanted

Thanks David.

Sven, does this resolve your issue?

#3 Updated by David van der Spoel over 6 years ago

  • Target version changed from 4.5.6 to 4.6

#4 Updated by Sven Steinhauer over 6 years ago

Thank you for working on this.
Currently, I can not tell you, if this solves the problem.
(When pulling in the development branch of the last few months, my compilation automatism has broken.)
I hope getting this fixed soon and give you feedback then.

#5 Updated by Erik Lindahl over 6 years ago

  • Status changed from Feedback wanted to Closed

Also available in: Atom PDF