Index: src/Actions/FragmentationAction/FragmentationAutomationAction.cpp
===================================================================
--- src/Actions/FragmentationAction/FragmentationAutomationAction.cpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Actions/FragmentationAction/FragmentationAutomationAction.cpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -163,5 +163,6 @@
         near_field_cells,
         interpolation_degree,
-        params.DoValenceOnly.get() ? MPQCData::DoSampleValenceOnly : MPQCData::DontSampleValenceOnly))
+        params.DoValenceOnly.get() ? MPQCData::DoSampleValenceOnly : MPQCData::DontSampleValenceOnly,
+        params.DoPrintDebug.get()))
       // can't use the summed up fragment due to saturation hydrogens
 //      summedChargeDensity.getFragment(),
@@ -178,17 +179,19 @@
   }
 
-  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();
+  if (params.DoPrintDebug.get()) {
+    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();
+      }
     }
   }
Index: src/Actions/FragmentationAction/FragmentationAutomationAction.def
===================================================================
--- src/Actions/FragmentationAction/FragmentationAutomationAction.def	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Actions/FragmentationAction/FragmentationAutomationAction.def	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -18,9 +18,9 @@
 // ValueStorage by the token "Z" -> first column: int, Z, "Z"
 // "undefine" if no parameters are required, use (NOPARAM_DEFAULT) for each (undefined) default value
-#define paramtypes (std::string)(std::string)(boost::filesystem::path)(unsigned int)(unsigned int)(unsigned int)(bool)(boost::filesystem::path)(bool)
-#define paramtokens ("server-address")("server-port")("fragment-executable")("grid-level")("near-field-cells")("interpolation-degree")("DoLongrange")("fragment-resultfile")("DoValenceOnly")
-#define paramdescriptions ("hostname of server")("controller port of server")("executable to launch on clients")("resolution of multigrid")("number of cells used in smearing out core charge")("interpolation degree for getting the nuclei potential from the grid")("whether to calculate long-range contributions")("parse fragment results from the given file")("whether the sampling uses only the valence electron and nuclei charge") 
-#define paramdefaults (PARAM_DEFAULT("127.0.0.1"))(NOPARAM_DEFAULT)(PARAM_DEFAULT("mpqc"))(PARAM_DEFAULT(5))(PARAM_DEFAULT(3))(PARAM_DEFAULT(3))(PARAM_DEFAULT("0"))(PARAM_DEFAULT(""))(PARAM_DEFAULT("0"))
-#define paramreferences (host)(port)(executable)(level)(near_field_cells)(interpolation_degree)(DoLongrange)(resultsfile)(DoValenceOnly)
+#define paramtypes (std::string)(std::string)(boost::filesystem::path)(unsigned int)(unsigned int)(unsigned int)(bool)(boost::filesystem::path)(bool)(bool)
+#define paramtokens ("server-address")("server-port")("fragment-executable")("grid-level")("near-field-cells")("interpolation-degree")("DoLongrange")("fragment-resultfile")("DoValenceOnly")("DoPrintDebug")
+#define paramdescriptions ("hostname of server")("controller port of server")("executable to launch on clients")("resolution of multigrid")("number of cells used in smearing out core charge")("interpolation degree for getting the nuclei potential from the grid")("whether to calculate long-range contributions")("parse fragment results from the given file")("whether the sampling uses only the valence electron and nuclei charge")("whether to print grids for debug visualization") 
+#define paramdefaults (PARAM_DEFAULT("127.0.0.1"))(NOPARAM_DEFAULT)(PARAM_DEFAULT("mpqc"))(PARAM_DEFAULT(5))(PARAM_DEFAULT(3))(PARAM_DEFAULT(3))(PARAM_DEFAULT("0"))(PARAM_DEFAULT(""))(PARAM_DEFAULT("0"))(PARAM_DEFAULT("0"))
+#define paramreferences (host)(port)(executable)(level)(near_field_cells)(interpolation_degree)(DoLongrange)(resultsfile)(DoValenceOnly)(DoPrintDebug)
 #define paramvalids \
 (DummyValidator< std::string >()) \
@@ -32,4 +32,5 @@
 (DummyValidator< bool >()) \
 (DummyValidator< boost::filesystem::path >()) \
+(DummyValidator< bool >()) \
 (DummyValidator< bool >())
 
Index: src/Fragmentation/Automation/VMGFragmentController.cpp
===================================================================
--- src/Fragmentation/Automation/VMGFragmentController.cpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Fragmentation/Automation/VMGFragmentController.cpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -71,9 +71,10 @@
     const size_t near_field_cells,
     const size_t interpolation_degree,
