C++ MatchVectType::size方法代码示例

void test6(){
  std::cout << " ----------------- Test 6 (Issue71 related)" << std::endl;
  MatchVectType matchV;
  std::vector< MatchVectType > matches;
  int n;

  RWMol *m,*q1;
  Atom *a6 = new Atom(6);

  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(0,1,Bond::SINGLE);
  m->addBond(1,2,Bond::SINGLE);
  m->addBond(0,2,Bond::SINGLE);

  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6),true);
  q1->addAtom(new QueryAtom(6),true);
  q1->addAtom(new QueryAtom(6),true);
  q1->addBond(0,1,Bond::UNSPECIFIED);
  q1->addBond(1,2,Bond::UNSPECIFIED);

  bool found = SubstructMatch(*m,*q1,matchV);
  CHECK_INVARIANT(found,"");
  CHECK_INVARIANT(matchV.size()==3,"");
  n = SubstructMatch(*m,*q1,matches,true);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches[0].size()==3,"");

  // close the loop and try again (we should still match)
  q1->addBond(0,2,Bond::UNSPECIFIED);
  found = SubstructMatch(*m,*q1,matchV);
  CHECK_INVARIANT(found,"");
  CHECK_INVARIANT(matchV.size()==3,"");
  n = SubstructMatch(*m,*q1,matches,true);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches[0].size()==3,"");

  std::cout << "Done\n" << std::endl;
}

void test4(){
  std::cout << " ----------------- Test 4" << std::endl;
  MatchVectType matchV;
  std::vector< MatchVectType > matches;
  int n;

  RWMol *m,*q1,*q2;
  Atom *a6 = new Atom(6);
  Atom *a8 = new Atom(8);
  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a8);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(1,0,Bond::SINGLE);
  m->addBond(1,2,Bond::SINGLE);
  m->addBond(1,3,Bond::SINGLE);
  m->addBond(2,4,Bond::SINGLE);

  // this will be the recursive query
  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6),true);
  q1->addAtom(new QueryAtom(8),true);
  q1->addBond(0,1,Bond::UNSPECIFIED);

  // here's the main query
  q2 = new RWMol();
  QueryAtom *qA = new QueryAtom(6);
  RecursiveStructureQuery *rsq = new RecursiveStructureQuery(q1);
  qA->expandQuery(rsq,Queries::COMPOSITE_AND);
  //std::cout << "post expand: " << qA->getQuery() << std::endl;
  q2->addAtom(qA,true,true);
  //std::cout << "mol: " << q2->getAtomWithIdx(0)->getQuery() << std::endl;
  q2->addAtom(new QueryAtom(6),true,true);
  q2->addBond(0,1,Bond::UNSPECIFIED);

  bool found = SubstructMatch(*m,*q2,matchV);
  CHECK_INVARIANT(found,"");
  CHECK_INVARIANT(matchV.size()==2,"");
  TEST_ASSERT(matchV[0].first==0);
  TEST_ASSERT(matchV[0].second==1);
  TEST_ASSERT(matchV[1].first==1);
  TEST_ASSERT(matchV[1].second==0||matchV[1].second==3);
  n = SubstructMatch(*m,*q2,matches,true);
  TEST_ASSERT(n==2);
  TEST_ASSERT(matches.size()==n);
  TEST_ASSERT(matches[0].size()==2);
  TEST_ASSERT(matches[1].size()==2);
  TEST_ASSERT(matches[0][0].second==matches[1][0].second);
  TEST_ASSERT(matches[0][1].second!=matches[1][1].second);

  std::cout << "Done\n" << std::endl;
}

void test5QueryRoot() {
  std::cout << " ----------------- Test 5 QueryRoot" << std::endl;
  MatchVectType matchV;
  std::vector<MatchVectType> matches;
  int n;
  bool updateLabel = true;
  bool takeOwnership = true;

  RWMol *m, *q1, *q2;
  auto *a6 = new Atom(6);
  auto *a8 = new Atom(8);
  // CC(OC)C
  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a8);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(0, 1, Bond::SINGLE);
  m->addBond(1, 2, Bond::SINGLE);
  m->addBond(1, 4, Bond::SINGLE);
  m->addBond(2, 3, Bond::SINGLE);

  // this will be the recursive query
  q1 = new RWMol();
  q1->addAtom(new QueryAtom(8), updateLabel, takeOwnership);
  q1->addAtom(new QueryAtom(6), updateLabel, takeOwnership);
  q1->addBond(0, 1, Bond::UNSPECIFIED);
  q1->setProp(common_properties::_queryRootAtom, 1);

  // here's the main query
  q2 = new RWMol();
  auto *qA = new QueryAtom();
  auto *rsq = new RecursiveStructureQuery(q1);
  qA->setQuery(rsq);
  q2->addAtom(qA, true, true);
  q2->addAtom(new QueryAtom(6), true, true);
  q2->addBond(0, 1, Bond::UNSPECIFIED);

  bool found = SubstructMatch(*m, *q2, matchV);
  CHECK_INVARIANT(found, "");
  CHECK_INVARIANT(matchV.size() == 2, "");
  n = SubstructMatch(*m, *q2, matches, true);
  CHECK_INVARIANT(n == 2, "");
  CHECK_INVARIANT(matches[0].size() == 2, "");
  delete m;
  delete a6;
  delete a8;
  delete q2;

  std::cout << "Done\n" << std::endl;
}

bool check(MatchVectType v, MatchVectType match) {
  dump("v", v);
  dump("match", match);
  for (size_t i = 0; i < v.size(); ++i) {
    if (v[i].first != match[i].first) {
      return false;
    }
    if (v[i].second != match[i].second) {
      return false;
    }
  }
  return true;
}

void testCrippenO3AConstraints() {
  ROMol *m = SmilesToMol("n1ccc(cc1)-c1ccccc1");
  TEST_ASSERT(m);
  ROMol *m1 = MolOps::addHs(*m);
  delete m;
  TEST_ASSERT(m1);
  DGeomHelpers::EmbedMolecule(*m1);
  MMFF::sanitizeMMFFMol((RWMol &)(*m1));
  MMFF::MMFFMolProperties mp(*m1);
  TEST_ASSERT(mp.isValid());
  ForceFields::ForceField *field = MMFF::constructForceField(*m1, &mp);
  field->initialize();
  field->minimize();
  RWMol *patt = SmartsToMol("nccc-cccc");
  MatchVectType matchVect;
  TEST_ASSERT(SubstructMatch(*m1, (ROMol &)*patt, matchVect));
  unsigned int nIdx = matchVect[0].second;
  unsigned int cIdx = matchVect[matchVect.size() - 1].second;
  MolTransforms::setDihedralDeg(m1->getConformer(), matchVect[2].second,
    matchVect[3].second, matchVect[4].second, matchVect[5].second, 0.0);
  ROMol m2(*m1);
  MolAlign::randomTransform(m2);
  ROMol m3(m2);
  unsigned int prbNAtoms = m2.getNumAtoms();
  std::vector<double> prbLogpContribs(prbNAtoms);
  std::vector<double> prbMRContribs(prbNAtoms);
  std::vector<unsigned int> prbAtomTypes(prbNAtoms);
  std::vector<std::string> prbAtomTypeLabels(prbNAtoms);
  Descriptors::getCrippenAtomContribs(m2, prbLogpContribs,
    prbMRContribs, true, &prbAtomTypes, &prbAtomTypeLabels);
  MolAlign::O3A *o3a = new MolAlign::O3A(m2, *m1,
    &prbLogpContribs, &prbLogpContribs, MolAlign::O3A::CRIPPEN);
  TEST_ASSERT(o3a);
  o3a->align();
  delete o3a;
  double d = (m2.getConformer().getAtomPos(cIdx)
    - m1->getConformer().getAtomPos(cIdx)).length();
  TEST_ASSERT(feq(d, 0.0, 1));
  MatchVectType constraintMap;
  constraintMap.push_back(std::make_pair(cIdx, nIdx));
  o3a = new MolAlign::O3A(m3, *m1, &prbLogpContribs, &prbLogpContribs,
    MolAlign::O3A::CRIPPEN, -1, -1, false, 50, 0, &constraintMap);
  TEST_ASSERT(o3a);
  o3a->align();
  delete o3a;
  d = (m3.getConformer().getAtomPos(cIdx)
    - m1->getConformer().getAtomPos(cIdx)).length();
  TEST_ASSERT(feq(d, 7.0, 1.0));
  delete m1;
}

void test5QueryRoot(){
  std::cout << " ----------------- Test 5 QueryRoot" << std::endl;
  MatchVectType matchV;
  std::vector< MatchVectType > matches;
  int n;

  RWMol *m,*q1,*q2;
  Atom *a6 = new Atom(6);
  Atom *a8 = new Atom(8);
  // CC(OC)C
  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a8);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(0,1,Bond::SINGLE);
  m->addBond(1,2,Bond::SINGLE);
  m->addBond(1,4,Bond::SINGLE);
  m->addBond(2,3,Bond::SINGLE);

  // this will be the recursive query
  q1 = new RWMol();
  q1->addAtom(new QueryAtom(8),true);
  q1->addAtom(new QueryAtom(6),true);
  q1->addBond(0,1,Bond::UNSPECIFIED);
  q1->setProp("_queryRootAtom",1);
  
  // here's the main query
  q2 = new RWMol();
  QueryAtom *qA = new QueryAtom();
  RecursiveStructureQuery *rsq = new RecursiveStructureQuery(q1);
  qA->setQuery(rsq);
  q2->addAtom(qA,true,true);
  q2->addAtom(new QueryAtom(6),true,true);
  q2->addBond(0,1,Bond::UNSPECIFIED);

  bool found = SubstructMatch(*m,*q2,matchV);
  CHECK_INVARIANT(found,"");
  CHECK_INVARIANT(matchV.size()==2,"");
  n = SubstructMatch(*m,*q2,matches,true);
  CHECK_INVARIANT(n==2,"");
  CHECK_INVARIANT(matches[0].size()==2,"");


  std::cout << "Done\n" << std::endl;
}

void test7() {
  std::cout << " ----------------- Test 7 (leak check)" << std::endl;
  MatchVectType matchV;
  int n;

  RWMol *m, *q1;
  auto *a6 = new Atom(6);
  bool updateLabel = true;
  bool takeOwnership = true;

  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(0, 1, Bond::SINGLE);
  m->addBond(1, 2, Bond::SINGLE);
  m->addBond(0, 2, Bond::SINGLE);

  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6), updateLabel, takeOwnership);
  q1->addAtom(new QueryAtom(6), updateLabel, takeOwnership);
  q1->addAtom(new QueryAtom(6), updateLabel, takeOwnership);
  q1->addBond(0, 1, Bond::UNSPECIFIED);
  q1->addBond(1, 2, Bond::UNSPECIFIED);

  bool found = SubstructMatch(*m, *q1, matchV);
  CHECK_INVARIANT(found, "");
  CHECK_INVARIANT(matchV.size() == 3, "");
  std::vector<MatchVectType> matches;
  for (int i = 0; i < 300000; i++) {
    n = SubstructMatch(*m, *q1, matches, true, true);
    CHECK_INVARIANT(n == 1, "");
    CHECK_INVARIANT(matches[0].size() == 3, "");
    if (!(i % 500)) std::cout << i << std::endl;
  }
  delete m;
  delete a6;
  delete q1;

  std::cout << "Done\n" << std::endl;
}

void test8(){
  std::cout << " ----------------- Test 8 (molgraph cache)" << std::endl;
  MatchVectType matchV;
  int n;

  RWMol *m,*q1;
  Atom *a6 = new Atom(6);

  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(0,1,Bond::SINGLE);
  m->addBond(1,2,Bond::SINGLE);
  m->addBond(0,2,Bond::SINGLE);

  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6),true);
  q1->addAtom(new QueryAtom(6),true);
  q1->addAtom(new QueryAtom(6),true);
  q1->addBond(0,1,Bond::UNSPECIFIED);
  q1->addBond(1,2,Bond::UNSPECIFIED);


  bool found = SubstructMatch(*m,*q1,matchV);
  CHECK_INVARIANT(found,"");
  CHECK_INVARIANT(matchV.size()==3,"");
  std::vector< MatchVectType > matches;
  for(int i=0;i<30000;i++){
    n = SubstructMatch(*m,*q1,matches,true,true);
    CHECK_INVARIANT(n==1,"");
    CHECK_INVARIANT(matches[0].size()==3,"");
    if(! (i%500) ) std::cout << i << std::endl;
  }
  std::cout << "Done\n" << std::endl;
}

