1
2
3
4
5
6
7
8
9
10
11 """ Basic EState definitions
12
13 """
14 from __future__ import print_function
15 import numpy
16 from rdkit import Chem
17
19 if atNum<=2: return 1
20 elif atNum <= 10: return 2
21 elif atNum <= 18: return 3
22 elif atNum <= 36: return 4
23 elif atNum <= 54: return 5
24 elif atNum <= 86: return 6
25 else: return 7
26
27
29 """ returns a tuple of EState indices for the molecule
30
31 Reference: Hall, Mohney and Kier. JCICS _31_ 76-81 (1991)
32
33 """
34 if not force and hasattr(mol,'_eStateIndices'):
35 return mol._eStateIndices
36
37 tbl = Chem.GetPeriodicTable()
38 nAtoms = mol.GetNumAtoms()
39 Is = numpy.zeros(nAtoms,numpy.float)
40 for i in range(nAtoms):
41 at = mol.GetAtomWithIdx(i)
42 atNum = at.GetAtomicNum()
43 d = at.GetDegree()
44 if d>0:
45 h = at.GetTotalNumHs()
46 dv = tbl.GetNOuterElecs(atNum)-h
47 N = GetPrincipleQuantumNumber(atNum)
48 Is[i] = (4./(N*N) * dv + 1)/d
49 dists = Chem.GetDistanceMatrix(mol,useBO=0,useAtomWts=0)
50 dists += 1
51 accum = numpy.zeros(nAtoms,numpy.float)
52 for i in range(nAtoms):
53 for j in range(i+1,nAtoms):
54 p = dists[i,j]
55 if p < 1e6:
56 tmp = (Is[i]-Is[j])/(p*p)
57 accum[i] += tmp
58 accum[j] -= tmp
59
60 res = accum+Is
61 mol._eStateIndices=res
62 return res
63 EStateIndices.version='1.0.0'
64
67 MaxEStateIndex.version="1.0.0"
68
71 MinEStateIndex.version="1.0.0"
72
75 MaxAbsEStateIndex.version="1.0.0"
76
79 MinAbsEStateIndex.version="1.0.0"
80
81 if __name__ =='__main__':
82 smis = ['CCCC','CCCCC','CCCCCC','CC(N)C(=O)O','CC(N)C(=O)[O-].[Na+]']
83 for smi in smis:
84 m = Chem.MolFromSmiles(smi)
85 print(smi)
86 inds = EStateIndices(m)
87 print('\t',inds)
88