diff --git a/src/feature_creation/feature_space/FeatureSpace.cpp b/src/feature_creation/feature_space/FeatureSpace.cpp
index 54fccb2989c44b7ae461f55d308165063014ea33..5877355c7d994da3a594f1eafa69580e8d92d770 100644
--- a/src/feature_creation/feature_space/FeatureSpace.cpp
+++ b/src/feature_creation/feature_space/FeatureSpace.cpp
@@ -60,15 +60,13 @@ void FeatureSpace::generate_feature_space()
     for(int nn = 1; nn <= _max_phi; ++nn)
     {
         std::vector<node_ptr> next_phi;
-
         _n_feat = _phi.size();
         // std::array<int, 2> start_end = _mpi_comm->get_start_end_for_iterator(_phi.size() - _start_gen[_start_gen.size()-1], _start_gen[_start_gen.size()-1]);
 
-
-        // next_phi.reserve(node_value_arrs::get_max_number_features(_allowed_ops, 1, start_end[1] - start_end[0]));
         int feat_ind = _phi.size();
 
-        for(auto feat_1 = _phi.begin() + _start_gen.back() + _mpi_comm->rank(); feat_1 != _phi.end(); feat_1 += _mpi_comm->size())
+        // for(auto feat_1 = _phi.begin() + start_end[0]; feat_1 != _phi.begin() + start_end[1]; ++feat_1)
+        for(auto feat_1 = _phi.begin() + _mpi_comm->rank(); feat_1 < _phi.end(); feat_1 += _mpi_comm->size())
         {
             for(auto& op : _un_operators)
             {
@@ -126,7 +124,7 @@ void FeatureSpace::generate_feature_space()
         }
         _mpi_comm->barrier();
         _start_gen.push_back(_phi.size());
-        if(nn < _max_phi)
+        if((nn < _max_phi) || (_mpi_comm->size() == 1))
         {
             std::vector<std::vector<node_ptr>> next_phi_gathered;
             mpi::all_gather(*_mpi_comm, next_phi, next_phi_gathered);
@@ -153,7 +151,6 @@ void FeatureSpace::generate_feature_space()
         }
         else
         {
-            std::cout << "HERE" << std::endl;
             std::vector<size_t> next_phi_sizes;
             mpi::all_gather(*_mpi_comm, next_phi.size(), next_phi_sizes);
             size_t n_feat = std::accumulate(next_phi_sizes.begin(), next_phi_sizes.end(), _phi.size());
@@ -197,9 +194,6 @@ void FeatureSpace::generate_feature_space()
             std::vector<size_t> recv_sizes(next_phi_sizes.size(), 0);
 
             // Is this rank sending or receiving?
-            std::cout << "HERE 2" << std::endl;
-            std::cout << std::accumulate(next_phi_excess.begin(), next_phi_excess.end(), 0);
-            std::cout << std::accumulate(next_phi_needed.begin(), next_phi_needed.end(), 0);
             if(next_phi_excess[_mpi_comm->rank()] > 0)
             {
                 size_t total_sent = std::accumulate(next_phi_excess.begin(), next_phi_excess.begin() + _mpi_comm->rank(), 0);
@@ -211,16 +205,16 @@ void FeatureSpace::generate_feature_space()
                     prev_sent_recv += next_phi_needed[ind];
                     ++ind;
                 }
-                send_size = std::min(next_phi.size(), next_phi_needed[ind-1]);
-                send_sizes[ind-1] = send_size;
+                send_size = std::min(next_phi.size(), next_phi_needed[ind]);
+                send_sizes[ind] = send_size;
                 total_sent = send_size;
                 while((total_sent < next_phi.size()) && (ind < _mpi_comm->size()))
                 {
                     send_size = std::min(next_phi.size() - total_sent, next_phi_needed[ind]);
                     send_sizes[ind] = send_size;
                     total_sent += send_size;
+                    ++ind;
                 }
-
                 total_sent = 0;
                 for(int pp = 0; pp < send_sizes.size(); ++pp)
                 {
@@ -244,14 +238,15 @@ void FeatureSpace::generate_feature_space()
                     prev_recv_sent += next_phi_excess[ind];
                     ++ind;
                 }
-                recv_size = std::min(next_phi.size(), next_phi_excess[ind-1]);
-                recv_sizes[ind-1] = recv_size;
+                recv_size = std::min(next_phi.size(), next_phi_excess[ind]);
+                recv_sizes[ind] = recv_size;
                 total_recv = recv_size;
                 while((_phi.size() < n_feat_rank) && (ind < _mpi_comm->size()))
                 {
                     recv_size = std::min(_phi.size() + total_recv, next_phi_excess[ind]);
                     recv_sizes[ind] = recv_size;
                     total_recv += recv_size;
+                    ++ind;
                 }
 
                 total_recv = 0;
@@ -269,7 +264,6 @@ void FeatureSpace::generate_feature_space()
                     }
                 }
             }
-
         }
     }
     _n_feat = _phi.size();
@@ -295,13 +289,13 @@ void FeatureSpace::sis(std::vector<double>& prop)
     project_r(prop.data());
 
     std::vector<int> inds = util_funcs::argsort(_scores);
-    std::vector<double> D_temp(_n_sis_select * prop.size());
+    std::vector<double> D_selected(_n_sis_select * prop.size());
     std::vector<double> scores_selected(_n_sis_select, 0.0);
-    std::vector<node_ptr> phi_slected(_n_sis_select, nullptr);
+    std::vector<node_ptr> phi_selected(_n_sis_select, nullptr);
 
     int ii = 0;
     int cur_feat_local = 0;