void addReactantAtomsAndBonds(const ChemicalReaction& rxn,
                              RWMOL_SPTR product,const ROMOL_SPTR reactantSptr,
                              const MatchVectType &match,
                              const ROMOL_SPTR reactantTemplate,
                              Conformer *productConf)
{
    // start by looping over all matches and marking the reactant atoms that
    // have already been "added" by virtue of being in the product. We'll also
    // mark "skipped" atoms: those that are in the match, but not in this
    // particular product (or, perhaps, not in any product)
    // At the same time we'll set up a map between the indices of those
    // atoms and their index in the product.
    ReactantProductAtomMapping* mapping =
        getAtomMappingsReactantProduct(match,*reactantTemplate,product,reactantSptr->getNumAtoms());

    boost::dynamic_bitset<> visitedAtoms(reactantSptr->getNumAtoms());

    const ROMol *reactant=reactantSptr.get();

    // ---------- ---------- ---------- ---------- ---------- ----------
    // Loop over the bonds in the product and look for those that have
    // the NullBond property set. These are bonds for which no information
    // (other than their existance) was provided in the template:
    setReactantBondPropertiesToProduct(product, *reactant, mapping);

    // ---------- ---------- ---------- ---------- ---------- ----------
    // Loop over the atoms in the match that were added to the product
    // From the corresponding atom in the reactant, do a graph traversal
    // to find other connected atoms that should be added:

    std::vector<const Atom *> chiralAtomsToCheck;
    for(unsigned matchIdx=0; matchIdx<match.size(); matchIdx++) {

        int reactantAtomIdx = match[matchIdx].second;
        if(mapping->mappedAtoms[reactantAtomIdx]) {
            CHECK_INVARIANT(mapping->reactProdAtomMap.find(reactantAtomIdx)!=mapping->reactProdAtomMap.end(),
                            "mapped reactant atom not present in product.");

            const Atom *reactantAtom=reactant->getAtomWithIdx(reactantAtomIdx);
            for(unsigned i = 0; i < mapping->reactProdAtomMap[reactantAtomIdx].size(); i++) {
                // here's a pointer to the atom in the product:
                unsigned productAtomIdx = mapping->reactProdAtomMap[reactantAtomIdx][i];
                Atom *productAtom = product->getAtomWithIdx(productAtomIdx);
                setReactantAtomPropertiesToProduct(productAtom, *reactantAtom,rxn.getImplicitPropertiesFlag());
            }
            // now traverse:
            addReactantNeighborsToProduct(*reactant, *reactantAtom, product, visitedAtoms,
                                          chiralAtomsToCheck, mapping);

            // now that we've added all the reactant's neighbors, check to see if
            // it is chiral in the reactant but is not in the reaction. If so
            // we need to worry about its chirality
            checkAndCorrectChiralityOfMatchingAtomsInProduct(*reactant, reactantAtomIdx, *reactantAtom, product, mapping);
        }
    } // end of loop over matched atoms

    // ---------- ---------- ---------- ---------- ---------- ----------
    // now we need to loop over atoms from the reactants that were chiral but not
    // directly involved in the reaction in order to make sure their chirality hasn't
    // been disturbed
    checkAndCorrectChiralityOfProduct(chiralAtomsToCheck, product, mapping);

    // ---------- ---------- ---------- ---------- ---------- ----------
    // finally we may need to set the coordinates in the product conformer:
    if(productConf) {
        productConf->resize(product->getNumAtoms());
        generateProductConformers(productConf, *reactant, mapping);
    }
    delete(mapping);
} // end of addReactantAtomsAndBonds

void testRecursiveSerialNumbers() {
  std::cout << " ----------------- Testing serial numbers on recursive queries"
            << std::endl;
  MatchVectType matchV;
  std::vector<MatchVectType> matches;
  int n;
  bool updateLabel = true;
  bool takeOwnership = true;

  RWMol *m, *q1, *q2;
  auto *a6 = new Atom(6);
  auto *a8 = new Atom(8);
  m = new RWMol();
  m->addAtom(a6);
  m->addAtom(a6);
  m->addAtom(a8);
  m->addAtom(a6);
  m->addAtom(a6);
  m->addBond(1, 0, Bond::SINGLE);
  m->addBond(1, 2, Bond::SINGLE);
  m->addBond(1, 3, Bond::SINGLE);
  m->addBond(2, 4, Bond::SINGLE);

  {
    // this will be the recursive query
    q1 = new RWMol();
    q1->addAtom(new QueryAtom(6), updateLabel, takeOwnership);
    q1->addAtom(new QueryAtom(8), updateLabel, takeOwnership);
    q1->addBond(0, 1, Bond::UNSPECIFIED);

    // here's the main query
    q2 = new RWMol();
    auto *qA = new QueryAtom(6);
    auto *rsq = new RecursiveStructureQuery(new RWMol(*q1), 1);
    qA->expandQuery(rsq, Queries::COMPOSITE_AND);
    // std::cout << "post expand: " << qA->getQuery() << std::endl;
    q2->addAtom(qA, true, true);
    // std::cout << "mol: " << q2->getAtomWithIdx(0)->getQuery() << std::endl;
    q2->addAtom(new QueryAtom(8), true, true);
    q2->addBond(0, 1, Bond::UNSPECIFIED);

    qA = new QueryAtom(6);
    rsq = new RecursiveStructureQuery(new RWMol(*q1), 1);
    qA->expandQuery(rsq, Queries::COMPOSITE_AND);
    q2->addAtom(qA, true, true);
    q2->addBond(1, 2, Bond::UNSPECIFIED);

    bool found = SubstructMatch(*m, *q2, matchV);
    CHECK_INVARIANT(found, "");
    CHECK_INVARIANT(matchV.size() == 3, "");
    n = SubstructMatch(*m, *q2, matches, true);
    TEST_ASSERT(n == 1);
    TEST_ASSERT(matches.size() == 1);
    TEST_ASSERT(matches[0].size() == 3);
    delete q1;
    delete q2;
  }
  delete m;
  delete a6;
  delete a8;
  std::cout << "Done\n" << std::endl;
}

void test1(){
  std::cout << " ----------------- Test 1" << std::endl;
  MatchVectType matchV;
  std::vector< MatchVectType > matches;
  unsigned int n;

  RWMol *m,*q1;
  m = new RWMol();
  m->addAtom(new Atom(8));
  m->addAtom(new Atom(6));
  m->addAtom(new Atom(6));
  m->addBond(0,1,Bond::SINGLE);
  m->addBond(1,2,Bond::SINGLE);

  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6));
  q1->addAtom(new QueryAtom(6));
  q1->addBond(0,1,Bond::SINGLE);
  n = SubstructMatch(*m,*q1,matches,false);
  CHECK_INVARIANT(n==2,"");
  CHECK_INVARIANT(matches.size()==n,"");
  CHECK_INVARIANT(matches[0].size()==2,"");
  TEST_ASSERT(matches[0][0].first==0);
  TEST_ASSERT(matches[0][0].second==1||matches[0][0].second==2);
  TEST_ASSERT(matches[0][1].first==1);
  TEST_ASSERT(matches[0][1].second!=matches[0][0].second);

  TEST_ASSERT(matches[1][1].second==1||matches[0][1].second==2);
  TEST_ASSERT(matches[1][0].first==0);
  TEST_ASSERT(matches[1][0].second==1||matches[1][0].second==2);
  TEST_ASSERT(matches[1][0].second!=matches[0][0].second);
  TEST_ASSERT(matches[1][0].second==matches[0][1].second);
  TEST_ASSERT(matches[1][1].first==1);
  TEST_ASSERT(matches[1][1].second!=matches[1][0].second);
  TEST_ASSERT(matches[1][1].second==matches[0][0].second);

  n = SubstructMatch(*m,*q1,matches,true);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches.size()==n,"");
  CHECK_INVARIANT(matches[0].size()==2,"");
  TEST_ASSERT(matches[0][0].first==0);
  TEST_ASSERT(matches[0][0].second==1||matches[0][0].second==2);
  TEST_ASSERT(matches[0][1].first==1);
  TEST_ASSERT(matches[0][1].second!=matches[0][0].second);
  TEST_ASSERT(matches[0][1].second==1||matches[0][1].second==2);
  
  CHECK_INVARIANT(SubstructMatch(*m,*q1,matchV),"");
  CHECK_INVARIANT(matchV.size()==2,"");

  // make sure we reset the match vectors.
  // build a query we won't match:
  q1->addAtom(new QueryAtom(6));
  q1->addBond(1,2,Bond::SINGLE);
  q1->addAtom(new QueryAtom(6));
  q1->addBond(2,3,Bond::SINGLE);

  TEST_ASSERT(!SubstructMatch(*m,*q1,matchV));
  TEST_ASSERT(matchV.size()==0);

  n = SubstructMatch(*m,*q1,matches,false);
  TEST_ASSERT(n==0);
  TEST_ASSERT(matches.size()==0);
  
  std::cout << "Done\n" << std::endl;
}

void test5() {
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "FeatureFactory testing." << std::endl;

  ROMol *testMol;
  MatchVectType mv;
  std::string inText;
  FeatSPtrList featSPtrs;
  boost::shared_ptr<MolChemicalFeature> featSPtr;

  MolChemicalFeatureFactory *factory;

  MolChemicalFeatureDef::CollectionType::value_type featDef;

  inText =
      "DefineFeature HDonor1 [N,O;!H0]\n"
      "  Family HBondDonor\n"
      "  Weights 1.0\n"
      "EndFeature\n"
      "DefineFeature HAcceptor1 [N,O]\n"
      "  Family HBondAcceptor\n"
      "  Weights 1.0\n"
      "EndFeature\n";

  factory = buildFeatureFactory(inText);
  TEST_ASSERT(factory);
  TEST_ASSERT(factory->getNumFeatureDefs() == 2);

  testMol = SmilesToMol("COCN");
  TEST_ASSERT(testMol);
  featDef = *factory->beginFeatureDefs();
  TEST_ASSERT(SubstructMatch(*testMol, *featDef->getPattern(), mv));
  BOOST_LOG(rdErrorLog) << "1" << std::endl;
  TEST_ASSERT(mv.size() == 1);
  TEST_ASSERT(mv[0].first == 0);
  TEST_ASSERT(mv[0].second == 3);
  featDef = *(++factory->beginFeatureDefs());
  TEST_ASSERT(SubstructMatch(*testMol, *featDef->getPattern(), mv));
  BOOST_LOG(rdErrorLog) << "2" << std::endl;
  TEST_ASSERT(mv.size() == 1);
  TEST_ASSERT(mv[0].first == 0);
  TEST_ASSERT(mv[0].second == 1 || mv[0].second == 3);

  // Test using the factory to find features:
  featSPtrs = factory->getFeaturesForMol(*testMol);
  BOOST_LOG(rdErrorLog) << "3" << std::endl;
  TEST_ASSERT(featSPtrs.size() == 3);
  featSPtr = *featSPtrs.begin();
  TEST_ASSERT(featSPtr->getFamily() == "HBondDonor");
  TEST_ASSERT(featSPtr->getType() == "HDonor1");
  featSPtr = *(++featSPtrs.begin());
  TEST_ASSERT(featSPtr->getFamily() == "HBondAcceptor");
  TEST_ASSERT(featSPtr->getType() == "HAcceptor1");
  featSPtr = *(++featSPtrs.begin());
  TEST_ASSERT(featSPtr->getFamily() == "HBondAcceptor");
  TEST_ASSERT(featSPtr->getType() == "HAcceptor1");

  // Test limiting stuff with includeOnly
  featSPtrs = factory->getFeaturesForMol(*testMol, "HBondAcceptor");
  BOOST_LOG(rdErrorLog) << "4" << std::endl;
  TEST_ASSERT(featSPtrs.size() == 2);
  featSPtr = *featSPtrs.begin();
  TEST_ASSERT(featSPtr->getFamily() == "HBondAcceptor");
  TEST_ASSERT(featSPtr->getType() == "HAcceptor1");
  featSPtr = *(++featSPtrs.begin());
  TEST_ASSERT(featSPtr->getFamily() == "HBondAcceptor");
  TEST_ASSERT(featSPtr->getType() == "HAcceptor1");

  featSPtrs = factory->getFeaturesForMol(*testMol, "HBondDonor");
  BOOST_LOG(rdErrorLog) << "5" << std::endl;
  TEST_ASSERT(featSPtrs.size() == 1);
  featSPtr = *featSPtrs.begin();
  TEST_ASSERT(featSPtr->getFamily() == "HBondDonor");
  TEST_ASSERT(featSPtr->getType() == "HDonor1");

  featSPtrs = factory->getFeaturesForMol(*testMol, "NotPresent");
  BOOST_LOG(rdErrorLog) << "6" << std::endl;
  TEST_ASSERT(featSPtrs.size() == 0);

  delete testMol;
  delete factory;

  BOOST_LOG(rdErrorLog) << "  done" << std::endl;
}

void testGitHubIssue688() {
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Test GitHub issue 688: partially specified "
                           "chiral substructure queries don't work properly"
                        << std::endl;
  {
    std::string smi = "C1CC[[email protected]](Cl)(N)O1";
    ROMol *mol = SmilesToMol(smi);
    TEST_ASSERT(mol);
    // mol->debugMol(std::cerr);
    std::string sma = "C1CC[[email protected]](N)O1";
    ROMol *qmol = SmartsToMol(sma);
    TEST_ASSERT(qmol);
    // qmol->updatePropertyCache();
    // qmol->debugMol(std::cerr);

    MatchVectType match;
    TEST_ASSERT(SubstructMatch(*mol, *qmol, match, true, false));
    TEST_ASSERT(match.size() == qmol->getNumAtoms());

    TEST_ASSERT(SubstructMatch(*mol, *qmol, match, true, true));
    TEST_ASSERT(match.size() == qmol->getNumAtoms());

    delete mol;
    delete qmol;
  }
  {
    std::string smi = "C1CC[[email protected]](Cl)(N)O1";
    ROMol *mol = SmilesToMol(smi);
    TEST_ASSERT(mol);
    // mol->debugMol(std::cerr);
    std::string sma = "C1CC[[email protected]@](N)O1";
    ROMol *qmol = SmartsToMol(sma);
    TEST_ASSERT(qmol);
    // qmol->updatePropertyCache();
    // qmol->debugMol(std::cerr);

    MatchVectType match;
    TEST_ASSERT(SubstructMatch(*mol, *qmol, match, true, false));
    TEST_ASSERT(match.size() == qmol->getNumAtoms());

    TEST_ASSERT(!SubstructMatch(*mol, *qmol, match, true, true));

    delete mol;
    delete qmol;
  }
  {
    std::string smi = "N[[email protected]]1(Cl)CCCO1";
    ROMol *mol = SmilesToMol(smi);
    TEST_ASSERT(mol);
    // mol->debugMol(std::cerr);
    std::string sma = "N[[email protected]]1CCCO1";
    ROMol *qmol = SmartsToMol(sma);
    TEST_ASSERT(qmol);
    // qmol->updatePropertyCache();
    // qmol->debugMol(std::cerr);
    // std::cerr << MolToSmiles(*qmol, true) << std::endl;

    MatchVectType match;
    TEST_ASSERT(SubstructMatch(*mol, *qmol, match, true, false));
    TEST_ASSERT(match.size() == qmol->getNumAtoms());

    TEST_ASSERT(SubstructMatch(*mol, *qmol, match, true, true));
    TEST_ASSERT(match.size() == qmol->getNumAtoms());

    delete mol;
    delete qmol;
  }
  {
    std::string smi = "N[[email protected]]1(Cl)CCCO1";
    ROMol *mol = SmilesToMol(smi);
    TEST_ASSERT(mol);
    // mol->debugMol(std::cerr);
    std::string sma = "N[[email protected]@]1CCCO1";
    ROMol *qmol = SmartsToMol(sma);
    TEST_ASSERT(qmol);
    // qmol->updatePropertyCache();
    // qmol->debugMol(std::cerr);
    // std::cerr << MolToSmiles(*qmol, true) << std::endl;

    MatchVectType match;
    TEST_ASSERT(SubstructMatch(*mol, *qmol, match, true, false));
    TEST_ASSERT(match.size() == qmol->getNumAtoms());

    TEST_ASSERT(!SubstructMatch(*mol, *qmol, match, true, true));

    delete mol;
    delete qmol;
  }
  BOOST_LOG(rdErrorLog) << "  done" << std::endl;
}

