Skip to content

Commit

Permalink
added aromaticity recognition for symmetry calcs
Browse files Browse the repository at this point in the history
Aromaticity is now perceived in calculating cyclic symmetry numbers.
This commit copies a molecule inside `calculateCyclicSymmetryNumber`,
turns it into an rdkit object, finds all the aromatic bonds and atoms.
The same bonds and atoms in the corresponding rmg molecule object are
then redifined so that the alogorithm is using the correct information.

There is one strange thing about this fix that I dont yet understand.
The unittest thinks cyclic symmetry number of dimethylbenzene is 1
(should be 2). On my local machine, the same
`calculateCyclicSymmetryNumber` function gets it correct, at 2. Also,
i'm inadvertantly printing a bit to stdout, and not sure why it's
happening?

relevant to issues ReactionMechanismGenerator#12, ReactionMechanismGenerator#24, ReactionMechanismGenerator#119, ReactionMechanismGenerator#70, ReactionMechanismGenerator#141
  • Loading branch information
enochd committed Jul 17, 2014
1 parent 4599552 commit d14d7fc
Show file tree
Hide file tree
Showing 3 changed files with 19 additions and 4 deletions.
2 changes: 2 additions & 0 deletions rmgpy/molecule/molecule.py
Original file line number Diff line number Diff line change
Expand Up @@ -567,6 +567,8 @@ def applyAction(self, action):
self.incrementOrder()
elif action[2] == -1:
self.decrementOrder()
elif action[2] == 'B':
self.order = 'B'
else:
raise ActionError('Unable to update Bond due to CHANGE_BOND action: Invalid order "{0}".'.format(action[2]))
else:
Expand Down
20 changes: 17 additions & 3 deletions rmgpy/molecule/symmetry.py
Original file line number Diff line number Diff line change
Expand Up @@ -343,17 +343,31 @@ def calculateCyclicSymmetryNumber(molecule):
Get the symmetry number correction for cyclic regions of a molecule.
For complicated fused rings the smallest set of smallest rings is used.
"""

from rdkit.Chem.rdmolops import SanitizeMol
from rdkit.Chem.rdchem import Mol
mcopy = molecule.toRDKitMol(removeHs=True, returnMapping=False)
SanitizeMol(mcopy)
symmetryNumber = 1

# Get symmetry number for each ring in structure
rings = molecule.getSmallestSetOfSmallestRings()
for ring0 in rings:

# Make copy of structure
structure = molecule.copy(True)
ring = [structure.atoms[molecule.atoms.index(atom)] for atom in ring0]

# Figure out which atoms and bonds are aromatic and reassign appropriately:
for i, atom1 in enumerate(ring0):
for atom2 in ring0[i+1:]:
if molecule.hasBond(atom1, atom2):
if str(mcopy.GetBondBetweenAtoms(i,i+1)) != 'None':
if str(mcopy.GetBondBetweenAtoms(i,i+1).GetBondType()) == 'AROMATIC':
bond = molecule.getBond(atom1, atom2)
bond.applyAction(['CHANGE_BOND', atom1, 'B', atom2])
atom1.atomType.label = 'Cb'
atom2.atomType.label = 'Cb'
else:
pass
# Remove bonds of ring from structure
for i, atom1 in enumerate(ring):
for atom2 in ring[i+1:]:
Expand Down
1 change: 0 additions & 1 deletion rmgpy/molecule/symmetryTest.py
Original file line number Diff line number Diff line change
Expand Up @@ -244,7 +244,6 @@ def testCyclicSymmetryNumberCyclohexane(self):
symmetryNumber = calculateCyclicSymmetryNumber(molecule)
self.assertEqual(symmetryNumber, 12)

@work_in_progress
def testCyclicSymmetryNumberBenzene(self):
"""
Test the Molecule.calculateCyclicSymmetryNumber() method.
Expand Down

0 comments on commit d14d7fc

Please sign in to comment.