Changes in src/Actions/FragmentationAction/FragmentationAutomationAction.cpp [94d5ac6:db6f7d]

File:

• src/Actions/FragmentationAction/FragmentationAutomationAction.cpp (modified) (6 diffs)

src/Actions/FragmentationAction/FragmentationAutomationAction.cpp

@@ -39 +39 @@
 #include "CodePatterns/MemDebug.hpp"
 
+#include <boost/mpl/remove.hpp>
+#include <boost/lambda/lambda.hpp>
+
+#include "CodePatterns/Assert.hpp"
 #include "CodePatterns/Info.hpp"
 #include "CodePatterns/Log.hpp"
-#include "JobMarket/Controller/FragmentController.hpp"
 #include "JobMarket/Jobs/FragmentJob.hpp"
 
-#include "Atom/atom.hpp"
+#include "Fragmentation/Automation/createMatrixNrLookup.hpp"
+#include "Fragmentation/Automation/extractJobIds.hpp"
+#include "Fragmentation/Automation/FragmentationChargeDensity.hpp"
+#include "Fragmentation/Automation/FragmentationResults.hpp"
+#include "Fragmentation/Automation/MPQCFragmentController.hpp"
+#include "Fragmentation/Automation/VMGDebugGridFragmentController.hpp"
+#include "Fragmentation/Automation/VMGFragmentController.hpp"
 #include "Fragmentation/EnergyMatrix.hpp"
 #include "Fragmentation/ForceMatrix.hpp"
 #include "Fragmentation/Fragmentation.hpp"
+#include "Fragmentation/SetValues/Fragment.hpp"
+#include "Fragmentation/SetValues/Histogram.hpp"
+#include "Fragmentation/SetValues/IndexedVectors.hpp"
 #include "Fragmentation/HydrogenSaturation_enum.hpp"
+#include "Fragmentation/KeySet.hpp"
 #include "Fragmentation/KeySetsContainer.hpp"
+#include "Fragmentation/Summation/writeTable.hpp"
 #include "Graph/DepthFirstSearchAnalysis.hpp"
-#include "Jobs/MPQCCommandJob.hpp"
-#include "molecule.hpp"
+#include "Helpers/defs.hpp"
+#include "Jobs/MPQCJob.hpp"
+#include "Jobs/MPQCData.hpp"
+#include "Jobs/MPQCData_printKeyNames.hpp"
+#ifdef HAVE_VMG
+#include "Jobs/VMGDebugGridJob.hpp"
+#include "Jobs/VMGJob.hpp"
+#include "Jobs/VMGData.hpp"
+#include "Jobs/VMGDataFused.hpp"
+#include "Jobs/VMGDataMap.hpp"
+#include "Jobs/VMGData_printKeyNames.hpp"
+#endif
 #include "World.hpp"
 
+#include <fstream>
 #include <iostream>
 #include <string>
 #include <vector>
+
+#include <boost/mpl/for_each.hpp>
 
 #include "Actions/FragmentationAction/FragmentationAutomationAction.hpp"
@@ -74 +101 @@
 {
   return NULL;
-}
-
-/** Creates a MPQCCommandJob with argument \a filename.
- *
- * @param jobs created job is added to this vector
- * @param command mpqc command to execute
- * @param filename filename being argument to job
- * @param nextid id for this job
- */
-void parsejob(
-    std::vector<FragmentJob::ptr> &jobs,
-    const std::string &command,
-    const std::string &filename,
-    const JobId_t nextid)
-{
-  std::ifstream file;
-  file.open(filename.c_str());
-  ASSERT( file.good(), "parsejob() - file "+filename+" does not exist.");
-  std::string output((std::istreambuf_iterator<char>(file)),
-      std::istreambuf_iterator<char>());
-  FragmentJob::ptr testJob( new MPQCCommandJob(output, nextid, command) );
-  jobs.push_back(testJob);
-  file.close();
-  LOG(1, "INFO: Added MPQCCommandJob from file "+filename+".");
 }
 