-    const MPQCData::DoValenceOnly_t _DoValenceOnly)
+    const MPQCData::DoValenceOnly_t _DoValenceOnly,
+    const bool _DoPrintDebug)
 {
   std::vector<FragmentJob::ptr> jobs;
-  // add one job for each fragment as the short-range correction which we need
-  // to subtract from the obtained full potential to get the long-range part only
+  /// add one job for each fragment as the short-range correction which we need
+  /// to subtract from the obtained full potential to get the long-range part only
   for (std::map<JobId_t, MPQCData>::const_iterator iter = fragmentData.begin();
       iter != fragmentData.end(); ++iter) {
@@ -83,9 +84,18 @@
         << " gridpoints and " << data.charges.size() << " particle charges.");
     FragmentJob::ptr testJob(
-        new VMGJob(next_id, data.sampled_grid, data.positions, data.charges, near_field_cells, interpolation_degree) );
+        new VMGJob(
+            next_id,
+            data.sampled_grid,
+            data.positions,
+            data.charges,
+            near_field_cells,
+            interpolation_degree,
+            _DoPrintDebug) );
     jobs.push_back(testJob);
   }
 
-  // prepare positions and charges of full system
+  /// prepare positions and charges of full system
+  /// \note we cannot use the summed up Fragment here, as the saturation hydrogens
+  /// are in the way and cannot be sorted out properly/in a simple fashion.
   std::vector< std::vector<double> > positions;
   std::vector<double> charges;
@@ -108,4 +118,5 @@
     }
   }
+  /// and submit full job
   for(std::vector<SamplingGrid>::const_iterator iter = full_sampled_grid.begin();
       iter != full_sampled_grid.end();
@@ -115,9 +126,16 @@
         << " gridpoints and " << charges.size() << " particle charges.");
     FragmentJob::ptr testJob(
-        new VMGJob(next_id, *iter, positions, charges, near_field_cells, interpolation_degree) );
+        new VMGJob(
+            next_id,
+            *iter,
+            positions,
+            charges,
+            near_field_cells,
+            interpolation_degree,
+            _DoPrintDebug) );
     jobs.push_back(testJob);
   }
 
-  // then send jobs to controller
+  /// then send jobs to controller
   addJobs(jobs);
   sendJobs(host, port);
Index: src/Fragmentation/Automation/VMGFragmentController.hpp
===================================================================
--- src/Fragmentation/Automation/VMGFragmentController.hpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Fragmentation/Automation/VMGFragmentController.hpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -46,4 +46,5 @@
    * \param interpolation_degree interpolation degree for newton interpolation from grid
    * \param _DoValenceOnly whether to neglect core electron and nuclei charges
+   * \param _DoPrintDebug whether we do print grid for debug visualization or not
    */
   bool createLongRangeJobs(
@@ -52,5 +53,6 @@
       const size_t near_field_cells,
       const size_t interpolation_degree,
-      const MPQCData::DoValenceOnly_t _DoValenceOnly);
+      const MPQCData::DoValenceOnly_t _DoValenceOnly,
+      const bool _DoPrintDebug);
 
   /** Get results of finished long-range jobs.
Index: src/Jobs/InterfaceVMGJob.cpp
===================================================================
--- src/Jobs/InterfaceVMGJob.cpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Jobs/InterfaceVMGJob.cpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -71,4 +71,5 @@
     vmg_float _box_end,
     const int& near_field_cells,
+    const bool _DoPrintDebug,
     int coarseningSteps,
     double alpha) :
@@ -78,5 +79,6 @@
   sampled_input(_sampled_input),
   returndata(_returndata),
-  level(levelMax)
+  level(levelMax),
+  DoPrintDebug(_DoPrintDebug)
 {
   for (size_t i=0;i<3;++i) {
@@ -144,6 +146,8 @@
   comm.CommFromGhosts(particle_grid);
 
-  // print nuclei grid to vtk
-  comm.PrintGrid(particle_grid, "Sampled Nuclei Density");
+  if (DoPrintDebug) {
+    // print nuclei grid to vtk
+    comm.PrintGrid(particle_grid, "Sampled Nuclei Density");
+  }
 
   // add sampled electron charge density onto grid
@@ -153,6 +157,8 @@
       -1.);
 
-  // print electron grid to vtk
-  comm.PrintGrid(grid, "Sampled Electron Density");
+  if (DoPrintDebug) {
+    // print electron grid to vtk
+    comm.PrintGrid(grid, "Sampled Electron Density");
+  }
 
   // add particle_grid onto grid
@@ -180,6 +186,8 @@
   comm.PrintStringOnce("Grid charge integral: %e", charge_sum/(4.0 * VMG::Math::pi));
 
-  // print total grid to vtk
-  comm.PrintGrid(grid, "Total Charge Density");
+  if (DoPrintDebug) {
+    // print total grid to vtk
+    comm.PrintGrid(grid, "Total Charge Density");
+  }
 
 //  delete temp_grid;
@@ -197,6 +205,8 @@
   Index i;
 
-  // print output grid to vtk
-  comm.PrintGrid(grid, "Potential Solution");
+  if (DoPrintDebug) {
+    // print output grid to vtk
+    comm.PrintGrid(grid, "Potential Solution");
+  }
 
   // obtain sampled potential from grid
Index: src/Jobs/InterfaceVMGJob.hpp
===================================================================
--- src/Jobs/InterfaceVMGJob.hpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Jobs/InterfaceVMGJob.hpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -46,6 +46,8 @@
       vmg_float _box_end,
       const int& near_field_cells,
+      const bool _DoPrintDebug=false,
       int coarseningSteps=0,
-      double alpha=1.6);
+      double alpha=1.6
+      );
 
   ~InterfaceVMGJob() {}
@@ -71,4 +73,7 @@
   //!> contains end of grid intervals
   double box_end[3];
+
+  //!> whether we do print grid for debug visualization or not
+  const bool DoPrintDebug;
 };
 
Index: src/Jobs/VMGJob.cpp
===================================================================
--- src/Jobs/VMGJob.cpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Jobs/VMGJob.cpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -82,5 +82,6 @@
     const std::vector< double > &_particle_charges,
     const size_t _near_field_cells,
-    const size_t _interpolation_degree) :
+    const size_t _interpolation_degree,
+    const bool _DoPrintDebug) :
   FragmentJob(_JobId),
   density_grid(_density_grid),
@@ -90,5 +91,6 @@
   interpolation_degree(_interpolation_degree),
   returndata(static_cast<const SamplingGridProperties &>(_density_grid)),
-  particles()
+  particles(),
+  DoPrintDebug(_DoPrintDebug)
 {}
 
@@ -97,5 +99,6 @@
   near_field_cells(3),
   interpolation_degree(3),
-  particles()
+  particles(),
+  DoPrintDebug(false)
 {}
 
@@ -216,5 +219,6 @@
       Vector(density_grid.begin),
       density_grid.end[0]-density_grid.begin[0],
-      near_field_cells);
+      near_field_cells,
+      DoPrintDebug);
   new LevelOperatorCS(Stencils::RestrictionFullWeight, Stencils::InterpolationTrilinear);
   new Givens<SolverSingular>();
Index: src/Jobs/VMGJob.hpp
===================================================================
--- src/Jobs/VMGJob.hpp	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ src/Jobs/VMGJob.hpp	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -42,4 +42,5 @@
    * @param _interpolation_degree degree of interpolation polynomial for getting nuclei
    *        potential from grid
+   * @param _DoPrintDebug whether we do print grid for debug visualization or not
    */
   VMGJob(const JobId_t _JobId,
@@ -48,5 +49,6 @@
       const std::vector< double > &_particle_charges,
       const size_t _near_field_cells,
-      const size_t _interpolation_degree
+      const size_t _interpolation_degree,
+      const bool _DoPrintDebug=false
       );
   virtual ~VMGJob();
@@ -70,4 +72,6 @@
   //!> interpolation degree used in sampling the potential of the nuclei
   const size_t interpolation_degree;
+  //!> whether we do print grid for debug visualization or not
+  const bool DoPrintDebug;
 
 private:
@@ -120,4 +124,5 @@
     ar & const_cast< size_t &>(near_field_cells);
     ar & const_cast< size_t &>(interpolation_degree);
+    ar & const_cast< bool &>(DoPrintDebug);
     ar & returndata;
   }
Index: tests/Python/AllActions/options.dat
===================================================================
--- tests/Python/AllActions/options.dat	(revision 6ff62c8bb90891ea66dfc4882a5f7298de44ba0c)
+++ tests/Python/AllActions/options.dat	(revision b6b21a5af3ef8085052c11ae4c00c67da135353a)
@@ -60,4 +60,5 @@
 domain-position	"10. 10. 10."
 DoLongrange	"0"
+DoPrintDebug	"0"
 DoRotate	"0"
 DoSaturate	"0"
