break // break the outer loop as well
if end_branch:
if len(branch_queue) > 0:
pdbqt_lines.append(branch_queue.pop())
if old_roots:
current_root = old_roots.pop()
else: // go to next disconnected fragment
next_frag = (frag for frag_num, frag in enumerate(frags)
After Change
mol = Chem.Mol(mol)
// if flexible molecule contains multiple fragments write them separately
if flexible and len(Chem.GetMolFrags(mol)) > 1:
return "".join(MolToPDBQTBlock(frag, flexible=flexible, addHs=addHs,
computeCharges=computeCharges)
for frag in Chem.GetMolFrags(mol, asMols=True))
// Identify donors and acceptors for atom typing
// Acceptors
patt = Chem.MolFromSmarts("[$([O;H1;v2]),"
"$([O;H0;v2;!$(O=N-*),"
"$([O;-;!$(*-N=O)]),"
"$([o;+0])]),"
"$([n;+0;!X3;!$([n;H1](cc)cc),"
"$([$([N;H0]//[C&v4])]),"
"$([N&v3;H0;$(Nc)])]),"
"$([F;$(F-[/ǜ]);!$(FC[F,Cl,Br,I])])]")
acceptors = list(map(lambda x: x[0],
mol.GetSubstructMatches(patt, maxMatches=mol.GetNumAtoms())))
// Donors
patt = Chem.MolFromSmarts("[$([N&!H0&v3,N&!H0&+1&v4,n&H1&+0,$([$([Nv3](-C)(-C)-C)]),"
"$([$(n[n;H1]),"
"$(nc[n;H1])])]),"
// Guanidine can be tautormeic - e.g. Arginine
"$([NX3,NX2]([!O,!S])!@C(!@[NX3,NX2]([!O,!S]))!@[NX3,NX2]([!O,!S])),"
"$([O,S;H1;+0])]")
donors = list(map(lambda x: x[0],
mol.GetSubstructMatches(patt, maxMatches=mol.GetNumAtoms())))
if addHs:
mol = Chem.AddHs(mol, addCoords=True, onlyOnAtoms=donors, )
if addHs or computeCharges:
AllChem.ComputeGasteigerCharges(mol)
atom_lines = PDBQTAtomLines(mol, donors, acceptors)
assert len(atom_lines) == mol.GetNumAtoms()
pdbqt_lines = []
// vina scores
if (mol.HasProp("vina_affinity") and mol.HasProp("vina_rmsd_lb") and
mol.HasProp("vina_rmsd_lb")):
pdbqt_lines.append("REMARK VINA RESULT: " +
("%.1f" % float(mol.GetProp("vina_affinity"))).rjust(8) +
("%.3f" % float(mol.GetProp("vina_rmsd_lb"))).rjust(11) +
("%.3f" % float(mol.GetProp("vina_rmsd_ub"))).rjust(11))
pdbqt_lines.append("REMARK Name = " +
(mol.GetProp("_Name") if mol.HasProp("_Name") else ""))
if flexible:
// Find rotatable bonds
rot_bond = Chem.MolFromSmarts("[!$(*//*)&!D1&!$(C(F)(F)F)&"
"!$(C(Cl)(Cl)Cl)&"
"!$(C(Br)(Br)Br)&"
"!$(C([CH3])([CH3])[CH3])&"
"!$([CD3](=[N,O,S])-!@[/ǝ,O,S!D1])&"
"!$([/ǝ,O,S!D1]-!@[CD3]=[N,O,S])&"
"!$([CD3](=[N+])-!@[/ǝ!D1])&"
"!$([/ǝ!D1]-!@[CD3]=[N+])]-!@[!$(*//*)&"
"!D1&!$(C(F)(F)F)&"
"!$(C(Cl)(Cl)Cl)&"
"!$(C(Br)(Br)Br)&"
"!$(C([CH3])([CH3])[CH3])]")
bond_atoms = list(mol.GetSubstructMatches(rot_bond))
num_torsions = len(bond_atoms)
// Active torsions header
pdbqt_lines.append("REMARK %i active torsions:" % num_torsions)
pdbqt_lines.append("REMARK status: (\"A\" for Active; \"I\" for Inactive)")
for i, (a1, a2) in enumerate(bond_atoms):
pdbqt_lines.append("REMARK%5.0i A between atoms: _%i and _%i"
% (i + 1, a1 + 1, a2 + 1))
// Fragment molecule on bonds to ge rigid fragments
bond_ids = [mol.GetBondBetweenAtoms(a1, a2).GetIdx()
for a1, a2 in bond_atoms]
if bond_ids:
mol_rigid_frags = Chem.FragmentOnBonds(mol, bond_ids, addDummies=False)
else:
mol_rigid_frags = mol
frags = list(Chem.GetMolFrags(mol_rigid_frags))
def weigh_frags(frag):
sort by the fragment size and the number of bonds (secondary)