@@ -130 +133 @@
 
 
-/** Print MPQCData from received results.
- *
- * @param results received results to extract MPQCData from
- * @param KeySetFilename filename with keysets to associate forces correctly
+
+/** 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 printReceivedMPQCResults(
-    const std::vector<FragmentResult::ptr> &results,
-    const std::string &KeySetFilename,
-    size_t NoAtoms)
+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)
 {
-  EnergyMatrix Energy;
-  EnergyMatrix EnergyFragments;
-  ForceMatrix Force;
-  ForceMatrix ForceFragments;
-  KeySetsContainer KeySet;
-
-  // align fragments
-  std::map< JobId_t, size_t > MatrixNrLookup;
-  size_t FragmentCounter = 0;
-  {
-    // bring ids in order ...
-    typedef std::map< JobId_t, FragmentResult::ptr> IdResultMap_t;
-    IdResultMap_t IdResultMap;
-    for (std::vector<FragmentResult::ptr>::const_iterator iter = results.begin();
-        iter != results.end(); ++iter) {
-  #ifndef NDEBUG
-      std::pair< IdResultMap_t::iterator, bool> inserter =
-  #endif
-      IdResultMap.insert( make_pair((*iter)->getId(), *iter) );
-      ASSERT( inserter.second,
-          "printReceivedMPQCResults() - two results have same id "
-          +toString((*iter)->getId())+".");
-    }
-    // ... and fill lookup
-    for(IdResultMap_t::const_iterator iter = IdResultMap.begin();
-        iter != IdResultMap.end(); ++iter)
-      MatrixNrLookup.insert( make_pair(iter->first, FragmentCounter++) );
-  }
-  LOG(1, "INFO: There are " << FragmentCounter << " fragments.");
-
-  // extract results
-  std::vector<MPQCData> fragmentData(results.size());
-  MPQCData combinedData;
-
-  LOG(2, "DEBUG: Parsing now through " << results.size() << " results.");
-  for (std::vector<FragmentResult::ptr>::const_iterator iter = results.begin();
-      iter != results.end(); ++iter) {
-    LOG(1, "RESULT: job #"+toString((*iter)->getId())+": "+toString((*iter)->result));
-    MPQCData extractedData;
-    std::stringstream inputstream((*iter)->result);
-    LOG(2, "DEBUG: First 50 characters FragmentResult's string: "+(*iter)->result.substr(0, 50));
-    boost::archive::text_iarchive ia(inputstream);
-    ia >> extractedData;
-    LOG(1, "INFO: extracted data is " << extractedData << ".");
-
+  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.energy);
+      matrix[0].resize(1, extractedData.energies.total);
       if (!Energy.AddMatrix(
-          std::string("MPQCJob ")+toString((*iter)->getId()),
+          std::string("MPQCJob ")+toString(jobid),
           matrix,
-          MatrixNrLookup[(*iter)->getId()])) {
+          nriter->second)) {
         ELOG(1, "Adding energy matrix failed.");
         return false;
@@ -213 +187 @@
       }
       if (!Force.AddMatrix(
-          std::string("MPQCJob ")+toString((*iter)->getId()),
+          std::string("MPQCJob ")+toString(jobid),
           matrix,
-          MatrixNrLookup[(*iter)->getId()])) {
+          nriter->second)) {
         ELOG(1, "Adding force matrix failed.");
         return false;
@@ -234 +208 @@
     return false;
 
+  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
+ * @param full_sample summed up charge from fragments on return
+ */
+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);
+
+  // place results into maps
+  EnergyMatrix Energy;
+  ForceMatrix Force;
+  if (!putResultsintoMatrices(fragmentData, MatrixNrLookup, FragmentCounter, NoAtoms, Energy, Force))
+    return false;
+
+  // initialise keysets
+  KeySetsContainer KeySet;
+  KeySetsContainer ForceKeySet;
+  if (!Energy.InitialiseIndices()) return false;
+
+  if (!Force.ParseIndices(KeySetFilename.c_str())) return false;
+
+  {
+    // else needs keysets without hydrogens
+    std::stringstream filename;
+    filename << FRAGMENTPREFIX << KEYSETFILE;
+    if (!KeySet.ParseKeySets(KeySetFilename, filename.str(), FragmentCounter)) return false;
+  }
+
+  {
+    // forces need keysets including hydrogens
+    std::stringstream filename;
+    filename << FRAGMENTPREFIX << FORCESFILE;
+    if (!ForceKeySet.ParseKeySets(KeySetFilename, filename.str(), FragmentCounter)) return false;
+  }
 
   // combine all found data
-  if (!Energy.InitialiseIndices()) return false;
-
-  if (!Force.ParseIndices(KeySetFilename.c_str())) return false;
-
-  if (!KeySet.ParseKeySets(KeySetFilename.c_str(), Force.RowCounter, Force.MatrixCounter)) return false;
-
   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;
@@ -275 +298 @@
 }
 
