Hi Lewis, at the moment the Python ForceField class (either MMFF or UFF) does not have accessors for force-field contributions. There are C++ assessors:
https://www.rdkit.org/docs/cppapi/classForceFields_1_1ForceField.html#afbf46190b614f8906aae883d1259e026 but they are not exposed to Python. This means that you can access the force-field parameters, but you will have to enumerate torsions yourself. That's easily done by doing the same enumeration that the C++ code does: from rdkit import Chem, ForceField from rdkit.Chem import ChemicalForceFields, rdDistGeom mol = Chem.MolFromSmiles("O=C(O)C(c1ccc(cc1)CCN4CCC(c2nc3ccccc3n2CCOCC)CC4)(C)C") mol_h = Chem.AddHs(mol) rdDistGeom.EmbedMolecule(mol_h) mp = ChemicalForceFields.MMFFGetMoleculeProperties(mol_h) ff = ChemicalForceFields.MMFFGetMoleculeForceField(mol_h, mp) torsion_smarts = "[!$(*#*)&!D1]~[!$(*#*)&!D1]" torsion_query = Chem.MolFromSmarts(torsion_smarts) torsion_bonds = mol_h.GetSubstructMatches(torsion_query) torsion_params = [] allowed_hyb = (Chem.HybridizationType.SP2, Chem.HybridizationType.SP3) for j, k in torsion_bonds: b = mol_h.GetBondBetweenAtoms(j, k) assert b j_atom = mol_h.GetAtomWithIdx(j) k_atom = mol_h.GetAtomWithIdx(k) if (j_atom.GetHybridization() not in allowed_hyb or k_atom.GetHybridization() not in allowed_hyb): continue for jb in j_atom.GetBonds(): if jb.GetIdx() == b.GetIdx(): continue i = jb.GetOtherAtomIdx(j) for kb in k_atom.GetBonds(): if kb.GetIdx() == b.GetIdx() or kb.GetIdx() == jb.GetIdx(): continue l = kb.GetOtherAtomIdx(k) if (i == l): continue params = mp.GetMMFFTorsionParams(mol_h, i, j, k, l) if params: torsion_params.append((i, j, k, l, *params)) len(torsion_params) 191 torsion_params[0] (0, 1, 2, 34, 0, 1.662, 6.152, -0.058) Likewise you may access the other force-field parameters by using the relevant Get*Params function. This will also work for UFF terms. Cheers, p. On Sat, Nov 20, 2021 at 6:41 AM Lewis Martin <lewis.marti...@gmail.com> wrote: > Hi all, > Does anyone have a way to access the atom indices and the parameters for > all of the dihedrals (aka torsions) in an MMFF force field object? I > noticed that one can set mmffVerbosity=2 to print them to stdout, but I'd > like to access them programmatically. I believe Paolo Tosco must have done > this as part of an effort to port MMFF to OpenMM via rdkit ( > https://github.com/ptosco/rdkit/tree/openmm ) , but I can't find the > portion of code that does it! > > Thanks a lot > Lewis > _______________________________________________ > Rdkit-discuss mailing list > Rdkit-discuss@lists.sourceforge.net > https://lists.sourceforge.net/lists/listinfo/rdkit-discuss >
_______________________________________________ Rdkit-discuss mailing list Rdkit-discuss@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