-    while((cur_feat_local != D_temp.size() / prop.size()) && (ii < _scores.size()))
+    while((cur_feat_local != D_selected.size() / prop.size()) && (ii < _scores.size()))
     {
         bool is_valid = true;
         for(int dd = 0; dd < cur_feat; ++dd)
@@ -314,7 +308,7 @@ void FeatureSpace::sis(std::vector<double>& prop)
         }
         for(int dd = 0; dd < cur_feat_local; ++dd)
         {
-            if(1.0 - std::abs(util_funcs::r(&D_temp[dd*prop.size()], _phi[inds[ii]]->value_ptr(), prop.size())) < 1e-13)
+            if(1.0 - std::abs(util_funcs::r(&D_selected[dd*prop.size()], _phi[inds[ii]]->value_ptr(), prop.size())) < 1e-13)
             {
                 is_valid = false;
                 break;
@@ -322,48 +316,57 @@ void FeatureSpace::sis(std::vector<double>& prop)
         }
         if(is_valid)
         {
-            std::copy_n(_phi[inds[ii]]->value_ptr(), prop.size(), &D_temp[cur_feat_local * prop.size()]);
+            std::copy_n(_phi[inds[ii]]->value_ptr(), prop.size(), &D_selected[cur_feat_local * prop.size()]);
             scores_selected[cur_feat_local] = _scores[inds[ii]];
-            phi_slected[cur_feat_local] = _phi[inds[ii]];
+            phi_selected[cur_feat_local] = _phi[inds[ii]];
             ++cur_feat_local;
         }
         ++ii;
     }
-    phi_slected.resize(cur_feat_local);
+
+    phi_selected.resize(cur_feat_local);
     scores_selected.resize(cur_feat_local);
-    D_temp.resize(cur_feat_local * prop.size());
+    D_selected.resize(cur_feat_local * prop.size());
 
-    std::vector<std::vector<double>> all_scores;
-    std::vector<std::vector<double>> all_D_temp;
-    std::vector<std::vector<node_ptr>> all_phi_selected;
+    if(_mpi_comm->size() > 1)
+    {
+        std::vector<std::vector<double>> all_scores;
+        std::vector<std::vector<double>> all_D;
+        std::vector<std::vector<node_ptr>> all_phi;
 
-    mpi::all_gather(*_mpi_comm, scores_selected, all_scores);
-    mpi::all_gather(*_mpi_comm, D_temp, all_scores);
-    mpi::all_gather(*_mpi_comm, phi_slected, all_phi_selected);
+        mpi::all_gather(*_mpi_comm, scores_selected, all_scores);
+        mpi::all_gather(*_mpi_comm, phi_selected, all_phi);
+        mpi::all_gather(*_mpi_comm, D_selected, all_D);
 
-    int iter_start = 0;
+        int iter_start = 0;
 
-    scores_selected = {};
-    D_temp = {};
-    phi_slected = {};
+        scores_selected = {};
+        D_selected = {};
+        phi_selected = {};
 
-    for(int sv = 0; sv < _mpi_comm->size(); ++sv)
-    {
-        scores_selected.resize(all_scores[sv].size() + scores_selected.size());
-        D_temp.resize(all_D_temp[sv].size() + D_temp.size());
-        std::copy_n(all_scores[sv].begin(), all_scores[sv].size(), scores_selected.end() - all_scores[sv].size());
-        std::copy_n(all_D_temp[sv].begin(), all_D_temp[sv].size(), D_temp.end() - all_D_temp[sv].size());
-        std::copy_n(all_phi_selected[sv].begin(), all_phi_selected[sv].size(), phi_slected.end() - all_phi_selected[sv].size());
+        for(int sv = 0; sv < _mpi_comm->size(); ++sv)
+        {
+            int cur_ind = scores_selected.size();
+            scores_selected.resize(all_scores[sv].size() + scores_selected.size());
+            phi_selected.resize(all_phi[sv].size() + phi_selected.size());
+            D_selected.resize(all_D[sv].size() * prop.size() + D_selected.size());
+
+            std::copy_n(all_scores[sv].begin(), all_scores[sv].size(), &scores_selected[cur_ind]);
+            std::copy_n(all_phi[sv].begin(), all_phi[sv].size(), &phi_selected[cur_ind]);
+            std::copy_n(all_D[sv].begin(), all_D[sv].size(), &D_selected[cur_ind * prop.size()]);
+
+        }
     }
 
     inds = util_funcs::argsort(scores_selected);
     cur_feat_local = 0;
+    ii = 0;
     while((cur_feat != _D.size() / prop.size()) && (ii < scores_selected.size()))
     {
         bool is_valid = true;
         for(int dd = 0; dd < cur_feat_local; ++dd)
         {
-            if(1.0 - std::abs(util_funcs::r(&_D[previous_size + dd * prop.size()], &D_temp[inds[ii] * prop.size()], prop.size())) < 1e-13)
+            if(1.0 - std::abs(util_funcs::r(&_D[previous_size + dd * prop.size()], &D_selected[inds[ii] * prop.size()], prop.size())) < 1e-13)
             {
                 is_valid = false;
                 break;
@@ -372,8 +375,8 @@ void FeatureSpace::sis(std::vector<double>& prop)
 
         if(is_valid)
         {
-            std::copy_n(&D_temp[inds[ii] * prop.size()], prop.size(), &_D[cur_feat * prop.size()]);
-            _phi_selected.push_back(phi_slected[inds[ii]]);
+            std::copy_n(&D_selected[inds[ii] * prop.size()], prop.size(), &_D[cur_feat * prop.size()]);
+            _phi_selected.push_back(phi_selected[inds[ii]]);
             _phi_selected.back()->set_value();
             ++cur_feat_local;
             ++cur_feat;