From ec366c1ab78b461afc5cdf4b8352de8679f0338e Mon Sep 17 00:00:00 2001 From: Greg Landrum Date: Wed, 22 Jan 2020 15:10:58 +0100 Subject: [PATCH] Dev/pvs studio cleanups2 (#2877) * a round of cleanups courtesy of PVS studio * add a test to make sure that a warning is a false alarm * bug fix * Fix a UFF bug * more PVS studio cleanups * next round of PVS studio cleanups * completely remove the chances for that bug * changes in response to review * add an additional test + a bit of reformatting that snuck in --- Code/DataStructs/SparseIntVect.h | 7 +- Code/ForceField/Wrap/PyForceField.h | 2 +- Code/Geometry/UniformGrid3D.cpp | 5 +- .../Enumerate/EvenSamplePairs.cpp | 28 +- .../GraphMol/ChemReactions/ReactionRunner.cpp | 3 - .../ChemTransforms/ChemTransforms.cpp | 24 +- .../GraphMol/ChemTransforms/MolFragmenter.cpp | 4 +- Code/GraphMol/DistGeomHelpers/Embedder.cpp | 5 +- Code/GraphMol/FMCS/FMCS.cpp | 115 +- Code/GraphMol/FMCS/FMCS.h | 23 +- Code/GraphMol/FMCS/MaximumCommonSubgraph.cpp | 55 +- Code/GraphMol/FMCS/MaximumCommonSubgraph.h | 21 +- Code/GraphMol/FileParsers/MaeMolSupplier.cpp | 3 +- .../FileParsers/SmilesMolSupplier.cpp | 2 +- Code/GraphMol/FileParsers/TplFileParser.cpp | 6 +- Code/GraphMol/Fingerprints/test1.cpp | 33 + .../ForceFieldHelpers/MMFF/AtomTyper.cpp | 138 +- .../ForceFieldHelpers/UFF/Builder.cpp | 27 +- .../ForceFieldHelpers/Wrap/rdForceFields.cpp | 20 +- .../GraphMol/FragCatalog/FragCatGenerator.cpp | 5 +- .../GraphMol/FragCatalog/FragCatalogUtils.cpp | 4 +- Code/GraphMol/FragCatalog/FragFPGenerator.cpp | 7 +- Code/GraphMol/MolAlign/O3AAlignMolecules.cpp | 94 +- Code/GraphMol/MolAlign/O3AAlignMolecules.h | 6 +- Code/GraphMol/MolCatalog/MolCatalogEntry.cpp | 4 +- Code/GraphMol/MolCatalog/MolCatalogEntry.h | 2 +- .../MolChemicalFeatures/FeatureParser.cpp | 20 +- Code/GraphMol/MolDraw2D/MolDraw2D.cpp | 12 - Code/GraphMol/MolHash/MolHash.cpp | 26 +- Code/GraphMol/MolHash/catch_tests.cpp | 45 +- Code/GraphMol/MolHash/hashfunctions.cpp | 39 +- Code/GraphMol/MolInterchange/Parser.cpp | 7 +- Code/GraphMol/MolPickler.cpp | 9 +- .../AcidBaseCatalog/AcidBaseCatalogEntry.h | 2 - .../AcidBaseCatalog/AcidBaseCatalogUtils.cpp | 2 +- Code/GraphMol/MolStandardize/Charge.cpp | 4 +- .../FragmentCatalog/FragmentCatalogUtils.cpp | 3 +- Code/GraphMol/MolStandardize/Metal.cpp | 23 +- .../TautomerCatalog/TautomerCatalogUtils.cpp | 2 +- .../TransformCatalogUtils.cpp | 2 +- Code/GraphMol/Resonance.cpp | 32 +- Code/GraphMol/SmilesParse/SmilesWrite.cpp | 5 +- Code/GraphMol/Substruct/SubstructUtils.cpp | 4 - Code/GraphMol/Substruct/vf2.hpp | 1169 ++++++++--------- Code/GraphMol/Wrap/EditableMol.cpp | 5 +- Code/GraphMol/Wrap/ForwardSDMolSupplier.cpp | 5 +- Code/GraphMol/Wrap/MaeMolSupplier.cpp | 3 +- Code/GraphMol/Wrap/MolOps.cpp | 2 +- Code/ML/Cluster/Murtagh/Clustering.cpp | 1 + Code/ML/Data/cQuantize.cpp | 8 +- Code/SimDivPickers/LeaderPicker.h | 10 +- Code/SimDivPickers/MaxMinPicker.h | 4 - 52 files changed, 1063 insertions(+), 1024 deletions(-) diff --git a/Code/DataStructs/SparseIntVect.h b/Code/DataStructs/SparseIntVect.h index 6188118e5..ef5965d61 100644 --- a/Code/DataStructs/SparseIntVect.h +++ b/Code/DataStructs/SparseIntVect.h @@ -480,11 +480,8 @@ double DiceSimilarity(const SparseIntVect &v1, v2Sum = v2.getTotalVal(true); double denom = v1Sum + v2Sum; if (fabs(denom) < 1e-6) { - if (returnDistance) { - return 1.0; - } else { - return 0.0; - } + // no need to worry about returnDistance here + return 0.0; } double minV = v1Sum < v2Sum ? v1Sum : v2Sum; if (2. * minV / denom < bounds) { diff --git a/Code/ForceField/Wrap/PyForceField.h b/Code/ForceField/Wrap/PyForceField.h index 852c43fe8..63de3c84f 100644 --- a/Code/ForceField/Wrap/PyForceField.h +++ b/Code/ForceField/Wrap/PyForceField.h @@ -34,8 +34,8 @@ class PyForceField { } int addExtraPoint(double x, double y, double z, bool fixed = true) { - RDGeom::Point3D *pt = new RDGeom::Point3D(x, y, z); PRECONDITION(this->field, "no force field"); + RDGeom::Point3D *pt = new RDGeom::Point3D(x, y, z); this->extraPoints.push_back(boost::shared_ptr(pt)); unsigned int ptIdx = this->extraPoints.size() - 1; RDGeom::Point3D *ptr = this->extraPoints[ptIdx].get(); diff --git a/Code/Geometry/UniformGrid3D.cpp b/Code/Geometry/UniformGrid3D.cpp index 707d471e4..e19664506 100644 --- a/Code/Geometry/UniformGrid3D.cpp +++ b/Code/Geometry/UniformGrid3D.cpp @@ -156,8 +156,9 @@ Point3D UniformGrid3D::getGridPointLoc(unsigned int pointId) const { Point3D res; res.x = (pointId % d_numX) * d_spacing; // the rounding here is intentional, we want the coordinates of a grid point - res.y = ((pointId % (d_numX * d_numY)) / d_numX) * d_spacing; - res.z = (pointId / (d_numX * d_numY)) * d_spacing; + res.y = + static_cast((pointId % (d_numX * d_numY)) / d_numX) * d_spacing; + res.z = static_cast(pointId / (d_numX * d_numY)) * d_spacing; res += d_offSet; // d_origin; return res; } diff --git a/Code/GraphMol/ChemReactions/Enumerate/EvenSamplePairs.cpp b/Code/GraphMol/ChemReactions/Enumerate/EvenSamplePairs.cpp index 4f489e9ef..85bff6a15 100644 --- a/Code/GraphMol/ChemReactions/Enumerate/EvenSamplePairs.cpp +++ b/Code/GraphMol/ChemReactions/Enumerate/EvenSamplePairs.cpp @@ -40,9 +40,10 @@ using namespace EnumerationTypes; void EvenSamplePairsStrategy::initializeStrategy(const ChemicalReaction &, const BBS &bbs) { // If we fail here, someone has a ridiculous amount of memory - PRECONDITION(m_numPermutations != EnumerationStrategyBase::EnumerationOverflow, - "Cannot represent all permutations for the even sampler"); - + PRECONDITION( + m_numPermutations != EnumerationStrategyBase::EnumerationOverflow, + "Cannot represent all permutations for the even sampler"); + boost::uint64_t npos = bbs.size(); used_count.resize(npos); std::fill(used_count.begin(), used_count.end(), 0); @@ -97,8 +98,8 @@ void EvenSamplePairsStrategy::initializeStrategy(const ChemicalReaction &, // This is fairly suboptimal for large collections // of building blocks and may take a while to // terminate... -bool EvenSamplePairsStrategy::try_add(boost::uint64_t seed) { - const RGROUPS &digits = decode(seed); +bool EvenSamplePairsStrategy::try_add(boost::uint64_t iseed) { + const RGROUPS &digits = decode(iseed); const RGROUPS &rgroups = m_permutationSizes; boost::uint64_t islack = 0; boost::uint64_t num_rgroups = m_permutationSizes.size(); @@ -145,7 +146,7 @@ bool EvenSamplePairsStrategy::try_add(boost::uint64_t seed) { if (used_count[i] == rdcast(rgroups[i])) { // complete variable scan => initialize if (nslack > min_nslack && rgroups[i] > 1) // cleared slack on i - nslack = min_nslack; + nslack = min_nslack; used_count[i] = 0; for (boost::uint64_t j = 0; j < rgroups[i]; ++j) { @@ -186,14 +187,14 @@ bool EvenSamplePairsStrategy::try_add(boost::uint64_t seed) { } } - selected.insert(seed); + selected.insert(iseed); return true; } const RGROUPS &EvenSamplePairsStrategy::next() { nslack = 0; - while (m_numPermutationsProcessed < rdcast(m_numPermutations)) { - bool added = false; + while (m_numPermutationsProcessed < + rdcast(m_numPermutations)) { for (boost::uint64_t l = 0; l < M; ++l) { seed = ((seed * a + b) % M); if (seed > rdcast(m_numPermutations)) { @@ -204,16 +205,13 @@ const RGROUPS &EvenSamplePairsStrategy::next() { continue; } else if (try_add(seed)) { m_numPermutationsProcessed++; - added = true; return decode(seed); } } - if (!added) { - // loosen heuristic - nslack += 1; - min_nslack += 1; - } + // loosen heuristic + nslack += 1; + min_nslack += 1; } throw EnumerationStrategyException("Ran out of molecules"); diff --git a/Code/GraphMol/ChemReactions/ReactionRunner.cpp b/Code/GraphMol/ChemReactions/ReactionRunner.cpp index 92b1ffe48..862352b7c 100644 --- a/Code/GraphMol/ChemReactions/ReactionRunner.cpp +++ b/Code/GraphMol/ChemReactions/ReactionRunner.cpp @@ -1600,9 +1600,6 @@ ROMol *reduceProductToSideChains(const ROMOL_SPTR &product, if (atommapno) { Atom *at = mol->getAtomWithIdx(idx); - PRECONDITION(at, - "Internal error in reduceProductToSideChains, bad " - "atom retrieval"); at->setProp(common_properties::molAtomMapNumber, atommapno); } } diff --git a/Code/GraphMol/ChemTransforms/ChemTransforms.cpp b/Code/GraphMol/ChemTransforms/ChemTransforms.cpp index 99b1bf9c8..a0c81053c 100644 --- a/Code/GraphMol/ChemTransforms/ChemTransforms.cpp +++ b/Code/GraphMol/ChemTransforms/ChemTransforms.cpp @@ -51,14 +51,15 @@ void updateSubMolConfs(const ROMol &mol, RWMol &res, res.addConformer(newConf, false); } } -} +} // namespace ROMol *deleteSubstructs(const ROMol &mol, const ROMol &query, bool onlyFrags, bool useChirality) { - RWMol *res = static_cast(new ROMol(mol, false)); + RWMol *res = new RWMol(mol, false); std::vector fgpMatches; std::vector::const_iterator mati; - VECT_INT_VECT matches; // all matches on the molecule - list of list of atom ids + VECT_INT_VECT + matches; // all matches on the molecule - list of list of atom ids MatchVectType::const_iterator mi; // do the substructure matching and get the atoms that match the query const bool uniquify = true; @@ -259,7 +260,7 @@ ROMol *replaceSidechains(const ROMol &mol, const ROMol &coreQuery, } std::vector keepList; - RWMol *newMol = static_cast(new ROMol(mol)); + RWMol *newMol = new RWMol(mol); unsigned int nDummies = 0; for (MatchVectType::const_iterator mvit = matchV.begin(); mvit != matchV.end(); mvit++) { @@ -349,7 +350,7 @@ ROMol *replaceCore(const ROMol &mol, const ROMol &core, allIndices[mvit.second] = mvit.first; } - RWMol *newMol = static_cast(new ROMol(mol)); + RWMol *newMol = new RWMol(mol); std::vector keepList; std::map dummyAtomMap; @@ -374,7 +375,7 @@ ROMol *replaceCore(const ROMol &mol, const ROMol &core, int mapping = -1; // if we were not in the matching list, still keep // the original indices (replaceDummies=False) - if (matchingIndices[i] == -1 && allIndices[i] != -1) { + if (allIndices[i] != -1) { mapping = allIndices[i]; } // loop over our neighbors and see if any are in the match: @@ -650,7 +651,8 @@ void addRecursiveQueries( if (!at->hasProp(propName)) continue; std::string pval; at->getProp(propName, pval); - std::string maybeSmarts = pval; // keep unmodified in case we are a smarts string + std::string maybeSmarts = + pval; // keep unmodified in case we are a smarts string boost::algorithm::to_lower(pval); if (reactantLabels != nullptr) { std::pair label(at->getIdx(), pval); @@ -678,8 +680,7 @@ void addRecursiveQueries( auto iter = queries.find(pval); if (iter == queries.end()) { notFound = true; - } - else { + } else { qToAdd = new RecursiveStructureQuery(new ROMol(*(iter->second))); } } @@ -689,8 +690,7 @@ void addRecursiveQueries( RWMol *m = 0; try { m = SmartsToMol(maybeSmarts); - if (!m) - throw KeyErrorException(pval); + if (!m) throw KeyErrorException(pval); qToAdd = new RecursiveStructureQuery(m); } catch (...) { throw KeyErrorException(pval); @@ -718,7 +718,7 @@ void parseQueryDefFile(std::istream *inStream, boost::char_separator sep(delimiter.c_str()); unsigned int line = 0; std::string tempStr; - while (!inStream->eof()) { + while (!inStream->eof() && !inStream->fail()) { line++; tempStr = getLine(inStream); if (tempStr == "" || tempStr.find(comment) == 0) continue; diff --git a/Code/GraphMol/ChemTransforms/MolFragmenter.cpp b/Code/GraphMol/ChemTransforms/MolFragmenter.cpp index e5e423a3d..2ac085774 100644 --- a/Code/GraphMol/ChemTransforms/MolFragmenter.cpp +++ b/Code/GraphMol/ChemTransforms/MolFragmenter.cpp @@ -52,7 +52,7 @@ void constructFragmenterAtomTypes( PRECONDITION(inStream, "no stream"); defs.clear(); unsigned int line = 0; - while (!inStream->eof()) { + while (!inStream->eof() && !inStream->fail()) { ++line; std::string tempStr = getLine(inStream); if (tempStr == "" || tempStr.find(comment) == 0) continue; @@ -141,7 +141,7 @@ void constructFragmenterBondTypes( defs.resize(0); unsigned int line = 0; - while (!inStream->eof()) { + while (!inStream->eof() && !inStream->fail()) { ++line; std::string tempStr = getLine(inStream); if (tempStr == "" || tempStr.find(comment) == 0) continue; diff --git a/Code/GraphMol/DistGeomHelpers/Embedder.cpp b/Code/GraphMol/DistGeomHelpers/Embedder.cpp index 72ade6e8f..35cb817a5 100644 --- a/Code/GraphMol/DistGeomHelpers/Embedder.cpp +++ b/Code/GraphMol/DistGeomHelpers/Embedder.cpp @@ -647,7 +647,7 @@ bool embedPoints(RDGeom::PointPtrVect *positions, detail::EmbedArgs eargs, } } // if(gotCoords) } // while - if (seed > -1 && rng) { + if (seed > -1) { delete rng; delete generator; delete distrib; @@ -1009,7 +1009,8 @@ void EmbedMultipleConfs(ROMol &mol, INT_VECT &res, unsigned int numConfs, "size of boundsMat provided does not match the number of atoms in " "the molecule."); } - collectBondsAndAngles((*piece.get()),etkdgDetails.bonds, etkdgDetails.angles); + collectBondsAndAngles((*piece.get()), etkdgDetails.bonds, + etkdgDetails.angles); mmat.reset(new DistGeom::BoundsMatrix(*params.boundsMat)); } diff --git a/Code/GraphMol/FMCS/FMCS.cpp b/Code/GraphMol/FMCS/FMCS.cpp index 0c963e45b..a98169cd5 100644 --- a/Code/GraphMol/FMCS/FMCS.cpp +++ b/Code/GraphMol/FMCS/FMCS.cpp @@ -92,13 +92,13 @@ MCSResult findMCS_P(const std::vector& mols, } MCSResult findMCS(const std::vector& mols, bool maximizeBonds, - double threshold, unsigned timeout, bool verbose, - bool matchValences, bool ringMatchesRingOnly, - bool completeRingsOnly, bool matchChiralTag, - AtomComparator atomComp, BondComparator bondComp) { + double threshold, unsigned timeout, bool verbose, + bool matchValences, bool ringMatchesRingOnly, + bool completeRingsOnly, bool matchChiralTag, + AtomComparator atomComp, BondComparator bondComp) { return findMCS(mols, maximizeBonds, threshold, timeout, verbose, - matchValences, ringMatchesRingOnly, completeRingsOnly, - matchChiralTag, atomComp, bondComp, IgnoreRingFusion); + matchValences, ringMatchesRingOnly, completeRingsOnly, + matchChiralTag, atomComp, bondComp, IgnoreRingFusion); } MCSResult findMCS(const std::vector& mols, bool maximizeBonds, @@ -143,7 +143,8 @@ MCSResult findMCS(const std::vector& mols, bool maximizeBonds, ps->BondCompareParameters.RingMatchesRingOnly = ringMatchesRingOnly; ps->BondCompareParameters.CompleteRingsOnly = completeRingsOnly; ps->BondCompareParameters.MatchFusedRings = (ringComp != IgnoreRingFusion); - ps->BondCompareParameters.MatchFusedRingsStrict = (ringComp == StrictRingFusion); + ps->BondCompareParameters.MatchFusedRingsStrict = + (ringComp == StrictRingFusion); MCSResult res = findMCS(mols, ps); delete ps; return res; @@ -242,13 +243,15 @@ bool MCSAtomCompareIsotopes(const MCSAtomCompareParameters& p, } bool MCSAtomCompareAnyHeavyAtom(const MCSAtomCompareParameters& p, - const ROMol& mol1, unsigned int atom1, - const ROMol& mol2, unsigned int atom2, void*) { + const ROMol& mol1, unsigned int atom1, + const ROMol& mol2, unsigned int atom2, void*) { const Atom& a1 = *mol1.getAtomWithIdx(atom1); const Atom& a2 = *mol2.getAtomWithIdx(atom2); - //Any atom, including H, matches another atom of the same type, according to the other flags - if (a1.getAtomicNum() == a2.getAtomicNum() || (a1.getAtomicNum() > 1 && a2.getAtomicNum() > 1)){ - return MCSAtomCompareAny(p,mol1,atom1,mol2,atom2,nullptr); + // Any atom, including H, matches another atom of the same type, according to + // the other flags + if (a1.getAtomicNum() == a2.getAtomicNum() || + (a1.getAtomicNum() > 1 && a2.getAtomicNum() > 1)) { + return MCSAtomCompareAny(p, mol1, atom1, mol2, atom2, nullptr); } return false; } @@ -367,12 +370,13 @@ bool MCSBondCompareOrderExact(const MCSBondCompareParameters& p, } //=== RING COMPARE ======================================================== -inline static bool ringFusionCheck(const short unsigned c1[], - const short unsigned c2[], const ROMol& mol1, - const FMCS::Graph& query, const ROMol& mol2, - const FMCS::Graph& target, const MCSParameters* p) { - const RingInfo *ri2 = mol2.getRingInfo(); - const VECT_INT_VECT &br2 = ri2->bondRings(); +inline static bool ringFusionCheck(const std::uint32_t c1[], + const std::uint32_t c2[], const ROMol& mol1, + const FMCS::Graph& query, const ROMol& mol2, + const FMCS::Graph& target, + const MCSParameters* p) { + const RingInfo* ri2 = mol2.getRingInfo(); + const VECT_INT_VECT& br2 = ri2->bondRings(); std::vector nonFusedBonds(br2.size(), 0); std::vector numMcsBondRings(br2.size(), 0); RDKit::FMCS::Graph::BOND_ITER_PAIR bpIter = boost::edges(query); @@ -380,8 +384,9 @@ inline static bool ringFusionCheck(const short unsigned c1[], // numMcsBondRings stores the number of bonds which // are part of the MCS for each ring for (auto it = bpIter.first; it != bpIter.second; ++it) { - const Bond *b = mol2.getBondBetweenAtoms( - target[c2[boost::source(*it, query)]], target[c2[boost::target(*it, query)]]); + const Bond* b = + mol2.getBondBetweenAtoms(target[c2[boost::source(*it, query)]], + target[c2[boost::target(*it, query)]]); if (!b) continue; unsigned int bi = b->getIdx(); if (!ri2->numBondRings(bi)) continue; @@ -393,9 +398,8 @@ inline static bool ringFusionCheck(const short unsigned c1[], // nonFusedBonds stores the number of non-fused bonds // (i.e., which belong to a single ring) for each ring for (i = 0; i < br2.size(); ++i) { - for (auto bi: br2[i]) { - if (ri2->numBondRings(bi) == 1) - ++nonFusedBonds[i]; + for (auto bi : br2[i]) { + if (ri2->numBondRings(bi) == 1) ++nonFusedBonds[i]; } } /* if a ring has at least one bond which is part of the MCS, @@ -410,16 +414,16 @@ inline static bool ringFusionCheck(const short unsigned c1[], indeed part of the MCS. */ bool missingFusedBond = false; for (i = 0; i < br2.size(); ++i) { - if (numMcsBondRings[i] - && numMcsBondRings[i] > (br2[i].size() - nonFusedBonds[i]) - && numMcsBondRings[i] < nonFusedBonds[i]) { + if (numMcsBondRings[i] && + numMcsBondRings[i] > (br2[i].size() - nonFusedBonds[i]) && + numMcsBondRings[i] < nonFusedBonds[i]) { if (p->BondCompareParameters.CompleteRingsOnly) return true; else continue; } - if (numMcsBondRings[i] == nonFusedBonds[i] - && numMcsBondRings[i] < br2[i].size()) + if (numMcsBondRings[i] == nonFusedBonds[i] && + numMcsBondRings[i] < br2[i].size()) missingFusedBond = true; } /* If we found that the MCS is missing a fused bond, we may need to @@ -451,13 +455,14 @@ inline static bool ringFusionCheck(const short unsigned c1[], // fused bonds, but not both. In strict mode we allow neither. if (!p->BondCompareParameters.MatchFusedRingsStrict) missingFusedBond = false; - const RingInfo *ri1 = mol1.getRingInfo(); - const VECT_INT_VECT &br1 = ri1->bondRings(); + const RingInfo* ri1 = mol1.getRingInfo(); + const VECT_INT_VECT& br1 = ri1->bondRings(); std::set mcsRingBondIdxSet; // put all MCS bond indices which belong to one or more rings in a set for (auto it = bpIter.first; it != bpIter.second; ++it) { - const Bond *b = mol1.getBondBetweenAtoms( - query[c1[boost::source(*it, query)]], query[c1[boost::target(*it, query)]]); + const Bond* b = + mol1.getBondBetweenAtoms(query[c1[boost::source(*it, query)]], + query[c1[boost::target(*it, query)]]); if (b && ri1->numBondRings(b->getIdx())) mcsRingBondIdxSet.insert(b->getIdx()); } @@ -468,23 +473,22 @@ inline static bool ringFusionCheck(const short unsigned c1[], // for each ring, check how many bonds belong to MCS (mcsBonds), // how many bonds in the MCS belong to multiple rings (mcsFusedBonds). // and how many bonds in each ring belong to a single ring (nonFusedBonds) - for (auto bi: br1[i]) { + for (auto bi : br1[i]) { bool isFused = (ri1->numBondRings(bi) > 1); - if (!isFused) - ++nonFusedBonds; - if (std::find(mcsRingBondIdxSet.begin(), - mcsRingBondIdxSet.end(), bi) != mcsRingBondIdxSet.end()) { + if (!isFused) ++nonFusedBonds; + if (std::find(mcsRingBondIdxSet.begin(), mcsRingBondIdxSet.end(), bi) != + mcsRingBondIdxSet.end()) { ++mcsBonds; - if (isFused) - ++mcsFusedBonds; + if (isFused) ++mcsFusedBonds; } } - if (!p->BondCompareParameters.CompleteRingsOnly - && mcsBonds > mcsFusedBonds && mcsBonds - mcsFusedBonds < nonFusedBonds) + if (!p->BondCompareParameters.CompleteRingsOnly && + mcsBonds > mcsFusedBonds && mcsBonds - mcsFusedBonds < nonFusedBonds) continue; - // if the ring is part of the MCS, and we found more bonds than the fused ones - // (i.e., the ring is not simply adjacent to an MCS ring) but less than - // the bonds which are part of this ring, then some fused bond is missing. + // if the ring is part of the MCS, and we found more bonds than the fused + // ones (i.e., the ring is not simply adjacent to an MCS ring) but less + // than the bonds which are part of this ring, then some fused bond is + // missing. if (mcsBonds && mcsBonds > mcsFusedBonds && mcsBonds < br1[i].size()) { missingFusedBond2 = true; break; @@ -495,23 +499,22 @@ inline static bool ringFusionCheck(const short unsigned c1[], return (!missingFusedBond && !missingFusedBond2); } -bool FinalMatchCheckFunction(const short unsigned c1[], - const short unsigned c2[], const ROMol& mol1, - const FMCS::Graph& query, const ROMol& mol2, - const FMCS::Graph& target, - const MCSParameters *p) { - if ((p->BondCompareParameters.MatchFusedRings - || p->BondCompareParameters.MatchFusedRingsStrict) - && !ringFusionCheck(c1, c2, mol1, query, mol2, target, p)) +bool FinalMatchCheckFunction(const std::uint32_t c1[], const std::uint32_t c2[], + const ROMol& mol1, const FMCS::Graph& query, + const ROMol& mol2, const FMCS::Graph& target, + const MCSParameters* p) { + if ((p->BondCompareParameters.MatchFusedRings || + p->BondCompareParameters.MatchFusedRingsStrict) && + !ringFusionCheck(c1, c2, mol1, query, mol2, target, p)) return false; - if (p->AtomCompareParameters.MatchChiralTag - && !FinalChiralityCheckFunction(c1, c2, mol1, query, mol2, target, p)) + if (p->AtomCompareParameters.MatchChiralTag && + !FinalChiralityCheckFunction(c1, c2, mol1, query, mol2, target, p)) return false; return true; } -bool FinalChiralityCheckFunction(const short unsigned c1[], - const short unsigned c2[], const ROMol& mol1, +bool FinalChiralityCheckFunction(const std::uint32_t c1[], + const std::uint32_t c2[], const ROMol& mol1, const FMCS::Graph& query, const ROMol& mol2, const FMCS::Graph& target, const MCSParameters* /*unused*/) { diff --git a/Code/GraphMol/FMCS/FMCS.h b/Code/GraphMol/FMCS/FMCS.h index 24cfd9044..6fa25f1fb 100644 --- a/Code/GraphMol/FMCS/FMCS.h +++ b/Code/GraphMol/FMCS/FMCS.h @@ -34,7 +34,7 @@ struct RDKIT_FMCS_EXPORT MCSBondCompareParameters { }; typedef bool (*MCSFinalMatchCheckFunction)( - const short unsigned c1[], const short unsigned c2[], const ROMol& mol1, + const std::uint32_t c1[], const std::uint32_t c2[], const ROMol& mol1, const FMCS::Graph& query, const ROMol& mol2, const FMCS::Graph& target, const MCSParameters* p); typedef bool (*MCSAtomCompareFunction)(const MCSAtomCompareParameters& p, @@ -51,10 +51,9 @@ RDKIT_FMCS_EXPORT bool MCSAtomCompareAny(const MCSAtomCompareParameters& p, const ROMol& mol1, unsigned int atom1, const ROMol& mol2, unsigned int atom2, void* userData); -RDKIT_FMCS_EXPORT bool MCSAtomCompareAnyHeavyAtom(const MCSAtomCompareParameters& p, - const ROMol& mol1, unsigned int atom1, - const ROMol& mol2, unsigned int atom2, - void* userData); +RDKIT_FMCS_EXPORT bool MCSAtomCompareAnyHeavyAtom( + const MCSAtomCompareParameters& p, const ROMol& mol1, unsigned int atom1, + const ROMol& mol2, unsigned int atom2, void* userData); RDKIT_FMCS_EXPORT bool MCSAtomCompareElements( const MCSAtomCompareParameters& p, const ROMol& mol1, unsigned int atom1, @@ -116,6 +115,7 @@ struct RDKIT_FMCS_EXPORT MCSResult { bool Canceled; // interrupted by timeout or user defined progress callback. // Contains valid current MCS ! RWMOL_SPTR QueryMol; + public: MCSResult() : NumAtoms(0), NumBonds(0), Canceled(false) {} bool isCompleted() const { return !Canceled; } @@ -145,14 +145,11 @@ typedef enum { PermissiveRingFusion, StrictRingFusion } RingComparator; -RDKIT_FMCS_EXPORT MCSResult -findMCS(const std::vector& mols, bool maximizeBonds, - double threshold, unsigned timeout, bool verbose, - bool matchValences, bool ringMatchesRingOnly, - bool completeRingsOnly, bool matchChiralTag, - AtomComparator atomComp, - BondComparator bondComp, - RingComparator ringComp); +RDKIT_FMCS_EXPORT MCSResult findMCS( + const std::vector& mols, bool maximizeBonds, double threshold, + unsigned timeout, bool verbose, bool matchValences, + bool ringMatchesRingOnly, bool completeRingsOnly, bool matchChiralTag, + AtomComparator atomComp, BondComparator bondComp, RingComparator ringComp); RDKIT_FMCS_EXPORT MCSResult findMCS(const std::vector& mols, bool maximizeBonds, double threshold = 1.0, unsigned timeout = 3600, bool verbose = false, diff --git a/Code/GraphMol/FMCS/MaximumCommonSubgraph.cpp b/Code/GraphMol/FMCS/MaximumCommonSubgraph.cpp index 90247af37..6ae51d940 100755 --- a/Code/GraphMol/FMCS/MaximumCommonSubgraph.cpp +++ b/Code/GraphMol/FMCS/MaximumCommonSubgraph.cpp @@ -37,8 +37,8 @@ MaximumCommonSubgraph::MaximumCommonSubgraph(const MCSParameters* params) { Parameters.ProgressCallbackUserData = &To; } if ((Parameters.AtomCompareParameters.MatchChiralTag || - Parameters.BondCompareParameters.MatchFusedRings || - Parameters.BondCompareParameters.MatchFusedRingsStrict) && + Parameters.BondCompareParameters.MatchFusedRings || + Parameters.BondCompareParameters.MatchFusedRingsStrict) && nullptr == Parameters.FinalMatchChecker) { Parameters.FinalMatchChecker = FinalMatchCheckFunction; if (Parameters.AtomCompareParameters.MatchChiralTag) @@ -600,7 +600,8 @@ struct AtomMatch { // for each seed atom (matched) }; typedef std::vector AtomMatchSet; -std::pair MaximumCommonSubgraph::generateResultSMARTSAndQueryMol( +std::pair +MaximumCommonSubgraph::generateResultSMARTSAndQueryMol( const MCS& mcsIdx) const { // match the result MCS with all targets to check if it is exact match or // template @@ -656,8 +657,8 @@ std::pair MaximumCommonSubgraph::generateResultSMARTSAndQue // Generate result's SMARTS // create molecule from MCS for MolToSmarts() - RWMol *mol = new RWMol(); - const RingInfo *ri = QueryMolecule->getRingInfo(); + RWMol* mol = new RWMol(); + const RingInfo* ri = QueryMolecule->getRingInfo(); unsigned ai = 0; // SeedAtomIdx for (auto atom = mcsIdx.Atoms.begin(); atom != mcsIdx.Atoms.end(); atom++, ai++) { @@ -682,7 +683,7 @@ std::pair MaximumCommonSubgraph::generateResultSMARTSAndQue } } if (Parameters.AtomCompareParameters.RingMatchesRingOnly) { - ATOM_EQUALS_QUERY *q = makeAtomInRingQuery(); + ATOM_EQUALS_QUERY* q = makeAtomInRingQuery(); q->setNegation(!ri->numAtomRings((*atom)->getIdx())); a.expandQuery(q, Queries::COMPOSITE_AND, true); } @@ -708,7 +709,7 @@ std::pair MaximumCommonSubgraph::generateResultSMARTSAndQue b.setStereo(bm->second->getStereo()); } if (Parameters.BondCompareParameters.RingMatchesRingOnly) { - BOND_EQUALS_QUERY *q = makeBondIsInRingQuery(); + BOND_EQUALS_QUERY* q = makeBondIsInRingQuery(); q->setNegation(!ri->numBondRings((*bond)->getIdx())); b.expandQuery(q, Queries::COMPOSITE_AND, true); } @@ -718,32 +719,32 @@ std::pair MaximumCommonSubgraph::generateResultSMARTSAndQue return std::make_pair(MolToSmarts(*mol, true), mol); } -bool MaximumCommonSubgraph::addFusedBondQueries(const MCS& mcsIdx, RWMol *rwMol) const { - const RingInfo *ri = QueryMolecule->getRingInfo(); +bool MaximumCommonSubgraph::addFusedBondQueries(const MCS& mcsIdx, + RWMol* rwMol) const { + const RingInfo* ri = QueryMolecule->getRingInfo(); unsigned bi = 0; // Seed Idx bool haveFusedBondQuery = false; std::map> bondRingMembership; - const VECT_INT_VECT &br = ri->bondRings(); + const VECT_INT_VECT& br = ri->bondRings(); std::vector mcsRingSize(br.size(), 0); for (size_t ringIdx = 0; ringIdx < br.size(); ++ringIdx) { - for (int bondIdx: br[ringIdx]) - bondRingMembership[bondIdx].insert(ringIdx); + for (int bondIdx : br[ringIdx]) bondRingMembership[bondIdx].insert(ringIdx); } for (auto bond = mcsIdx.Bonds.begin(); bond != mcsIdx.Bonds.end(); ++bond) { unsigned int bondIdx = (*bond)->getIdx(); if (!ri->numBondRings(bondIdx)) continue; - for (size_t ringIdx: bondRingMembership[bondIdx]) - ++mcsRingSize[ringIdx]; + for (size_t ringIdx : bondRingMembership[bondIdx]) ++mcsRingSize[ringIdx]; } for (auto bond = mcsIdx.Bonds.begin(); bond != mcsIdx.Bonds.end(); ++bond, ++bi) { unsigned int bondIdx = (*bond)->getIdx(); if (!ri->numBondRings(bondIdx)) continue; haveFusedBondQuery = true; - for (int ringIdx: bondRingMembership[bondIdx]) { + for (int ringIdx : bondRingMembership[bondIdx]) { if (mcsRingSize[ringIdx] < br[ringIdx].size()) continue; rwMol->getBondWithIdx(bi)->expandQuery( - makeBondInRingOfSizeQuery(br[ringIdx].size()), Queries::COMPOSITE_AND, true); + makeBondInRingOfSizeQuery(br[ringIdx].size()), Queries::COMPOSITE_AND, + true); } } return haveFusedBondQuery; @@ -855,15 +856,16 @@ MCSResult MaximumCommonSubgraph::find(const std::vector& src_mols) { res.Canceled = growSeeds() ? false : true; // verify what MCS is equal to one of initial seed for chirality match if ((FinalMatchCheckFunction == Parameters.FinalMatchChecker && - 1 == getMaxNumberBonds()) || 0 == getMaxNumberBonds()) { + 1 == getMaxNumberBonds()) || + 0 == getMaxNumberBonds()) { McsIdx = MCS(); // clear makeInitialSeeds(); // check all possible initial seeds if (!Seeds.empty()) { const Seed& fs = Seeds.front(); - if (1 == getMaxNumberBonds() - || !(Parameters.BondCompareParameters.CompleteRingsOnly - && fs.MoleculeFragment.Bonds.size() == 1 - && queryIsBondInRing(fs.MoleculeFragment.Bonds.front()))) { + if (1 == getMaxNumberBonds() || + !(Parameters.BondCompareParameters.CompleteRingsOnly && + fs.MoleculeFragment.Bonds.size() == 1 && + queryIsBondInRing(fs.MoleculeFragment.Bonds.front()))) { McsIdx.QueryMolecule = QueryMolecule; McsIdx.Atoms = fs.MoleculeFragment.Atoms; McsIdx.Bonds = fs.MoleculeFragment.Bonds; @@ -885,8 +887,8 @@ MCSResult MaximumCommonSubgraph::find(const std::vector& src_mols) { res.NumBonds = getMaxNumberBonds(); ; if (res.NumBonds > 0) { - std::pair smartsQueryMolPair = - generateResultSMARTSAndQueryMol(McsIdx); + std::pair smartsQueryMolPair = + generateResultSMARTSAndQueryMol(McsIdx); res.SmartsString = smartsQueryMolPair.first; res.QueryMol = RWMOL_SPTR(smartsQueryMolPair.second); if (Parameters.BondCompareParameters.MatchFusedRingsStrict) @@ -900,7 +902,8 @@ MCSResult MaximumCommonSubgraph::find(const std::vector& src_mols) { res.NumAtoms > 0 && tag != Targets.end(); tag++, itarget++) { MatchVectType match; - bool target_matched = SubstructMatch(*tag->Molecule, *res.QueryMol, match); + bool target_matched = + SubstructMatch(*tag->Molecule, *res.QueryMol, match); if (!target_matched) std::cout << "Target " << itarget + 1 << (target_matched ? " matched " : " MISMATCHED ") @@ -1257,9 +1260,9 @@ bool MaximumCommonSubgraph::matchIncrementalFast(Seed& seed, unsigned itarget) { } // CHIRALITY: FinalMatchCheck if (matched && Parameters.FinalMatchChecker) { - std::vector c1; + std::vector c1; c1.reserve(4096); - std::vector c2; + std::vector c2; c2.reserve(4096); for (unsigned si = 0; si < seed.getNumAtoms(); si++) { // index in the seed topology diff --git a/Code/GraphMol/FMCS/MaximumCommonSubgraph.h b/Code/GraphMol/FMCS/MaximumCommonSubgraph.h index 72c2a5ba3..474c7325e 100755 --- a/Code/GraphMol/FMCS/MaximumCommonSubgraph.h +++ b/Code/GraphMol/FMCS/MaximumCommonSubgraph.h @@ -25,16 +25,14 @@ namespace RDKit { -inline bool FinalChiralityCheckFunction(const short unsigned c1[], - const short unsigned c2[], const ROMol& mol1, - const FMCS::Graph& query, const ROMol& mol2, - const FMCS::Graph& target, - const MCSParameters* p); +inline bool FinalChiralityCheckFunction( + const std::uint32_t c1[], const std::uint32_t c2[], const ROMol& mol1, + const FMCS::Graph& query, const ROMol& mol2, const FMCS::Graph& target, + const MCSParameters* p); -bool FinalMatchCheckFunction(const short unsigned c1[], - const short unsigned c2[], const ROMol& mol1, - const FMCS::Graph& query, const ROMol& mol2, - const FMCS::Graph& target, +bool FinalMatchCheckFunction(const std::uint32_t c1[], const std::uint32_t c2[], + const ROMol& mol1, const FMCS::Graph& query, + const ROMol& mol2, const FMCS::Graph& target, const MCSParameters* p); namespace FMCS { @@ -98,8 +96,9 @@ class RDKIT_FMCS_EXPORT MaximumCommonSubgraph { void makeInitialSeeds(); bool createSeedFromMCS(size_t newQueryTarget, Seed& seed); bool growSeeds(); // returns false if canceled - std::pair generateResultSMARTSAndQueryMol(const MCS& mcsIdx) const; - bool addFusedBondQueries(const MCS& McsIdx, RWMol *rwMol) const; + std::pair generateResultSMARTSAndQueryMol( + const MCS& mcsIdx) const; + bool addFusedBondQueries(const MCS& McsIdx, RWMol* rwMol) const; bool match(Seed& seed); bool matchIncrementalFast(Seed& seed, unsigned itarget); diff --git a/Code/GraphMol/FileParsers/MaeMolSupplier.cpp b/Code/GraphMol/FileParsers/MaeMolSupplier.cpp index 2461bf6e1..f28c4356a 100644 --- a/Code/GraphMol/FileParsers/MaeMolSupplier.cpp +++ b/Code/GraphMol/FileParsers/MaeMolSupplier.cpp @@ -391,7 +391,8 @@ MaeMolSupplier::MaeMolSupplier(const std::string &fileName, bool sanitize, bool removeHs) { df_owner = true; auto *ifs = new std::ifstream(fileName.c_str(), std::ios_base::binary); - if (!ifs || !(*ifs) || ifs->bad()) { + if (!(*ifs) || ifs->bad()) { + delete ifs; std::ostringstream errout; errout << "Bad input file " << fileName; throw BadFileException(errout.str()); diff --git a/Code/GraphMol/FileParsers/SmilesMolSupplier.cpp b/Code/GraphMol/FileParsers/SmilesMolSupplier.cpp index 5b5300a27..c63d86dd6 100644 --- a/Code/GraphMol/FileParsers/SmilesMolSupplier.cpp +++ b/Code/GraphMol/FileParsers/SmilesMolSupplier.cpp @@ -38,7 +38,7 @@ SmilesMolSupplier::SmilesMolSupplier(const std::string &fileName, // Need to check if this has been fixed in VC++ 7.0 auto *tmpStream = new std::ifstream(fileName.c_str(), std::ios_base::binary); - if (!tmpStream || (!(*tmpStream)) || (tmpStream->bad())) { + if (!(*tmpStream) || tmpStream->bad()) { std::ostringstream errout; errout << "Bad input file " << fileName; delete tmpStream; diff --git a/Code/GraphMol/FileParsers/TplFileParser.cpp b/Code/GraphMol/FileParsers/TplFileParser.cpp index f1c89cf4a..85b1c1af3 100644 --- a/Code/GraphMol/FileParsers/TplFileParser.cpp +++ b/Code/GraphMol/FileParsers/TplFileParser.cpp @@ -1,4 +1,3 @@ -// $Id$ // // Copyright (C) 2007-2010 Greg Landrum // @@ -146,6 +145,7 @@ Conformer *ParseConfData(std::istream *inStream, unsigned int &line, RWMol *mol, boost::split(splitLine, tempStr, boost::is_any_of(" \t"), boost::token_compress_on); if (splitLine.size() < 3) { + delete conf; std::ostringstream errout; errout << "Did not find enough fields on line " << line << " while reading conformer " << confId << std::endl; @@ -266,7 +266,7 @@ RWMol *TPLDataStreamToMol(std::istream *inStream, unsigned int &line, throw FileParseException("Found a non-blank line between conformers."); } } - if (sanitize && res) { + if (sanitize) { MolOps::sanitizeMol(*res); } @@ -307,7 +307,7 @@ RWMol *TPLFileToMol(const std::string &fName, bool sanitize, //------------------------------------------------ RWMol *MolBlockToMol(const std::string &molBlock, bool sanitize){ std::istringstream inStream(molBlock); - RWMol *res=NULL; + RWMol *res = nullptr; unsigned int line = 0; return MolDataStreamToMol(inStream, line, sanitize); } diff --git a/Code/GraphMol/Fingerprints/test1.cpp b/Code/GraphMol/Fingerprints/test1.cpp index 36d39caf5..6658c8cd4 100644 --- a/Code/GraphMol/Fingerprints/test1.cpp +++ b/Code/GraphMol/Fingerprints/test1.cpp @@ -2321,6 +2321,39 @@ void testMACCS() { delete m1; delete fp1; } + { + // check bits 124 and 130 part 1 + std::string smi = "CNNCOOC"; + RWMol *m1 = SmilesToMol(smi); + TEST_ASSERT(m1); + ExplicitBitVect *fp1 = MACCSFingerprints::getFingerprintAsBitVect(*m1); + TEST_ASSERT((*fp1)[124]); + TEST_ASSERT((*fp1)[130]); + delete m1; + delete fp1; + } + { + // check bits 124 and 130 part 2 + std::string smi = "CNNCC"; + RWMol *m1 = SmilesToMol(smi); + TEST_ASSERT(m1); + ExplicitBitVect *fp1 = MACCSFingerprints::getFingerprintAsBitVect(*m1); + TEST_ASSERT((*fp1)[124]); + TEST_ASSERT(!(*fp1)[130]); + delete m1; + delete fp1; + } + { + // check bits 124 and 130 part 3 + std::string smi = "CNCCC"; + RWMol *m1 = SmilesToMol(smi); + TEST_ASSERT(m1); + ExplicitBitVect *fp1 = MACCSFingerprints::getFingerprintAsBitVect(*m1); + TEST_ASSERT(!(*fp1)[124]); + TEST_ASSERT(!(*fp1)[130]); + delete m1; + delete fp1; + } BOOST_LOG(rdInfoLog) << "done" << std::endl; } diff --git a/Code/GraphMol/ForceFieldHelpers/MMFF/AtomTyper.cpp b/Code/GraphMol/ForceFieldHelpers/MMFF/AtomTyper.cpp index f145c46bb..6347b71a8 100644 --- a/Code/GraphMol/ForceFieldHelpers/MMFF/AtomTyper.cpp +++ b/Code/GraphMol/ForceFieldHelpers/MMFF/AtomTyper.cpp @@ -45,18 +45,18 @@ const MMFFBondCollection *getMMFFBond() { static MMFFBondCollection ds_instance; return &ds_instance; } - + const MMFFChgCollection *getMMFFChg() { static MMFFChgCollection ds_instance; return &ds_instance; } - + const MMFFDefCollection *getMMFFDef() { static MMFFDefCollection ds_instance; return &ds_instance; } -const MMFFHerschbachLaurieCollection * getMMFFHerschbachLaurie() { +const MMFFHerschbachLaurieCollection *getMMFFHerschbachLaurie() { static MMFFHerschbachLaurieCollection ds_instance; return &ds_instance; } @@ -71,7 +71,6 @@ const MMFFAngleCollection *getMMFFAngle() { return &ds_instance; } - const MMFFStbnCollection *getMMFFStbn() { static MMFFStbnCollection ds_instance; return &ds_instance; @@ -90,13 +89,13 @@ const MMFFCovRadPauEleCollection *getMMFFCovRadPauEle() { const MMFFTorCollection *getMMFFTor(const bool isMMFFs) { static MMFFTorCollection MMFF94(false, ""); static MMFFTorCollection MMFF94s(true, ""); - return (isMMFFs) ? &MMFF94s : &MMFF94; + return (isMMFFs) ? &MMFF94s : &MMFF94; } const MMFFOopCollection *getMMFFOop(const bool isMMFFs) { static MMFFOopCollection MMFF94(false, ""); static MMFFOopCollection MMFF94s(true, ""); - return (isMMFFs) ? &MMFF94s : &MMFF94; + return (isMMFFs) ? &MMFF94s : &MMFF94; } const MMFFVdWCollection *getMMFFVdW() { @@ -104,37 +103,35 @@ const MMFFVdWCollection *getMMFFVdW() { return &ds_instance; } -} +} // namespace DefaultParameters class RingMembership { public: - RingMembership() : d_isInAromaticRing(false) {}; - bool getIsInAromaticRing() const { - return d_isInAromaticRing; - } + RingMembership() : d_isInAromaticRing(false){}; + bool getIsInAromaticRing() const { return d_isInAromaticRing; } void setIsInAromaticRing(bool isInAromaticRing) { d_isInAromaticRing = isInAromaticRing; } - const std::set& getRingIdxSet() const { - return d_ringIdxSet; - } - std::set& getRingIdxSet() { - return d_ringIdxSet; - } + const std::set &getRingIdxSet() const { return d_ringIdxSet; } + std::set &getRingIdxSet() { return d_ringIdxSet; } + private: bool d_isInAromaticRing; std::set d_ringIdxSet; }; class RingMembershipSize { - typedef std::map RingMembershipMap; - typedef std::map RingSizeMembershipMap; + typedef std::map RingMembershipMap; + typedef std::map RingSizeMembershipMap; public: static const std::uint32_t IS_AROMATIC_BIT; RingMembershipSize(const ROMol &mol); - bool isAtomInAromaticRingOfSize(const Atom *atom, const unsigned int ringSize) const; + bool isAtomInAromaticRingOfSize(const Atom *atom, + const unsigned int ringSize) const; bool areAtomsInSameAromaticRing(const Atom *atom1, const Atom *atom2) const; - bool areAtomsInSameRingOfSize(const unsigned int ringSize, const unsigned int numAtoms, ...) const; + bool areAtomsInSameRingOfSize(const unsigned int ringSize, + const unsigned int numAtoms, ...) const; + private: RingSizeMembershipMap d_ringSizeMembershipMap; }; @@ -150,19 +147,19 @@ RingMembershipSize::RingMembershipSize(const ROMol &mol) { unsigned int ringSize = atomRings[ringIdx].size(); std::uint32_t ringIdxWithAromaticFlag = ringIdx; bool ringIsAromatic = isRingAromatic(mol, atomRings[ringIdx]); - if (ringIsAromatic) - ringIdxWithAromaticFlag |= IS_AROMATIC_BIT; + if (ringIsAromatic) ringIdxWithAromaticFlag |= IS_AROMATIC_BIT; auto it = d_ringSizeMembershipMap.find(ringSize); if (it == d_ringSizeMembershipMap.end()) - it = d_ringSizeMembershipMap.insert(std::make_pair(ringSize, RingMembershipMap())).first; + it = d_ringSizeMembershipMap + .insert(std::make_pair(ringSize, RingMembershipMap())) + .first; for (int atomIdxIt : atomRings[ringIdx]) { auto it2 = it->second.find(atomIdxIt); if (it2 == it->second.end()) it2 = it->second.insert(std::make_pair(atomIdxIt, RingMembership())) .first; it2->second.getRingIdxSet().insert(ringIdxWithAromaticFlag); - if (ringIsAromatic) - it2->second.setIsInAromaticRing(true); + if (ringIsAromatic) it2->second.setIsInAromaticRing(true); } } } @@ -175,15 +172,14 @@ bool RingMembershipSize::isAtomInAromaticRingOfSize( if (isAromatic) { auto it2 = it->second.find(atom->getIdx()); isAromatic = (it2 != it->second.end()); - if (isAromatic) - isAromatic = it2->second.getIsInAromaticRing(); + if (isAromatic) isAromatic = it2->second.getIsInAromaticRing(); } - + return isAromatic; } -bool RingMembershipSize::areAtomsInSameAromaticRing( - const Atom *atom1, const Atom *atom2) const { +bool RingMembershipSize::areAtomsInSameAromaticRing(const Atom *atom1, + const Atom *atom2) const { bool areInSameAromaticRing = false; for (auto it = d_ringSizeMembershipMap.begin(); @@ -193,10 +189,13 @@ bool RingMembershipSize::areAtomsInSameAromaticRing( auto it2 = it->second.find(atom2->getIdx()); if ((it1 != it->second.end()) && (it2 != it->second.end())) { std::set_intersection(it1->second.getRingIdxSet().begin(), - it1->second.getRingIdxSet().end(), it2->second.getRingIdxSet().begin(), - it2->second.getRingIdxSet().end(), std::back_inserter(intersectVect)); - for (std::vector::const_iterator ivIt = intersectVect.begin(); - !areInSameAromaticRing && (ivIt != intersectVect.end()); ++ ivIt) + it1->second.getRingIdxSet().end(), + it2->second.getRingIdxSet().begin(), + it2->second.getRingIdxSet().end(), + std::back_inserter(intersectVect)); + for (std::vector::const_iterator ivIt = + intersectVect.begin(); + !areInSameAromaticRing && (ivIt != intersectVect.end()); ++ivIt) areInSameAromaticRing = *ivIt & IS_AROMATIC_BIT; } } @@ -204,8 +203,9 @@ bool RingMembershipSize::areAtomsInSameAromaticRing( return areInSameAromaticRing; } -bool RingMembershipSize::areAtomsInSameRingOfSize( - const unsigned int ringSize, const unsigned int numAtoms, ...) const { +bool RingMembershipSize::areAtomsInSameRingOfSize(const unsigned int ringSize, + const unsigned int numAtoms, + ...) const { va_list atoms; bool areInSameRingOfSize = false; @@ -223,8 +223,9 @@ bool RingMembershipSize::areAtomsInSameRingOfSize( if (it2 == it->second.end()) break; std::set intersect; std::set_intersection(commonSet.begin(), commonSet.end(), - it2->second.getRingIdxSet().begin(), it2->second.getRingIdxSet().end(), - std::inserter(intersect, intersect.end())); + it2->second.getRingIdxSet().begin(), + it2->second.getRingIdxSet().end(), + std::inserter(intersect, intersect.end())); commonSet = intersect; areInSameRingOfSize = !commonSet.empty(); } @@ -282,8 +283,8 @@ unsigned int getPeriodicTableRowHL(const int atomicNum) { // given the MMFF atom type, this function returns true // if it is aromatic bool isAromaticAtomType(const unsigned int atomType) { - static const unsigned int aromatic_array[] = {37, 38, 39, 44, 58, 59, 63, 64, 65, - 66, 69, 76, 78, 79, 80, 81, 82}; + static const unsigned int aromatic_array[] = { + 37, 38, 39, 44, 58, 59, 63, 64, 65, 66, 69, 76, 78, 79, 80, 81, 82}; const std::set aromaticTypes( aromatic_array, aromatic_array + sizeof(aromatic_array) / sizeof(aromatic_array[0])); @@ -294,8 +295,8 @@ bool isAromaticAtomType(const unsigned int atomType) { bool isRingAromatic(const ROMol &mol, const INT_VECT &ringIndxVect) { bool isAromatic = true; for (unsigned int i = 0; isAromatic && (i < ringIndxVect.size() - 1); ++i) - isAromatic = (mol.getBondBetweenAtoms( - ringIndxVect[i], ringIndxVect[i + 1])->getBondType() == Bond::AROMATIC); + isAromatic = (mol.getBondBetweenAtoms(ringIndxVect[i], ringIndxVect[i + 1]) + ->getBondType() == Bond::AROMATIC); return isAromatic; } @@ -663,8 +664,8 @@ void setMMFFAromaticity(RWMol &mol) { } // sets the MMFF atomType for a heavy atom -void MMFFMolProperties::setMMFFHeavyAtomType( - const RingMembershipSize &rmSize, const Atom *atom) { +void MMFFMolProperties::setMMFFHeavyAtomType(const RingMembershipSize &rmSize, + const Atom *atom) { unsigned int atomType = 0; unsigned int i; unsigned int j; @@ -800,8 +801,7 @@ void MMFFMolProperties::setMMFFHeavyAtomType( } // if there are neither alpha nor beta heteroatoms // or if there are both, but they belong to different rings - if (((!(alphaHet.size())) && (!(betaHet.size()))) || - (alphaHet.size() && betaHet.size())) { + if (alphaHet.size() == betaHet.size()) { bool surroundedByBenzeneC = true; bool surroundedByArom = true; // loop over neighbors @@ -879,8 +879,7 @@ void MMFFMolProperties::setMMFFHeavyAtomType( break; } if ((atom->getTotalDegree() == 3) && - ((alphaHet.size() && (!(betaHet.size()))) || - (betaHet.size() && (!(alphaHet.size()))))) { + alphaHet.size() != betaHet.size()) { // NIM+ // Aromatic nitrogen in imidazolium // N5A+ @@ -1713,7 +1712,7 @@ void MMFFMolProperties::setMMFFHeavyAtomType( unsigned int nObondedToCorNorS = 0; unsigned int nSbondedToCorNorS = 0; bool isOxideOBondedToH = - atom->getNumExplicitHs() + atom->getNumImplicitHs(); + atom->getNumExplicitHs() + atom->getNumImplicitHs() > 0; bool isCarboxylateO = false; bool isCarbonylO = false; bool isOxideOBondedToC = false; @@ -1728,11 +1727,11 @@ void MMFFMolProperties::setMMFFHeavyAtomType( // loop over neighbors boost::tie(nbrIdx, endNbrs) = mol.getAtomNeighbors(atom); for (; (nbrIdx != endNbrs) && (!isOxideOBondedToC) && - (!isOxideOBondedToN) && (!isOxideOBondedToH) && - (!isCarboxylateO) && (!isNitroO) && (!isNOxideO) && - (!isThioSulfinateO) && (!isSulfateO) && - (!isPhosphateOrPerchlorateO) && (!isCarbonylO) && - (!isNitrosoO) && (!isSulfoxideO); + (!isOxideOBondedToN) && (!isOxideOBondedToH) && + (!isCarboxylateO) && (!isNitroO) && (!isNOxideO) && + (!isThioSulfinateO) && (!isSulfateO) && + (!isPhosphateOrPerchlorateO) && (!isCarbonylO) && + (!isNitrosoO) && (!isSulfoxideO); ++nbrIdx) { const Atom *nbrAtom = mol[*nbrIdx]; const Bond *bond = @@ -2542,7 +2541,8 @@ MMFFMolProperties::MMFFMolProperties(ROMol &mol, const std::string &mmffVariant, "A T O M T Y P E S A N D C H A R G E S\n\n" " ATOM FORMAL PARTIAL\n" " ATOM TYPE CHARGE CHARGE\n" - "-----------------------------------" << std::endl; + "-----------------------------------" + << std::endl; for (idx = 0; idx < mol.getNumAtoms(); ++idx) { oStream << std::left << std::setw(2) << mol.getAtomWithIdx(idx)->getSymbol() << std::left << " #" @@ -2737,9 +2737,9 @@ MMFFMolProperties::getMMFFBondStretchEmpiricalRuleParams(const ROMol &mol, const MMFFCovRadPauEle *mmffAtomCovRadPauEleParams[2]; const MMFFBndkCollection *mmffBndk = DefaultParameters::getMMFFBndk(); const MMFFHerschbachLaurieCollection *mmffHerschbachLaurie = - DefaultParameters::getMMFFHerschbachLaurie(); + DefaultParameters::getMMFFHerschbachLaurie(); const MMFFCovRadPauEleCollection *mmffCovRadPauEle = - DefaultParameters::getMMFFCovRadPauEle(); + DefaultParameters::getMMFFCovRadPauEle(); const MMFFPropCollection *mmffProp = DefaultParameters::getMMFFProp(); unsigned int atomicNum1 = bond->getBeginAtom()->getAtomicNum(); @@ -3692,7 +3692,8 @@ bool MMFFMolProperties::getMMFFAngleBendParams(const ROMol &mol, atomType[i] = getMMFFAtomType(idx[i]); } const MMFFAngle *mmffAngleParams = - (*mmffAngle)(DefaultParameters::getMMFFDef(), angleType, atomType[0], atomType[1], atomType[2]); + (*mmffAngle)(DefaultParameters::getMMFFDef(), angleType, atomType[0], + atomType[1], atomType[2]); const MMFFProp *mmffPropParamsCentralAtom = (*mmffProp)(atomType[1]); if ((!mmffAngleParams) || (isDoubleZero(mmffAngleParams->ka))) { areMMFFAngleParamsEmpirical = true; @@ -3798,9 +3799,9 @@ bool MMFFMolProperties::getMMFFTorsionParams( getMMFFTorsionType(mol, idx1, idx2, idx3, idx4); bool areMMFFTorParamsEmpirical = false; const std::pair mmffTorPair = - mmffTor->getMMFFTorParams(DefaultParameters::getMMFFDef(), - torTypePair, atomType[0], atomType[1], - atomType[2], atomType[3]); + mmffTor->getMMFFTorParams(DefaultParameters::getMMFFDef(), torTypePair, + atomType[0], atomType[1], atomType[2], + atomType[3]); torsionType = (mmffTorPair.first ? mmffTorPair.first : torTypePair.first); const MMFFTor *mmffTorParams = mmffTorPair.second; if (!mmffTorParams) { @@ -3841,8 +3842,8 @@ bool MMFFMolProperties::getMMFFOopBendParams(const ROMol &mol, atomType[i] = getMMFFAtomType(idx[i]); } const MMFFOop *mmffOopParams = - (*mmffOop)(DefaultParameters::getMMFFDef(), - atomType[0], atomType[1], atomType[2], atomType[3]); + (*mmffOop)(DefaultParameters::getMMFFDef(), atomType[0], atomType[1], + atomType[2], atomType[3]); // if no parameters could be found, we exclude this term (SURDOX02) if (mmffOopParams) { mmffOopBendParams = *mmffOopParams; @@ -3870,12 +3871,13 @@ bool MMFFMolProperties::getMMFFVdWParams(const unsigned int idx1, mmffVdWParamsJAtom); mmffVdWParams.R_ij_star = mmffVdWParams.R_ij_starUnscaled; mmffVdWParams.epsilon = mmffVdWParams.epsilonUnscaled; - MMFF::Utils::scaleVdWParams(mmffVdWParams.R_ij_star, mmffVdWParams.epsilon, - mmffVdW, mmffVdWParamsIAtom, mmffVdWParamsJAtom); + MMFF::Utils::scaleVdWParams(mmffVdWParams.R_ij_star, + mmffVdWParams.epsilon, mmffVdW, + mmffVdWParamsIAtom, mmffVdWParamsJAtom); res = true; } } return res; } -} -} +} // namespace MMFF +} // namespace RDKit diff --git a/Code/GraphMol/ForceFieldHelpers/UFF/Builder.cpp b/Code/GraphMol/ForceFieldHelpers/UFF/Builder.cpp index 12e328219..487dc585d 100644 --- a/Code/GraphMol/ForceFieldHelpers/UFF/Builder.cpp +++ b/Code/GraphMol/ForceFieldHelpers/UFF/Builder.cpp @@ -67,7 +67,7 @@ void setTwoBitCell(boost::shared_array &res, unsigned int pos, } std::uint8_t getTwoBitCell(boost::shared_array &res, - unsigned int pos) { + unsigned int pos) { unsigned int twoBitPos = pos / 4; unsigned int shift = 2 * (pos % 4); std::uint8_t twoBitMask = 3 << shift; @@ -89,16 +89,17 @@ std::uint8_t getTwoBitCell(boost::shared_array &res, // ------------------------------------------------------------------------ boost::shared_array buildNeighborMatrix(const ROMol &mol) { const std::uint8_t RELATION_1_X_INIT = RELATION_1_X | (RELATION_1_X << 2) | - (RELATION_1_X << 4) | - (RELATION_1_X << 6); + (RELATION_1_X << 4) | + (RELATION_1_X << 6); unsigned int nAtoms = mol.getNumAtoms(); unsigned nTwoBitCells = (nAtoms * (nAtoms + 1) - 1) / 8 + 1; boost::shared_array res(new std::uint8_t[nTwoBitCells]); std::memset(res.get(), RELATION_1_X_INIT, nTwoBitCells); for (ROMol::ConstBondIterator bondi = mol.beginBonds(); bondi != mol.endBonds(); ++bondi) { - setTwoBitCell(res, twoBitCellPos(nAtoms, (*bondi)->getBeginAtomIdx(), - (*bondi)->getEndAtomIdx()), + setTwoBitCell(res, + twoBitCellPos(nAtoms, (*bondi)->getBeginAtomIdx(), + (*bondi)->getEndAtomIdx()), RELATION_1_2); unsigned int bondiBeginAtomIdx = (*bondi)->getBeginAtomIdx(); unsigned int bondiEndAtomIdx = (*bondi)->getEndAtomIdx(); @@ -182,10 +183,8 @@ void addAngles(const ROMol &mol, const AtomicParamVect ¶ms, // if the central atom and one of the bonded atoms, but not the // other one are inside a ring, then this angle is between a // ring substituent and a ring edge - if ((rings->isAtomInRingOfSize(i, 3) && - !rings->isAtomInRingOfSize(k, 3)) || - (!rings->isAtomInRingOfSize(i, 3) && - rings->isAtomInRingOfSize(k, 3))) { + if (rings->isAtomInRingOfSize(i, 3) != + rings->isAtomInRingOfSize(k, 3)) { order = 30; } // if all atoms are inside the ring, then this is one of ring @@ -200,10 +199,8 @@ void addAngles(const ROMol &mol, const AtomicParamVect ¶ms, // if the central atom and one of the bonded atoms, but not the // other one are inside a ring, then this angle is between a // ring substituent and a ring edge - if ((rings->isAtomInRingOfSize(i, 4) && - !rings->isAtomInRingOfSize(k, 4)) || - (!rings->isAtomInRingOfSize(i, 4) && - rings->isAtomInRingOfSize(k, 4))) { + if (rings->isAtomInRingOfSize(i, 4) != + rings->isAtomInRingOfSize(k, 4)) { order = 40; } // if all atoms are inside the ring, then this is one of ring @@ -718,5 +715,5 @@ ForceFields::ForceField *constructForceField(ROMol &mol, double vdwThresh, return constructForceField(mol, params, vdwThresh, confId, ignoreInterfragInteractions); } -} -} +} // namespace UFF +} // namespace RDKit diff --git a/Code/GraphMol/ForceFieldHelpers/Wrap/rdForceFields.cpp b/Code/GraphMol/ForceFieldHelpers/Wrap/rdForceFields.cpp index bc6d5f309..96461badc 100644 --- a/Code/GraphMol/ForceFieldHelpers/Wrap/rdForceFields.cpp +++ b/Code/GraphMol/ForceFieldHelpers/Wrap/rdForceFields.cpp @@ -71,11 +71,13 @@ int FFHelper(ForceFields::ForceField &ff, int maxIters) { return ForceFieldsHelper::OptimizeMolecule(ff, maxIters).first; } -python::object FFConfsHelper(ROMol &mol, ForceFields::ForceField &ff, int numThreads, int maxIters) { +python::object FFConfsHelper(ROMol &mol, ForceFields::ForceField &ff, + int numThreads, int maxIters) { std::vector> res; { NOGIL gil; - ForceFieldsHelper::OptimizeMoleculeConfs(mol, ff, res, numThreads, maxIters); + ForceFieldsHelper::OptimizeMoleculeConfs(mol, ff, res, numThreads, + maxIters); } python::list pyres; for (auto &itm : res) { @@ -132,8 +134,10 @@ ForceFields::PyMMFFMolProperties *GetMMFFMolProperties( new MMFF::MMFFMolProperties(mol, mmffVariant, mmffVerbosity); ForceFields::PyMMFFMolProperties *pyMP = nullptr; - if (mmffMolProperties && mmffMolProperties->isValid()) { + if (mmffMolProperties->isValid()) { pyMP = new ForceFields::PyMMFFMolProperties(mmffMolProperties); + } else { + delete mmffMolProperties; } return pyMP; @@ -153,9 +157,7 @@ ForceFields::PyForceField *MMFFGetMoleculeForceField( MMFF::constructForceField(mol, mmffMolProperties, nonBondedThresh, confId, ignoreInterfragInteractions); pyFF = new ForceFields::PyForceField(ff); - if (pyFF) { - pyFF->initialize(); - } + pyFF->initialize(); } return pyFF; @@ -167,7 +169,7 @@ bool MMFFHasAllMoleculeParams(const ROMol &mol) { return mmffMolProperties.isValid(); } -}; +}; // namespace RDKit namespace ForceFields { PyObject *getUFFBondStretchParams(const RDKit::ROMol &mol, @@ -236,7 +238,7 @@ PyObject *getUFFVdWParams(const RDKit::ROMol &mol, const unsigned int idx1, } return res; }; -} +} // namespace ForceFields BOOST_PYTHON_MODULE(rdForceFieldHelpers) { python::scope().attr("__doc__") = @@ -477,6 +479,6 @@ RETURNS: a list of (not_converged, energy) 2-tuples. \n\ \n"; python::def("OptimizeMoleculeConfs", RDKit::FFConfsHelper, (python::arg("mol"), python::arg("ff"), - python::arg("numThreads") = 1, python::arg("maxIters") = 200), + python::arg("numThreads") = 1, python::arg("maxIters") = 200), docString.c_str()); } diff --git a/Code/GraphMol/FragCatalog/FragCatGenerator.cpp b/Code/GraphMol/FragCatalog/FragCatGenerator.cpp index c8eec788c..d9b95f85f 100644 --- a/Code/GraphMol/FragCatalog/FragCatGenerator.cpp +++ b/Code/GraphMol/FragCatalog/FragCatGenerator.cpp @@ -1,4 +1,3 @@ -// $Id$ // // Copyright (C) 2003-2006 Rational Discovery LLC // @@ -54,6 +53,7 @@ unsigned int addOrder1Paths(PATH_LIST &paths, const ROMol &mol, invar = computeIntVectPrimesProduct(*pi); mapkm1[invar] = (*eti); delete nent; + nent = nullptr; break; } } @@ -178,6 +178,7 @@ unsigned int addHigherOrderPaths(const INT_PATH_LIST_MAP &allPaths, found = true; mEntId = (*iti); delete nent; + nent = nullptr; break; } } @@ -251,4 +252,4 @@ unsigned int FragCatGenerator::addFragsFromMol(const ROMol &mol, delete coreMol; return (nO1Pths + nremPths); } -} +} // namespace RDKit diff --git a/Code/GraphMol/FragCatalog/FragCatalogUtils.cpp b/Code/GraphMol/FragCatalog/FragCatalogUtils.cpp index 33c2ca05c..8f82b085c 100644 --- a/Code/GraphMol/FragCatalog/FragCatalogUtils.cpp +++ b/Code/GraphMol/FragCatalog/FragCatalogUtils.cpp @@ -1,4 +1,3 @@ -// $Id$ // // Copyright (C) 2003-2006 Rational Discovery LLC // @@ -73,7 +72,8 @@ MOL_SPTR_VECT readFuncGroups(std::istream &inStream, int nToRead) { int nRead = 0; MOL_SPTR_VECT funcGroups; - while (!inStream.eof() && (nToRead < 0 || nRead < nToRead)) { + while (!inStream.eof() && !inStream.fail() && + (nToRead < 0 || nRead < nToRead)) { std::string tmpstr; std::getline(inStream, tmpstr); // parse the molecule on this line (if there is one) diff --git a/Code/GraphMol/FragCatalog/FragFPGenerator.cpp b/Code/GraphMol/FragCatalog/FragFPGenerator.cpp index ed0f894ed..d9a36cc34 100644 --- a/Code/GraphMol/FragCatalog/FragFPGenerator.cpp +++ b/Code/GraphMol/FragCatalog/FragFPGenerator.cpp @@ -1,4 +1,3 @@ -// $Id$ // // Copyright (C) 2003-2006 Rational Discovery LLC // @@ -82,10 +81,10 @@ void FragFPGenerator::computeFP(const ROMol &mol, const FragCatalog &fcat, fp->setBit(bitId); } mapkm1[invar] = (*eti); - delete nent; break; } } + delete nent; } // now deal with the higher order stuff. @@ -170,10 +169,10 @@ void FragFPGenerator::computeFP(const ROMol &mol, const FragCatalog &fcat, if (bitId >= 0) { fp->setBit(bitId); } - delete nent; break; } } + delete nent; } // overwrite mapkm1 with mapk before we move on to order k+1 @@ -235,4 +234,4 @@ void FragFPGenerator::computeFP(const ROMol &mol, const FragCatalog &fcat, **************************/ #endif } -} +} // namespace RDKit diff --git a/Code/GraphMol/MolAlign/O3AAlignMolecules.cpp b/Code/GraphMol/MolAlign/O3AAlignMolecules.cpp index c1b79a4c1..b1390c8fb 100644 --- a/Code/GraphMol/MolAlign/O3AAlignMolecules.cpp +++ b/Code/GraphMol/MolAlign/O3AAlignMolecules.cpp @@ -1095,14 +1095,13 @@ double o3aMMFFWeightFunc(const unsigned int prbIdx, const unsigned int refIdx, return static_cast(O3_MAX_WEIGHT_COEFF - (static_cast(data))->weight) * - ((1.0 + - O3_CHARGE_COEFF * - fabs((static_cast( - (static_cast(data))->refProp)) - ->getMMFFPartialCharge(refIdx) + - (static_cast( - (static_cast(data))->prbProp)) - ->getMMFFPartialCharge(prbIdx))) / + ((1.0 + O3_CHARGE_COEFF * + fabs((static_cast( + (static_cast(data))->refProp)) + ->getMMFFPartialCharge(refIdx) + + (static_cast( + (static_cast(data))->prbProp)) + ->getMMFFPartialCharge(prbIdx))) / (1.0 + fabs((static_cast( (static_cast(data))->refProp)) ->getMMFFPartialCharge(refIdx) - @@ -1221,8 +1220,8 @@ O3A::O3A(int (*costFunc)(const unsigned int, const unsigned int, double, double (*weightFunc)(const unsigned int, const unsigned int, void *), double (*scoringFunc)(const unsigned int, const unsigned int, void *), void *data, ROMol &prbMol, const ROMol &refMol, const int prbCid, - const int refCid, boost::dynamic_bitset<> *prbHvyAtoms, - boost::dynamic_bitset<> *refHvyAtoms, const bool reflect, + const int refCid, const boost::dynamic_bitset<> &prbHvyAtoms, + const boost::dynamic_bitset<> &refHvyAtoms, const bool reflect, const unsigned int maxIters, unsigned int options, O3AConstraintVect *o3aConstraintVect, ROMol *extWorkPrbMol, LAP *extLAP, MolHistogram *extPrbHist, MolHistogram *extRefHist) @@ -1243,21 +1242,6 @@ O3A::O3A(int (*costFunc)(const unsigned int, const unsigned int, double, unsigned int prbNHeavyAtoms = prbMol.getNumHeavyAtoms(); unsigned int largestNHeavyAtoms = std::max(refNHeavyAtoms, prbNHeavyAtoms); unsigned int i; - if (!refHvyAtoms) { - refHvyAtoms = new boost::dynamic_bitset<>(refNAtoms); - for (i = 0; i < refNAtoms; ++i) { - if (refMol[i]->getAtomicNum() != 1) { - refHvyAtoms->set(i); - } - } - } - if (!prbHvyAtoms) { - for (i = 0; i < prbNAtoms; ++i) { - if (prbMol[i]->getAtomicNum() != 1) { - prbHvyAtoms->set(i); - } - } - } std::vector pairs(4, 0); std::vector score(3, 0.0); std::vector pairsRMSD(2, 0.0); @@ -1267,18 +1251,20 @@ O3A::O3A(int (*costFunc)(const unsigned int, const unsigned int, double, MolHistogram *prbHist = nullptr; LAP *lap = nullptr; if (local) { - startSDM.fillFromDist(O3_SDM_THRESHOLD_START, *refHvyAtoms, *prbHvyAtoms); + startSDM.fillFromDist(O3_SDM_THRESHOLD_START, refHvyAtoms, prbHvyAtoms); } else { - refHist = (extRefHist ? extRefHist - : new MolHistogram( - refMol, MolOps::get3DDistanceMat( - refMol, refCid, false, false, ""), - true)); - prbHist = (extPrbHist ? extPrbHist - : new MolHistogram( - prbMol, MolOps::get3DDistanceMat( - prbMol, prbCid, false, false, ""), - true)); + refHist = + (extRefHist ? extRefHist + : new MolHistogram(refMol, + MolOps::get3DDistanceMat( + refMol, refCid, false, false, ""), + true)); + prbHist = + (extPrbHist ? extPrbHist + : new MolHistogram(prbMol, + MolOps::get3DDistanceMat( + prbMol, prbCid, false, false, ""), + true)); lap = (extLAP ? extLAP : new LAP(largestNHeavyAtoms)); lap->computeCostMatrix(prbMol, *prbHist, refMol, *refHist, o3aConstraintVect, costFunc, data); @@ -1303,7 +1289,7 @@ O3A::O3A(int (*costFunc)(const unsigned int, const unsigned int, double, (double)sdmThresholdIt * O3_SDM_THRESHOLD_STEP; while (flag && (iter < O3_MAX_SDM_ITERATIONS)) { SDM progressSDM(&prbConf, &refConf, o3aConstraintVect); - progressSDM.fillFromDist(sdmThresholdDist, *refHvyAtoms, *prbHvyAtoms); + progressSDM.fillFromDist(sdmThresholdDist, refHvyAtoms, prbHvyAtoms); pairs[3] = progressSDM.size(); if (pairs[3] < 3) { break; @@ -1380,16 +1366,16 @@ O3A::O3A(int (*costFunc)(const unsigned int, const unsigned int, double, delete bestSDM[i]; } } - if ((!extLAP) && lap) { + if (!extLAP) { delete lap; } - if ((!extRefHist) && refHist) { + if (!extRefHist) { delete refHist; } - if ((!extPrbHist) && prbHist) { + if (!extPrbHist) { delete prbHist; } - if ((!extWorkPrbMol) && workPrbMol) { + if (!extWorkPrbMol) { delete workPrbMol; } } @@ -1476,16 +1462,18 @@ O3A::O3A(ROMol &prbMol, const ROMol &refMol, void *prbProp, void *refProp, MolHistogram *prbHist = nullptr; LAP *lap = nullptr; if (!local) { - refHist = (extRefHist ? extRefHist - : new MolHistogram( - refMol, MolOps::get3DDistanceMat( - refMol, refCid, false, false, ""), - true)); - prbHist = (extPrbHist ? extPrbHist - : new MolHistogram( - prbMol, MolOps::get3DDistanceMat( - prbMol, prbCid, false, false, ""), - true)); + refHist = + (extRefHist ? extRefHist + : new MolHistogram(refMol, + MolOps::get3DDistanceMat( + refMol, refCid, false, false, ""), + true)); + prbHist = + (extPrbHist ? extPrbHist + : new MolHistogram(prbMol, + MolOps::get3DDistanceMat( + prbMol, prbCid, false, false, ""), + true)); lap = (extLAP ? extLAP : new LAP(largestNHeavyAtoms)); } for (l = 0, score[0] = 0.0; @@ -1496,7 +1484,7 @@ O3A::O3A(ROMol &prbMol, const ROMol &refMol, void *prbProp, void *refProp, data.weight = (l ? c : 0); data.coeff = c; auto *o3a = new O3A(costFunc, weightFunc, scoringFunc, &data, prbMol, - refMol, prbCid, refCid, &prbHvyAtoms, &refHvyAtoms, + refMol, prbCid, refCid, prbHvyAtoms, refHvyAtoms, reflect, maxIters, options, &o3aConstraintVect, &extWorkPrbMol, lap, prbHist, refHist); score[1] = o3a->score(); @@ -1681,7 +1669,7 @@ void O3AHelper_(ROMol *prbMol, const ROMol *refMol, void *prbProp, (*res)[i].reset(lres); } } -} // end of detail namespace +} // namespace detail #endif void getO3AForProbeConfs(ROMol &prbMol, const ROMol &refMol, void *prbProp, diff --git a/Code/GraphMol/MolAlign/O3AAlignMolecules.h b/Code/GraphMol/MolAlign/O3AAlignMolecules.h index 1faf90ea2..5b429980d 100644 --- a/Code/GraphMol/MolAlign/O3AAlignMolecules.h +++ b/Code/GraphMol/MolAlign/O3AAlignMolecules.h @@ -1,4 +1,3 @@ -// $Id$ // // Copyright (C) 2013-2014 Paolo Tosco // @@ -208,6 +207,7 @@ class RDKIT_MOLALIGN_EXPORT SDM { }; // assignment operator SDM &operator=(const SDM &other) { + if (this == &other) return *this; d_prbConf = other.d_prbConf; d_refConf = other.d_refConf; d_o3aConstraintVect = other.d_o3aConstraintVect; @@ -286,8 +286,8 @@ class RDKIT_MOLALIGN_EXPORT O3A { double (*weightFunc)(const unsigned int, const unsigned int, void *), double (*scoringFunc)(const unsigned int, const unsigned int, void *), void *data, ROMol &prbMol, const ROMol &refMol, const int prbCid, - const int refCid, boost::dynamic_bitset<> *prbHvyAtoms = NULL, - boost::dynamic_bitset<> *refHvyAtoms = NULL, const bool reflect = false, + const int refCid, const boost::dynamic_bitset<> &prbHvyAtoms, + const boost::dynamic_bitset<> &refHvyAtoms, const bool reflect = false, const unsigned int maxIters = 50, unsigned int options = 0, O3AConstraintVect *o3aConstraintVect = NULL, ROMol *extWorkPrbMol = NULL, LAP *extLAP = NULL, MolHistogram *extPrbHist = NULL, diff --git a/Code/GraphMol/MolCatalog/MolCatalogEntry.cpp b/Code/GraphMol/MolCatalog/MolCatalogEntry.cpp index 22712dafd..48fc486db 100644 --- a/Code/GraphMol/MolCatalog/MolCatalogEntry.cpp +++ b/Code/GraphMol/MolCatalog/MolCatalogEntry.cpp @@ -26,6 +26,7 @@ MolCatalogEntry::MolCatalogEntry(const ROMol *omol) { dp_props = new Dict(); d_descrip = ""; dp_mol = omol; + d_order = 0; } MolCatalogEntry::MolCatalogEntry(const MolCatalogEntry &other) { @@ -39,6 +40,7 @@ MolCatalogEntry::MolCatalogEntry(const MolCatalogEntry &other) { if (other.dp_mol) { dp_mol = new ROMol(*other.dp_mol); } + d_order = other.d_order; } MolCatalogEntry::~MolCatalogEntry() { @@ -119,4 +121,4 @@ void MolCatalogEntry::initFromString(const std::string &text) { // now start reading out values: initFromStream(ss); } -} +} // namespace RDKit diff --git a/Code/GraphMol/MolCatalog/MolCatalogEntry.h b/Code/GraphMol/MolCatalog/MolCatalogEntry.h index 1c7cfe62d..6439eb659 100644 --- a/Code/GraphMol/MolCatalog/MolCatalogEntry.h +++ b/Code/GraphMol/MolCatalog/MolCatalogEntry.h @@ -16,7 +16,7 @@ class ROMol; //! This class is used to store ROMol objects in a MolCatalog class RDKIT_MOLCATALOG_EXPORT MolCatalogEntry : public RDCatalog::CatalogEntry { public: - MolCatalogEntry() : dp_mol(0), d_descrip("") { + MolCatalogEntry() : dp_mol(0), d_order(0), d_descrip("") { dp_props = new Dict(); setBitId(-1); } diff --git a/Code/GraphMol/MolChemicalFeatures/FeatureParser.cpp b/Code/GraphMol/MolChemicalFeatures/FeatureParser.cpp index ef5258ddf..c399ce24a 100644 --- a/Code/GraphMol/MolChemicalFeatures/FeatureParser.cpp +++ b/Code/GraphMol/MolChemicalFeatures/FeatureParser.cpp @@ -19,7 +19,7 @@ #include #include #include -typedef boost::tokenizer > tokenizer; +typedef boost::tokenizer> tokenizer; #include #include @@ -27,14 +27,14 @@ typedef boost::tokenizer > tokenizer; namespace RDKit { namespace Local { -typedef boost::tokenizer > CommaTokenizer; +typedef boost::tokenizer> CommaTokenizer; void getNextLine(std::istream &inStream, std::string &line, unsigned int &lineNo) { - if (inStream.eof()) return; + if (inStream.eof() || inStream.fail()) return; line = ""; bool continuationLine = false; - while (!inStream.eof()) { + while (!inStream.eof() && !inStream.fail()) { std::string tmpLine; std::getline(inStream, tmpLine); lineNo++; @@ -90,8 +90,8 @@ void parseAtomType(const std::string &inLine, std::map &atomTypeDefs, const unsigned int &lineNo) { boost::char_separator sep(" \t"); - boost::tokenizer > tok(inLine, sep); - boost::tokenizer >::iterator tokIt = tok.begin(); + boost::tokenizer> tok(inLine, sep); + boost::tokenizer>::iterator tokIt = tok.begin(); if (tokIt == tok.end()) { throw FeatureFileParseException(lineNo, inLine, "empty input line for AtomType"); @@ -167,12 +167,12 @@ MolChemicalFeatureDef *parseFeatureDef( } boost::char_separator sep(" \t"); - boost::tokenizer > tok(nextLine, sep); + boost::tokenizer> tok(nextLine, sep); if (tok.begin() == tok.end()) { throw FeatureFileParseException(lineNo, inLine, "bad DefineFeature line, no tokens found"); } - boost::tokenizer >::iterator tokIt = tok.begin(); + boost::tokenizer>::iterator tokIt = tok.begin(); tokIt++; if (tokIt == tok.end()) { throw FeatureFileParseException(lineNo, inLine, @@ -281,7 +281,7 @@ int parseFeatureData(std::istream &inStream, std::string inLine; Local::getNextLine(inStream, inLine, lineNo); std::map atomTypeDefs; - while (!inStream.eof()) { + while (!inStream.eof() && !inStream.fail()) { // clean any whitespace off the line: boost::trim_if(inLine, boost::is_any_of(" \t\r\n")); if (!inLine.empty() && inLine[0] != '#' && inLine[0] != '\n') { @@ -314,4 +314,4 @@ int parseFeatureFile(const std::string &fileName, } return parseFeatureData(inStream, res); } -} +} // namespace RDKit diff --git a/Code/GraphMol/MolDraw2D/MolDraw2D.cpp b/Code/GraphMol/MolDraw2D/MolDraw2D.cpp index aa9a0ce68..af7f036ea 100644 --- a/Code/GraphMol/MolDraw2D/MolDraw2D.cpp +++ b/Code/GraphMol/MolDraw2D/MolDraw2D.cpp @@ -1635,18 +1635,6 @@ pair MolDraw2D::getAtomSymbolAndOrientation( orient = MolDraw2D::W; } } - // last check: degree zero atoms from the last three periods should have - // the Hs first - if (!atom.getDegree()) { - static int HsListedFirstSrc[] = {8, 9, 16, 17, 34, 35, 52, 53, 84, 85}; - std::vector HsListedFirst( - HsListedFirstSrc, - HsListedFirstSrc + sizeof(HsListedFirstSrc) / sizeof(int)); - if (std::find(HsListedFirst.begin(), HsListedFirst.end(), - atom.getAtomicNum()) != HsListedFirst.end()) { - orient = MolDraw2D::W; - } - } if (orient == MolDraw2D::W) { preText.push_back(h); } else { diff --git a/Code/GraphMol/MolHash/MolHash.cpp b/Code/GraphMol/MolHash/MolHash.cpp index b4046ef5a..49d7befb6 100644 --- a/Code/GraphMol/MolHash/MolHash.cpp +++ b/Code/GraphMol/MolHash/MolHash.cpp @@ -22,15 +22,15 @@ namespace RDKit { namespace MolHash { struct MolFragment // Reference to a fragment of source molecule - { +{ std::vector Atoms; std::vector Bonds; std::vector AtomsIdx; std::vector BondsIdx; std::map MolAtomIdxMap; // Full Molecule to - // fragment indices - // backward - // conversion map + // fragment indices + // backward + // conversion map public: std::uint32_t getNumAtoms() const { return AtomsIdx.size(); } std::uint32_t getNumBonds() const { return BondsIdx.size(); } @@ -251,7 +251,8 @@ std::string generateMoleculeHashSet(const ROMol &mol, // snprintf(buf, sizeof(buf),"%u-%u-%u-", res.Version, // res.NumAtoms,res.NumBonds); str = (boost::format("%u-%u-%u-") % (res.Version) % (res.NumAtoms) % - (res.NumBonds)).str(); + (res.NumBonds)) + .str(); str += encode(&res.FormulaCRC32, sizeof(res.FormulaCRC32)); str += "-"; str += encode(&res.NonChiralAtomsHash, sizeof(res.NonChiralAtomsHash)); @@ -320,9 +321,8 @@ static void prepareMolFragment(MolFragment &m, const ROMol &mol, n = mol.getNumBonds(); m.BondsIdx.resize(n); for (unsigned i = 0; i < n; i++) m.BondsIdx[i] = i; - } else if (nullptr != - atomsToUse) // selected atoms only and all/selected bonds - // between them + } else if (nullptr != atomsToUse) // selected atoms only and all/selected + // bonds between them { std::map addedBonds; unsigned n = atomsToUse->size(); @@ -375,11 +375,11 @@ static void prepareMolFragment(MolFragment &m, const ROMol &mol, unsigned endAtoms[2]; endAtoms[0] = bond->getBeginAtomIdx(); endAtoms[1] = bond->getEndAtomIdx(); - for (unsigned ai = 0; ai < 2 && nullptr != bond; - ai++) // both ending bonds of the atom + for (unsigned ai = 0; ai < 2; ai++) // both ending bonds of the atom { if (addedAtoms.end() == addedAtoms.find(endAtoms[ai])) { - addedAtoms[endAtoms[ai]] = addedAtoms.size(); + auto val = addedAtoms.size(); + addedAtoms[endAtoms[ai]] = val; m.AtomsIdx.push_back( endAtoms[ai]); // the atom has NOT been already added } @@ -422,5 +422,5 @@ static void prepareLabels(std::vector &atomLabels, } } //============================================================================= -} -} +} // namespace MolHash +} // namespace RDKit diff --git a/Code/GraphMol/MolHash/catch_tests.cpp b/Code/GraphMol/MolHash/catch_tests.cpp index db45cdb76..6da28eb07 100644 --- a/Code/GraphMol/MolHash/catch_tests.cpp +++ b/Code/GraphMol/MolHash/catch_tests.cpp @@ -185,26 +185,37 @@ TEST_CASE("Basic MolHash", "[molhash]") { } } - TEST_CASE("Tautomers and chirality", "[molhash]") { SECTION("basics") { auto om = "C[C@H](C(=O)O)C(=O)[O-]"_smiles; REQUIRE(om); -{ std::unique_ptr m(new RWMol(*om)); - auto hsh = - MolHash::MolHash(m.get(), MolHash::HashFunction::CanonicalSmiles); - CHECK(hsh == "C[C@@H](C(=O)[O-])C(=O)O"); -} -{ std::unique_ptr m(new RWMol(*om)); - auto hsh = - MolHash::MolHash(m.get(), MolHash::HashFunction::HetAtomTautomer); - CHECK(hsh == "CC([C]([O])[O])[C]([O])[O]_1_-1"); -} -{ std::unique_ptr m(new RWMol(*om)); - auto hsh = - MolHash::MolHash(m.get(), MolHash::HashFunction::HetAtomProtomer); - CHECK(hsh == "CC([C]([O])[O])[C]([O])[O]_2"); -} + { + std::unique_ptr m(new RWMol(*om)); + auto hsh = + MolHash::MolHash(m.get(), MolHash::HashFunction::CanonicalSmiles); + CHECK(hsh == "C[C@@H](C(=O)[O-])C(=O)O"); + } + { + std::unique_ptr m(new RWMol(*om)); + auto hsh = + MolHash::MolHash(m.get(), MolHash::HashFunction::HetAtomTautomer); + CHECK(hsh == "CC([C]([O])[O])[C]([O])[O]_1_-1"); + } + { + std::unique_ptr m(new RWMol(*om)); + auto hsh = + MolHash::MolHash(m.get(), MolHash::HashFunction::HetAtomProtomer); + CHECK(hsh == "CC([C]([O])[O])[C]([O])[O]_2"); + } + } +} + +TEST_CASE("Molecular formula with fragments", "[molhash]") { + SECTION("basics") { + auto om = "CC(=O)[O-].C[N+](C)(C)C"_smiles; + REQUIRE(om); + auto hsh = MolHash::MolHash(om.get(), MolHash::HashFunction::MolFormula); + CHECK(hsh == "C6H15NO2"); + } } -} \ No newline at end of file diff --git a/Code/GraphMol/MolHash/hashfunctions.cpp b/Code/GraphMol/MolHash/hashfunctions.cpp index 3a7550683..4a56fd891 100644 --- a/Code/GraphMol/MolHash/hashfunctions.cpp +++ b/Code/GraphMol/MolHash/hashfunctions.cpp @@ -24,6 +24,7 @@ namespace { unsigned int NMRDKitBondGetOrder(const RDKit::Bond *bnd) { + PRECONDITION(bnd, "bad bond"); switch (bnd->getBondType()) { case RDKit::Bond::AROMATIC: case RDKit::Bond::SINGLE: @@ -46,6 +47,9 @@ unsigned int NMRDKitBondGetOrder(const RDKit::Bond *bnd) { RDKit::Bond *NMRDKitMolNewBond(RDKit::RWMol *mol, RDKit::Atom *src, RDKit::Atom *dst, unsigned int order, bool arom) { + PRECONDITION(mol, "bad molecule"); + PRECONDITION(src, "bad src atom"); + PRECONDITION(dst, "bad dest atom"); RDKit::Bond *result; result = mol->getBondBetweenAtoms(src->getIdx(), dst->getIdx()); if (result) { @@ -92,6 +96,7 @@ RDKit::Bond *NMRDKitMolNewBond(RDKit::RWMol *mol, RDKit::Atom *src, } void NMRDKitSanitizeHydrogens(RDKit::RWMol *mol) { + PRECONDITION(mol, "bad molecule"); // Move all of the implicit Hs into one box for (auto aptr : mol->atoms()) { unsigned int hcount = aptr->getTotalNumHs(); @@ -107,6 +112,8 @@ namespace RDKit { namespace MolHash { static unsigned int NMDetermineComponents(RWMol *mol, unsigned int *parts, unsigned int acount) { + PRECONDITION(mol, "bad molecule"); + PRECONDITION(parts, "bad parts pointer"); memset(parts, 0, acount * sizeof(unsigned int)); std::vector todo; @@ -137,6 +144,8 @@ static unsigned int NMDetermineComponents(RWMol *mol, unsigned int *parts, static std::string NMMolecularFormula(RWMol *mol, const unsigned int *parts, unsigned int part) { + PRECONDITION(mol, "bad molecule"); + PRECONDITION((!part || parts), "bad parts pointer"); unsigned int hist[256]; int charge = 0; @@ -187,6 +196,7 @@ static std::string NMMolecularFormula(RWMol *mol, const unsigned int *parts, } static std::string NMMolecularFormula(RWMol *mol, bool sep = false) { + PRECONDITION(mol, "bad molecule"); if (!sep) return NMMolecularFormula(mol, 0, 0); unsigned int acount = mol->getNumAtoms(); @@ -215,6 +225,7 @@ static std::string NMMolecularFormula(RWMol *mol, bool sep = false) { } static void NormalizeHCount(Atom *aptr) { + PRECONDITION(aptr, "bad atom pointer"); unsigned int hcount; switch (aptr->getAtomicNum()) { @@ -256,6 +267,7 @@ static void NormalizeHCount(Atom *aptr) { } static std::string AnonymousGraph(RWMol *mol, bool elem) { + PRECONDITION(mol, "bad molecule"); std::string result; int charge = 0; @@ -281,6 +293,7 @@ static std::string AnonymousGraph(RWMol *mol, bool elem) { } static std::string MesomerHash(RWMol *mol, bool netq) { + PRECONDITION(mol, "bad molecule"); std::string result; char buffer[32]; int charge = 0; @@ -305,6 +318,7 @@ static std::string MesomerHash(RWMol *mol, bool netq) { } static std::string TautomerHash(RWMol *mol, bool proto) { + PRECONDITION(mol, "bad molecule"); std::string result; char buffer[32]; int hcount = 0; @@ -334,9 +348,9 @@ static std::string TautomerHash(RWMol *mol, bool proto) { MolOps::assignRadicals(*mol); // we may have just destroyed some stereocenters/bonds // clean that up: - bool cleanIt=true; - bool force=true; - MolOps::assignStereochemistry(*mol,cleanIt,force); + bool cleanIt = true; + bool force = true; + MolOps::assignStereochemistry(*mol, cleanIt, force); result = MolToSmiles(*mol); if (!proto) { sprintf(buffer, "_%d_%d", hcount, charge); @@ -347,6 +361,8 @@ static std::string TautomerHash(RWMol *mol, bool proto) { } static bool TraverseForRing(Atom *atom, unsigned char *visit) { + PRECONDITION(atom, "bad atom pointer"); + PRECONDITION(visit, "bad pointer"); visit[atom->getIdx()] = 1; for (auto nbri : boost::make_iterator_range( atom->getOwningMol().getAtomNeighbors(atom))) { @@ -362,6 +378,9 @@ static bool TraverseForRing(Atom *atom, unsigned char *visit) { static bool DepthFirstSearchForRing(Atom *root, Atom *nbor, unsigned int maxatomidx) { + PRECONDITION(root, "bad atom pointer"); + PRECONDITION(nbor, "bad atom pointer"); + unsigned int natoms = maxatomidx; unsigned char *visit = (unsigned char *)alloca(natoms); memset(visit, 0, natoms); @@ -371,6 +390,7 @@ static bool DepthFirstSearchForRing(Atom *root, Atom *nbor, } bool IsInScaffold(Atom *atom, unsigned int maxatomidx) { + PRECONDITION(atom, "bad atom pointer"); if (RDKit::queryIsAtomInRing(atom)) return true; unsigned int count = 0; @@ -383,6 +403,8 @@ bool IsInScaffold(Atom *atom, unsigned int maxatomidx) { } static bool HasNbrInScaffold(Atom *aptr, unsigned char *is_in_scaffold) { + PRECONDITION(aptr, "bad atom pointer"); + PRECONDITION(is_in_scaffold, "bad pointer"); for (auto nbri : boost::make_iterator_range( aptr->getOwningMol().getAtomNeighbors(aptr))) { auto nptr = aptr->getOwningMol()[nbri]; @@ -392,6 +414,7 @@ static bool HasNbrInScaffold(Atom *aptr, unsigned char *is_in_scaffold) { } static std::string ExtendedMurckoScaffold(RWMol *mol) { + PRECONDITION(mol, "bad molecule"); RDKit::MolOps::fastFindRings(*mol); unsigned int maxatomidx = mol->getNumAtoms(); @@ -424,6 +447,7 @@ static std::string ExtendedMurckoScaffold(RWMol *mol) { } static std::string MurckoScaffoldHash(RWMol *mol) { + PRECONDITION(mol, "bad molecule"); std::vector for_deletion; do { for_deletion.clear(); @@ -454,6 +478,7 @@ static std::string MurckoScaffoldHash(RWMol *mol) { } static std::string NetChargeHash(RWMol *mol) { + PRECONDITION(mol, "bad molecule"); int totalq = 0; for (auto aptr : mol->atoms()) { @@ -466,6 +491,7 @@ static std::string NetChargeHash(RWMol *mol) { } static std::string SmallWorldHash(RWMol *mol, bool brl) { + PRECONDITION(mol, "bad molecule"); char buffer[64]; unsigned int acount = mol->getNumAtoms(); @@ -505,6 +531,7 @@ static void DegreeVector(RWMol *mol, unsigned int *v) { } static bool HasDoubleBond(Atom *atom) { + PRECONDITION(atom, "bad atom"); for (const auto &nbri : boost::make_iterator_range(atom->getOwningMol().getAtomBonds(atom))) { auto bptr = (atom->getOwningMol())[nbri]; @@ -521,6 +548,7 @@ static bool HasDoubleBond(Atom *atom) { // 3 means break, with asterisks on both beg and end static int RegioisomerBond(Bond *bnd) { + PRECONDITION(bnd, "bad bond"); if (NMRDKitBondGetOrder(bnd) != 1) return -1; if (RDKit::queryIsBondInRing(bnd)) return -1; @@ -544,6 +572,7 @@ static int RegioisomerBond(Bond *bnd) { } static void ClearEZStereo(Atom *atm) { + PRECONDITION(atm, "bad atom"); for (const auto &nbri : boost::make_iterator_range(atm->getOwningMol().getAtomBonds(atm))) { auto bptr = (atm->getOwningMol())[nbri]; @@ -553,6 +582,8 @@ static void ClearEZStereo(Atom *atm) { } static std::string RegioisomerHash(RWMol *mol) { + PRECONDITION(mol, "bad molecule"); + // we need a copy of the molecule so that we can loop over the bonds of // something while modifying something else RDKit::MolOps::fastFindRings(*mol); @@ -597,6 +628,7 @@ static std::string RegioisomerHash(RWMol *mol) { } static std::string ArthorSubOrderHash(RWMol *mol) { + PRECONDITION(mol, "bad molecule"); char buffer[256]; unsigned int acount = mol->getNumAtoms(); @@ -681,6 +713,7 @@ static std::string ArthorSubOrderHash(RWMol *mol) { } std::string MolHash(RWMol *mol, HashFunction func) { + PRECONDITION(mol, "bad molecule"); std::string result; char buffer[32]; NMRDKitSanitizeHydrogens(mol); diff --git a/Code/GraphMol/MolInterchange/Parser.cpp b/Code/GraphMol/MolInterchange/Parser.cpp index aea4e6e6e..d30930f58 100644 --- a/Code/GraphMol/MolInterchange/Parser.cpp +++ b/Code/GraphMol/MolInterchange/Parser.cpp @@ -145,8 +145,10 @@ void readAtom(RWMol *mol, const rj::Value &atomVal, at->setNumRadicalElectrons(getIntDefaultValue("nRad", atomVal, atomDefaults)); at->setIsotope(getIntDefaultValue("isotope", atomVal, atomDefaults)); std::string stereo = getStringDefaultValue("stereo", atomVal, atomDefaults); - if (chilookup.find(stereo) == chilookup.end()) + if (chilookup.find(stereo) == chilookup.end()) { + delete at; throw FileParseException("Bad Format: bad stereo value for atom"); + } at->setChiralTag(chilookup.find(stereo)->second); bool updateLabel = false, takeOwnership = true; mol->addAtom(at, updateLabel, takeOwnership); @@ -235,7 +237,8 @@ void readPartialCharges(RWMol *mol, const rj::Value &repVal, throw FileParseException( "Bad Format: size of values array != num atoms"); for (unsigned int idx = 0; idx != mol->getNumAtoms(); ++idx) { - const auto &val = miter->value.GetArray()[idx]; + const auto &aval = miter->value.GetArray(); + const auto &val = aval[idx]; if (!val.IsDouble()) throw FileParseException("Bad Format: partial charge not double"); mol->getAtomWithIdx(idx)->setProp(common_properties::_GasteigerCharge, diff --git a/Code/GraphMol/MolPickler.cpp b/Code/GraphMol/MolPickler.cpp index c7c4ee17d..7c08b99a1 100644 --- a/Code/GraphMol/MolPickler.cpp +++ b/Code/GraphMol/MolPickler.cpp @@ -511,12 +511,12 @@ Query *unpickleQuery(std::istream &ss, ROMol *tmpMol; switch (tag) { case MolPickler::QUERY_ATOMRING: - res = new AtomRingQuery(); streamRead(ss, tag, version); if (tag != MolPickler::QUERY_VALUE) { throw MolPicklerException( "Bad pickle format: QUERY_VALUE tag not found."); } + res = new AtomRingQuery(); streamRead(ss, val, version); static_cast *>(res)->setVal(val); streamRead(ss, val, version); @@ -723,12 +723,12 @@ AtomMonomerInfo *unpickleAtomMonomerInfo(std::istream &ss, int version) { switch (typ) { case AtomMonomerInfo::UNKNOWN: case AtomMonomerInfo::OTHER: - res = - new AtomMonomerInfo(RDKit::AtomMonomerInfo::AtomMonomerType(typ), nm); streamRead(ss, tag, version); if (tag != MolPickler::END_ATOM_MONOMER) throw MolPicklerException( "did not find expected end of atom monomer info"); + res = + new AtomMonomerInfo(RDKit::AtomMonomerInfo::AtomMonomerType(typ), nm); break; case AtomMonomerInfo::PDBRESIDUE: res = static_cast(new AtomPDBResidueInfo(nm)); @@ -945,9 +945,6 @@ void MolPickler::_pickle(const ROMol *mol, std::ostream &ss, const Conformer *conf = ci->get(); _pickleProperties(ss, *conf, propertyFlags); } - } - - if (propertyFlags & PicklerOps::MolProps) { streamWrite(ss, BEGINPROPS); _pickleProperties(ss, *mol, propertyFlags); streamWrite(ss, ENDPROPS); diff --git a/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogEntry.h b/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogEntry.h index 41ea9a8ac..e58374657 100644 --- a/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogEntry.h +++ b/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogEntry.h @@ -25,8 +25,6 @@ class RDKIT_MOLSTANDARDIZE_EXPORT AcidBaseCatalogEntry : public RDCatalog::CatalogEntry { public: AcidBaseCatalogEntry() { - dp_pair->first = nullptr; - dp_pair->second = nullptr; d_descrip = ""; dp_props = new Dict(); setBitId(-1); diff --git a/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogUtils.cpp b/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogUtils.cpp index 6ce9ce457..c2628b058 100644 --- a/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogUtils.cpp +++ b/Code/GraphMol/MolStandardize/AcidBaseCatalog/AcidBaseCatalogUtils.cpp @@ -87,7 +87,7 @@ std::vector> readPairs(std::istream& inStream, char inLine[MAX_LINE_LEN]; std::string tmpstr; int nRead = 0; - while (!inStream.eof() && (nToRead < 0 || nRead < nToRead)) { + while (!inStream.eof() && !inStream.fail() && (nToRead < 0 || nRead < nToRead)) { inStream.getline(inLine, MAX_LINE_LEN, '\n'); tmpstr = inLine; // parse the molpair on this line (if there is one) diff --git a/Code/GraphMol/MolStandardize/Charge.cpp b/Code/GraphMol/MolStandardize/Charge.cpp index b74d2c28a..c8b0370b3 100644 --- a/Code/GraphMol/MolStandardize/Charge.cpp +++ b/Code/GraphMol/MolStandardize/Charge.cpp @@ -330,7 +330,7 @@ ROMol *Uncharger::uncharge(const ROMol &mol) { for (const auto pr : a_atoms) nonAcids.reset(pr.second); unsigned int midx = 0; // zwitterion with more negative charges than quaternary positive centres - while (neg_surplus > 0 && n_matched > 0 && midx < n_atoms.size()) { + while (neg_surplus > 0 && midx < n_atoms.size()) { unsigned int idx = n_atoms[midx++].second; if (!nonAcids[idx]) continue; Atom *atom = omol->getAtomWithIdx(idx); @@ -361,7 +361,7 @@ ROMol *Uncharger::uncharge(const ROMol &mol) { if (a_matched > 0 && neg_surplus > 0) { unsigned int midx = 0; // zwitterion with more negative charges than quaternary positive centres - while (neg_surplus > 0 && a_matched > 0 && midx < a_atoms.size()) { + while (neg_surplus > 0 && midx < a_atoms.size()) { // Add hydrogen to first negative acidic atom, increase formal charge // Until quaternary positive == negative total or no more negative atoms Atom *atom = omol->getAtomWithIdx(a_atoms[midx++].second); diff --git a/Code/GraphMol/MolStandardize/FragmentCatalog/FragmentCatalogUtils.cpp b/Code/GraphMol/MolStandardize/FragmentCatalog/FragmentCatalogUtils.cpp index 6850e5a55..d58b29290 100644 --- a/Code/GraphMol/MolStandardize/FragmentCatalog/FragmentCatalogUtils.cpp +++ b/Code/GraphMol/MolStandardize/FragmentCatalog/FragmentCatalogUtils.cpp @@ -83,7 +83,8 @@ std::vector> readFuncGroups(std::istream &inStream, int nRead = 0; std::vector> funcGroups; - while (!inStream.eof() && (nToRead < 0 || nRead < nToRead)) { + while (!inStream.eof() && !inStream.fail() && + (nToRead < 0 || nRead < nToRead)) { inStream.getline(inLine, MAX_LINE_LEN, '\n'); std::string tmpstr(inLine); // parse the molecule on this line (if there is one) diff --git a/Code/GraphMol/MolStandardize/Metal.cpp b/Code/GraphMol/MolStandardize/Metal.cpp index fb0b9961f..128547b71 100644 --- a/Code/GraphMol/MolStandardize/Metal.cpp +++ b/Code/GraphMol/MolStandardize/Metal.cpp @@ -32,9 +32,9 @@ MetalDisconnector::MetalDisconnector() "W,Re,Os,Ir,Pt,Au,Hg,Tl,Pb,Bi]~[N,O,F]")), metal_non(SmartsToMol( "[Al,Sc,Ti,V,Cr,Mn,Fe,Co,Ni,Cu,Zn,Y,Zr,Nb,Mo,Tc,Ru,Rh,Pd,Ag,Cd,Hf,Ta," - "W,Re,Os,Ir,Pt,Au]~[B,C,Si,P,As,Sb,S,Se,Te,Cl,Br,I,At]")){ - BOOST_LOG(rdInfoLog) << "Initializing MetalDisconnector\n"; - }; + "W,Re,Os,Ir,Pt,Au]~[B,C,Si,P,As,Sb,S,Se,Te,Cl,Br,I,At]")) { + BOOST_LOG(rdInfoLog) << "Initializing MetalDisconnector\n"; +}; MetalDisconnector::MetalDisconnector(const MetalDisconnector &other) { metal_nof = other.metal_nof; @@ -43,16 +43,16 @@ MetalDisconnector::MetalDisconnector(const MetalDisconnector &other) { MetalDisconnector::~MetalDisconnector(){}; -ROMol *MetalDisconnector::getMetalNof() {return metal_nof.get();} +ROMol *MetalDisconnector::getMetalNof() { return metal_nof.get(); } -ROMol *MetalDisconnector::getMetalNon() {return metal_non.get();} +ROMol *MetalDisconnector::getMetalNon() { return metal_non.get(); } void MetalDisconnector::setMetalNof(const ROMol &mol) { - this->metal_nof = static_cast(new ROMol(mol)); + this->metal_nof.reset(new ROMol(mol)); } void MetalDisconnector::setMetalNon(const ROMol &mol) { - this->metal_non = static_cast(new ROMol(mol)); + this->metal_non.reset(new ROMol(mol)); } ROMol *MetalDisconnector::disconnect(const ROMol &mol) { @@ -62,7 +62,7 @@ ROMol *MetalDisconnector::disconnect(const ROMol &mol) { } void MetalDisconnector::disconnect(RWMol &mol) { - BOOST_LOG(rdInfoLog) << "Running MetalDisconnector\n"; + BOOST_LOG(rdInfoLog) << "Running MetalDisconnector\n"; std::list metalList = {metal_nof, metal_non}; for (auto &query : metalList) { std::vector matches; @@ -83,11 +83,12 @@ void MetalDisconnector::disconnect(RWMol &mol) { a1->setFormalCharge(a1->getFormalCharge() + int(order)); Atom *a2 = mol.getAtomWithIdx(non_idx); a2->setFormalCharge(a2->getFormalCharge() - int(order)); - BOOST_LOG(rdInfoLog) << "Removed covalent bond between " - << a1->getSymbol() << " and " << a2->getSymbol() << "\n"; + BOOST_LOG(rdInfoLog) << "Removed covalent bond between " + << a1->getSymbol() << " and " << a2->getSymbol() + << "\n"; } // std::cout << "After removing bond and charge adjustment: " << - //MolToSmiles(mol) << std::endl; + // MolToSmiles(mol) << std::endl; } } } // namespace MolStandardize diff --git a/Code/GraphMol/MolStandardize/TautomerCatalog/TautomerCatalogUtils.cpp b/Code/GraphMol/MolStandardize/TautomerCatalog/TautomerCatalogUtils.cpp index 10371acfa..74facb8bd 100644 --- a/Code/GraphMol/MolStandardize/TautomerCatalog/TautomerCatalogUtils.cpp +++ b/Code/GraphMol/MolStandardize/TautomerCatalog/TautomerCatalogUtils.cpp @@ -151,7 +151,7 @@ std::vector readTautomers(std::istream& inStream, char inLine[MAX_LINE_LEN]; std::string tmpstr; int nRead = 0; - while (!inStream.eof() && (nToRead < 0 || nRead < nToRead)) { + while (!inStream.eof() && !inStream.fail() && (nToRead < 0 || nRead < nToRead)) { inStream.getline(inLine, MAX_LINE_LEN, '\n'); tmpstr = inLine; // parse the molpair on this line (if there is one) diff --git a/Code/GraphMol/MolStandardize/TransformCatalog/TransformCatalogUtils.cpp b/Code/GraphMol/MolStandardize/TransformCatalog/TransformCatalogUtils.cpp index 242a0d715..080ec448f 100644 --- a/Code/GraphMol/MolStandardize/TransformCatalog/TransformCatalogUtils.cpp +++ b/Code/GraphMol/MolStandardize/TransformCatalog/TransformCatalogUtils.cpp @@ -80,7 +80,7 @@ std::vector> readTransformations( char inLine[MAX_LINE_LEN]; std::string tmpstr; int nRead = 0; - while (!inStream.eof() && (nToRead < 0 || nRead < nToRead)) { + while (!inStream.eof() && !inStream.fail() && (nToRead < 0 || nRead < nToRead)) { inStream.getline(inLine, MAX_LINE_LEN, '\n'); tmpstr = inLine; // parse the reaction on this line (if there is one) diff --git a/Code/GraphMol/Resonance.cpp b/Code/GraphMol/Resonance.cpp index 8f37d59e4..2f0391b05 100644 --- a/Code/GraphMol/Resonance.cpp +++ b/Code/GraphMol/Resonance.cpp @@ -790,22 +790,21 @@ void ConjElectrons::enumerateNonBonded(CEMap &ceMap) { void ConjElectrons::computeMetrics() { // 1000 * Electronegativity according to the Allen scale // (Allen, L.C. J. Am. Chem. Soc. 1989, 111, 9003-9014) - static const unsigned int en[] = {1000, - 2300, 4160, 912, 1576, 2051, 2544, 3066, 3610, 4193, 4789, 869, - 1293, 1613, 1916, 2253, 2589, 2869, 3242, 734, 1034, 1190, 1380, - 1530, 1650, 1750, 1800, 1840, 1880, 1850, 1590, 1756, 1994, 2211, - 2434, 2685, 2966, 706, 963, 1120, 1320, 1410, 1470, 1510, 1540, - 1560, 1590, 1870, 1520, 1656, 1824, 1984, 2158, 2359, 2582, 659, - 881, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, - 1000, 1000, 1000, 1000, 1090, 1160, 1340, 1470, 1600, 1650, 1680, - 1720, 1920, 1760, 1789, 1854, 2010, 2190, 2390, 2600, 670, 890}; - const size_t enSize = sizeof(en) / sizeof(double); + static const std::vector en{ + 1000, 2300, 4160, 912, 1576, 2051, 2544, 3066, 3610, 4193, 4789, 869, + 1293, 1613, 1916, 2253, 2589, 2869, 3242, 734, 1034, 1190, 1380, 1530, + 1650, 1750, 1800, 1840, 1880, 1850, 1590, 1756, 1994, 2211, 2434, 2685, + 2966, 706, 963, 1120, 1320, 1410, 1470, 1510, 1540, 1560, 1590, 1870, + 1520, 1656, 1824, 1984, 2158, 2359, 2582, 659, 881, 1000, 1000, 1000, + 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1000, 1090, + 1160, 1340, 1470, 1600, 1650, 1680, 1720, 1920, 1760, 1789, 1854, 2010, + 2190, 2390, 2600, 670, 890}; for (ConjAtomMap::const_iterator it = d_conjAtomMap.begin(); it != d_conjAtomMap.end(); ++it) { d_ceMetrics.d_absFormalCharges += abs(it->second->fc()); size_t anIdx = it->second->atom()->getAtomicNum(); d_ceMetrics.d_wtdFormalCharges += - (it->second->fc() * ((anIdx >= enSize) ? 1000 : en[anIdx])); + (it->second->fc() * ((anIdx >= en.size()) ? 1000 : en[anIdx])); d_ceMetrics.d_nbMissing += it->second->neededNbForOctet(); } computeDistFormalCharges(); @@ -1196,10 +1195,9 @@ void ResonanceMolSupplier::assignConjGrpIdx() { unsigned int na = d_mol->getNumAtoms(); d_atomConjGrpIdx.resize(na, -1); for (d_nConjGrp = 0; true; ++d_nConjGrp) { - for (const auto b: d_mol->bonds()) { + for (const auto b : d_mol->bonds()) { unsigned int i = b->getIdx(); - if (b->getIsConjugated() && (d_bondConjGrpIdx[i] == -1)) - { + if (b->getIsConjugated() && (d_bondConjGrpIdx[i] == -1)) { bondIndexStack.push(i); break; } @@ -1209,10 +1207,10 @@ void ResonanceMolSupplier::assignConjGrpIdx() { unsigned int i = bondIndexStack.top(); bondIndexStack.pop(); const Bond *bi = d_mol->getBondWithIdx(i); - for (const Atom *bondAtom: { bi->getBeginAtom(), bi->getEndAtom() }) { + for (const Atom *bondAtom : {bi->getBeginAtom(), bi->getEndAtom()}) { // loop over bonds sprouting from bondAtom - for (const auto &bNbri: boost::make_iterator_range( - d_mol->getAtomBonds(bondAtom))) { + for (const auto &bNbri : + boost::make_iterator_range(d_mol->getAtomBonds(bondAtom))) { const auto &bNbr = (*d_mol)[bNbri]; unsigned int biNbr = bNbr->getIdx(); if (bNbr->getIsConjugated() && (d_bondConjGrpIdx[biNbr] == -1)) { diff --git a/Code/GraphMol/SmilesParse/SmilesWrite.cpp b/Code/GraphMol/SmilesParse/SmilesWrite.cpp index d1c1d0561..36315e131 100644 --- a/Code/GraphMol/SmilesParse/SmilesWrite.cpp +++ b/Code/GraphMol/SmilesParse/SmilesWrite.cpp @@ -440,11 +440,8 @@ std::string MolToSmiles(const ROMol &mol, bool doIsomericSmiles, bool doKekule, // adding randomness without setting the rootedAtAtom if (doRandom) { if (rootedAtAtom == -1) { - rootedAtAtom = std::rand() % mol.getNumAtoms(); // need to find an atom id between 0 and mol.getNumAtoms() exclusively - PRECONDITION(rootedAtAtom < 0 || static_cast( - rootedAtAtom) < mol.getNumAtoms(), - "rootedAtomAtom must be less than the number of atoms"); + rootedAtAtom = std::rand() % mol.getNumAtoms(); } } diff --git a/Code/GraphMol/Substruct/SubstructUtils.cpp b/Code/GraphMol/Substruct/SubstructUtils.cpp index 9e733c3ad..0d7b04bce 100644 --- a/Code/GraphMol/Substruct/SubstructUtils.cpp +++ b/Code/GraphMol/Substruct/SubstructUtils.cpp @@ -29,10 +29,6 @@ bool atomCompat(const Atom *a1, const Atom *a2, res = a1->Match(a2); } return res; - std::cerr << "\t\tatomCompat: " << a1 << " " << a1->getIdx() << "-" << a2 - << " " << a2->getIdx() << std::endl; - std::cerr << "\t\t " << res << std::endl; - return res; } bool chiralAtomCompat(const Atom *&a1, const Atom *&a2) { diff --git a/Code/GraphMol/Substruct/vf2.hpp b/Code/GraphMol/Substruct/vf2.hpp index 16b91bb3d..5a8b9216b 100644 --- a/Code/GraphMol/Substruct/vf2.hpp +++ b/Code/GraphMol/Substruct/vf2.hpp @@ -1,4 +1,3 @@ -// $Id$ /* * This is an extensive modification by Greg Landrum of * pieces from several files in the vflib-2.0 distribution @@ -21,226 +20,231 @@ //#define RDK_VF2_PRUNING #define RDK_ADJ_ITER typename Graph::adjacency_iterator -namespace boost{ - namespace detail { - typedef unsigned short node_id; - const node_id NULL_NODE=0xFFFF; - struct NodeInfo { - node_id id; - node_id in; - node_id out; - }; - - template - struct Pair { - node_id n1, n2; - bool hasiter; - RDK_ADJ_ITER nbrbeg, nbrend; - - Pair() : n1(NULL_NODE), n2(NULL_NODE), hasiter(false) { - } - }; +namespace boost { +namespace detail { +typedef std::uint32_t node_id; +const node_id NULL_NODE = 0xFFFF; +struct NodeInfo { + node_id id; + node_id in; + node_id out; +}; - /** - * The ordering by in/out degree - */ - static bool nodeInfoComp1(const NodeInfo &a, const NodeInfo &b) { - if(a.out < b.out) return true; - if(a.out > b.out) return false; - if(a.in < b.in) return true; - if(a.in > b.in) return false; - return false; +template +struct Pair { + node_id n1, n2; + bool hasiter; + RDK_ADJ_ITER nbrbeg, nbrend; + + Pair() : n1(NULL_NODE), n2(NULL_NODE), hasiter(false) {} +}; + +/** + * The ordering by in/out degree + */ +static bool nodeInfoComp1(const NodeInfo &a, const NodeInfo &b) { + if (a.out < b.out) return true; + if (a.out > b.out) return false; + if (a.in < b.in) return true; + if (a.in > b.in) return false; + return false; +} + +/** + * The ordering by frequency/valence. + * The frequency is in the out field, the valence in `in'. + */ +static int nodeInfoComp2(const NodeInfo &a, const NodeInfo &b) { + if (!a.in && b.in) return 1; + if (a.in && !b.in) return -1; + if (a.out < b.out) return -1; + if (a.out > b.out) return 1; + if (a.in < b.in) return -1; + if (a.in > b.in) return 1; + return 0; +} + +template +VertexDescr getOtherIdx(const Graph &g, const EdgeDescr &edge, + const VertexDescr &vertex) { + VertexDescr tmp = boost::source(edge, g); + if (tmp == vertex) { + tmp = boost::target(edge, g); + } + return tmp; +} + +/*---------------------------------------------------- + * Sorts the nodes of a graphs, returning a + * heap-allocated vector (using new) with the node ids + * in the proper orders. + * The sorting criterion takes into account: + * 1 - The number of nodes with the same in/out + * degree. + * 2 - The valence of the nodes. + * The nodes at the beginning of the vector are + * the most singular, from which the matching should + * start. + *--------------------------------------------------*/ +template +node_id *SortNodesByFrequency(const Graph *g) { + std::vector vect; + vect.reserve(boost::num_vertices(*g)); + typename Graph::vertex_iterator bNode, eNode; + boost::tie(bNode, eNode) = boost::vertices(*g); + while (bNode != eNode) { + NodeInfo t; + t.id = vect.size(); + t.in = boost::out_degree(*bNode, *g); // <- assuming undirected graph + t.out = boost::out_degree(*bNode, *g); + vect.push_back(t); + ++bNode; + } + std::sort(vect.begin(), vect.end(), nodeInfoComp1); + + unsigned int run = 1; + for (unsigned int i = 0; i < vect.size(); i += run) { + for (run = 1; i + run < vect.size() && vect[i + run].in == vect[i].in && + vect[i + run].out == vect[i].out; + ++run) + ; + for (unsigned int j = 0; j < run; ++j) { + vect[i + j].in += vect[i + j].out; + vect[i + j].out = run; + } + } + std::sort(vect.begin(), vect.end(), nodeInfoComp2); + + node_id *nodes = new node_id[vect.size()]; + for (unsigned int i = 0; i < vect.size(); ++i) { + nodes[i] = vect[i].id; + } + + return nodes; +} + +/*---------------------------------------------------------- + * class VF2SubState + * A representation of the SSS current state + ---------------------------------------------------------*/ +template +class VF2SubState { + private: + Graph *g1, *g2; + VertexCompatible &vc; + EdgeCompatible &ec; + MatchChecking &mc; + unsigned int n1, n2; + + unsigned int core_len; + unsigned int t1_len; + unsigned int t2_len; // Core nodes are also counted by these... + node_id *core_1; + node_id *core_2; + node_id *term_1; + node_id *term_2; + + node_id *order; + + long *share_count; + int *vs_compared; + + public: + VF2SubState(Graph *ag1, Graph *ag2, VertexCompatible &avc, + EdgeCompatible &aec, MatchChecking &amc, bool sortNodes = false) + : g1(ag1), + g2(ag2), + vc(avc), + ec(aec), + mc(amc), + n1(num_vertices(*ag1)), + n2(num_vertices(*ag2)) { + if (sortNodes) { + order = SortNodesByFrequency(ag1); + } else { + order = NULL; } - /** - * The ordering by frequency/valence. - * The frequency is in the out field, the valence in `in'. - */ - static int nodeInfoComp2(const NodeInfo &a, const NodeInfo &b) { - if (!a.in && b.in ) return 1; - if (a.in && !b.in) return -1; - if (a.out < b.out) return -1; - if (a.out > b.out) return 1; - if( a.in < b.in ) return -1; - if (a.in > b.in) return 1; - return 0; + core_len = 0; + t1_len = 0; + t2_len = 0; + + core_1 = new node_id[n1]; + core_2 = new node_id[n2]; + term_1 = new node_id[n1]; + term_2 = new node_id[n2]; + share_count = new long; + + for (unsigned int i = 0; i < n1; i++) { + core_1[i] = NULL_NODE; + term_1[i] = 0; } - - template - VertexDescr getOtherIdx(const Graph &g,const EdgeDescr &edge,const VertexDescr &vertex) { - VertexDescr tmp=boost::source(edge,g); - if(tmp==vertex){ - tmp=boost::target(edge,g); - } - return tmp; + for (unsigned int i = 0; i < n2; i++) { + core_2[i] = NULL_NODE; + term_2[i] = 0; } - - /*---------------------------------------------------- - * Sorts the nodes of a graphs, returning a - * heap-allocated vector (using new) with the node ids - * in the proper orders. - * The sorting criterion takes into account: - * 1 - The number of nodes with the same in/out - * degree. - * 2 - The valence of the nodes. - * The nodes at the beginning of the vector are - * the most singular, from which the matching should - * start. - *--------------------------------------------------*/ - template - node_id* SortNodesByFrequency(const Graph *g) { - std::vector vect; - vect.reserve(boost::num_vertices(*g)); - typename Graph::vertex_iterator bNode,eNode; - boost::tie(bNode,eNode) = boost::vertices(*g); - while(bNode!=eNode){ - NodeInfo t; - t.id=vect.size(); - t.in=boost::out_degree(*bNode,*g);// <- assuming undirected graph - t.out=boost::out_degree(*bNode,*g); - vect.push_back(t); - ++bNode; - } - std::sort(vect.begin(),vect.end(),nodeInfoComp1); - - unsigned int run=1; - for(unsigned int i=0; i(); + *share_count = 1; + }; + + VF2SubState(const VF2SubState &state) + : g1(state.g1), + g2(state.g2), + vc(state.vc), + ec(state.ec), + mc(state.mc), + n1(state.n1), + n2(state.n2), + order(state.order), + vs_compared(state.vs_compared) + // es_compared(state.es_compared) + { + core_len = state.core_len; + t1_len = state.t1_len; + t2_len = state.t2_len; + + core_1 = state.core_1; + core_2 = state.core_2; + term_1 = state.term_1; + term_2 = state.term_2; + share_count = state.share_count; + + ++(*share_count); + }; + + ~VF2SubState() { + if (--*share_count == 0) { + delete[] core_1; + delete[] core_2; + delete[] term_1; + delete[] term_2; + delete share_count; + delete[] order; + // delete [] vs_compared; + // delete es_compared; } + }; - /*---------------------------------------------------------- - * class VF2SubState - * A representation of the SSS current state - ---------------------------------------------------------*/ - template - class VF2SubState - { - private: - Graph *g1, *g2; - VertexCompatible &vc; - EdgeCompatible &ec; - MatchChecking &mc; - unsigned int n1, n2; + bool IsGoal() { return core_len == n1; }; + bool MatchChecks(const node_id c1[], const node_id c2[]) { + return mc(c1, c2); + }; + bool IsDead() { return n1 > n2 || t1_len > t2_len; }; + unsigned int CoreLen() { return core_len; } + Graph *GetGraph1() { return g1; } + Graph *GetGraph2() { return g2; } - unsigned int core_len; - unsigned int t1_len; - unsigned int t2_len; // Core nodes are also counted by these... - node_id *core_1; - node_id *core_2; - node_id *term_1; - node_id *term_2; - - node_id *order; - - long *share_count; - int *vs_compared; - - public: - VF2SubState(Graph *ag1, Graph *ag2, - VertexCompatible &avc, - EdgeCompatible &aec, - MatchChecking &amc, - bool sortNodes=false) : g1(ag1), g2(ag2), vc(avc), ec(aec), mc(amc), - n1(num_vertices(*ag1)),n2(num_vertices(*ag2)) { - if (sortNodes){ - order = SortNodesByFrequency(ag1); - } else { - order = NULL; - } - - core_len=0; - t1_len=0; - t2_len=0; - - core_1=new node_id[n1]; - core_2=new node_id[n2]; - term_1=new node_id[n1]; - term_2=new node_id[n2]; - share_count = new long; - - for(unsigned int i=0; i(); - *share_count = 1; - }; - - VF2SubState(const VF2SubState &state) : - g1(state.g1), g2(state.g2), vc(state.vc), ec(state.ec), mc(state.mc), - n1(state.n1),n2(state.n2), order(state.order),vs_compared(state.vs_compared) - //es_compared(state.es_compared) - { - - core_len=state.core_len; - t1_len=state.t1_len; - t2_len=state.t2_len; - - core_1=state.core_1; - core_2=state.core_2; - term_1=state.term_1; - term_2=state.term_2; - share_count=state.share_count; - - ++(*share_count); - }; - - ~VF2SubState(){ - if (-- *share_count == 0) { - delete [] core_1; - delete [] core_2; - delete [] term_1; - delete [] term_2; - delete share_count; - delete [] order; - //delete [] vs_compared; - //delete es_compared; - } - }; - - bool IsGoal() { return core_len==n1 ; }; - bool MatchChecks(const node_id c1[],const node_id c2[]){ - return mc(c1,c2); - }; - bool IsDead() { return n1>n2 || - t1_len>t2_len; - }; - unsigned int CoreLen() { return core_len; } - Graph *GetGraph1() { return g1; } - Graph *GetGraph2() { return g2; } - - bool NextPair(Pair &pair){ - if (pair.n1==NULL_NODE) - pair.n1=0; - if (pair.n2==NULL_NODE) - pair.n2=0; - else - pair.n2++; + bool NextPair(Pair &pair) { + if (pair.n1 == NULL_NODE) pair.n1 = 0; + if (pair.n2 == NULL_NODE) + pair.n2 = 0; + else + pair.n2++; #if 0 std::cerr<<" **** np: "<< prev_n1<<","<core_len && t2_len>core_len) { - while (pair.n1core_len && t2_len>core_len) { - while (pair.n2size()<boost::out_degree(node2,*g2)) - return false; - if(!vc(node1,node2)) return false; - - unsigned int other1, other2; -#ifdef RDK_VF2_PRUNING - unsigned int term1 = 0, term2 = 0; - unsigned int new1 = 0, new2 = 0; -#endif - - // Check the out edges of node1 - typename Graph::out_edge_iterator bNbrs,eNbrs; - boost::tie(bNbrs,eNbrs) = boost::out_edges(node1,*g1); - while(bNbrs!=eNbrs){ - other1=getOtherIdx(*g1,*bNbrs,node1); - if (core_1[other1] != NULL_NODE) { - other2 = core_1[other1]; - typename Graph::edge_descriptor oEdge; - bool found; - boost::tie(oEdge,found) = boost::edge(node2,other2,*g2); - if(!found || !ec(*bNbrs,oEdge) ){ - //std::cerr<<" short2"< pair; - while (NextPair(pair)) { - if (IsFeasiblePair(pair.n1, pair.n2)){ - AddPair(pair.n1, pair.n2); - if (Match(c1, c2)) // recurse - return true; - BackTrack(pair.n1, pair.n2); - } - } - return false; + if (t1_len > core_len && t2_len > core_len) { + while (pair.n1 < n1 && + (core_1[pair.n1] != NULL_NODE || term_1[pair.n1] == 0)) { + pair.n1++; + pair.n2 = 0; } - - template - bool MatchAll(node_id c1[], node_id c2[], - DoubleBackInsertionSequence &res, unsigned int lim=0) - { - if (IsGoal()) { - GetCoreSet(c1, c2); - if(MatchChecks(c1,c2)) { - typename DoubleBackInsertionSequence::value_type newSeq; - for(unsigned int i=0;i(c1[i],c2[i])); - } - res.push_back(newSeq); - return lim && res.size() >= lim; - } - } - - if (IsDead()) - return false; - - Pair pair; - while (NextPair(pair)) { - if (IsFeasiblePair(pair.n1, pair.n2)){ - AddPair(pair.n1, pair.n2); - if (MatchAll(c1, c2, res, lim)) // recurse - return true; - BackTrack(pair.n1, pair.n2); - } - } - return false; - } - }; - /*------------------------------------------------------------- - * static bool match(pn, c1, c2, s) - * Finds a matching between two graphs, if it exists, starting - * from state s. - * Returns true a match has been found. - * *pn is assigned the numbero of matched nodes, and - * c1 and c2 will contain the ids of the corresponding nodes - * in the two graphs. - ------------------------------------------------------------*/ - template - bool match(int *pn, node_id c1[], node_id c2[], SubState &s) - { - if (s.Match(c1, c2)) { - // not needed, pn = num query atoms (n1)... - *pn=s.CoreLen(); - return true; + /* Initialize VF2 Plus neighbor iterator. + * The next query node (pair.n1) has been selected from the terminal + * set and is therefore adjacent to an already mapped atom (in + * core_1). Rather than select pair.n2 from all atoms (0...n2) we can + * select it from the neighbors of this mapped atom (0...deg(nbor)) + * since it must also be adajcent to this mapped atom! + */ + if (!pair.hasiter) { + RDK_ADJ_ITER n1iter_beg, n1iter_end; + boost::tie(n1iter_beg, n1iter_end) = + boost::adjacent_vertices(pair.n1, *g1); + + while (n1iter_beg != n1iter_end && core_1[*n1iter_beg] == NULL_NODE) + ++n1iter_beg; + + assert(n1iter_beg != n1iter_end); + + boost::tie(pair.nbrbeg, pair.nbrend) = + boost::adjacent_vertices(core_1[*n1iter_beg], *g2); + pair.hasiter = true; } + } else if (pair.n1 == 0 && order != NULL) { + // Optimisation: if the order vector is laid out in a DFS/BFS then this + // loop can be replaced with: + // pair.n1=order[core_len]; + // :) + unsigned int i = 0; + while (i < n1 && core_1[pair.n1 = order[i]] != NULL_NODE) i++; + if (i == n1) pair.n1 = n1; + } else { + while (pair.n1 < n1 && core_1[pair.n1] != NULL_NODE) { + pair.n1++; + pair.n2 = 0; + } + } + + /* VF2 Plus iterator available? */ + if (pair.hasiter) { + while (pair.nbrbeg < pair.nbrend && core_2[*pair.nbrbeg] != NULL_NODE) { + ++pair.nbrbeg; + } + + if (pair.nbrbeg < pair.nbrend) { + pair.n2 = *pair.nbrbeg; + ++pair.nbrbeg; + } else { + pair.n2 = n2; + } + } else if (t1_len > core_len && t2_len > core_len) { + while (pair.n2 < n2 && + (core_2[pair.n2] != NULL_NODE || term_2[pair.n2] == 0)) { + pair.n2++; + } + } else { + while (pair.n2 < n2 && core_2[pair.n2] != NULL_NODE) { + pair.n2++; + } + } + return pair.n1 < n1 && pair.n2 < n2; + }; + bool IsFeasiblePair(node_id node1, node_id node2) { + assert(node1 < n1); + assert(node2 < n2); + assert(core_1[node1] == NULL_NODE); + assert(core_2[node2] == NULL_NODE); + + // std::cerr<<" ifp:"<size()< boost::out_degree(node2, *g2)) return false; + if (!vc(node1, node2)) return false; + + unsigned int other1, other2; +#ifdef RDK_VF2_PRUNING + unsigned int term1 = 0, term2 = 0; + unsigned int new1 = 0, new2 = 0; +#endif + + // Check the out edges of node1 + typename Graph::out_edge_iterator bNbrs, eNbrs; + boost::tie(bNbrs, eNbrs) = boost::out_edges(node1, *g1); + while (bNbrs != eNbrs) { + other1 = getOtherIdx(*g1, *bNbrs, node1); + if (core_1[other1] != NULL_NODE) { + other2 = core_1[other1]; + typename Graph::edge_descriptor oEdge; + bool found; + boost::tie(oEdge, found) = boost::edge(node2, other2, *g2); + if (!found || !ec(*bNbrs, oEdge)) { + // std::cerr<<" short2"< - bool match(node_id c1[], node_id c2[], SubState &s, DoubleBackInsertionSequence &res, - unsigned int max_results) { - s.MatchAll(c1, c2, res, max_results); - return !res.empty(); +#ifdef RDK_VF2_PRUNING + // Check the out edges of node2 + boost::tie(bNbrs, eNbrs) = boost::out_edges(node2, *g2); + while (bNbrs != eNbrs) { + other2 = getOtherIdx(*g2, *bNbrs, node2); + if (core_2[other2] != NULL_NODE) { + // do nothing + } else { + if (term_2[other2]) ++term2; + if (!term_2[other2]) ++new2; + } + ++bNbrs; } - }; //end of namespace detail + // std::cerr<<(termin1 <= termin2 && termout1 <= termout2 && + // (termin1+termout1+new1)<=(termin2+termout2+new2))< - > - bool vf2(const Graph &g1,const Graph &g2, - VertexLabeling& vertex_labeling, - EdgeLabeling& edge_labeling, - MatchChecking& match_checking, - BackInsertionSequence& F){ - detail::VF2SubState s0(&g1,&g2,vertex_labeling, - edge_labeling,match_checking,false); - detail::node_id *ni1 = new detail::node_id[num_vertices(g1)]; - detail::node_id *ni2 = new detail::node_id[num_vertices(g2)]; - int n=0; - - F.clear(); - F.resize(0); - if(match(&n,ni1,ni2,s0)){ - for(unsigned int i=0;i(ni1[i],ni2[i])); + // n.b. term1+new1 == boost::out_degree(node1) and + // term2+new2 == boost::out_degree(node2) + return term1 <= term2 && (term1 + new1) <= (term2 + new2); +#else + return true; +#endif + }; + void AddPair(node_id node1, node_id node2) { + assert(node1 < n1); + assert(node2 < n2); + assert(core_len < n1); + assert(core_len < n2); + + ++core_len; + if (!term_1[node1]) { + term_1[node1] = core_len; + ++t1_len; + } + + if (!term_2[node2]) { + term_2[node2] = core_len; + ++t2_len; + } + + core_1[node1] = node2; + core_2[node2] = node1; + + typename Graph::out_edge_iterator bNbrs, eNbrs; + // FIX: this is explicitly ignoring directionality + boost::tie(bNbrs, eNbrs) = boost::out_edges(node1, *g1); + while (bNbrs != eNbrs) { + unsigned int other = getOtherIdx(*g1, *bNbrs, node1); + if (!term_1[other]) { + term_1[other] = core_len; + ++t1_len; + } + ++bNbrs; + } + + // FIX: this is explicitly ignoring directionality + boost::tie(bNbrs, eNbrs) = boost::out_edges(node2, *g2); + while (bNbrs != eNbrs) { + unsigned int other = getOtherIdx(*g2, *bNbrs, node2); + if (!term_2[other]) { + term_2[other] = core_len; + ++t2_len; + } + ++bNbrs; + } + }; + void GetCoreSet(node_id c1[], node_id c2[]) { + unsigned int i, j; + for (i = 0, j = 0; i < n1; ++i) { + if (core_1[i] != NULL_NODE) { + c1[j] = i; + c2[j] = core_1[i]; + ++j; } } - delete [] ni1; - delete [] ni2; - - return !F.empty(); }; - template < class Graph - , class VertexLabeling // binary predicate - , class EdgeLabeling // binary predicate - , class MatchChecking // binary predicate - , class DoubleBackInsertionSequence // contains a back insertion sequence - > - bool vf2_all(const Graph& g1, const Graph& g2, - VertexLabeling& vertex_labeling, - EdgeLabeling& edge_labeling, - MatchChecking& match_checking, - DoubleBackInsertionSequence& F, - unsigned int max_results=1000) { - detail::VF2SubState s0(&g1,&g2,vertex_labeling, - edge_labeling,match_checking,false); - detail::node_id *ni1 = new detail::node_id[num_vertices(g1)]; - detail::node_id *ni2 = new detail::node_id[num_vertices(g2)]; - - F.clear(); - F.resize(0); + VF2SubState *Clone() { return new VF2SubState(*this); }; + void BackTrack(node_id node1, node_id node2) { + if (term_1[node1] == core_len) { + term_1[node1] = 0; + --t1_len; + } - match(ni1,ni2,s0,F,max_results); + typename Graph::out_edge_iterator bNbrs, eNbrs; + boost::tie(bNbrs, eNbrs) = boost::out_edges(node1, *g1); + while (bNbrs != eNbrs) { + unsigned int other = getOtherIdx(*g1, *bNbrs, node1); + if (term_1[other] == core_len) { + term_1[other] = 0; + --t1_len; + } + ++bNbrs; + } - delete [] ni1; - delete [] ni2; - - return !F.empty(); + if (term_2[node2] == core_len) { + term_2[node2] = 0; + --t2_len; + } + + boost::tie(bNbrs, eNbrs) = boost::out_edges(node2, *g2); + while (bNbrs != eNbrs) { + unsigned int other = getOtherIdx(*g2, *bNbrs, node2); + if (term_2[other] == core_len) { + term_2[other] = 0; + --t2_len; + } + ++bNbrs; + } + + core_1[node1] = NULL_NODE; + core_2[node2] = NULL_NODE; + --core_len; }; -} // end of namespace boost + bool Match(node_id c1[], node_id c2[]) { + if (IsGoal()) { + GetCoreSet(c1, c2); + if (MatchChecks(c1, c2)) return true; + } + + if (IsDead()) return false; + + Pair pair; + while (NextPair(pair)) { + if (IsFeasiblePair(pair.n1, pair.n2)) { + AddPair(pair.n1, pair.n2); + if (Match(c1, c2)) // recurse + return true; + BackTrack(pair.n1, pair.n2); + } + } + return false; + } + + template + bool MatchAll(node_id c1[], node_id c2[], DoubleBackInsertionSequence &res, + unsigned int lim = 0) { + if (IsGoal()) { + GetCoreSet(c1, c2); + if (MatchChecks(c1, c2)) { + typename DoubleBackInsertionSequence::value_type newSeq; + for (unsigned int i = 0; i < core_len; ++i) { + newSeq.push_back(std::pair(c1[i], c2[i])); + } + res.push_back(newSeq); + return lim && res.size() >= lim; + } + } + + if (IsDead()) return false; + + Pair pair; + while (NextPair(pair)) { + if (IsFeasiblePair(pair.n1, pair.n2)) { + AddPair(pair.n1, pair.n2); + if (MatchAll(c1, c2, res, lim)) // recurse + return true; + BackTrack(pair.n1, pair.n2); + } + } + return false; + } +}; + +/*------------------------------------------------------------- + * static bool match(pn, c1, c2, s) + * Finds a matching between two graphs, if it exists, starting + * from state s. + * Returns true a match has been found. + * *pn is assigned the numbero of matched nodes, and + * c1 and c2 will contain the ids of the corresponding nodes + * in the two graphs. + ------------------------------------------------------------*/ +template +bool match(int *pn, node_id c1[], node_id c2[], SubState &s) { + if (s.Match(c1, c2)) { + // not needed, pn = num query atoms (n1)... + *pn = s.CoreLen(); + return true; + } + return false; +} + +/*------------------------------------------------------------- + * static bool match(c1, c2, vis, usr_data, pcount) + * Visits all the matchings between two graphs, starting + * from state s. + * Returns true if the caller must stop the visit. + * Stops when max_results is reached, set max_results to 0 to + * keep going until there are no more matches + * + ------------------------------------------------------------*/ +template +bool match(node_id c1[], node_id c2[], SubState &s, + DoubleBackInsertionSequence &res, unsigned int max_results) { + s.MatchAll(c1, c2, res, max_results); + return !res.empty(); +} +}; // end of namespace detail + +template < + class Graph, class VertexLabeling // binary predicate + , + class EdgeLabeling // binary predicate + , + class MatchChecking // binary predicate + , + class + BackInsertionSequence // contains + // std::pair + > +bool vf2(const Graph &g1, const Graph &g2, VertexLabeling &vertex_labeling, + EdgeLabeling &edge_labeling, MatchChecking &match_checking, + BackInsertionSequence &F) { + detail::VF2SubState + s0(&g1, &g2, vertex_labeling, edge_labeling, match_checking, false); + detail::node_id *ni1 = new detail::node_id[num_vertices(g1)]; + detail::node_id *ni2 = new detail::node_id[num_vertices(g2)]; + int n = 0; + + F.clear(); + F.resize(0); + if (match(&n, ni1, ni2, s0)) { + for (unsigned int i = 0; i < num_vertices(g1); i++) { + F.push_back(std::pair(ni1[i], ni2[i])); + } + } + delete[] ni1; + delete[] ni2; + + return !F.empty(); +}; +template +bool vf2_all(const Graph &g1, const Graph &g2, VertexLabeling &vertex_labeling, + EdgeLabeling &edge_labeling, MatchChecking &match_checking, + DoubleBackInsertionSequence &F, unsigned int max_results = 1000) { + detail::VF2SubState + s0(&g1, &g2, vertex_labeling, edge_labeling, match_checking, false); + detail::node_id *ni1 = new detail::node_id[num_vertices(g1)]; + detail::node_id *ni2 = new detail::node_id[num_vertices(g2)]; + + F.clear(); + F.resize(0); + + match(ni1, ni2, s0, F, max_results); + + delete[] ni1; + delete[] ni2; + + return !F.empty(); +}; +} // end of namespace boost #endif #undef RDK_VF2_PRUNING diff --git a/Code/GraphMol/Wrap/EditableMol.cpp b/Code/GraphMol/Wrap/EditableMol.cpp index f13f9ef1a..8dd057919 100644 --- a/Code/GraphMol/Wrap/EditableMol.cpp +++ b/Code/GraphMol/Wrap/EditableMol.cpp @@ -29,10 +29,7 @@ namespace { class EditableMol : boost::noncopyable { public: EditableMol(const ROMol &m) { dp_mol = new RWMol(m); }; - ~EditableMol() noexcept { - PRECONDITION(dp_mol, "no molecule"); - delete dp_mol; - }; + ~EditableMol() noexcept { delete dp_mol; }; void RemoveAtom(unsigned int idx) { PRECONDITION(dp_mol, "no molecule"); diff --git a/Code/GraphMol/Wrap/ForwardSDMolSupplier.cpp b/Code/GraphMol/Wrap/ForwardSDMolSupplier.cpp index 9bfcb821e..902d3baf8 100644 --- a/Code/GraphMol/Wrap/ForwardSDMolSupplier.cpp +++ b/Code/GraphMol/Wrap/ForwardSDMolSupplier.cpp @@ -34,7 +34,7 @@ class LocalForwardSDMolSupplier : public RDKit::ForwardSDMolSupplier { LocalForwardSDMolSupplier(python::object &input, bool sanitize, bool removeHs, bool strictParsing) { // FIX: minor leak here - auto *sb = new streambuf(input,'b'); + auto *sb = new streambuf(input, 'b'); dp_inStream = new streambuf::istream(*sb); df_owner = true; df_sanitize = sanitize; @@ -56,7 +56,8 @@ class LocalForwardSDMolSupplier : public RDKit::ForwardSDMolSupplier { std::istream *tmpStream = nullptr; tmpStream = static_cast( new std::ifstream(filename.c_str(), std::ios_base::binary)); - if (!tmpStream || (!(*tmpStream)) || (tmpStream->bad())) { + if (!(*tmpStream) || tmpStream->bad()) { + delete tmpStream; std::ostringstream errout; errout << "Bad input file " << filename; throw RDKit::BadFileException(errout.str()); diff --git a/Code/GraphMol/Wrap/MaeMolSupplier.cpp b/Code/GraphMol/Wrap/MaeMolSupplier.cpp index a5dcbe578..e1749e790 100644 --- a/Code/GraphMol/Wrap/MaeMolSupplier.cpp +++ b/Code/GraphMol/Wrap/MaeMolSupplier.cpp @@ -77,7 +77,8 @@ class LocalMaeMolSupplier : public RDKit::MaeMolSupplier { bool removeHs = true) { df_owner = true; auto *ifs = new std::ifstream(fname.c_str(), std::ios_base::binary); - if (!ifs || !(*ifs) || ifs->bad()) { + if (!(*ifs) || ifs->bad()) { + delete ifs; std::ostringstream errout; errout << "Bad input file " << fname; throw RDKit::BadFileException(errout.str()); diff --git a/Code/GraphMol/Wrap/MolOps.cpp b/Code/GraphMol/Wrap/MolOps.cpp index e8cbabed0..907fa5c3f 100644 --- a/Code/GraphMol/Wrap/MolOps.cpp +++ b/Code/GraphMol/Wrap/MolOps.cpp @@ -372,7 +372,7 @@ MolOps::SanitizeFlags sanitizeMol(ROMol &mol, boost::uint64_t sanitizeOps, } RWMol *getEditable(const ROMol &mol) { - RWMol *res = static_cast(new ROMol(mol, false)); + RWMol *res = new RWMol(mol, false); return res; } diff --git a/Code/ML/Cluster/Murtagh/Clustering.cpp b/Code/ML/Cluster/Murtagh/Clustering.cpp index 3e6b6b34e..5f35d2a53 100644 --- a/Code/ML/Cluster/Murtagh/Clustering.cpp +++ b/Code/ML/Cluster/Murtagh/Clustering.cpp @@ -38,6 +38,7 @@ void clusterit(real *dataP, boost::int64_t n, boost::int64_t m, double tmp; len = (n * (n - 1)) / 2; dists = (real *)calloc(len, sizeof(real)); + CHECK_INVARIANT(dists, "failed to allocate memory"); for (i = 1; i < n; i++) { iTab = i * m; for (j = 0; j < i; j++) { diff --git a/Code/ML/Data/cQuantize.cpp b/Code/ML/Data/cQuantize.cpp index a82e3ba38..2faab6f71 100644 --- a/Code/ML/Data/cQuantize.cpp +++ b/Code/ML/Data/cQuantize.cpp @@ -1,4 +1,3 @@ -// $Id$ // // Copyright 2003-2008 Rational Discovery LLC and Greg Landrum // @@ All Rights Reserved @@ @@ -124,6 +123,8 @@ long int *GenVarTable(double *vals, int nVals, long int *cuts, int nCuts, double RecurseHelper(double *vals, int nVals, long int *cuts, int nCuts, int which, long int *starts, int nStarts, long int *results, int nPossibleRes) { + PRECONDITION(vals, "bad vals pointer"); + double maxGain = -1e6, gainHere; long int *bestCuts, *tCuts; long int *varTable = nullptr; @@ -133,6 +134,9 @@ double RecurseHelper(double *vals, int nVals, long int *cuts, int nCuts, varTable = (long int *)calloc((nCuts + 1) * nPossibleRes, sizeof(long int)); bestCuts = (long int *)calloc(nCuts, sizeof(long int)); tCuts = (long int *)calloc(nCuts, sizeof(long int)); + CHECK_INVARIANT(varTable, "failed to allocate memory"); + CHECK_INVARIANT(bestCuts, "failed to allocate memory"); + CHECK_INVARIANT(tCuts, "failed to allocate memory"); GenVarTable(vals, nVals, cuts, nCuts, starts, results, nPossibleRes, varTable); while (cuts[which] <= highestCutHere) { @@ -244,6 +248,7 @@ static python::tuple cQuantize_RecurseOnBounds(python::object vals, python::ssize_t nCuts = python::len(pyCuts); cuts = (long int *)calloc(nCuts, sizeof(long int)); + CHECK_INVARIANT(cuts, "failed to allocate memory"); for (python::ssize_t i = 0; i < nCuts; i++) { python::object elem = pyCuts[i]; cuts[i] = python::extract(elem); @@ -251,6 +256,7 @@ static python::tuple cQuantize_RecurseOnBounds(python::object vals, python::ssize_t nStarts = python::len(pyStarts); starts = (long int *)calloc(nStarts, sizeof(long int)); + CHECK_INVARIANT(starts, "failed to allocate memory"); for (python::ssize_t i = 0; i < nStarts; i++) { python::object elem = pyStarts[i]; starts[i] = python::extract(elem); diff --git a/Code/SimDivPickers/LeaderPicker.h b/Code/SimDivPickers/LeaderPicker.h index f0e71691d..18a7495aa 100644 --- a/Code/SimDivPickers/LeaderPicker.h +++ b/Code/SimDivPickers/LeaderPicker.h @@ -116,9 +116,10 @@ class LeaderPicker : public DistPicker { }; #ifdef USE_THREADED_LEADERPICKER -// Note that this block of code currently only works on linux (which is why it's disabled by default) -// We will revisit this during the 2020.03 release cycle in order to get a multi-threaded version of -// the LeaderPicker that works on all supported platforms +// Note that this block of code currently only works on linux (which is why it's +// disabled by default) We will revisit this during the 2020.03 release cycle in +// order to get a multi-threaded version of the LeaderPicker that works on all +// supported platforms template void *LeaderPickerWork(void *arg); @@ -325,10 +326,9 @@ struct LeaderPickerState { int query; T *func; - LeaderPickerState(unsigned int count, int) { + LeaderPickerState(unsigned int count, int) : left(count), threshold(0.0), query(0), func(nullptr) { v.resize(count); for (unsigned int i = 0; i < count; i++) v[i] = i; - left = count; } bool empty() { return left == 0; } diff --git a/Code/SimDivPickers/MaxMinPicker.h b/Code/SimDivPickers/MaxMinPicker.h index f9ffa28f5..21acb671d 100644 --- a/Code/SimDivPickers/MaxMinPicker.h +++ b/Code/SimDivPickers/MaxMinPicker.h @@ -151,10 +151,6 @@ RDKit::INT_VECT MaxMinPicker::lazyPick(T &func, unsigned int poolSize, unsigned int memsize = (unsigned int)(poolSize * sizeof(MaxMinPickInfo)); MaxMinPickInfo *pinfo = new MaxMinPickInfo[memsize]; - if (!pinfo) { - threshold = -1.0; - return picks; - } memset(pinfo, 0, memsize); picks.reserve(pickSize);