Context Navigation

-              r503acc1
+              r0cd8cf
 #include "Fragmentation/Automation/FragmentationChargeDensity.hpp"
 #include "Fragmentation/Automation/FragmentJobQueue.hpp"
+#include "Fragmentation/Automation/FragmentationResults.hpp"
+#include "Fragmentation/Automation/FragmentationLongRangeResults.hpp"
+#include "Fragmentation/Automation/FragmentationShortRangeResults.hpp"
 #include "Fragmentation/Automation/MPQCFragmentController.hpp"
 #include "Fragmentation/Automation/VMGDebugGridFragmentController.hpp"
 …
+}
-/** Helper function to get number of atoms somehow.
+ *
- * Here, we just parse the number of lines in the adjacency file as
- * it should correspond to the number of atoms, except when some atoms
- * are not bonded, but then fragmentation makes no sense.
+ *
- * @param path path to the adjacency file
- */
-size_t getNoAtomsFromAdjacencyFile(const std::string &path)
+{
-  size_t NoAtoms = 0;
-  // parse in special file to get atom count (from line count)
-  std::string filename(path);
-  filename += FRAGMENTPREFIX;
-  filename += ADJACENCYFILE;
-  std::ifstream adjacency(filename.c_str());
-  if (adjacency.fail()) {
-    LOG(0, endl << "getNoAtomsFromAdjacencyFile() - Unable to open " << filename << ", is the directory correct?");
-    return false;
+  }
-  std::string buffer;
-  while (getline(adjacency, buffer))
-    NoAtoms++;
-  LOG(1, "INFO: There are " << NoAtoms << " atoms.");
-  return NoAtoms;
+}
-/** Place results from FragmentResult into EnergyMatrix and ForceMatrix.
+ *
- * @param fragmentData MPQCData resulting from the jobs
- * @param MatrixNrLookup Lookup up-map from job id to fragment number
- * @param FragmentCounter total number of fragments
- * @param NoAtoms total number of atoms
- * @param Energy energy matrix to be filled on return
- * @param Force force matrix to be filled on return
- * @return true - everything ok, false - else
- */
-bool putResultsintoMatrices(
-    const std::map<JobId_t, MPQCData> &fragmentData,
-    const std::map< JobId_t, size_t > &MatrixNrLookup,
-    const size_t FragmentCounter,
-    const size_t NoAtoms,
-    EnergyMatrix &Energy,
-    ForceMatrix &Force)
+{
-  for (std::map<JobId_t, MPQCData>::const_iterator dataiter = fragmentData.begin();
-      dataiter != fragmentData.end(); ++dataiter) {
-    const MPQCData &extractedData = dataiter->second;
-    const JobId_t &jobid = dataiter->first;
-    std::map< JobId_t, size_t >::const_iterator nriter = MatrixNrLookup.find(jobid);
-    ASSERT( nriter != MatrixNrLookup.end(),
-        "putResultsintoMatrices() - MatrixNrLookup does not contain id "
-        +toString(jobid)+".");
-    // place results into EnergyMatrix ...
+    {
-      MatrixContainer::MatrixArray matrix;
-      matrix.resize(1);
-      matrix[0].resize(1, extractedData.energies.total);
-      if (!Energy.AddMatrix(
-          std::string("MPQCJob ")+toString(jobid),
-          matrix,
-          nriter->second)) {
-        ELOG(1, "Adding energy matrix failed.");
-        return false;
+      }
+    }
-    // ... and ForceMatrix (with two empty columns in front)
+    {
-      MatrixContainer::MatrixArray matrix;
-      const size_t rows = extractedData.forces.size();
-      matrix.resize(rows);
-      for (size_t i=0;i<rows;++i) {
-        const size_t columns = 2+extractedData.forces[i].size();
-        matrix[i].resize(columns, 0.);
-  //      for (size_t j=0;j<2;++j)
-  //        matrix[i][j] = 0.;
-        for (size_t j=2;j<columns;++j)
-          matrix[i][j] = extractedData.forces[i][j-2];
+      }
-      if (!Force.AddMatrix(
-          std::string("MPQCJob ")+toString(jobid),
-          matrix,
-          nriter->second)) {
-        ELOG(1, "Adding force matrix failed.");
-        return false;
+      }
+    }
+  }
-  // add one more matrix (not required for energy)
-  MatrixContainer::MatrixArray matrix;
-  matrix.resize(1);
-  matrix[0].resize(1, 0.);
-  if (!Energy.AddMatrix(std::string("MPQCJob total"), matrix, FragmentCounter))
-    return false;
-  // but for energy because we need to know total number of atoms
-  matrix.resize(NoAtoms);
-  for (size_t i = 0; i< NoAtoms; ++i)
-    matrix[i].resize(2+NDIM, 0.);
-  if (!Force.AddMatrix(std::string("MPQCJob total"), matrix, FragmentCounter))
-    return false;
-  return true;
+}
-/** Print MPQCData from received results.
+ *
- * @param fragmentData MPQCData resulting from the jobs, associated to job id
- * @param MatrixNrLookup map with enumerated jobids
- * @param KeySet KeySets of all (non-hydrogen) atoms
- * @param ForceKeySet KeySets of all atoms except those added by saturation
- * @param NoAtoms total number of atoms
- */
-bool printReceivedMPQCResults(
-    const std::map<JobId_t, MPQCData> &fragmentData,
-    const std::map< JobId_t, size_t > MatrixNrLookup,
-    KeySetsContainer KeySet,
-    const KeySetsContainer &ForceKeySet,
-    size_t NoAtoms)
+{
-  size_t FragmentCounter = MatrixNrLookup.size();
-  // place results into maps
-  EnergyMatrix Energy;
-  ForceMatrix Force;
-  if (!putResultsintoMatrices(fragmentData, MatrixNrLookup, FragmentCounter, NoAtoms, Energy, Force))
-    return false;
-  if (!Energy.InitialiseIndices()) return false;
-  if (!Force.ParseIndices(ForceKeySet)) return false;
-  // combine all found data
-  if (!KeySet.ParseManyBodyTerms()) return false;
-  EnergyMatrix EnergyFragments;
-  ForceMatrix ForceFragments;
-  if (!EnergyFragments.AllocateMatrix(Energy.Header, Energy.MatrixCounter, Energy.RowCounter, Energy.ColumnCounter)) return false;
-  if (!ForceFragments.AllocateMatrix(Force.Header, Force.MatrixCounter, Force.RowCounter, Force.ColumnCounter)) return false;
-  if(!Energy.SetLastMatrix(0., 0)) return false;
-  if(!Force.SetLastMatrix(0., 2)) return false;
-  for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
-    // --------- sum up energy --------------------
-    LOG(1, "INFO: Summing energy of order " << BondOrder+1 << " ...");
-    if (!EnergyFragments.SumSubManyBodyTerms(Energy, KeySet, BondOrder)) return false;
-    if (!Energy.SumSubEnergy(EnergyFragments, NULL, KeySet, BondOrder, 1.)) return false;
-    // --------- sum up Forces --------------------
-    LOG(1, "INFO: Summing forces of order " << BondOrder+1 << " ...");
-    if (!ForceFragments.SumSubManyBodyTerms(Force, KeySet, BondOrder)) return false;
-    if (!Force.SumSubForces(ForceFragments, KeySet, BondOrder, 1.)) return false;
+  }
-  // for debugging print resulting energy and forces
-  LOG(1, "INFO: Resulting energy is " << Energy.Matrix[ FragmentCounter ][0][0]);
-  std::stringstream output;
-  for (int i=0; i< Force.RowCounter[FragmentCounter]; ++i) {
-    for (int j=0; j< Force.ColumnCounter[FragmentCounter]; ++j)
-      output << Force.Matrix[ FragmentCounter ][i][j] << " ";
-    output << "\n";
+  }
-  LOG(1, "INFO: Resulting forces are " << std::endl << output.str());
-  return true;
+}
-/** Print MPQCData from received results.
+ *
- * @param fragmentData MPQCData resulting from the jobs, associated to job id
- * @param KeySetFilename filename with keysets to associate forces correctly
- * @param NoAtoms total number of atoms
- */
-bool printReceivedMPQCResults(
-    const std::map<JobId_t, MPQCData> &fragmentData,
-    const std::string &KeySetFilename,
-    size_t NoAtoms)
+{
-  // create a vector of all job ids
-  std::vector<JobId_t> jobids;
-  std::transform(fragmentData.begin(),fragmentData.end(),
-      std::back_inserter(jobids),
-      boost::bind( &std::map<JobId_t,MPQCData>::value_type::first, boost::lambda::_1 )
-  );
-  // create lookup from job nr to fragment number
-  size_t FragmentCounter = 0;
-  const std::map< JobId_t, size_t > MatrixNrLookup =
-      createMatrixNrLookup(jobids, FragmentCounter);
-  // initialise keysets
-  KeySetsContainer KeySet;
-  parseKeySetFile(KeySet, KeySetFilename, FragmentCounter, NonHydrogenKeySets);
-  KeySetsContainer ForceKeySet;
-  parseKeySetFile(ForceKeySet, KeySetFilename, FragmentCounter, HydrogenKeySets);
-  return printReceivedMPQCResults(
-      fragmentData,
-      MatrixNrLookup,
-      KeySet,
-      ForceKeySet,
-      NoAtoms);
+}
-/** Print MPQCData from received results.
+ *
- * @param fragmentData MPQCData resulting from the jobs, associated to job id
- * @param KeySet KeySets of all (non-hydrogen) atoms
- * @param ForceKeySet KeySets of all atoms except those added by saturation
- * @param NoAtoms total number of atoms
- */
-bool printReceivedMPQCResults(
-    const std::map<JobId_t, MPQCData> &fragmentData,
-    const KeySetsContainer &KeySet,
-    const KeySetsContainer &ForceKeySet,
-    size_t NoAtoms)
+{
-  // create a vector of all job ids
-  std::vector<JobId_t> jobids;
-  std::transform(fragmentData.begin(),fragmentData.end(),
-      std::back_inserter(jobids),
-      boost::bind( &std::map<JobId_t,MPQCData>::value_type::first, boost::lambda::_1 )
-  );
-  // create lookup from job nr to fragment number
-  size_t FragmentCounter = 0;
-  const std::map< JobId_t, size_t > MatrixNrLookup =
-      createMatrixNrLookup(jobids, FragmentCounter);
-  return printReceivedMPQCResults(
-      fragmentData,
-      MatrixNrLookup,
-      KeySet,
-      ForceKeySet,
-      NoAtoms);
+}
 void writeToFile(const std::string &filename, const std::string contents)