+void writeToFile(const std::string &filename, const std::string contents)
+{
+  std::ofstream tablefile(filename.c_str());
+  tablefile << contents;
+  tablefile.close();
+}
+
+/** Print MPQCData from received results.
+ *
+ * @param results summed up results container
+ */
+void printReceivedFullResults(
+    const FragmentationResults &results)
+{
+  // print tables (without eigenvalues, they go extra)
+  {
+    typedef boost::mpl::remove<MPQCDataEnergyVector_t, MPQCDataFused::energy_eigenvalues>::type
+      MPQCDataEnergyVector_noeigenvalues_t;
+    const std::string energyresult =
+        writeTable<MPQCDataEnergyMap_t, MPQCDataEnergyVector_noeigenvalues_t >()(
+            results.Result_Energy_fused, results.getMaxLevel());
+    LOG(0, "Energy table is \n" << energyresult);
+    std::string filename;
+    filename += FRAGMENTPREFIX + std::string("_Energy.dat");
+    writeToFile(filename, energyresult);
+  }
+
+  {
+    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);
+  }
+
+  {
+    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);
+  }
+
+  {
+    const std::string eigenvalueresult;
+    LOG(0, "Eigenvalue table is \n" << eigenvalueresult);
+    std::string filename;
+    filename += FRAGMENTPREFIX + std::string("_Eigenvalues.dat");
+    writeToFile(filename, eigenvalueresult);
+  }
+
+  {
+    const std::string forceresult =
+        writeTable<MPQCDataForceMap_t, MPQCDataForceVector_t>()(
+            results.Result_Force_fused, results.getMaxLevel());
+    LOG(0, "Force table is \n" << forceresult);
+    std::string filename;
+    filename += FRAGMENTPREFIX + std::string("_Forces.dat");
+    writeToFile(filename, forceresult);
+  }
+  // we don't want to print grid to a table
+  {
+    // print times (without flops for now)
+    typedef boost::mpl::remove<
+        boost::mpl::remove<MPQCDataTimeVector_t, MPQCDataFused::times_total_flops>::type,
+        MPQCDataFused::times_gather_flops>::type
+        MPQCDataTimeVector_noflops_t;
+    const std::string timesresult =
+        writeTable<MPQCDataTimeMap_t, MPQCDataTimeVector_noflops_t >()(
+            results.Result_Time_fused, results.getMaxLevel());
+    LOG(0, "Times table is \n" << timesresult);
+    std::string filename;
+    filename += FRAGMENTPREFIX + std::string("_Times.dat");
+    writeToFile(filename, timesresult);
+  }
+}
+
+
 Action::state_ptr FragmentationFragmentationAutomationAction::performCall() {
   boost::asio::io_service io_service;
-  FragmentController controller(io_service);
 
   // TODO: Have io_service run in second thread and merge with current again eventually
 
-  // Phase One: obtain ids
-  std::vector< boost::filesystem::path > jobfiles = params.jobfiles.get();
-  controller.requestIds(params.host.get(), params.port.get(), jobfiles.size());
-  {
-    io_service.reset();
-    Info info("io_service: Phase One");
-    io_service.run();
-  }
-  // Phase Two: create and add jobs
-  {
-    std::vector<FragmentJob::ptr> jobs;
-    for (std::vector< boost::filesystem::path >::const_iterator iter = jobfiles .begin();
-        iter != jobfiles .end(); ++iter) {
-      const std::string &filename = (*iter).string();
-      if (boost::filesystem::exists(filename)) {
-        const JobId_t next_id = controller.getAvailableId();
-        LOG(1, "INFO: Creating MPQCCommandJob with filename'"
-            +filename+"', and id "+toString(next_id)+".");
-        parsejob(jobs, params.executable.get().string(), filename, next_id);
-      } else {
-        ELOG(1, "Fragment job "+filename+" does not exist.");
+  size_t Exitflag = 0;
+  std::map<JobId_t, MPQCData> fragmentData;
+  {
+    MPQCFragmentController mpqccontroller(io_service);
+    mpqccontroller.setHost(params.host.get());
+    mpqccontroller.setPort(params.port.get());
+    mpqccontroller.setLevel(params.level.get());
+    // Phase One: obtain ids
+    std::vector< boost::filesystem::path > jobfiles = params.jobfiles.get();
+    mpqccontroller.requestIds(jobfiles.size());
+
+    // Phase Two: create and add MPQCJobs
+    if (!mpqccontroller.addJobsFromFiles(params.executable.get().string(), jobfiles))
+      return Action::failure;
+
+    // Phase Three: calculate result
+    mpqccontroller.waitforResults(jobfiles.size());
+    mpqccontroller.getResults(fragmentData);
+
+    Exitflag += mpqccontroller.getExitflag();
+  }
+
+#ifdef HAVE_VMG
+  if (params.DoLongrange.get()) {
+  if ( World::getInstance().getAllAtoms().size() == 0) {
+    ELOG(1, "Please load the full molecule into the world before starting this action.");
+    return Action::failure;
+  }
+
+  // obtain combined charge density
+  FragmentationChargeDensity summedChargeDensity(
+      fragmentData,
+      params.path.get());
+  const std::vector<SamplingGrid> full_sample = summedChargeDensity.getFullSampledGrid();
+
+  LOG(1, "INFO: There are " << fragmentData.size() << " short-range and "
+      << full_sample.size() << " level-wise long-range jobs.");
+
+  // Phase Four: obtain more ids
+  std::map<JobId_t, VMGData> longrangeData;
+  {
+    VMGFragmentController vmgcontroller(io_service);
+    vmgcontroller.setHost(params.host.get());
+    vmgcontroller.setPort(params.port.get());
+    const size_t NoJobs = fragmentData.size()+full_sample.size();
+    vmgcontroller.requestIds(NoJobs);
+
+    // Phase Five: create VMGJobs
+    const size_t near_field_cells = params.near_field_cells.get();
+    const size_t interpolation_degree = params.interpolation_degree.get();
+    if (!vmgcontroller.createLongRangeJobs(
+        fragmentData,
+        full_sample,
+        summedChargeDensity.getFragment(),
+        near_field_cells,
+        interpolation_degree))
+      return Action::failure;
+
+    // Phase Six: calculate result
+    vmgcontroller.waitforResults(NoJobs);
+    vmgcontroller.getResults(longrangeData);
+    ASSERT( NoJobs == longrangeData.size(),
+        "FragmentationFragmentationAutomationAction::performCall() - number of MPQCresults+"
+        +toString(full_sample.size())+"="+toString(NoJobs)
+        +" and VMGresults "+toString(longrangeData.size())+" don't match.");
+    Exitflag += vmgcontroller.getExitflag();
+  }
+
+  // remove full solution corresponding to full_sample from map (must be highest ids), has to be treated extra
+  std::map<JobId_t, VMGData>::iterator iter = longrangeData.end();
+  for (size_t i=0;i<full_sample.size();++i)
+    --iter;
+  std::map<JobId_t, VMGData>::iterator remove_iter = iter;
+  std::vector<VMGData> fullsolutionData;
+  for (; iter != longrangeData.end(); ++iter)
+    fullsolutionData.push_back(iter->second);
+  longrangeData.erase(remove_iter, longrangeData.end());
+
+  // Final phase: sum up and print result
+  {
+    FragmentationResults results(
+        fragmentData,
+        longrangeData,
+        fullsolutionData,
+        params.path.get(),
+        getNoAtomsFromAdjacencyFile(params.path.get()),
+        full_sample);
+
+    LOG(1, "INFO: Parsing fragment files from " << params.path.get() << ".");
+    printReceivedFullResults(results);
+  }
+
+  std::map<JobId_t, std::string> debugData;
+  {
+    if (!full_sample.empty()) {
+      // create debug jobs for each level to print the summed-up potential to vtk files
+      VMGDebugGridFragmentController debugcontroller(io_service);
+      debugcontroller.setHost(params.host.get());
+      debugcontroller.setPort(params.port.get());
+      debugcontroller.requestIds(full_sample.size());
+      if (!debugcontroller.createDebugJobs(full_sample))
         return Action::failure;
-      }
+      debugcontroller.waitforResults(full_sample.size());
+      debugcontroller.getResults(debugData);
+      Exitflag += debugcontroller.getExitflag();
     }
-    controller.addJobs(jobs);
-    controller.sendJobs(params.host.get(), params.port.get());
-  }
-  {
-    io_service.reset();
-    Info info("io_service: Phase Two");
-    io_service.run();
-  }
-  // Phase Three: calculate result
-  size_t NoCalculatedResults = 0;
-  while (NoCalculatedResults != jobfiles.size()) {
-    // wait a bit
-    boost::asio::deadline_timer timer(io_service);
-    timer.expires_from_now(boost::posix_time::milliseconds(500));
-    timer.wait();
-    // then request status
-    controller.checkResults(params.host.get(), params.port.get());
-    {
-      io_service.reset();
-      Info info("io_service: Phase Three");
-      io_service.run();
-    }
-    const std::pair<size_t, size_t> JobStatus = controller.getJobStatus();
-    LOG(1, "INFO: #" << JobStatus.first << " are waiting in the queue and #" << JobStatus.second << " jobs are calculated so far.");
-    NoCalculatedResults = JobStatus.second;
-  }
-  // Phase Three: get result
-  controller.receiveResults(params.host.get(), params.port.get());
-  {
-    io_service.reset();
-    Info info("io_service: Phase Four");
-    io_service.run();
-  }
+  }
+  }
+#else
   // Final phase: print result
   {
     LOG(1, "INFO: Parsing fragment files from " << params.path.get() << ".");
-    std::vector<FragmentResult::ptr> results = controller.getReceivedResults();
     printReceivedMPQCResults(
-        results,
+        fragmentData,
         params.path.get(),
         getNoAtomsFromAdjacencyFile(params.path.get()));
   }
-  size_t Exitflag = controller.getExitflag();
+#endif
 
   return (Exitflag == 0) ? Action::success : Action::failure;