void test3(){
  std::cout << " ----------------- Test 3" << std::endl;
  MatchVectType matchV;
  std::vector< MatchVectType > matches;
  unsigned int n;

  RWMol *m,*q1;
  m = new RWMol();
  m->addAtom(new Atom(6));
  m->addAtom(new Atom(6));
  m->addAtom(new Atom(8));
  m->addBond(0,1,Bond::SINGLE);
  m->addBond(1,2,Bond::SINGLE);

  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6));
  q1->addAtom(new QueryAtom(8));
  q1->addBond(0,1,Bond::UNSPECIFIED);
  n = SubstructMatch(*m,*q1,matches,false);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches.size()==n,"");
  CHECK_INVARIANT(matches[0].size()==2,"");
  n = SubstructMatch(*m,*q1,matches,true);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches.size()==n,"");
  CHECK_INVARIANT(matches[0].size()==2,"");

  CHECK_INVARIANT(SubstructMatch(*m,*q1,matchV),"");
  CHECK_INVARIANT(matchV.size()==2,"");

  delete m;
  m = new RWMol();
  m->addAtom(new Atom(6));
  m->addAtom(new Atom(6));
  m->addAtom(new Atom(8));
  m->addBond(0,1,Bond::SINGLE);
  m->addBond(1,2,Bond::DOUBLE);

  matches.clear();
  n = SubstructMatch(*m,*q1,matches,false);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches.size()==n,"");
  CHECK_INVARIANT(matches[0].size()==2,"");
  n = SubstructMatch(*m,*q1,matches,true);
  CHECK_INVARIANT(n==1,"");
  CHECK_INVARIANT(matches.size()==n,"");
  CHECK_INVARIANT(matches[0].size()==2,"");

  CHECK_INVARIANT(SubstructMatch(*m,*q1,matchV),"");
  CHECK_INVARIANT(matchV.size()==2,"");


  delete q1;
  q1 = new RWMol();
  q1->addAtom(new QueryAtom(6));
  q1->addAtom(new QueryAtom(6));
  q1->addBond(0,1,Bond::UNSPECIFIED);
  n = SubstructMatch(*m,*q1,matches,false);
  TEST_ASSERT(n==2);
  TEST_ASSERT(matches.size()==n);
  TEST_ASSERT(matches[0].size()==2);
  TEST_ASSERT(matches[1].size()==2);
  TEST_ASSERT(matches[0][0].second!=matches[1][0].second);
  TEST_ASSERT(matches[0][1].second!=matches[1][1].second);
  n = SubstructMatch(*m,*q1,matches,true);
  TEST_ASSERT(n==1);
  TEST_ASSERT(matches.size()==n);

  std::cout << "Done\n" << std::endl;
}