+{
 …
+}
 /** Print MPQCData from received results.
+/** Print (short range) energy, forces, and timings from received results.
+ *
  * @param results summed up results container
  */
 void printReceivedFullResults(
     const FragmentationResults &results)
+void printReceivedShortResults(
+    const FragmentationShortRangeResults &results)
+{
   // print tables (without eigenvalues, they go extra)
 …
     filename += FRAGMENTPREFIX + std::string("_Energy.dat");
     writeToFile(filename, energyresult);
+  }
-  if (results.Result_LongRange_fused.size() >= (results.getMaxLevel()-2))
+  {
-    const std::string gridresult =
-        writeTable<VMGDataMap_t, VMGDataVector_t >()(
-            results.Result_LongRange_fused, results.getMaxLevel(), 2);
-    LOG(0, "VMG table is \n" << gridresult);
-    std::string filename;
-    filename += FRAGMENTPREFIX + std::string("_VMGEnergy.dat");
-    writeToFile(filename, gridresult);
+  }
-  if (results.Result_LongRangeIntegrated_fused.size() >= (results.getMaxLevel()-2))
+  {
-    const std::string gridresult =
-        writeTable<VMGDataLongRangeMap_t, VMGDataLongRangeVector_t >()(
-            results.Result_LongRangeIntegrated_fused, results.getMaxLevel(), 2);
-    LOG(0, "LongRange table is \n" << gridresult);
-    std::string filename;
-    filename += FRAGMENTPREFIX + std::string("_LongRangeEnergy.dat");
-    writeToFile(filename, gridresult);
+  }
 …
