Index: src/Actions/FragmentationAction/FragmentationAutomationAction.cpp
===================================================================
--- src/Actions/FragmentationAction/FragmentationAutomationAction.cpp	(revision 94d530feefb1ab3f22599705d28b424606dab346)
+++ src/Actions/FragmentationAction/FragmentationAutomationAction.cpp	(revision 849cd8e107c07fa8b744cd0f691da136460b27f0)
@@ -65,4 +65,5 @@
 #include "Jobs/MPQCData.hpp"
 #include "Jobs/MPQCData_printKeyNames.hpp"
+#include "Jobs/Grid/SamplingGrid.hpp"
 #include "LinearAlgebra/RealSpaceMatrix.hpp"
 #ifdef HAVE_VMG
@@ -428,4 +429,57 @@
  * @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 sumUpChargeDensity(
+    const std::vector<FragmentResult::ptr> &results,
+    const std::vector<MPQCData> &fragmentData,
+    const std::string &KeySetFilename,
+    SamplingGrid &full_sample)
+{
+  // create lookup from job nr to fragment number
+  std::map< JobId_t, size_t > MatrixNrLookup;
+  size_t FragmentCounter = 0;
+  createMatrixNrLookup(results, MatrixNrLookup, FragmentCounter);
+
+  // initialise keysets
+  KeySetsContainer KeySet;
+  {
+    // else needs keysets without hydrogens
+    std::stringstream filename;
+    filename << FRAGMENTPREFIX << KEYSETFILE;
+    if (!KeySet.ParseKeySets(KeySetFilename, filename.str(), FragmentCounter)) return false;
+  }
+
+  /// prepare for OrthogonalSummation
+
+  // convert KeySetContainer to IndexSetContainer
+  IndexSetContainer::ptr container(new IndexSetContainer(KeySet));
+  // create the map of all keysets
+  SubsetMap::ptr subsetmap(new SubsetMap(*container));
+
+  // create a vector of all job ids
+  std::vector<JobId_t> jobids(results.size(), JobId::IllegalJob);
+  std::transform(results.begin(), results.end(), jobids.begin(),
+      boost::bind(&FragmentResult::getId,
+          boost::bind(&FragmentResult::ptr::operator->, _1)));
+
+  /// convert all MPQCData to MPQCDataMap_t
+  std::vector<MPQCDataGridMap_t> Result_Grid_fused(
+    SumUpPerLevel<MPQCDataGridMap_t, MPQCDataGridVector_t>(
+        fragmentData, jobids, MatrixNrLookup, container, subsetmap));
+  // obtain full grid
+  full_sample = boost::fusion::at_key<MPQCDataFused::sampled_grid>(Result_Grid_fused.back());
+
+  return true;
+}
+
+
+/** Print MPQCData from received results.
+ *
+ * @param results results with ids to associate with fragment number
+ * @param fragmentData MPQCData resulting from the jobs
+ * @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(
@@ -433,5 +487,6 @@
     const std::vector<MPQCData> &fragmentData,
     const std::string &KeySetFilename,
-    size_t NoAtoms)
+    size_t NoAtoms,
+    SamplingGrid &full_sample)
 {
   // create lookup from job nr to fragment number
@@ -507,4 +562,7 @@
                 Result_Force_fused);
     boost::mpl::for_each<MPQCDataForceVector_t>(boost::ref(forceSummer));
+
+    // obtain full grid
+    full_sample = boost::fusion::at_key<MPQCDataFused::sampled_grid>(Result_Grid_fused.back());
 
     // print tables (without eigenvalues, they go extra)
@@ -614,5 +672,6 @@
     FragmentController &controller,
     const FragmentationFragmentationAutomationAction::FragmentationFragmentationAutomationParameters &params,
-    const std::vector<MPQCData> &fragmentData)
+    const std::vector<MPQCData> &fragmentData,
+    const SamplingGrid &full_sampled_grid)
 {
   std::vector<FragmentJob::ptr> jobs;
@@ -629,20 +688,16 @@
   }
 
-  // add one more job for the full calculation
-  // TODO: Here, we actually have to combine all the other sampled_grids
   {
-    const int level = LEVEL;
-    const int GRID = pow(2, level);
-    std::vector<double> full_sample(GRID*GRID*GRID, 0.);
-    double begin[NDIM] = { 0., 0., 0. };
-    const RealSpaceMatrix& M = World::getInstance().getDomain().getM();
-    const double size = M.at(0,0);
-    ASSERT( M.determinant() == size*size*size,
-        "createLongRangeJobs() - current domain matrix "+toString(M)+" is not cubic.");
-    const SamplingGrid full_sampled_grid(begin, size, level, full_sample);
+//    const int level = LEVEL;
+//    const int GRID = pow(2, level);
+//    double begin[NDIM] = { 0., 0., 0. };
+//    const RealSpaceMatrix& M = World::getInstance().getDomain().getM();
+//    const double size = M.at(0,0);
+//    ASSERT( M.determinant() == size*size*size,
+//        "createLongRangeJobs() - current domain matrix "+toString(M)+" is not cubic.");
     const std::vector< std::vector<double> > positions;
     const std::vector<double> charges;
     const JobId_t next_id = controller.getAvailableId();
-    LOG(1, "INFO: Creating full VMGJob with " << full_sample.size()
+    LOG(1, "INFO: Creating full VMGJob with " << full_sampled_grid.sampled_grid.size()
         << " gridpoints and " << charges.size() << " particle charges.");
     FragmentJob::ptr testJob(
@@ -706,16 +761,15 @@
   ConvertFragmentResultToMPQCData(MPQCresults, fragmentData);
 
-  // print intermediate short-range results
-  {
-    LOG(1, "INFO: Parsing fragment files from " << params.path.get() << ".");
-    printReceivedMPQCResults(
-        MPQCresults,
-        fragmentData,
-        params.path.get(),
-        getNoAtomsFromAdjacencyFile(params.path.get()));
-  }
-
 #ifdef HAVE_VMG
   if (params.DoLongrange.get()) {
+  // obtain combined charge density
+  LOG(1, "INFO: Parsing fragment files from " << params.path.get() << ".");
+  SamplingGrid full_sample;
+  sumUpChargeDensity(
+      MPQCresults,
+      fragmentData,
+      params.path.get(),
+      full_sample);
+
   // Phase Four: obtain more ids
   requestIds(controller, params, fragmentData.size()+1);
@@ -723,5 +777,5 @@
 
   // Phase Five: create VMGJobs
-  if (!createLongRangeJobs(controller, params, fragmentData))
+  if (!createLongRangeJobs(controller, params, fragmentData, full_sample))
     return Action::failure;
   RunService(io_service, "Adding VMGJobs");
@@ -742,5 +796,6 @@
         fragmentData,
         params.path.get(),
-        getNoAtomsFromAdjacencyFile(params.path.get()));
+        getNoAtomsFromAdjacencyFile(params.path.get()),
+        full_sample);
   }
   }
