Merge pull request #2464 from peastman/element

Improved determination of elements from Gromacs top files

wrappers/python/simtk/openmm/app/gromacstopfile.py

@@ -563,20 +563,27 @@ class GromacsTopFile(object):
                     if atomName in atomReplacements:
                         atomName = atomReplacements[atomName]
 
-                    # Try to guess the element.
-
-                    upper = atomName.upper()
-                    if upper.startswith('CL'):
-                        element = elem.chlorine
-                    elif upper.startswith('NA'):
-                        element = elem.sodium
-                    elif upper.startswith('MG'):
-                        element = elem.magnesium
+                    # Try to determine the element.
+
+                    atomicNumber = self._atomTypes[fields[1]][2]
+                    if atomicNumber is None:
+                        # Try to guess the element from the name.
+                        upper = atomName.upper()
+                        if upper.startswith('CL'):
+                            element = elem.chlorine
+                        elif upper.startswith('NA'):
+                            element = elem.sodium
+                        elif upper.startswith('MG'):
+                            element = elem.magnesium
+                        else:
+                            try:
+                                element = elem.get_by_symbol(atomName[0])
+                            except KeyError:
+                                element = None
+                    elif atomicNumber == '0':
+                        element = None
                     else:
-                        try:
-                            element = elem.get_by_symbol(atomName[0])
-                        except KeyError:
-                            element = None
+                        element = elem.Element.getByAtomicNumber(int(atomicNumber))
                     atoms.append(top.addAtom(atomName, element, r))
 
                 # Add bonds to the topology

wrappers/python/tests/TestGromacsTopFile.py

@@ -172,6 +172,13 @@ class TestGromacsTopFile(unittest.TestCase):
 
         top = GromacsTopFile('systems/bnz.top')
         gro = GromacsGroFile('systems/bnz.gro')
+        for atom in top.topology.atoms():
+            if atom.name.startswith('C'):
+                self.assertEqual(elem.carbon, atom.element)
+            elif atom.name.startswith('H'):
+                self.assertEqual(elem.hydrogen, atom.element)
+            else:
+                self.assertIsNone(atom.element)
         system = top.createSystem()
 
         self.assertEqual(26, system.getNumParticles())