+}
+/** Print long range energy from received results.
+ *
+ * @param results summed up results container
+ */
+void printReceivedFullResults(
+    const FragmentationLongRangeResults &results)
+{
+  // print tables (without eigenvalues, they go extra)
+  if (results.Result_LongRange_fused.size() >= (results.getMaxLevel()-2))
+  {
+    const std::string gridresult =
+        writeTable<VMGDataMap_t, VMGDataVector_t >()(
+            results.Result_LongRange_fused, results.getMaxLevel(), 2);
+    LOG(0, "VMG table is \n" << gridresult);
+    std::string filename;
+    filename += FRAGMENTPREFIX + std::string("_VMGEnergy.dat");
+    writeToFile(filename, gridresult);
+  }
+  if (results.Result_LongRangeIntegrated_fused.size() >= (results.getMaxLevel()-2))
+  {
+    const std::string gridresult =
+        writeTable<VMGDataLongRangeMap_t, VMGDataLongRangeVector_t >()(
+            results.Result_LongRangeIntegrated_fused, results.getMaxLevel(), 2);
+    LOG(0, "LongRange table is \n" << gridresult);
+    std::string filename;
+    filename += FRAGMENTPREFIX + std::string("_LongRangeEnergy.dat");
+    writeToFile(filename, gridresult);
+  }
+}
 bool appendToHomologyFile(
     const boost::filesystem::path &homology_file,
+    const FragmentationResults &results)
+    const FragmentationShortRangeResults &shortrangeresults,
+    const FragmentationLongRangeResults &longrangeresults
+    )
+{
   /// read homology container (if present)
 …
   /// append all fragments to a HomologyContainer
   HomologyContainer::container_t values;
   const size_t FragmentCounter = results.Result_perIndexSet_Energy.size();
+  const size_t FragmentCounter = shortrangeresults.Result_perIndexSet_Energy.size();
   // convert KeySetContainer to IndexSetContainer
   IndexSetContainer::ptr ForceContainer(new IndexSetContainer(results.getForceKeySet()));
   const IndexSetContainer::Container_t &Indices = results.getContainer();
+  IndexSetContainer::ptr ForceContainer(new IndexSetContainer(shortrangeresults.getForceKeySet()));
+  const IndexSetContainer::Container_t &Indices = shortrangeresults.getContainer();
   const IndexSetContainer::Container_t &ForceIndices = ForceContainer->getContainer();
   IndexSetContainer::Container_t::const_iterator iter = Indices.begin();
 …
     // obtain fragment as key
     std::map<IndexSet::ptr, std::pair< MPQCDataFragmentMap_t, MPQCDataFragmentMap_t> >::const_iterator fragmentiter
         = results.Result_perIndexSet_Fragment.find(index);
     ASSERT( fragmentiter != results.Result_perIndexSet_Fragment.end(),
+        = longrangeresults.Result_perIndexSet_Fragment.find(index);
+    ASSERT( fragmentiter != longrangeresults.Result_perIndexSet_Fragment.end(),
         "appendToHomologyFile() - cannot find index "+toString(*index)
         +" in FragmentResults.");
 …
     // obtain energy as value
     std::map<IndexSet::ptr, std::pair< MPQCDataEnergyMap_t, MPQCDataEnergyMap_t> >::const_iterator energyiter
         = results.Result_perIndexSet_Energy.find(index);
     ASSERT( energyiter != results.Result_perIndexSet_Energy.end(),
+        = shortrangeresults.Result_perIndexSet_Energy.find(index);
+    ASSERT( energyiter != shortrangeresults.Result_perIndexSet_Energy.end(),
         "appendToHomologyFile() - cannot find index "+toString(*index)
         +" in FragmentResults.");
 …
+  }
+  FragmentationShortRangeResults shortrangeresults(
+      fragmentData,
+      FragmentJobQueue::getInstance().getKeySets(),
+      FragmentJobQueue::getInstance().getFullKeySets());
+  shortrangeresults(fragmentData);
+  printReceivedShortResults(shortrangeresults);
 #ifdef HAVE_VMG
   if (params.DoLongrange.get()) {
 …
   // Final phase: sum up and print result
   FragmentationResults results(
+  FragmentationLongRangeResults longrangeresults(
       fragmentData,
       longrangeData,
       FragmentJobQueue::getInstance().getKeySets(),
       FragmentJobQueue::getInstance().getFullKeySets());
   results(
+  longrangeresults(
       fragmentData,
       longrangeData,
       fullsolutionData,
       full_sample);
   printReceivedFullResults(results);
+  printReceivedFullResults(longrangeresults);
   // append all keysets to homology file
 …
     if (homology_file.string() != "") {
       LOG(1, "INFO: Appending HomologyGraphs to file " << homology_file.string() << ".");
       if (!appendToHomologyFile(homology_file, results))
+      if (!appendToHomologyFile(homology_file, shortrangeresults, longrangeresults))
         Exitflag = 1;
+    }
 …
+  }
+  }
-#else
-  // Final phase: print result
+  {
-    LOG(1, "INFO: Parsing fragment files from " << params.path.get() << ".");
-    printReceivedMPQCResults(
-        fragmentData,
-        FragmentJobQueue::getInstance().getKeySets(),
-        FragmentJobQueue::getInstance().getFullKeySets(),
-        World::getInstance().getAllAtoms().size()));
+  }
 #endif

Note: See TracChangeset for help on using the changeset viewer.

Context Navigation

Changeset 0cd8cf for src/Actions/FragmentationAction/FragmentationAutomationAction.cpp

Legend:

src/Actions/FragmentationAction/FragmentationAutomationAction.cpp

Download in other formats: