diff --git a/src/cluster/cluster.cpp b/src/cluster/cluster.cpp
index 46de0e4970fd32402cd1b9b96d1b9ae902131263..be7e59573bdb80b86427597697692325289763b9 100644
--- a/src/cluster/cluster.cpp
+++ b/src/cluster/cluster.cpp
@@ -440,8 +440,8 @@ void cluster(const string& config_file, const string& data_file, const string& o
 	    int jer = cluster_jxi488_cycle(myjxi488, sconf, gconf, p_scattering_angles, cid_2, p_output_2, output_path, vtppoanp_2);
 	  } else {
 	    if (myompthread > 0) {
-	      // If there is no input for this thread, set output pointer to NULL.
-	      p_outarray[myompthread] = NULL;
+	      // If there is no input for this thread, mark to skip.
+	      p_outarray[myompthread] = new ClusterOutputInfo(1);
 	    }
 	  }
 #pragma omp barrier
@@ -453,11 +453,11 @@ void cluster(const string& config_file, const string& data_file, const string& o
 	  // threads different from 0 append their virtual files to the one of thread 0, and delete them
 	  if (myompthread == 0) {
 	    for (int ti=1; ti<ompnumthreads; ti++) {
-	      if (p_outarray[ti] != NULL) {
+	      if (p_outarray[ti]->skip_flag != 1) {
 		p_outarray[0]->insert(*(p_outarray[ti]));
-		delete p_outarray[ti];
-		p_outarray[ti] = NULL;
 	      }
+	      delete p_outarray[ti];
+	      p_outarray[ti] = NULL;
 	      vtppoanarray[0]->append(*(vtppoanarray[ti]));
 	      delete vtppoanarray[ti];
 	    }
@@ -628,39 +628,27 @@ void cluster(const string& config_file, const string& data_file, const string& o
 	  }
 	  int jer = cluster_jxi488_cycle(myjxi488, sconf, gconf, p_scattering_angles, cid_2, p_output_2, output_path, vtppoanp_2);
 	} else {
-	  // if (myompthread > 0) {
-	    // If there is no input for this thread, set the output pointer to NULL.
-	    p_outarray[myompthread] = NULL;
-	  //}
+	  p_outarray[myompthread] = new ClusterOutputInfo(1);
 	}
 
 #pragma omp barrier
 	// threads different from 0 append their virtual files to the one of thread 0, and delete them
 	if (myompthread == 0) {
 	  for (int ti=1; ti<ompnumthreads; ti++) {
-	    if (p_outarray[ti] != NULL) {
+	    if (p_outarray[ti]->skip_flag != 1) {
 	      p_outarray[0]->insert(*(p_outarray[ti]));
-	      delete p_outarray[ti];
-	      p_outarray[ti] = NULL;
 	    }
+	    delete p_outarray[ti];
+	    p_outarray[ti] = NULL;
 	    vtppoanarray[0]->append(*(vtppoanarray[ti]));
 	    delete vtppoanarray[ti];
 	  }
 	  // thread 0 sends the collected virtualfiles to thread 0 of MPI process 0, then deletes them
 	  for (int rr=1; rr<mpidata->nprocs; rr++) {
 	    if (rr == mpidata->rank) {
-	      if (p_outarray[0] == NULL) {
-		// signal that we are not sending anything
-		int skip_flag = 1;
-		MPI_Send(&skip_flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-	      }
-	      else {
-		// signal that we are sending something
-		int skip_flag = 0;
-		MPI_Send(&skip_flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-		p_outarray[0]->mpisend(mpidata);
-		delete p_outarray[0];
-	      }
+	      p_outarray[0]->mpisend(mpidata);
+	      delete p_outarray[0];
+	      p_outarray[0] = NULL;
 	      vtppoanarray[0]->mpisend(mpidata);
 	      delete vtppoanarray[0];
 	    }
diff --git a/src/include/errors.h b/src/include/errors.h
index a986defa166fb9c7edbcdf5911f45084461e6590..1908ade318ecb90db059bde369971a27bcf85d16 100644
--- a/src/include/errors.h
+++ b/src/include/errors.h
@@ -22,6 +22,32 @@
 #ifndef INCLUDE_ERRORS_H_
 #define INCLUDE_ERRORS_H_
 
+/*! \brief Exception for wrong OutputInfo NULL calls.
+ */
+class UnrecognizedOutputInfo: public std::exception {
+protected:
+  //! Description of the problem.
+  std::string message;
+
+public:
+  /**
+   * \brief Exception instance constructor.
+   *
+   * \param requested: `int` The index that was requested.
+   */
+  UnrecognizedOutputInfo(int requested) {
+    message = "Error: passed parameter " + std::to_string(requested)
+      + ", but only 1 is allowed";
+  }
+  
+  /**
+   * \brief Exception message.
+   */
+  virtual const char* what() const throw() {
+    return message.c_str();
+  }
+};
+
 /*! \brief Exception for out of bounds List requests.
  */
 class ListOutOfBoundsException: public std::exception {
diff --git a/src/include/outputs.h b/src/include/outputs.h
index b4f9b02ed34cbb91a56cfbe9de2832e6fd28be52..82ee40ae0c2d19dbaa18135e680c1f001a94ecd6 100644
--- a/src/include/outputs.h
+++ b/src/include/outputs.h
@@ -33,6 +33,8 @@
  */
 class ClusterOutputInfo {
 protected:
+  //! \brief Flag for skipping mpisend() and mpireceive()
+  int _skip_flag;
   //! \brief Number of incident azimuth calculations.
   int _num_theta;
   //! \brief Number of scattered azimuth calculations.
@@ -63,6 +65,8 @@ protected:
   int write_legacy(const std::string &output);
   
 public:
+  //! \brief Read-only view on skip_flag
+  const int &skip_flag = _skip_flag;
   //! \brief Read-only view on the ID of the first scale
   const int &first_xi = _first_xi;
   //! \brief Number of spheres in the aggregate.
@@ -471,6 +475,12 @@ public:
    */   
   ClusterOutputInfo(const std::string &hdf5_name);
 
+  /*! \brief `ClusterOutputInfo` constructor for the dummy NULL case
+   *
+   * \param skip_flag: `const int` must be passed as the `1` constant.
+   */   
+  ClusterOutputInfo(const int skip_flag);
+
   /*! \brief `ClusterOutputInfo` instance destroyer.
    */
   ~ClusterOutputInfo();
@@ -551,6 +561,8 @@ public:
  */
 class InclusionOutputInfo {
 protected:
+  //! \brief Flag for skipping mpisend() and mpireceive()
+  int _skip_flag;
   //! \brief Number of incident azimuth calculations.
   int _num_theta;
   //! \brief Number of scattered azimuth calculations.
@@ -581,6 +593,8 @@ protected:
   int write_legacy(const std::string &output);
   
 public:
+  //! \brief Read-only view on skip_flag
+  const int &skip_flag = _skip_flag;
   //! \brief Read-only view on the ID of the first scale
   const int &first_xi = _first_xi;
   //! \brief Number of spheres in the aggregate.
@@ -901,6 +915,12 @@ public:
    */   
   InclusionOutputInfo(const std::string &hdf5_name);
 
+  /*! \brief `InclusionOutputInfo` constructor for the dummy NULL case
+   *
+   * \param skip_flag: `const int` must be passed as the `1` constant.
+   */   
+  InclusionOutputInfo(const int skip_flag);
+
   /*! \brief `InclusionOutputInfo` instance destroyer.
    */
   ~InclusionOutputInfo();
@@ -981,6 +1001,8 @@ public:
  */
 class SphereOutputInfo {
 protected:
+  //! \brief Flag for skipping mpisend() and mpireceive()
+  int _skip_flag;
   //! \brief Number of incident azimuth calculations.
   int _num_theta;
   //! \brief Number of scattered azimuth calculations.
@@ -1011,6 +1033,8 @@ protected:
   int write_legacy(const std::string &output);
   
 public:
+  //! \brief Read-only view on skip_flag
+  const int &skip_flag = _skip_flag;
   //! \brief Read-only view on the ID of the first scale
   const int &first_xi = _first_xi;
   //! \brief Number of spheres.
@@ -1181,6 +1205,12 @@ public:
    */   
   SphereOutputInfo(const std::string &hdf5_name);
 
+  /*! \brief `SphereOutputInfo` constructor for the dummy NULL case
+   *
+   * \param skip_flag: `const int` must be passed as the `1` constant.
+   */   
+  SphereOutputInfo(const int skip_flag);
+
   /*! \brief `InclusionOutputInfo` instance destroyer.
    */
   ~SphereOutputInfo();
diff --git a/src/inclusion/inclusion.cpp b/src/inclusion/inclusion.cpp
index 846cb953b780368598eb7c5a70ed9cfff337711a..dac0821198526b2e3adfa2ef70cfce57a03ca50f 100644
--- a/src/inclusion/inclusion.cpp
+++ b/src/inclusion/inclusion.cpp
@@ -437,7 +437,7 @@ void inclusion(const string& config_file, const string& data_file, const string&
 	  } else {
 	    if (myompthread > 0) {
 	      // If there is no input for this thread, set output pointer to NULL.
-	      p_outarray[myompthread] = NULL;
+	      p_outarray[myompthread] = new InclusionOutputInfo(1);
 	    }
 	  }
 #pragma omp barrier
@@ -449,11 +449,11 @@ void inclusion(const string& config_file, const string& data_file, const string&
 	  // threads different from 0 append their virtual files to the one of thread 0, and delete them
 	  if (myompthread == 0) {
 	    for (int ti=1; ti<ompnumthreads; ti++) {
-	      if (p_outarray[ti] != NULL) {
+	      if (p_outarray[ti]->skip_flag != 1) {
 		p_outarray[0]->insert(*(p_outarray[ti]));
-		delete p_outarray[ti];
-		p_outarray[ti] = NULL;
 	      }
+	      delete p_outarray[ti];
+	      p_outarray[ti] = NULL;
 	      vtppoanarray[0]->append(*(vtppoanarray[ti]));
 	      delete vtppoanarray[ti];
 	    }
@@ -624,39 +624,27 @@ void inclusion(const string& config_file, const string& data_file, const string&
 	  }
 	  int jer = inclusion_jxi488_cycle(myjxi488, sconf, gconf, p_scattering_angles, cid_2, p_output_2, output_path, vtppoanp_2);
 	} else {
-	  //if (myompthread > 0) {
-	    // If there is no input for this thread, set the output pointer to NULL.
-	    p_outarray[myompthread] = NULL;
-	    //}	  
+	  p_outarray[myompthread] = new InclusionOutputInfo(1);
 	}
 
 #pragma omp barrier
 	// threads different from 0 append their virtual files to the one of thread 0, and delete them
 	if (myompthread == 0) {
 	  for (int ti=1; ti<ompnumthreads; ti++) {
-	    if (p_outarray[ti] != NULL) {
+	    if (p_outarray[ti]->skip_flag != 1) {
 	      p_outarray[0]->insert(*(p_outarray[ti]));
-	      delete p_outarray[ti];
-	      p_outarray[ti] = NULL;
 	    }
+	    delete p_outarray[ti];
+	    p_outarray[ti] = NULL;
 	    vtppoanarray[0]->append(*(vtppoanarray[ti]));
 	    delete vtppoanarray[ti];
 	  }
 	  // thread 0 sends the collected virtualfiles to thread 0 of MPI process 0, then deletes them
 	  for (int rr=1; rr<mpidata->nprocs; rr++) {
 	    if (rr == mpidata->rank) {
-	      if (p_outarray[0] == NULL) {
-		// signal that we are not sending anything
-		int skip_flag = 1;
-		MPI_Send(&skip_flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-	      }
-	      else {
-		// signal that we are sending something
-		int skip_flag = 0;
-		MPI_Send(&skip_flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-		p_outarray[0]->mpisend(mpidata);
-		delete p_outarray[0];
-	      }
+	      p_outarray[0]->mpisend(mpidata);
+	      delete p_outarray[0];
+	      p_outarray[0] = NULL;
 	      vtppoanarray[0]->mpisend(mpidata);
 	      delete vtppoanarray[0];
 	    }
diff --git a/src/libnptm/outputs.cpp b/src/libnptm/outputs.cpp
index 3f29eb0638bb096bb4a5c24e531e738dba5227ff..1827626877a92121d446cea437ede638ad2ccf41 100644
--- a/src/libnptm/outputs.cpp
+++ b/src/libnptm/outputs.cpp
@@ -61,6 +61,7 @@ ClusterOutputInfo::ClusterOutputInfo(
   ScattererConfiguration *sc, GeometryConfiguration *gc,
   const mixMPI *mpidata, int first_xi, int xi_length
 ) {
+  _skip_flag = 0;
   nsph = gc->number_of_spheres;
   li = gc->li;
   le = gc->le;
@@ -293,6 +294,7 @@ ClusterOutputInfo::ClusterOutputInfo(const std::string &hdf5_name) {
   HDFFile *hdf_file = new HDFFile(hdf5_name, flags);
   herr_t status = hdf_file->get_status();
   string str_name, str_type;
+  _skip_flag = 0;
   if (status == 0) {
     status = hdf_file->read("NSPH", "INT32_(1)", &nsph);
     status = hdf_file->read("LI", "INT32_(1)", &li);
@@ -701,170 +703,184 @@ ClusterOutputInfo::ClusterOutputInfo(const std::string &hdf5_name) {
   }
 }
 
+ClusterOutputInfo::ClusterOutputInfo(const int flag) {
+  /*
+    create a dummy placeholder just to know I should skip MPI_Send and MPI_Recv
+  */
+  if (flag == 1) {
+    _skip_flag = 1;
+  } else {
+    UnrecognizedOutputInfo ex(flag);
+    throw ex;
+  }
+}
+
 ClusterOutputInfo::~ClusterOutputInfo() {
-  delete[] vec_x_coords;
-  delete[] vec_y_coords;
-  delete[] vec_z_coords;
-  delete[] vec_jxi;
-  delete[] vec_ier;
-  delete[] vec_vk;
-  delete[] vec_xi;
-  delete[] vec_sphere_sizes;
-  delete[] vec_sphere_ref_indices;
-  delete[] vec_sphere_scs;
-  delete[] vec_sphere_abs;
-  delete[] vec_sphere_exs;
-  delete[] vec_sphere_albs;
-  delete[] vec_sphere_sqscs;
-  delete[] vec_sphere_sqabs;
-  delete[] vec_sphere_sqexs;
-  delete[] vec_fsas;
-  delete[] vec_qschus;
-  delete[] vec_pschus;
-  delete[] vec_s0mags;
-  delete[] vec_cosavs;
-  delete[] vec_raprs;
-  delete[] vec_tqek1;
-  delete[] vec_tqsk1;
-  delete[] vec_tqek2;
-  delete[] vec_tqsk2;
-  delete[] vec_fsat;
-  delete[] vec_qschut;
-  delete[] vec_pschut;
-  delete[] vec_s0magt;
-  delete[] vec_scc1;
-  delete[] vec_scc2;
-  delete[] vec_abc1;
-  delete[] vec_abc2;
-  delete[] vec_exc1;
-  delete[] vec_exc2;
-  delete[] vec_albedc1;
-  delete[] vec_albedc2;
-  delete[] vec_sccrt1;
-  delete[] vec_sccrt2;
-  delete[] vec_abcrt1;
-  delete[] vec_abcrt2;
-  delete[] vec_excrt1;
-  delete[] vec_excrt2;
-  delete[] vec_fsac11;
-  delete[] vec_fsac21;
-  delete[] vec_fsac22;
-  delete[] vec_fsac12;
-  delete[] vec_qschuc1;
-  delete[] vec_qschuc2;
-  delete[] vec_pschuc1;
-  delete[] vec_pschuc2;
-  delete[] vec_s0magc1;
-  delete[] vec_s0magc2;
-  delete[] vec_cosavc1;
-  delete[] vec_cosavc2;
-  delete[] vec_raprc1;
-  delete[] vec_raprc2;
-  delete[] vec_fkc1;
-  delete[] vec_fkc2;
-  delete[] vec_dir_tidg;
-  delete[] vec_dir_pidg;
-  delete[] vec_dir_tsdg;
-  delete[] vec_dir_psdg;
-  delete[] vec_dir_scand;
-  delete[] vec_dir_cfmp;
-  delete[] vec_dir_sfmp;
-  delete[] vec_dir_cfsp;
-  delete[] vec_dir_sfsp;
-  delete[] vec_dir_un;
-  delete[] vec_dir_uns;
-  delete[] vec_dir_sas11;
-  delete[] vec_dir_sas21;
-  delete[] vec_dir_sas12;
-  delete[] vec_dir_sas22;
-  delete[] vec_dir_muls;
-  delete[] vec_dir_mulslr;
-  delete[] vec_dir_sat11;
-  delete[] vec_dir_sat21;
-  delete[] vec_dir_sat12;
-  delete[] vec_dir_sat22;
-  delete[] vec_dir_scc1;
-  delete[] vec_dir_scc2;
-  delete[] vec_dir_abc1;
-  delete[] vec_dir_abc2;
-  delete[] vec_dir_exc1;
-  delete[] vec_dir_exc2;
-  delete[] vec_dir_albedc1;
-  delete[] vec_dir_albedc2;
-  delete[] vec_dir_qscc1;
-  delete[] vec_dir_qscc2;
-  delete[] vec_dir_qabc1;
-  delete[] vec_dir_qabc2;
-  delete[] vec_dir_qexc1;
-  delete[] vec_dir_qexc2;
-  delete[] vec_qscamc1;
-  delete[] vec_qscamc2;
-  delete[] vec_qabsmc1;
-  delete[] vec_qabsmc2;
-  delete[] vec_qextmc1;
-  delete[] vec_qextmc2;
-  delete[] vec_dir_sccrt1;
-  delete[] vec_dir_sccrt2;
-  delete[] vec_dir_abcrt1;
-  delete[] vec_dir_abcrt2;
-  delete[] vec_dir_excrt1;
-  delete[] vec_dir_excrt2;
-  delete[] vec_dir_fsac11;
-  delete[] vec_dir_fsac21;
-  delete[] vec_dir_fsac12;
-  delete[] vec_dir_fsac22;
-  delete[] vec_dir_sac11;
-  delete[] vec_dir_sac21;
-  delete[] vec_dir_sac12;
-  delete[] vec_dir_sac22;
-  delete[] vec_dir_qschuc1;
-  delete[] vec_dir_qschuc2;
-  delete[] vec_dir_pschuc1;
-  delete[] vec_dir_pschuc2;
-  delete[] vec_dir_s0magc1;
-  delete[] vec_dir_s0magc2;
-  delete[] vec_dir_cosavc1;
-  delete[] vec_dir_cosavc2;
-  delete[] vec_dir_raprc1;
-  delete[] vec_dir_raprc2;
-  delete[] vec_dir_flc1;
-  delete[] vec_dir_flc2;
-  delete[] vec_dir_frc1;
-  delete[] vec_dir_frc2;
-  delete[] vec_dir_fkc1;
-  delete[] vec_dir_fkc2;
-  delete[] vec_dir_fxc1;
-  delete[] vec_dir_fxc2;
-  delete[] vec_dir_fyc1;
-  delete[] vec_dir_fyc2;
-  delete[] vec_dir_fzc1;
-  delete[] vec_dir_fzc2;
-  delete[] vec_dir_tqelc1;
-  delete[] vec_dir_tqelc2;
-  delete[] vec_dir_tqerc1;
-  delete[] vec_dir_tqerc2;
-  delete[] vec_dir_tqekc1;
-  delete[] vec_dir_tqekc2;
-  delete[] vec_dir_tqexc1;
-  delete[] vec_dir_tqexc2;
-  delete[] vec_dir_tqeyc1;
-  delete[] vec_dir_tqeyc2;
-  delete[] vec_dir_tqezc1;
-  delete[] vec_dir_tqezc2;
-  delete[] vec_dir_tqslc1;
-  delete[] vec_dir_tqslc2;
-  delete[] vec_dir_tqsrc1;
-  delete[] vec_dir_tqsrc2;
-  delete[] vec_dir_tqskc1;
-  delete[] vec_dir_tqskc2;
-  delete[] vec_dir_tqsxc1;
-  delete[] vec_dir_tqsxc2;
-  delete[] vec_dir_tqsyc1;
-  delete[] vec_dir_tqsyc2;
-  delete[] vec_dir_tqszc1;
-  delete[] vec_dir_tqszc2;
-  delete[] vec_dir_mulc;
-  delete[] vec_dir_mulclr;
+  if (_skip_flag != 1) {
+    delete[] vec_x_coords;
+    delete[] vec_y_coords;
+    delete[] vec_z_coords;
+    delete[] vec_jxi;
+    delete[] vec_ier;
+    delete[] vec_vk;
+    delete[] vec_xi;
+    delete[] vec_sphere_sizes;
+    delete[] vec_sphere_ref_indices;
+    delete[] vec_sphere_scs;
+    delete[] vec_sphere_abs;
+    delete[] vec_sphere_exs;
+    delete[] vec_sphere_albs;
+    delete[] vec_sphere_sqscs;
+    delete[] vec_sphere_sqabs;
+    delete[] vec_sphere_sqexs;
+    delete[] vec_fsas;
+    delete[] vec_qschus;
+    delete[] vec_pschus;
+    delete[] vec_s0mags;
+    delete[] vec_cosavs;
+    delete[] vec_raprs;
+    delete[] vec_tqek1;
+    delete[] vec_tqsk1;
+    delete[] vec_tqek2;
+    delete[] vec_tqsk2;
+    delete[] vec_fsat;
+    delete[] vec_qschut;
+    delete[] vec_pschut;
+    delete[] vec_s0magt;
+    delete[] vec_scc1;
+    delete[] vec_scc2;
+    delete[] vec_abc1;
+    delete[] vec_abc2;
+    delete[] vec_exc1;
+    delete[] vec_exc2;
+    delete[] vec_albedc1;
+    delete[] vec_albedc2;
+    delete[] vec_sccrt1;
+    delete[] vec_sccrt2;
+    delete[] vec_abcrt1;
+    delete[] vec_abcrt2;
+    delete[] vec_excrt1;
+    delete[] vec_excrt2;
+    delete[] vec_fsac11;
+    delete[] vec_fsac21;
+    delete[] vec_fsac22;
+    delete[] vec_fsac12;
+    delete[] vec_qschuc1;
+    delete[] vec_qschuc2;
+    delete[] vec_pschuc1;
+    delete[] vec_pschuc2;
+    delete[] vec_s0magc1;
+    delete[] vec_s0magc2;
+    delete[] vec_cosavc1;
+    delete[] vec_cosavc2;
+    delete[] vec_raprc1;
+    delete[] vec_raprc2;
+    delete[] vec_fkc1;
+    delete[] vec_fkc2;
+    delete[] vec_dir_tidg;
+    delete[] vec_dir_pidg;
+    delete[] vec_dir_tsdg;
+    delete[] vec_dir_psdg;
+    delete[] vec_dir_scand;
+    delete[] vec_dir_cfmp;
+    delete[] vec_dir_sfmp;
+    delete[] vec_dir_cfsp;
+    delete[] vec_dir_sfsp;
+    delete[] vec_dir_un;
+    delete[] vec_dir_uns;
+    delete[] vec_dir_sas11;
+    delete[] vec_dir_sas21;
+    delete[] vec_dir_sas12;
+    delete[] vec_dir_sas22;
+    delete[] vec_dir_muls;
+    delete[] vec_dir_mulslr;
+    delete[] vec_dir_sat11;
+    delete[] vec_dir_sat21;
+    delete[] vec_dir_sat12;
+    delete[] vec_dir_sat22;
+    delete[] vec_dir_scc1;
+    delete[] vec_dir_scc2;
+    delete[] vec_dir_abc1;
+    delete[] vec_dir_abc2;
+    delete[] vec_dir_exc1;
+    delete[] vec_dir_exc2;
+    delete[] vec_dir_albedc1;
+    delete[] vec_dir_albedc2;
+    delete[] vec_dir_qscc1;
+    delete[] vec_dir_qscc2;
+    delete[] vec_dir_qabc1;
+    delete[] vec_dir_qabc2;
+    delete[] vec_dir_qexc1;
+    delete[] vec_dir_qexc2;
+    delete[] vec_qscamc1;
+    delete[] vec_qscamc2;
+    delete[] vec_qabsmc1;
+    delete[] vec_qabsmc2;
+    delete[] vec_qextmc1;
+    delete[] vec_qextmc2;
+    delete[] vec_dir_sccrt1;
+    delete[] vec_dir_sccrt2;
+    delete[] vec_dir_abcrt1;
+    delete[] vec_dir_abcrt2;
+    delete[] vec_dir_excrt1;
+    delete[] vec_dir_excrt2;
+    delete[] vec_dir_fsac11;
+    delete[] vec_dir_fsac21;
+    delete[] vec_dir_fsac12;
+    delete[] vec_dir_fsac22;
+    delete[] vec_dir_sac11;
+    delete[] vec_dir_sac21;
+    delete[] vec_dir_sac12;
+    delete[] vec_dir_sac22;
+    delete[] vec_dir_qschuc1;
+    delete[] vec_dir_qschuc2;
+    delete[] vec_dir_pschuc1;
+    delete[] vec_dir_pschuc2;
+    delete[] vec_dir_s0magc1;
+    delete[] vec_dir_s0magc2;
+    delete[] vec_dir_cosavc1;
+    delete[] vec_dir_cosavc2;
+    delete[] vec_dir_raprc1;
+    delete[] vec_dir_raprc2;
+    delete[] vec_dir_flc1;
+    delete[] vec_dir_flc2;
+    delete[] vec_dir_frc1;
+    delete[] vec_dir_frc2;
+    delete[] vec_dir_fkc1;
+    delete[] vec_dir_fkc2;
+    delete[] vec_dir_fxc1;
+    delete[] vec_dir_fxc2;
+    delete[] vec_dir_fyc1;
+    delete[] vec_dir_fyc2;
+    delete[] vec_dir_fzc1;
+    delete[] vec_dir_fzc2;
+    delete[] vec_dir_tqelc1;
+    delete[] vec_dir_tqelc2;
+    delete[] vec_dir_tqerc1;
+    delete[] vec_dir_tqerc2;
+    delete[] vec_dir_tqekc1;
+    delete[] vec_dir_tqekc2;
+    delete[] vec_dir_tqexc1;
+    delete[] vec_dir_tqexc2;
+    delete[] vec_dir_tqeyc1;
+    delete[] vec_dir_tqeyc2;
+    delete[] vec_dir_tqezc1;
+    delete[] vec_dir_tqezc2;
+    delete[] vec_dir_tqslc1;
+    delete[] vec_dir_tqslc2;
+    delete[] vec_dir_tqsrc1;
+    delete[] vec_dir_tqsrc2;
+    delete[] vec_dir_tqskc1;
+    delete[] vec_dir_tqskc2;
+    delete[] vec_dir_tqsxc1;
+    delete[] vec_dir_tqsxc2;
+    delete[] vec_dir_tqsyc1;
+    delete[] vec_dir_tqsyc2;
+    delete[] vec_dir_tqszc1;
+    delete[] vec_dir_tqszc2;
+    delete[] vec_dir_mulc;
+    delete[] vec_dir_mulclr;
+  }
 }
 
 long ClusterOutputInfo::compute_size(
@@ -2356,13 +2372,13 @@ int ClusterOutputInfo::write_legacy(const std::string &output) {
 #ifdef MPI_VERSION
 int ClusterOutputInfo::mpireceive(const mixMPI *mpidata, int pid) {
   int result = 0;
-  int skip_flag;
+  int flag;
   int chk_nsph, chk_inpol, chk_iavm, chk_isam, chk_num_theta, chk_num_thetas;
   int chk_num_phi, chk_num_phis, chk_ndirs, chk_idfc, chk_configs;
   double chk_exri;
-  MPI_Recv(&skip_flag, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
-  // Proceed with the rest _only if__ skip_flag==0, else nothing is to be received
-  if (skip_flag == 0) {
+  MPI_Recv(&flag, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+  // Proceed with the rest _only if_ skip_flag==0, else nothing is to be received
+  if (flag == 0) {
     MPI_Recv(&chk_nsph, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     MPI_Recv(&chk_inpol, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     MPI_Recv(&chk_iavm, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
@@ -2562,191 +2578,200 @@ int ClusterOutputInfo::mpireceive(const mixMPI *mpidata, int pid) {
   return result;
 }
 
-
 int ClusterOutputInfo::mpisend(const mixMPI *mpidata) {
   int result = 0;
   int chunk_size;
-  // Send output metadata for configuration cross-check
-  MPI_Send(&nsph, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&inpol, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&iavm, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&isam, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_theta, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_thetas, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_phi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_phis, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&ndirs, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&exri, 1, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&idfc, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&configurations, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  // Wait for process 0 to cross-check the configuration
-  MPI_Recv(&result, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
-  if (result == 0) {
-    // Process 0 confirmed the consistency of configuration. Send the data.
-    // Send vectors of single values per scale
-    MPI_Send(&_first_xi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(&xi_block_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_jxi, xi_block_size, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_ier, xi_block_size, MPI_SHORT, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_vk, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_xi, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsat, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0magt, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_albedc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_albedc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qscamc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qscamc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qabsmc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qabsmc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qextmc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qextmc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sccrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sccrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abcrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abcrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_excrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_excrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsac11, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsac21, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsac22, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsac12, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschuc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschuc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschuc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschuc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0magc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0magc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_cosavc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_cosavc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_raprc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_raprc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fkc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fkc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+  if (_skip_flag == 1) {
+    // tell the receiver we are not sending anything
+    int flag = 1;
+    MPI_Send(&flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+  }
+  else {
+    // tell the receiver we are sending actual stuff
+    int flag = 0;
+    MPI_Send(&flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    // Send output metadata for configuration cross-check
+    MPI_Send(&nsph, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&inpol, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&iavm, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&isam, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_theta, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_thetas, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_phi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_phis, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&ndirs, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&exri, 1, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&idfc, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&configurations, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    // Wait for process 0 to cross-check the configuration
+    MPI_Recv(&result, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+    if (result == 0) {
+      // Process 0 confirmed the consistency of configuration. Send the data.
+      // Send vectors of single values per scale
+      MPI_Send(&_first_xi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(&xi_block_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_jxi, xi_block_size, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_ier, xi_block_size, MPI_SHORT, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_vk, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_xi, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsat, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0magt, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_albedc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_albedc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qscamc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qscamc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qabsmc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qabsmc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qextmc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qextmc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sccrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sccrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abcrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abcrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_excrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_excrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsac11, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsac21, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsac22, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsac12, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschuc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschuc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschuc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschuc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0magc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0magc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_cosavc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_cosavc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_raprc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_raprc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fkc1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fkc2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors of multiple configuration values per scale
-    chunk_size = xi_block_size * configurations;
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_sizes, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_ref_indices, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_scs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_abs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_exs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_albs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_sqscs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_sqabs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_sqexs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsas, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschus, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschus, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0mags, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_cosavs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_raprs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqek1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqsk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqek2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqsk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors of multiple configuration values per scale
+      chunk_size = xi_block_size * configurations;
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_sizes, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_ref_indices, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_scs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_abs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_exs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_albs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_sqscs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_sqabs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_sqexs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsas, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschus, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschus, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0mags, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_cosavs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_raprs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqek1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqsk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqek2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqsk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors whose sizes depend on directions and configurations.
-    chunk_size = ndirs * configurations * xi_block_size;
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_muls, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_mulslr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on directions and configurations.
+      chunk_size = ndirs * configurations * xi_block_size;
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_muls, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_mulslr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors whose sizes depend on directions and scales.
-    chunk_size = xi_block_size * ndirs;
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sat11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sat21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sat12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sat22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_scc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_scc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_abc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_abc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_exc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_exc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_albedc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_albedc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qscc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qscc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qabc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qabc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qexc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qexc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sccrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sccrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_abcrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_abcrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_excrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_excrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsac11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsac21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsac12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsac22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sac11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sac21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sac12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sac22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qschuc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qschuc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_pschuc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_pschuc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_s0magc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_s0magc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_cosavc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_cosavc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_raprc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_raprc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_flc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_flc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_frc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_frc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fkc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fkc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fxc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fxc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fyc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fyc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fzc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fzc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqelc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqelc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqerc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqerc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqekc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqekc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqexc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqexc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqeyc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqeyc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqezc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqezc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqslc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqslc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsrc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsrc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqskc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqskc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsxc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsxc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsyc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsyc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqszc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqszc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_mulc, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_mulclr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on directions and scales.
+      chunk_size = xi_block_size * ndirs;
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sat11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sat21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sat12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sat22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_scc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_scc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_abc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_abc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_exc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_exc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_albedc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_albedc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qscc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qscc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qabc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qabc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qexc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qexc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sccrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sccrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_abcrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_abcrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_excrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_excrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsac11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsac21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsac12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsac22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sac11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sac21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sac12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sac22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qschuc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qschuc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_pschuc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_pschuc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_s0magc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_s0magc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_cosavc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_cosavc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_raprc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_raprc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_flc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_flc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_frc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_frc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fkc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fkc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fxc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fxc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fyc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fyc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fzc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fzc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqelc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqelc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqerc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqerc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqekc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqekc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqexc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqexc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqeyc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqeyc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqezc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqezc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqslc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqslc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsrc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsrc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqskc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqskc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsxc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsxc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsyc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsyc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqszc1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqszc2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_mulc, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_mulclr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+    }
   }
   return result;
 }
@@ -2758,6 +2783,7 @@ InclusionOutputInfo::InclusionOutputInfo(
   ScattererConfiguration *sc, GeometryConfiguration *gc,
   const mixMPI *mpidata, int first_xi, int xi_length
 ) {
+  _skip_flag = 0;
   nsph = gc->number_of_spheres;
   li = gc->li;
   le = gc->le;
@@ -2946,6 +2972,7 @@ InclusionOutputInfo::InclusionOutputInfo(const std::string &hdf5_name) {
   HDFFile *hdf_file = new HDFFile(hdf5_name, flags);
   herr_t status = hdf_file->get_status();
   string str_name, str_type;
+  _skip_flag = 0;
   if (status == 0) {
     status = hdf_file->read("NSPH", "INT32_(1)", &nsph);
     status = hdf_file->read("LI", "INT32_(1)", &li);
@@ -3260,127 +3287,141 @@ InclusionOutputInfo::InclusionOutputInfo(const std::string &hdf5_name) {
   }
 }
 
+InclusionOutputInfo::InclusionOutputInfo(const int flag) {
+  /*
+    create a dummy placeholder just to know I should skip MPI_Send and MPI_Recv
+  */
+  if (flag == 1) {
+    _skip_flag = 1;
+  } else {
+    UnrecognizedOutputInfo ex(flag);
+    throw ex;
+  }
+}
+
 InclusionOutputInfo::~InclusionOutputInfo() {
-  delete[] vec_x_coords;
-  delete[] vec_y_coords;
-  delete[] vec_z_coords;
-  delete[] vec_jxi;
-  delete[] vec_ier;
-  delete[] vec_vk;
-  delete[] vec_xi;
-  delete[] vec_sphere_sizes;
-  delete[] vec_sphere_ref_indices;
-  delete[] vec_scs1;
-  delete[] vec_scs2;
-  delete[] vec_abs1;
-  delete[] vec_abs2;
-  delete[] vec_exs1;
-  delete[] vec_exs2;
-  delete[] vec_albeds1;
-  delete[] vec_albeds2;
-  delete[] vec_scsrt1;
-  delete[] vec_scsrt2;
-  delete[] vec_absrt1;
-  delete[] vec_absrt2;
-  delete[] vec_exsrt1;
-  delete[] vec_exsrt2;
-  delete[] vec_qschu1;
-  delete[] vec_qschu2;
-  delete[] vec_pschu1;
-  delete[] vec_pschu2;
-  delete[] vec_s0mag1;
-  delete[] vec_s0mag2;
-  delete[] vec_cosav1;
-  delete[] vec_cosav2;
-  delete[] vec_raprs1;
-  delete[] vec_raprs2;
-  delete[] vec_fk1;
-  delete[] vec_fk2;
-  delete[] vec_fsas11;
-  delete[] vec_fsas21;
-  delete[] vec_fsas22;
-  delete[] vec_fsas12;
-  delete[] vec_dir_tidg;
-  delete[] vec_dir_pidg;
-  delete[] vec_dir_tsdg;
-  delete[] vec_dir_psdg;
-  delete[] vec_dir_scand;
-  delete[] vec_dir_cfmp;
-  delete[] vec_dir_sfmp;
-  delete[] vec_dir_cfsp;
-  delete[] vec_dir_sfsp;
-  delete[] vec_dir_un;
-  delete[] vec_dir_uns;
-  delete[] vec_dir_scs1;
-  delete[] vec_dir_scs2;
-  delete[] vec_dir_abs1;
-  delete[] vec_dir_abs2;
-  delete[] vec_dir_exs1;
-  delete[] vec_dir_exs2;
-  delete[] vec_dir_albeds1;
-  delete[] vec_dir_albeds2;
-  delete[] vec_dir_scsrt1;
-  delete[] vec_dir_scsrt2;
-  delete[] vec_dir_absrt1;
-  delete[] vec_dir_absrt2;
-  delete[] vec_dir_exsrt1;
-  delete[] vec_dir_exsrt2;
-  delete[] vec_dir_fsas11;
-  delete[] vec_dir_fsas21;
-  delete[] vec_dir_fsas12;
-  delete[] vec_dir_fsas22;
-  delete[] vec_dir_sas11;
-  delete[] vec_dir_sas21;
-  delete[] vec_dir_sas12;
-  delete[] vec_dir_sas22;
-  delete[] vec_dir_qschu1;
-  delete[] vec_dir_qschu2;
-  delete[] vec_dir_pschu1;
-  delete[] vec_dir_pschu2;
-  delete[] vec_dir_s0mag1;
-  delete[] vec_dir_s0mag2;
-  delete[] vec_dir_cosav1;
-  delete[] vec_dir_cosav2;
-  delete[] vec_dir_rapr1;
-  delete[] vec_dir_rapr2;
-  delete[] vec_dir_fl1;
-  delete[] vec_dir_fl2;
-  delete[] vec_dir_fr1;
-  delete[] vec_dir_fr2;
-  delete[] vec_dir_fk1;
-  delete[] vec_dir_fk2;
-  delete[] vec_dir_fx1;
-  delete[] vec_dir_fx2;
-  delete[] vec_dir_fy1;
-  delete[] vec_dir_fy2;
-  delete[] vec_dir_fz1;
-  delete[] vec_dir_fz2;
-  delete[] vec_dir_tqel1;
-  delete[] vec_dir_tqel2;
-  delete[] vec_dir_tqer1;
-  delete[] vec_dir_tqer2;
-  delete[] vec_dir_tqek1;
-  delete[] vec_dir_tqek2;
-  delete[] vec_dir_tqex1;
-  delete[] vec_dir_tqex2;
-  delete[] vec_dir_tqey1;
-  delete[] vec_dir_tqey2;
-  delete[] vec_dir_tqez1;
-  delete[] vec_dir_tqez2;
-  delete[] vec_dir_tqsl1;
-  delete[] vec_dir_tqsl2;
-  delete[] vec_dir_tqsr1;
-  delete[] vec_dir_tqsr2;
-  delete[] vec_dir_tqsk1;
-  delete[] vec_dir_tqsk2;
-  delete[] vec_dir_tqsx1;
-  delete[] vec_dir_tqsx2;
-  delete[] vec_dir_tqsy1;
-  delete[] vec_dir_tqsy2;
-  delete[] vec_dir_tqsz1;
-  delete[] vec_dir_tqsz2;
-  delete[] vec_dir_mull;
-  delete[] vec_dir_mulllr;
+  if (_skip_flag != 1) {
+    delete[] vec_x_coords;
+    delete[] vec_y_coords;
+    delete[] vec_z_coords;
+    delete[] vec_jxi;
+    delete[] vec_ier;
+    delete[] vec_vk;
+    delete[] vec_xi;
+    delete[] vec_sphere_sizes;
+    delete[] vec_sphere_ref_indices;
+    delete[] vec_scs1;
+    delete[] vec_scs2;
+    delete[] vec_abs1;
+    delete[] vec_abs2;
+    delete[] vec_exs1;
+    delete[] vec_exs2;
+    delete[] vec_albeds1;
+    delete[] vec_albeds2;
+    delete[] vec_scsrt1;
+    delete[] vec_scsrt2;
+    delete[] vec_absrt1;
+    delete[] vec_absrt2;
+    delete[] vec_exsrt1;
+    delete[] vec_exsrt2;
+    delete[] vec_qschu1;
+    delete[] vec_qschu2;
+    delete[] vec_pschu1;
+    delete[] vec_pschu2;
+    delete[] vec_s0mag1;
+    delete[] vec_s0mag2;
+    delete[] vec_cosav1;
+    delete[] vec_cosav2;
+    delete[] vec_raprs1;
+    delete[] vec_raprs2;
+    delete[] vec_fk1;
+    delete[] vec_fk2;
+    delete[] vec_fsas11;
+    delete[] vec_fsas21;
+    delete[] vec_fsas22;
+    delete[] vec_fsas12;
+    delete[] vec_dir_tidg;
+    delete[] vec_dir_pidg;
+    delete[] vec_dir_tsdg;
+    delete[] vec_dir_psdg;
+    delete[] vec_dir_scand;
+    delete[] vec_dir_cfmp;
+    delete[] vec_dir_sfmp;
+    delete[] vec_dir_cfsp;
+    delete[] vec_dir_sfsp;
+    delete[] vec_dir_un;
+    delete[] vec_dir_uns;
+    delete[] vec_dir_scs1;
+    delete[] vec_dir_scs2;
+    delete[] vec_dir_abs1;
+    delete[] vec_dir_abs2;
+    delete[] vec_dir_exs1;
+    delete[] vec_dir_exs2;
+    delete[] vec_dir_albeds1;
+    delete[] vec_dir_albeds2;
+    delete[] vec_dir_scsrt1;
+    delete[] vec_dir_scsrt2;
+    delete[] vec_dir_absrt1;
+    delete[] vec_dir_absrt2;
+    delete[] vec_dir_exsrt1;
+    delete[] vec_dir_exsrt2;
+    delete[] vec_dir_fsas11;
+    delete[] vec_dir_fsas21;
+    delete[] vec_dir_fsas12;
+    delete[] vec_dir_fsas22;
+    delete[] vec_dir_sas11;
+    delete[] vec_dir_sas21;
+    delete[] vec_dir_sas12;
+    delete[] vec_dir_sas22;
+    delete[] vec_dir_qschu1;
+    delete[] vec_dir_qschu2;
+    delete[] vec_dir_pschu1;
+    delete[] vec_dir_pschu2;
+    delete[] vec_dir_s0mag1;
+    delete[] vec_dir_s0mag2;
+    delete[] vec_dir_cosav1;
+    delete[] vec_dir_cosav2;
+    delete[] vec_dir_rapr1;
+    delete[] vec_dir_rapr2;
+    delete[] vec_dir_fl1;
+    delete[] vec_dir_fl2;
+    delete[] vec_dir_fr1;
+    delete[] vec_dir_fr2;
+    delete[] vec_dir_fk1;
+    delete[] vec_dir_fk2;
+    delete[] vec_dir_fx1;
+    delete[] vec_dir_fx2;
+    delete[] vec_dir_fy1;
+    delete[] vec_dir_fy2;
+    delete[] vec_dir_fz1;
+    delete[] vec_dir_fz2;
+    delete[] vec_dir_tqel1;
+    delete[] vec_dir_tqel2;
+    delete[] vec_dir_tqer1;
+    delete[] vec_dir_tqer2;
+    delete[] vec_dir_tqek1;
+    delete[] vec_dir_tqek2;
+    delete[] vec_dir_tqex1;
+    delete[] vec_dir_tqex2;
+    delete[] vec_dir_tqey1;
+    delete[] vec_dir_tqey2;
+    delete[] vec_dir_tqez1;
+    delete[] vec_dir_tqez2;
+    delete[] vec_dir_tqsl1;
+    delete[] vec_dir_tqsl2;
+    delete[] vec_dir_tqsr1;
+    delete[] vec_dir_tqsr2;
+    delete[] vec_dir_tqsk1;
+    delete[] vec_dir_tqsk2;
+    delete[] vec_dir_tqsx1;
+    delete[] vec_dir_tqsx2;
+    delete[] vec_dir_tqsy1;
+    delete[] vec_dir_tqsy2;
+    delete[] vec_dir_tqsz1;
+    delete[] vec_dir_tqsz2;
+    delete[] vec_dir_mull;
+    delete[] vec_dir_mulllr;
+  }
 }
 
 long InclusionOutputInfo::compute_size() {
@@ -4485,13 +4526,13 @@ int InclusionOutputInfo::write_legacy(const std::string &output) {
 #ifdef MPI_VERSION
 int InclusionOutputInfo::mpireceive(const mixMPI* mpidata, int pid) {
   int result = 0;
-  int skip_flag;
+  int flag;
   int chk_nsph, chk_inpol, chk_iavm, chk_isam, chk_num_theta, chk_num_thetas;
   int chk_num_phi, chk_num_phis, chk_ndirs, chk_idfc, chk_configs;
   double chk_exri;
-  MPI_Recv(&skip_flag, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
-  // Proceed with the rest _only if__ skip_flag==0, else nothing is to be received
-  if (skip_flag == 0) {
+  MPI_Recv(&flag, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+  // Proceed with the rest _only if_ flag==0, else nothing is to be received
+  if (flag == 0) {
     MPI_Recv(&chk_nsph, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     MPI_Recv(&chk_inpol, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     MPI_Recv(&chk_iavm, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
@@ -4648,140 +4689,150 @@ int InclusionOutputInfo::mpireceive(const mixMPI* mpidata, int pid) {
 int InclusionOutputInfo::mpisend(const mixMPI *mpidata) {
   int result = 0;
   int chunk_size;
-  // Send output metadata for configuration cross-check
-  MPI_Send(&nsph, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&inpol, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&iavm, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&isam, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_theta, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_thetas, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_phi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_phis, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&ndirs, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&exri, 1, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&idfc, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&configurations, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  // Wait for process 0 to cross-check the configuration
-  MPI_Recv(&result, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
-  if (result == 0) {
-    // Process 0 confirmed the consistency of configuration. Send the data.
-    // Send vectors of single values per scale
-    MPI_Send(&_first_xi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(&xi_block_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_jxi, xi_block_size, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_ier, xi_block_size, MPI_SHORT, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_vk, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_xi, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_albeds1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_albeds2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scsrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scsrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_absrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_absrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exsrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exsrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschu1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschu2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschu1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschu2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0mag1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0mag2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_cosav1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_cosav2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_raprs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_raprs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fk1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fk2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsas11, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsas21, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsas22, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsas12, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+  if (_skip_flag == 1) {
+    // tell the receiver we are not sending anything
+    int flag = 1;
+    MPI_Send(&flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+  }
+  else {
+    // tell the receiver we are sending actual stuff
+    int flag = 0;
+    MPI_Send(&flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    // Send output metadata for configuration cross-check
+    MPI_Send(&nsph, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&inpol, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&iavm, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&isam, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_theta, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_thetas, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_phi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_phis, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&ndirs, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&exri, 1, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&idfc, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&configurations, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    // Wait for process 0 to cross-check the configuration
+    MPI_Recv(&result, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+    if (result == 0) {
+      // Process 0 confirmed the consistency of configuration. Send the data.
+      // Send vectors of single values per scale
+      MPI_Send(&_first_xi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(&xi_block_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_jxi, xi_block_size, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_ier, xi_block_size, MPI_SHORT, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_vk, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_xi, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_albeds1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_albeds2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scsrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scsrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_absrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_absrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exsrt1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exsrt2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschu1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschu2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschu1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschu2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0mag1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0mag2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_cosav1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_cosav2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_raprs1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_raprs2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fk1, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fk2, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsas11, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsas21, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsas22, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsas12, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors whose sizes depend on configurations
-    chunk_size = xi_block_size * (configurations + 1);
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_sizes, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_ref_indices, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on configurations
+      chunk_size = xi_block_size * (configurations + 1);
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_sizes, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_ref_indices, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors whose sizes depend on directions and wavelengths
-    chunk_size = xi_block_size * ndirs;
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_scs1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_scs2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_abs1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_abs2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_exs1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_exs2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_albeds1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_albeds2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_scsrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_scsrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_absrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_absrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_exsrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_exsrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fsas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qschu1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_qschu2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_pschu1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_pschu2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_s0mag1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_s0mag2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_cosav1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_cosav2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_rapr1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_rapr2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fl1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fl2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fr1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fr2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fx1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fx2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fy1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fy2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fz1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fz2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqel1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqel2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqer1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqer2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqek1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqek2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqex1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqex2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqey1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqey2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqez1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqez2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsl1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsl2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsr1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsr2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsx1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsx2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsy1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsy2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsz1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_tqsz2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_mull, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_mulllr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on directions and wavelengths
+      chunk_size = xi_block_size * ndirs;
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_scs1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_scs2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_abs1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_abs2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_exs1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_exs2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_albeds1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_albeds2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_scsrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_scsrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_absrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_absrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_exsrt1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_exsrt2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fsas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qschu1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_qschu2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_pschu1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_pschu2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_s0mag1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_s0mag2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_cosav1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_cosav2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_rapr1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_rapr2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fl1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fl2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fr1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fr2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fx1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fx2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fy1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fy2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fz1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fz2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqel1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqel2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqer1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqer2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqek1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqek2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqex1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqex2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqey1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqey2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqez1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqez2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsl1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsl2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsr1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsr2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsx1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsx2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsy1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsy2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsz1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_tqsz2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_mull, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_mulllr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+    }
   }
   return result;
 }
@@ -4793,6 +4844,7 @@ SphereOutputInfo::SphereOutputInfo(
   ScattererConfiguration *sc, GeometryConfiguration *gc,
   const mixMPI *mpidata, int first_xi, int xi_length
 ) {
+  _skip_flag = 0;
   _first_xi = first_xi;
   nsph = gc->number_of_spheres;
   lm = gc->l_max;
@@ -4905,6 +4957,7 @@ SphereOutputInfo::SphereOutputInfo(const std::string &hdf5_name) {
   HDFFile *hdf_file = new HDFFile(hdf5_name, flags);
   herr_t status = hdf_file->get_status();
   string str_name, str_type;
+  _skip_flag = 0;
   if (status == 0) {
     status = hdf_file->read("NSPH", "INT32_(1)", &nsph);
     status = hdf_file->read("LM", "INT32_(1)", &lm);
@@ -5077,56 +5130,70 @@ SphereOutputInfo::SphereOutputInfo(const std::string &hdf5_name) {
   }
 }
 
+SphereOutputInfo::SphereOutputInfo(const int flag) {
+  /*
+    create a dummy placeholder just to know I should skip MPI_Send and MPI_Recv
+  */
+  if (flag == 1) {
+    _skip_flag = 1;
+  } else {
+    UnrecognizedOutputInfo ex(flag);
+    throw ex;
+  }
+}
+
 SphereOutputInfo::~SphereOutputInfo() {
-  delete[] vec_jxi;
-  delete[] vec_ier;
-  delete[] vec_vk;
-  delete[] vec_xi;
-  delete[] vec_sphere_sizes;
-  delete[] vec_sphere_ref_indices;
-  delete[] vec_scs;
-  delete[] vec_abs;
-  delete[] vec_exs;
-  delete[] vec_albeds;
-  delete[] vec_scsrt;
-  delete[] vec_absrt;
-  delete[] vec_exsrt;
-  delete[] vec_fsas;
-  delete[] vec_qschu;
-  delete[] vec_pschu;
-  delete[] vec_s0mag;
-  delete[] vec_cosav;
-  delete[] vec_raprs;
-  delete[] vec_tqek1;
-  delete[] vec_tqek2;
-  delete[] vec_tqsk1;
-  delete[] vec_tqsk2;
-  if (nsph != 1) {
-    delete[] vec_fsat;
-    delete[] vec_qschut;
-    delete[] vec_pschut;
-    delete[] vec_s0magt;
+  if (_skip_flag != 1) {
+    delete[] vec_jxi;
+    delete[] vec_ier;
+    delete[] vec_vk;
+    delete[] vec_xi;
+    delete[] vec_sphere_sizes;
+    delete[] vec_sphere_ref_indices;
+    delete[] vec_scs;
+    delete[] vec_abs;
+    delete[] vec_exs;
+    delete[] vec_albeds;
+    delete[] vec_scsrt;
+    delete[] vec_absrt;
+    delete[] vec_exsrt;
+    delete[] vec_fsas;
+    delete[] vec_qschu;
+    delete[] vec_pschu;
+    delete[] vec_s0mag;
+    delete[] vec_cosav;
+    delete[] vec_raprs;
+    delete[] vec_tqek1;
+    delete[] vec_tqek2;
+    delete[] vec_tqsk1;
+    delete[] vec_tqsk2;
+    if (nsph != 1) {
+      delete[] vec_fsat;
+      delete[] vec_qschut;
+      delete[] vec_pschut;
+      delete[] vec_s0magt;
+    }
+    delete[] vec_dir_tidg;
+    delete[] vec_dir_pidg;
+    delete[] vec_dir_tsdg;
+    delete[] vec_dir_psdg;
+    delete[] vec_dir_scand;
+    delete[] vec_dir_cfmp;
+    delete[] vec_dir_cfsp;
+    delete[] vec_dir_sfmp;
+    delete[] vec_dir_sfsp;
+    delete[] vec_dir_un;
+    delete[] vec_dir_uns;
+    delete[] vec_dir_sas11;
+    delete[] vec_dir_sas21;
+    delete[] vec_dir_sas12;
+    delete[] vec_dir_sas22;
+    delete[] vec_dir_fx;
+    delete[] vec_dir_fy;
+    delete[] vec_dir_fz;
+    delete[] vec_dir_muls;
+    delete[] vec_dir_mulslr;
   }
-  delete[] vec_dir_tidg;
-  delete[] vec_dir_pidg;
-  delete[] vec_dir_tsdg;
-  delete[] vec_dir_psdg;
-  delete[] vec_dir_scand;
-  delete[] vec_dir_cfmp;
-  delete[] vec_dir_cfsp;
-  delete[] vec_dir_sfmp;
-  delete[] vec_dir_sfsp;
-  delete[] vec_dir_un;
-  delete[] vec_dir_uns;
-  delete[] vec_dir_sas11;
-  delete[] vec_dir_sas21;
-  delete[] vec_dir_sas12;
-  delete[] vec_dir_sas22;
-  delete[] vec_dir_fx;
-  delete[] vec_dir_fy;
-  delete[] vec_dir_fz;
-  delete[] vec_dir_muls;
-  delete[] vec_dir_mulslr;
 }
 
 long SphereOutputInfo::compute_size(
@@ -5746,13 +5813,13 @@ int SphereOutputInfo::write_legacy(const std::string &file_name) {
 #ifdef MPI_VERSION
 int SphereOutputInfo::mpireceive(const mixMPI *mpidata, int pid) {
   int result = 0;
-  int skip_flag;
+  int flag;
   int chk_nsph, chk_inpol, chk_isam, chk_num_theta, chk_num_thetas;
   int chk_num_phi, chk_num_phis, chk_ndirs, chk_idfc, chk_configs;
   double chk_exri;
-  MPI_Recv(&skip_flag, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
-  // Proceed with the rest _only if__ skip_flag==0, else nothing is to be received
-  if (skip_flag == 0) {
+  MPI_Recv(&flag, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+  // Proceed with the rest _only if_ flag==0, else nothing is to be received
+  if (flag == 0) {
     MPI_Recv(&chk_nsph, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     MPI_Recv(&chk_inpol, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     MPI_Recv(&chk_isam, 1, MPI_INT32_T, pid, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
@@ -5842,74 +5909,84 @@ int SphereOutputInfo::mpireceive(const mixMPI *mpidata, int pid) {
 int SphereOutputInfo::mpisend(const mixMPI *mpidata) {
   int result = 0;
   int chunk_size;
-  // Send output metadata for configuration cross-check
-  MPI_Send(&nsph, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&inpol, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&isam, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_theta, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_thetas, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_phi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&_num_phis, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&ndirs, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&exri, 1, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&idfc, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  MPI_Send(&configurations, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-  // Wait for process 0 to cross-check the configuration
-  MPI_Recv(&result, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
-  if (result == 0) {
-    // Process 0 confirmed the consistency of configuration. Send the data.
-    // Send vectors of single values per scale
-    MPI_Send(&_first_xi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(&xi_block_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_jxi, xi_block_size, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_ier, xi_block_size, MPI_SHORT, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_vk, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_xi, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    if (nsph != 1) {
-      MPI_Send(vec_fsat, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-      MPI_Send(vec_qschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-      MPI_Send(vec_pschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-      MPI_Send(vec_s0magt, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    }
+  if (_skip_flag == 1) {
+    // tell the receiver we are not sending anything
+    int flag = 1;
+    MPI_Send(&flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+  }
+  else {
+    // tell the receiver we are sending actual stuff
+    int flag = 0;
+    MPI_Send(&flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    // Send output metadata for configuration cross-check
+    MPI_Send(&nsph, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&inpol, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&isam, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_theta, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_thetas, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_phi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&_num_phis, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&ndirs, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&exri, 1, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&idfc, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    MPI_Send(&configurations, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+    // Wait for process 0 to cross-check the configuration
+    MPI_Recv(&result, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+    if (result == 0) {
+      // Process 0 confirmed the consistency of configuration. Send the data.
+      // Send vectors of single values per scale
+      MPI_Send(&_first_xi, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(&xi_block_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_jxi, xi_block_size, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_ier, xi_block_size, MPI_SHORT, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_vk, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_xi, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      if (nsph != 1) {
+	MPI_Send(vec_fsat, xi_block_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+	MPI_Send(vec_qschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+	MPI_Send(vec_pschut, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+	MPI_Send(vec_s0magt, xi_block_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      }
 
-    // Send vectors whose sizes depend on configurations and scales
-    chunk_size = xi_block_size * configurations;
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_sizes, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_sphere_ref_indices, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_abs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_albeds, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_scsrt, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_absrt, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_exsrt, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_fsas, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_qschu, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_pschu, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_s0mag, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_cosav, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_raprs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqek1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqek2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqsk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_tqsk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on configurations and scales
+      chunk_size = xi_block_size * configurations;
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_sizes, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_sphere_ref_indices, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_abs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_albeds, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_scsrt, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_absrt, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_exsrt, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_fsas, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_qschu, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_pschu, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_s0mag, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_cosav, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_raprs, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqek1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqek2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqsk1, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_tqsk2, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors whose sizes depend on NSPH, directions and wavelengths
-    chunk_size = xi_block_size * nsph * ndirs;
-    MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_sas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_muls, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_mulslr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on NSPH, directions and wavelengths
+      chunk_size = xi_block_size * nsph * ndirs;
+      MPI_Send(&chunk_size, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas11, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas21, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas12, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_sas22, chunk_size, MPI_C_DOUBLE_COMPLEX, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_muls, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_mulslr, 16 * chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
 
-    // Send vectors whose sizes depend on NSPH, incidence directions and wavelengths
-    chunk_size = xi_block_size * nsph * _num_theta * _num_phi;
-    MPI_Send(vec_dir_fx, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fy, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
-    MPI_Send(vec_dir_fz, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      // Send vectors whose sizes depend on NSPH, incidence directions and wavelengths
+      chunk_size = xi_block_size * nsph * _num_theta * _num_phi;
+      MPI_Send(vec_dir_fx, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fy, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+      MPI_Send(vec_dir_fz, chunk_size, MPI_DOUBLE, 0, 10, MPI_COMM_WORLD);
+    }
   }
   return result;
 }
diff --git a/src/sphere/sphere.cpp b/src/sphere/sphere.cpp
index 18d19b80bc716356e40f39449d21bde53c2af666..44c2a24c6c255719cb3a0a5b9d236680f07c7c25 100644
--- a/src/sphere/sphere.cpp
+++ b/src/sphere/sphere.cpp
@@ -330,19 +330,19 @@ void sphere(const string& config_file, const string& data_file, const string& ou
 	    int jer = sphere_jxi488_cycle(myjxi488 - 1, sconf, gconf, p_sa, sid_2, p_output_2, output_path, vtppoanp_2);
 	  } else {
 	    if (myompthread > 0) {
-	      // If there is no input for this thread, set output pointer to NULL.
-	      p_outarray[myompthread] = NULL;
+	      // If there is no input for this thread, mark to skip.
+	      p_outarray[myompthread] = new SphereOutputInfo(1);
 	    }
 	  }
 #pragma omp barrier
 	  // threads different from 0 append their virtual files to the one of thread 0, and delete them
 	  if (myompthread == 0) {
 	    for (int ti=1; ti<ompnumthreads; ti++) {
-	      if (p_outarray[ti] != NULL) {
+	      if (p_outarray[ti]->skip_flag != 1) {
 		p_outarray[0]->insert(*(p_outarray[ti]));
-		delete p_outarray[ti];
-		p_outarray[ti] = NULL;
 	      }
+	      delete p_outarray[ti];
+	      p_outarray[ti] = NULL;
 	      vtppoanarray[0]->append(*(vtppoanarray[ti]));
 	      delete vtppoanarray[ti];
 	    }
@@ -492,39 +492,27 @@ void sphere(const string& config_file, const string& data_file, const string& ou
 	  }
 	  int jer = sphere_jxi488_cycle(myjxi488 - 1, sconf, gconf, p_sa, sid_2, p_output_2, output_path, vtppoanp_2);
 	} else {
-	  // if (myompthread > 0) {
-	    // If there is no input for this thread, set the output pointer to NULL.
-	    p_outarray[myompthread] = NULL;
-	    //}	  
+	  p_outarray[myompthread] = new SphereOutputInfo(1);
 	}
 
 #pragma omp barrier
 	// threads different from 0 append their virtual files to the one of thread 0, and delete them
 	if (myompthread == 0) {
 	  for (int ti=1; ti<ompnumthreads; ti++) {
-	    if (p_outarray[ti] != NULL) {
+	    if (p_outarray[ti]->skip_flag != 1) {
 	      p_outarray[0]->insert(*(p_outarray[ti]));
-	      delete p_outarray[ti];
-	      p_outarray[ti] = NULL;
 	    }
+	    delete p_outarray[ti];
+	    p_outarray[ti] = NULL;
 	    vtppoanarray[0]->append(*(vtppoanarray[ti]));
 	    delete vtppoanarray[ti];
 	  }
 	  // thread 0 sends the collected virtualfiles to thread 0 of MPI process 0, then deletes them
 	  for (int rr=1; rr<mpidata->nprocs; rr++) {
 	    if (rr == mpidata->rank) {
-	      if (p_outarray[0] == NULL) {
-		// signal that we are not sending anything
-		int skip_flag = 1;
-		MPI_Send(&skip_flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-	      }
-	      else {
-		// signal that we are sending something
-		int skip_flag = 0;
-		MPI_Send(&skip_flag, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD);
-		p_outarray[0]->mpisend(mpidata);
-		delete p_outarray[0];
-	      }
+	      p_outarray[0]->mpisend(mpidata);
+	      delete p_outarray[0];
+	      p_outarray[0] = NULL;
 	      vtppoanarray[0]->mpisend(mpidata);
 	      delete vtppoanarray[0];
 	    }