C++ forcefields::ForceField类代码示例

void testUFFBuilderSpecialCases(){
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing UFF special cases." << std::endl;

  RWMol *mol;
  std::string key;
  int needMore;
  RDGeom::Point3D v1,v2;
  ForceFields::ForceField *field;

  std::string pathName=getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/UFF/test_data";
  // ----------
  //  Trigonal bipyramid
  // ----------
  mol = MolFileToMol(pathName+"/tbp.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  const Conformer &conf = mol->getConformer();
  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize(200,1e-8,1e-4);
  TEST_ASSERT(!needMore);
  v1 = conf.getAtomPos(0).directionVector(conf.getAtomPos(1));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(2));
  TEST_ASSERT(feq(v1.dotProduct(v2),-1.0,1e-3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(3));
  TEST_ASSERT(feq(v1.dotProduct(v2),0.0,1e-3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(4));
  TEST_ASSERT(feq(v1.dotProduct(v2),0.0,1e-3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(5));
  TEST_ASSERT(feq(v1.dotProduct(v2),0.0,1e-3));

  v1 = conf.getAtomPos(0).directionVector(conf.getAtomPos(2));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(3));
  TEST_ASSERT(feq(v1.dotProduct(v2),0.0,1e-3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(4));
  TEST_ASSERT(feq(v1.dotProduct(v2),0.0,1e-3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(5));
  TEST_ASSERT(feq(v1.dotProduct(v2),0.0,1e-3));

  v1 = conf.getAtomPos(0).directionVector(conf.getAtomPos(3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(4));
  TEST_ASSERT(feq(v1.dotProduct(v2),-0.5,1e-3));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(5));
  TEST_ASSERT(feq(v1.dotProduct(v2),-0.5,1e-3));

  v1 = conf.getAtomPos(0).directionVector(conf.getAtomPos(4));
  v2 = conf.getAtomPos(0).directionVector(conf.getAtomPos(5));
  TEST_ASSERT(feq(v1.dotProduct(v2),-0.5,1e-3));

  
  delete mol;
  delete field;
  

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

void runMol(ROMol *mol,int checkEvery=10,bool verbose=true){
  ForceFields::ForceField *field;

  std::cout << MolToMolBlock(*mol) << "$$$$" << std::endl;

  try{
    field=UFF::constructForceField(*mol,2.5);
  } catch (...) {
    field=0;
  }
  if(field){
    field->initialize();
    int needMore=1;
    int nPasses=0;
    while(needMore){
#if 1
      needMore = field->minimize(checkEvery);
      if(verbose) std::cerr << "\t" << ++nPasses << std::endl;
#else
      needMore = field->minimize(1);
      std::cout << MolToMolBlock(mol) << "$$$$" << std::endl;
#endif
    }
    std::cout << MolToMolBlock(*mol) << "$$$$" << std::endl;
    delete field;
  } else {
    std::cerr << "failed";
  }
  
}

void testIssue242()
{
// FIX: Changes to the forcefield (connected to Issue 408) have
// made it so that this particular problem no longer manifests
// in this molecule/embedding. A new test case is needed.
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing Issue242." << std::endl;

  RWMol *mol, *mol2;
  int needMore;
  ForceFields::ForceField *field = 0, *field2 = 0;
  std::string mb1,mb2;
  double e1,e2;

  std::string pathName = getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/MMFF/test_data";

  mol = MolFileToMol(pathName + "/Issue242-2.mol");
  TEST_ASSERT(mol);

  mol2 = MolFileToMol(pathName + "/Issue242-2.mol");
  TEST_ASSERT(mol2);

  TEST_ASSERT(DGeomHelpers::EmbedMolecule(*mol2, 30, 2370) >= 0);
  mb1 = MolToMolBlock(*mol);
  mb2 = MolToMolBlock(*mol2);

  //std::cout << mb1 << "------\n";
  //std::cout << mb2 << "------\n";

  field = MMFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  field2 = MMFF::constructForceField(*mol2);
  TEST_ASSERT(field2);
  field2->initialize();
  e1 = field->calcEnergy();
  e2 = field2->calcEnergy();
  BOOST_LOG(rdInfoLog) << "E1: " << e1 << std::endl;
  BOOST_LOG(rdInfoLog) << "E2: " << e2 << std::endl;
  //TEST_ASSERT(feq(e2,e1,0.1));

  needMore = field->minimize(200,1.0e-4);
  needMore = field2->minimize(200,1.0e-4);
  e1 = field->calcEnergy();
  e2 = field2->calcEnergy();
  BOOST_LOG(rdInfoLog) << "E1: " << e1 << std::endl;
  BOOST_LOG(rdInfoLog) << "E2: " << e2 << std::endl;
  TEST_ASSERT(feq(e2, e1, 0.1));

  delete mol;
  delete mol2;
  delete field;
  delete field2;
  

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

void testUFFMultiThread(){
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Test UFF multithreading" << std::endl;

  ForceFields::ForceField *field;

  std::string pathName = getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/UFF/test_data";
  SDMolSupplier suppl(pathName + "/bulk.sdf");
  std::vector<ROMol *> mols;
  while(!suppl.atEnd()&&mols.size()<100){
    ROMol *mol=0;
    try{
      mol=suppl.next();
    } catch(...){
      continue;
    }
    if(!mol) continue;
    mols.push_back(mol);
  }

  std::cerr<<"generating reference data"<<std::endl;
  std::vector<double> energies(mols.size(),0.0);
  for(unsigned int i=0;i<mols.size();++i){
    ROMol mol(*mols[i]);
    ForceFields::ForceField *field = 0;
    try {
      field = UFF::constructForceField(mol);
    } catch (...) {
      field = 0;
    }
    TEST_ASSERT(field);
    field->initialize();
    int failed = field->minimize(500);
    TEST_ASSERT(!failed);
    energies[i]=field->calcEnergy();
    delete field;
  }
  
  boost::thread_group tg;

  std::cerr<<"processing"<<std::endl;
  unsigned int count=4;
  for(unsigned int i=0;i<count;++i){
    std::cerr<<" launch :"<<i<<std::endl;std::cerr.flush();
    tg.add_thread(new boost::thread(runblock_uff,mols,energies,count,i));
  }
  tg.join_all();
  
  BOOST_FOREACH(ROMol *mol,mols){
    delete mol;
  }
  BOOST_LOG(rdErrorLog) << "  done" << std::endl;
}

ForceFields::ForceField *construct3DImproperForceField(
    const BoundsMatrix &mmat, RDGeom::Point3DPtrVect &positions,
    const std::vector<std::vector<int> > &improperAtoms,
    const std::vector<int> &atomNums) {
  (void)atomNums;
  unsigned int N = mmat.numRows();
  CHECK_INVARIANT(N == positions.size(), "");
  ForceFields::ForceField *field =
      new ForceFields::ForceField(positions[0]->dimension());
  for (unsigned int i = 0; i < N; ++i) {
    field->positions().push_back(positions[i]);
  }

  // improper torsions / out-of-plane bend / inversion
  double oobForceScalingFactor = 10.0;
  for (unsigned int t = 0; t < improperAtoms.size(); ++t) {
    std::vector<int> n(4);
    for (unsigned int i = 0; i < 3; ++i) {
      n[1] = 1;
      switch (i) {
        case 0:
          n[0] = 0;
          n[2] = 2;
          n[3] = 3;
          break;

        case 1:
          n[0] = 0;
          n[2] = 3;
          n[3] = 2;
          break;

        case 2:
          n[0] = 2;
          n[2] = 3;
          n[3] = 0;
          break;
      }
      ForceFields::UFF::InversionContrib *contrib =
          new ForceFields::UFF::InversionContrib(
              field, improperAtoms[t][n[0]], improperAtoms[t][n[1]],
              improperAtoms[t][n[2]], improperAtoms[t][n[3]],
              improperAtoms[t][4], improperAtoms[t][5], oobForceScalingFactor);
      field->contribs().push_back(ForceFields::ContribPtr(contrib));
    }
  }

  return field;

} // construct3DImproperForceField

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 testGithub308()
{
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing github 308: crash during MMFF parameterization ." << std::endl;
  ROMol *mol = SmilesToMol("FF");
  TEST_ASSERT(DGeomHelpers::EmbedMolecule(*mol) >= 0);
  int needMore;
  ForceFields::ForceField *field = 0;
  TEST_ASSERT(mol);
  MMFF::MMFFMolProperties mmffMolProperties(*mol);
  TEST_ASSERT(mmffMolProperties.isValid());
  field = MMFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize(200, 1.0e-6, 1.0e-3);
  TEST_ASSERT(!needMore);
}

 int MMFFOptimizeMolecule(ROMol &mol, std::string mmffVariant = "MMFF94",
   int maxIters = 200, double nonBondedThresh = 100.0, int confId = -1,
   bool ignoreInterfragInteractions = true)
 {
   int res = -1;
   
   MMFF::MMFFMolProperties mmffMolProperties(mol, mmffVariant);
   if (mmffMolProperties.isValid()) {
     ForceFields::ForceField *ff = MMFF::constructForceField(mol,
       &mmffMolProperties, nonBondedThresh, confId, ignoreInterfragInteractions);
     ff->initialize();
     res = ff->minimize(maxIters);
     delete ff;
   }
   
   return res;
 }

void testMMFFBatch() {
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog)
      << "    Testing bulk MMFF (regression to check that things run)."
      << std::endl;

  ROMol *mol;

  ForceFields::ForceField *field;

  std::string pathName = getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/MMFF/test_data";
  SDMolSupplier suppl(pathName + "/bulk.sdf");

  int count = 0;
  mol = suppl.next();
  while (mol && (!(suppl.atEnd()))) {
    ++count;
    std::string origMolBlock = MolToMolBlock(*mol);

    BOOST_LOG(rdErrorLog) << "Mol:" << count << std::endl;
    try {
      field = MMFF::constructForceField(*mol);
    } catch (...) {
      field = 0;
    }
    if (field) {
      field->initialize();
      int failed = field->minimize(500);
      if (failed) {
        BOOST_LOG(rdErrorLog) << " not converged (code = " << failed << ")"
                              << std::endl;
        std::cout << origMolBlock << "$$$$" << std::endl;
        std::cout << MolToMolBlock(*mol) << "$$$$" << std::endl;
      }
      delete field;
    }
    delete mol;
    mol = suppl.next();
  }

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

void testIssue239(){
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing Issue239." << std::endl;

  RWMol *mol;
  int needMore;
  ForceFields::ForceField *field;
  double e1,e2;

  std::string pathName=getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/UFF/test_data";
  mol = MolFileToMol(pathName+"/Issue239.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize(200,1e-6,1e-3);
  e1 = field->calcEnergy();
  needMore = field->minimize(200,1e-6,1e-3);
  e2 = field->calcEnergy();
  TEST_ASSERT(feq(e2,e1,0.1));
  
  delete mol;
  delete field;
  

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

void testGitHubIssue62() {
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing GitHubIssue62." << std::endl;

  std::string pathName=getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/UFF/test_data";
  {
    double energyValues[] =
      { 38.687, 174.698, 337.986, 115.248,
        2.482, 1.918, 10.165, 98.469, 39.078,
        267.236, 15.747, 202.121, 205.539,
        20.044, 218.986, 79.627 };
    SmilesMolSupplier smiSupplier(pathName + "/Issue62.smi");
    SDWriter *sdfWriter = new SDWriter(pathName + "/Issue62.sdf");
    for (unsigned int i = 0; i < smiSupplier.length(); ++i) {
      ROMol *mol = MolOps::addHs(*(smiSupplier[i]));
      TEST_ASSERT(mol);
      std::string molName = "";
      if (mol->hasProp(common_properties::_Name)) {
        mol->getProp(common_properties::_Name, molName);
      }
      DGeomHelpers::EmbedMolecule(*mol);
      ForceFields::ForceField *field = UFF::constructForceField(*mol);
      TEST_ASSERT(field);
      field->initialize();
      int needMore = field->minimize(200, 1.e-6, 1.e-3);
      TEST_ASSERT(!needMore);
      sdfWriter->write(*mol);
      double e = field->calcEnergy();
      BOOST_LOG(rdErrorLog) << molName << " " << e << std::endl;
      TEST_ASSERT(fabs(e - energyValues[i]) < 1.);
    }
    sdfWriter->close();
    BOOST_LOG(rdErrorLog) << "  done" << std::endl;
  }
}

void testSFIssue1653802(){
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing SFIssue1653802." << std::endl;

  RWMol *mol;
  int needMore;
  ForceFields::ForceField *field;

  std::string pathName=getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/UFF/test_data";

  mol = MolFileToMol(pathName+"/cyclobutadiene.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);


  UFF::AtomicParamVect types;
  bool foundAll;
  boost::shared_array<boost::uint8_t> nbrMat;
  boost::tie(types,foundAll)=UFF::getAtomTypes(*mol);
  TEST_ASSERT(foundAll);
  TEST_ASSERT(types.size()==mol->getNumAtoms());
  field=new ForceFields::ForceField();
  UFF::Tools::addBonds(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==8);

  nbrMat = UFF::Tools::buildNeighborMatrix(*mol);
  UFF::Tools::addAngles(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==20);
  UFF::Tools::addTorsions(*mol,types,field);
  //std::cout << field->contribs().size() << std::endl;
  TEST_ASSERT(field->contribs().size()==36);
  UFF::Tools::addNonbonded(*mol,0,types,field,nbrMat);
  delete field;

  field = UFF::constructForceField(*mol);
  field->initialize();
  needMore = field->minimize(200,1e-6,1e-3);
  TEST_ASSERT(!needMore);
  
  delete mol;
  delete field;
  BOOST_LOG(rdErrorLog) << "  done" << std::endl;
}

void testSFIssue1653802() {
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing SFIssue1653802." << std::endl;

  RWMol *mol;
  int needMore;
  ForceFields::ForceField *field;

  std::string pathName = getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/MMFF/test_data";

  mol = MolFileToMol(pathName + "/cyclobutadiene.mol", false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);
  MMFF::MMFFMolProperties *mmffMolProperties =
      new MMFF::MMFFMolProperties(*mol);
  TEST_ASSERT(mmffMolProperties);

  boost::shared_array<boost::uint8_t> nbrMat;
  field = new ForceFields::ForceField();
  MMFF::Tools::addBonds(*mol, mmffMolProperties, field);
  TEST_ASSERT(field->contribs().size() == 8);

  nbrMat = MMFF::Tools::buildNeighborMatrix(*mol);
  MMFF::Tools::addAngles(*mol, mmffMolProperties, field);
  TEST_ASSERT(field->contribs().size() == 20);
  MMFF::Tools::addTorsions(*mol, mmffMolProperties, field);
  // std::cout << field->contribs().size() << std::endl;
  TEST_ASSERT(field->contribs().size() == 36);
  MMFF::Tools::addVdW(*mol, 0, mmffMolProperties, field, nbrMat);
  delete field;

  field = MMFF::constructForceField(*mol);
  field->initialize();
  needMore = field->minimize(200, 1.0e-6, 1.0e-3);
  TEST_ASSERT(!needMore);

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

void testUFFBuilder2(){
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing UFF builder+minimization." << std::endl;

  RWMol *mol;
  std::string key;
  int needMore;

  ForceFields::ForceField *field;

  std::string pathName=getenv("RDBASE");
  pathName += "/Code/GraphMol/ForceFieldHelpers/UFF/test_data";
  mol = MolFileToMol(pathName+"/small1.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize();
  TEST_ASSERT(!needMore);
  //std::cout << MolToMolBlock(mol);
  delete mol;
  delete field;

  mol = MolFileToMol(pathName+"/small2.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize(150);
  TEST_ASSERT(!needMore);
  //std::cout << MolToMolBlock(mol);
  delete mol;
  delete field;

  mol = MolFileToMol(pathName+"/benzene.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize();
  TEST_ASSERT(!needMore);
  //std::cout << MolToMolBlock(mol);
  delete mol;
  delete field;
  
  mol = MolFileToMol(pathName+"/toluene.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize();
  TEST_ASSERT(!needMore);
  //std::cout << MolToMolBlock(mol);
  delete mol;
  delete field;

  mol = MolFileToMol(pathName+"/complex1.mol",false);
  TEST_ASSERT(mol);
  MolOps::sanitizeMol(*mol);

  field=UFF::constructForceField(*mol);
  TEST_ASSERT(field);
  field->initialize();
  needMore = field->minimize();
  TEST_ASSERT(!needMore);
  //std::cout << MolToMolBlock(mol);
  delete mol;
  delete field;

  

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

void testUFFBuilder1(){
  BOOST_LOG(rdErrorLog) << "-------------------------------------" << std::endl;
  BOOST_LOG(rdErrorLog) << "    Testing UFF builder tools." << std::endl;

  ROMol *mol,*mol2;
  std::string key;

  UFF::AtomicParamVect types;
  bool foundAll;
  ForceFields::ForceField *field;
  boost::shared_array<boost::uint8_t> nbrMat;

  mol = SmilesToMol("CC(O)C");
  Conformer *conf = new Conformer(mol->getNumAtoms());
  int cid = static_cast<int>(mol->addConformer(conf, true));
  TEST_ASSERT(mol);
  boost::tie(types,foundAll)=UFF::getAtomTypes(*mol);
  TEST_ASSERT(foundAll);
  TEST_ASSERT(types.size()==mol->getNumAtoms());
  field=new ForceFields::ForceField();
  UFF::Tools::addBonds(*mol,types,field);

  TEST_ASSERT(field->contribs().size()==3);

  nbrMat = UFF::Tools::buildNeighborMatrix(*mol);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,0)==UFF::Tools::RELATION_1_X);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,1)==UFF::Tools::RELATION_1_2);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,2)==UFF::Tools::RELATION_1_3);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,3)==UFF::Tools::RELATION_1_3);

  UFF::Tools::addAngles(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==6);

  // there are no non-bonded terms here:
  UFF::Tools::addNonbonded(*mol,cid,types,field,nbrMat);
  TEST_ASSERT(field->contribs().size()==6);

  // and no torsions either, until we add Hs:
  UFF::Tools::addTorsions(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==6);

  delete mol;
  delete field;
  mol = SmilesToMol("CCOC");
  Conformer *conf2 = new Conformer(mol->getNumAtoms());
  cid = static_cast<int>(mol->addConformer(conf2, true));
  TEST_ASSERT(mol);
  boost::tie(types,foundAll)=UFF::getAtomTypes(*mol);
  TEST_ASSERT(foundAll);
  TEST_ASSERT(types.size()==mol->getNumAtoms());
  field=new ForceFields::ForceField();
  UFF::Tools::addBonds(*mol,types,field);

  TEST_ASSERT(field->contribs().size()==3);

  nbrMat = UFF::Tools::buildNeighborMatrix(*mol);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,0)==UFF::Tools::RELATION_1_X);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,1)==UFF::Tools::RELATION_1_2);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,2)==UFF::Tools::RELATION_1_3);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,3)==UFF::Tools::RELATION_1_X);

  UFF::Tools::addAngles(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==5);
  UFF::Tools::addNonbonded(*mol,cid,types,field,nbrMat);
  TEST_ASSERT(field->contribs().size()==6);
  UFF::Tools::addTorsions(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==7);



  delete mol;
  delete field;
  mol = SmilesToMol("CO");
  Conformer *conf3 = new Conformer(mol->getNumAtoms());
  cid = static_cast<int>(mol->addConformer(conf3, true));
  TEST_ASSERT(mol);
  boost::tie(types,foundAll)=UFF::getAtomTypes(*mol);
  TEST_ASSERT(foundAll);
  TEST_ASSERT(types.size()==mol->getNumAtoms());

  field=new ForceFields::ForceField();
  UFF::Tools::addBonds(*mol,types,field);

  TEST_ASSERT(field->contribs().size()==1);

  nbrMat = UFF::Tools::buildNeighborMatrix(*mol);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,0)==UFF::Tools::RELATION_1_X);
  TEST_ASSERT(UFF::Tools::getTwoBitCell(nbrMat,1)==UFF::Tools::RELATION_1_2);

  UFF::Tools::addAngles(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==1);
  UFF::Tools::addNonbonded(*mol,cid,types,field,nbrMat);
  TEST_ASSERT(field->contribs().size()==1);
  UFF::Tools::addTorsions(*mol,types,field);
  TEST_ASSERT(field->contribs().size()==1);

  
  mol2 = MolOps::addHs(*mol);
  TEST_ASSERT(mol2->getNumAtoms()==6);
  delete field;
//.........这里部分代码省略.........