rdkit/Code/GraphMol/DetermineBonds/Wrap/testDetermineBonds.py at Release_2022_09 · rdkit/rdkit · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
#
#  Copyright (C) 2022 Greg Landrum
#   @@ All Rights Reserved @@
#
#  This file is part of the RDKit.
#  The contents are covered by the terms of the BSD license
#  which is included in the file license.txt, found at the root
#  of the RDKit source tree.

from rdkit import RDConfig
from rdkit import Chem
from rdkit.Chem import rdDetermineBonds
import unittest
import os
import glob


class TestCase(unittest.TestCase):

  def testVdWConnectivity(self):
    testDir = os.path.join(RDConfig.RDBaseDir, 'Code', 'GraphMol', 'DetermineBonds', 'test_data',
                           'connectivity')
    for fn in glob.glob(os.path.join(testDir, 'test*.xyz')):
      mol = Chem.MolFromXYZFile(fn)
      self.assertIsNotNone(mol)
      smi = mol.GetProp('_FileComments')
      omol = Chem.MolFromSmiles(smi)
      self.assertIsNotNone(omol)

      rdDetermineBonds.DetermineConnectivity(mol, useHueckel=False)
      mol = Chem.RemoveAllHs(mol, sanitize=False)
      self.assertEqual(mol.GetNumAtoms(), omol.GetNumAtoms())
      self.assertEqual(mol.GetNumBonds(), omol.GetNumBonds())
      for aid1 in range(mol.GetNumAtoms()):
        for aid2 in range(aid1 + 1, mol.GetNumAtoms()):
          if omol.GetBondBetweenAtoms(aid1, aid2):
            self.assertIsNotNone(mol.GetBondBetweenAtoms(aid1, aid2))

  def testHueckelConnectivity(self):
    testDir = os.path.join(RDConfig.RDBaseDir, 'Code', 'GraphMol', 'DetermineBonds', 'test_data',
                           'connectivity')
    for fn in glob.glob(os.path.join(testDir, 'test*.xyz')):
      mol = Chem.MolFromXYZFile(fn)
      self.assertIsNotNone(mol)
      smi = mol.GetProp('_FileComments')
      omol = Chem.MolFromSmiles(smi)
      self.assertIsNotNone(omol)
      charge = Chem.GetFormalCharge(omol)

      rdDetermineBonds.DetermineConnectivity(mol, useHueckel=True, charge=charge)
      mol = Chem.RemoveAllHs(mol, sanitize=False)
      self.assertEqual(mol.GetNumAtoms(), omol.GetNumAtoms())
      self.assertEqual(mol.GetNumBonds(), omol.GetNumBonds())
      for aid1 in range(mol.GetNumAtoms()):
        for aid2 in range(aid1 + 1, mol.GetNumAtoms()):
          if omol.GetBondBetweenAtoms(aid1, aid2):
            self.assertIsNotNone(mol.GetBondBetweenAtoms(aid1, aid2))

  def testVdWBonds(self):
    testDir = os.path.join(RDConfig.RDBaseDir, 'Code', 'GraphMol', 'DetermineBonds', 'test_data',
                           'connectivity')
    for fn in glob.glob(os.path.join(testDir, 'test*.xyz')):
      mol = Chem.MolFromXYZFile(fn)
      self.assertIsNotNone(mol)
      smi = mol.GetProp('_FileComments')
      omol = Chem.MolFromSmiles(smi)
      self.assertIsNotNone(omol)
      charge = Chem.GetFormalCharge(omol)

      rdDetermineBonds.DetermineBonds(mol, useHueckel=False, charge=charge)
      mol = Chem.RemoveAllHs(mol, sanitize=False)
      self.assertEqual(mol.GetNumAtoms(), omol.GetNumAtoms())
      self.assertEqual(mol.GetNumBonds(), omol.GetNumBonds())
      for aid1 in range(mol.GetNumAtoms()):
        for aid2 in range(aid1 + 1, mol.GetNumAtoms()):
          if omol.GetBondBetweenAtoms(aid1, aid2):
            self.assertIsNotNone(mol.GetBondBetweenAtoms(aid1, aid2))

  # FIX: this is problematic...
  # def testHueckelBonds(self):
  #   testDir = os.path.join(RDConfig.RDBaseDir, 'Code', 'GraphMol', 'DetermineBonds', 'test_data',
  #                          'connectivity')
  #   for fn in glob.glob(os.path.join(testDir, 'test*.xyz')):
  #     mol = Chem.MolFromXYZFile(fn)
  #     self.assertIsNotNone(mol)
  #     smi = mol.GetProp('_FileComments')
  #     omol = Chem.MolFromSmiles(smi)
  #     self.assertIsNotNone(omol)
  #     charge = Chem.GetFormalCharge(omol)

  #     try:
  #       rdDetermineBonds.DetermineBonds(mol, useHueckel=True, charge=charge)
  #     except ValueError:
  #       print(fn)
  #       print(' ', smi)
  #       print(' ', Chem.MolToSmiles(mol))

  #     mol = Chem.RemoveAllHs(mol, sanitize=False)
  #     self.assertEqual(mol.GetNumAtoms(), omol.GetNumAtoms())
  #     self.assertEqual(mol.GetNumBonds(), omol.GetNumBonds())
  #     for aid1 in range(mol.GetNumAtoms()):
  #       for aid2 in range(aid1 + 1, mol.GetNumAtoms()):
  #         if omol.GetBondBetweenAtoms(aid1, aid2):
  #           self.assertIsNotNone(mol.GetBondBetweenAtoms(aid1, aid2))


if __name__ == '__main__':  # pragma: nocover
  unittest.main()