void addTorsions(const ROMol &mol, MMFFMolProperties *mmffMolProperties,
                 ForceFields::ForceField *field,
                 const std::string &torsionBondSmarts) {
  PRECONDITION(field, "bad ForceField");
  PRECONDITION(mmffMolProperties, "bad MMFFMolProperties");
  PRECONDITION(mmffMolProperties->isValid(),
               "missing atom types - invalid force-field");

  std::ostream &oStream = mmffMolProperties->getMMFFOStream();
  ROMol::ADJ_ITER nbr1Idx;
  ROMol::ADJ_ITER end1Nbrs;
  ROMol::ADJ_ITER nbr2Idx;
  ROMol::ADJ_ITER end2Nbrs;
  double totalTorsionEnergy = 0.0;
  RDGeom::PointPtrVect points;
  if (mmffMolProperties->getMMFFVerbosity()) {
    if (mmffMolProperties->getMMFFVerbosity() == MMFF_VERBOSITY_HIGH) {
      oStream << "\n"
                 "T O R S I O N A L\n\n"
                 "--------------ATOMS---------------      ---ATOM TYPES---     "
                 "FF     TORSION              -----FORCE Params-----\n"
                 "  I        J        K        L          I    J    K    L   "
                 "CLASS     ANGLE    ENERGY       V1        V2        V3\n"
                 "-------------------------------------------------------------"
                 "-----------------------------------------------------"
              << std::endl;
    }
    points = field->positions();
  }
  std::vector<MatchVectType> matchVect;
  const ROMol *defaultQuery = DefaultTorsionBondSmarts::query();
  const ROMol *query = (torsionBondSmarts == DefaultTorsionBondSmarts::string())
                       ? defaultQuery : SmartsToMol(torsionBondSmarts);
  TEST_ASSERT(query);
  unsigned int nHits = SubstructMatch(mol, *query, matchVect);
  if (query != defaultQuery) delete query;

  for (unsigned int i = 0; i < nHits; ++i) {
    MatchVectType match = matchVect[i];
    TEST_ASSERT(match.size() == 2);
    int idx2 = match[0].second;
    int idx3 = match[1].second;
    const Bond *bond = mol.getBondBetweenAtoms(idx2, idx3);
    TEST_ASSERT(bond);
    const Atom *jAtom = mol.getAtomWithIdx(idx2);
    const Atom *kAtom = mol.getAtomWithIdx(idx3);
    if (((jAtom->getHybridization() == Atom::SP2) ||
         (jAtom->getHybridization() == Atom::SP3)) &&
        ((kAtom->getHybridization() == Atom::SP2) ||
         (kAtom->getHybridization() == Atom::SP3))) {
      ROMol::OEDGE_ITER beg1, end1;
      boost::tie(beg1, end1) = mol.getAtomBonds(jAtom);
      while (beg1 != end1) {
        const Bond *tBond1 = mol[*beg1].get();
        if (tBond1 != bond) {
          int idx1 = tBond1->getOtherAtomIdx(idx2);
          ROMol::OEDGE_ITER beg2, end2;
          boost::tie(beg2, end2) = mol.getAtomBonds(kAtom);
          while (beg2 != end2) {
            const Bond *tBond2 = mol[*beg2].get();
            if ((tBond2 != bond) && (tBond2 != tBond1)) {
              int idx4 = tBond2->getOtherAtomIdx(idx3);
              // make sure this isn't a three-membered ring:
              if (idx4 != idx1) {
                // we now have a torsion involving atoms (bonds):
                //  bIdx - (tBond1) - idx1 - (bond) - idx2 - (tBond2) - eIdx
                unsigned int torType;
                MMFFTor mmffTorParams;
                if (mmffMolProperties->getMMFFTorsionParams(
                        mol, idx1, idx2, idx3, idx4, torType, mmffTorParams)) {
                  TorsionAngleContrib *contrib = new TorsionAngleContrib(
                      field, idx1, idx2, idx3, idx4, &mmffTorParams);
                  field->contribs().push_back(ForceFields::ContribPtr(contrib));
                  if (mmffMolProperties->getMMFFVerbosity()) {
                    const Atom *iAtom = mol.getAtomWithIdx(idx1);
                    const Atom *lAtom = mol.getAtomWithIdx(idx4);
                    unsigned int iAtomType =
                        mmffMolProperties->getMMFFAtomType(idx1);
                    unsigned int jAtomType =
                        mmffMolProperties->getMMFFAtomType(idx2);
                    unsigned int kAtomType =
                        mmffMolProperties->getMMFFAtomType(idx3);
                    unsigned int lAtomType =
                        mmffMolProperties->getMMFFAtomType(idx4);
                    const RDGeom::Point3D p1((*(points[idx1]))[0],
                                             (*(points[idx1]))[1],
                                             (*(points[idx1]))[2]);
                    const RDGeom::Point3D p2((*(points[idx2]))[0],
                                             (*(points[idx2]))[1],
                                             (*(points[idx2]))[2]);
                    const RDGeom::Point3D p3((*(points[idx3]))[0],
                                             (*(points[idx3]))[1],
                                             (*(points[idx3]))[2]);
                    const RDGeom::Point3D p4((*(points[idx4]))[0],
                                             (*(points[idx4]))[1],
                                             (*(points[idx4]))[2]);
                    const double cosPhi =
                        MMFF::Utils::calcTorsionCosPhi(p1, p2, p3, p4);
                    const double torsionEnergy = MMFF::Utils::calcTorsionEnergy(
                        mmffTorParams.V1, mmffTorParams.V2, mmffTorParams.V3,
//.........这里部分代码省略.........