From d8e402832afa146987a6dbd8a36db129e95ad453 Mon Sep 17 00:00:00 2001
From: Thomas Purcell <purcell@fhi-berlin.mpg.de>
Date: Sun, 10 May 2020 16:57:10 +0200
Subject: [PATCH] Parallelized Code

Feature Creation, projection scores, and l0 normalization all parallelized
---
 src/Makefile.am                               |   3 +-
 src/descriptor_identifier/SISSORegressor.cpp  |  43 +++-
 src/descriptor_identifier/SISSORegressor.hpp  |   8 +-
 .../feature_space/FeatureSpace.cpp            | 206 ++++++++++--------
 .../feature_space/FeatureSpace.hpp            |  52 +++--
 src/feature_creation/node/FeatureNode.cpp     |   7 +-
 src/feature_creation/node/FeatureNode.hpp     |   4 +
 src/feature_creation/node/Node.cpp            |  12 +-
 src/feature_creation/node/Node.hpp            |   9 +-
 .../node/operator_nodes/OperatorNode.cpp      |  30 +--
 .../node/operator_nodes/OperatorNode.hpp      |  16 +-
 .../absolute_difference.cpp                   |  20 +-
 .../absolute_difference.hpp                   |   7 +
 .../allowed_operator_nodes/absolute_value.cpp |  24 +-
 .../allowed_operator_nodes/absolute_value.hpp |  10 +-
 .../allowed_operator_nodes/add.cpp            |  19 +-
 .../allowed_operator_nodes/add.hpp            |   8 +
 .../allowed_operator_nodes/cos.cpp            |  17 +-
 .../allowed_operator_nodes/cos.hpp            |   9 +
 .../allowed_operator_nodes/cube.cpp           |  21 +-
 .../allowed_operator_nodes/cube.hpp           |   9 +
 .../allowed_operator_nodes/cube_root.cpp      |  21 +-
 .../allowed_operator_nodes/cube_root.hpp      |   9 +
 .../allowed_operator_nodes/divide.cpp         |  21 +-
 .../allowed_operator_nodes/divide.hpp         |   8 +
 .../allowed_operator_nodes/exponential.cpp    |  17 +-
 .../allowed_operator_nodes/exponential.hpp    |   9 +
 .../allowed_operator_nodes/inverse.cpp        |  21 +-
 .../allowed_operator_nodes/inverse.hpp        |   8 +
 .../allowed_operator_nodes/log.cpp            |  17 +-
 .../allowed_operator_nodes/log.hpp            |   9 +
 .../allowed_operator_nodes/multiply.cpp       |  21 +-
 .../allowed_operator_nodes/multiply.hpp       |   9 +
 .../negative_exponential.cpp                  |  17 +-
 .../negative_exponential.hpp                  |   9 +
 .../allowed_operator_nodes/sin.cpp            |  17 +-
 .../allowed_operator_nodes/sin.hpp            |   9 +
 .../allowed_operator_nodes/sixth_power.cpp    |  21 +-
 .../allowed_operator_nodes/sixth_power.hpp    |   9 +
 .../allowed_operator_nodes/square.cpp         |  21 +-
 .../allowed_operator_nodes/square.hpp         |   9 +
 .../allowed_operator_nodes/square_root.cpp    |  21 +-
 .../allowed_operator_nodes/square_root.hpp    |   9 +
 .../allowed_operator_nodes/subtract.cpp       |  17 +-
 .../allowed_operator_nodes/subtract.hpp       |   8 +
 .../node/operator_nodes/functions.hpp         |  70 +++---
 src/inputs/InputParser.cpp                    |  22 +-
 src/inputs/InputParser.hpp                    |   2 +-
 src/main.cpp                                  |  33 +--
 src/mpi_interface/MPI_Interface.cpp           |   8 +-
 src/mpi_interface/MPI_Interface.hpp           |   2 +-
 src/utils/math_funcs.hpp                      |   5 +
 52 files changed, 734 insertions(+), 279 deletions(-)

diff --git a/src/Makefile.am b/src/Makefile.am
index 4e1dc086..5f650a54 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -4,7 +4,8 @@ AM_CXXFLAGS=-I$(top_srcdir)
 bin_PROGRAMS = $(top_builddir)/sisso_cpp
 
 __top_builddir__sisso_cpp_SOURCES = \
-feature_creation/node/value_storage/nodes_value_containers.cpp \
+    mpi_interface/MPI_Interface.cpp \
+    feature_creation/node/value_storage/nodes_value_containers.cpp \
     feature_creation/units/Unit.cpp \
     feature_creation/node/Node.cpp \
     feature_creation/node/FeatureNode.cpp \
diff --git a/src/descriptor_identifier/SISSORegressor.cpp b/src/descriptor_identifier/SISSORegressor.cpp
index c653dda8..da92ce2c 100644
--- a/src/descriptor_identifier/SISSORegressor.cpp
+++ b/src/descriptor_identifier/SISSORegressor.cpp
@@ -1,6 +1,8 @@
 #include <descriptor_identifier/SISSORegressor.hpp>
 
-SISSORegressor::SISSORegressor(FeatureSpace feat_space, std::vector<double> prop, int n_dim):
+SISSORegressor::SISSORegressor(std::shared_ptr<FeatureSpace> feat_space, std::vector<double> prop, int n_dim):
+    _feat_space(feat_space),
+    _mpi_comm(feat_space->mpi_comm()),
     _n_samp(prop.size()),
     _n_dim(n_dim),
     _lwork(-1),
@@ -12,8 +14,7 @@ SISSORegressor::SISSORegressor(FeatureSpace feat_space, std::vector<double> prop
     _work(nullptr),
     _s(new double[n_dim + 1]),
     _models(),
-    _prop(prop),
-    _feat_space(feat_space)
+    _prop(prop)
 {
 
     // Initialize a, b, ones, s, and _error arrays
@@ -35,7 +36,7 @@ SISSORegressor::SISSORegressor(FeatureSpace feat_space, std::vector<double> prop
 void  SISSORegressor::set_a(std::vector<int>& inds)
 {
     for(int ii = 0; ii < inds.size(); ++ii)
-        std::copy_n(_feat_space.D(inds[ii]), _n_samp, _a.get() + ii * _n_samp);
+        std::copy_n(_feat_space->D(inds[ii]), _n_samp, _a.get() + ii * _n_samp);
     std::copy_n(_ones.get(), _n_samp, _a.get() + inds.size() * _n_samp);
 }
 
@@ -79,14 +80,14 @@ void SISSORegressor::set_error(std::vector<int>& inds, double* coeffs)
 void SISSORegressor::fit()
 {
     std::vector<double> residual(_n_samp);
-    _feat_space.sis(_prop);
-    std::vector<node_ptr> min_node(1, _feat_space.phi_selected()[0]);
+    _feat_space->sis(_prop);
+    std::vector<node_ptr> min_node(1, _feat_space->phi_selected()[0]);
     _models.push_back(Model(_prop, min_node));
     _models[_models.size() - 1].copy_error(residual.data());
 
     for(int dd = 2; dd <= _n_dim; ++dd)
     {
-        _feat_space.sis(residual);
+        _feat_space->sis(residual);
         l0_norm(residual, dd);
         _models[_models.size() - 1].copy_error(residual.data());
     }
@@ -101,10 +102,21 @@ void SISSORegressor::l0_norm(std::vector<double>& prop, int n_dim)
 
     double min_error = util_funcs::norm(prop.data(), _n_samp);
 
-    int n = _feat_space.phi_selected().size();
+    int n = _feat_space->phi_selected().size();
 
+    int length = n;
+    for(int nn = n - 1; nn > n - n_dim; --nn)
+        length *= nn;
+    length /= n_dim;
+
+    std::array<int, 2> start_end = _mpi_comm->get_start_end_for_iterator(length, 0);
     std::vector<bool> v(n);
     std::fill(v.end() - n_dim, v.end(), true);
+    int rr = 0;
+    do{
+        ++rr;
+    }while(std::next_permutation(v.begin(), v.end()) && (rr < start_end[0]));
+
     do {
         int nn = 0;
         int ii = 0;
@@ -124,11 +136,22 @@ void SISSORegressor::l0_norm(std::vector<double>& prop, int n_dim)
             std::copy_n(inds.begin(), inds.size(), inds_min.begin());
             min_error = util_funcs::norm(_error.get(), _n_samp);
         }
-    } while (std::next_permutation(v.begin(), v.end()));
+
+        ++rr;
+    } while(std::next_permutation(v.begin(), v.end()) && (rr < start_end[1]));
+
+
+    std::vector<double> all_min_error;
+    std::vector<std::vector<int>> all_inds_min;
+
+    mpi::all_gather(*_mpi_comm, min_error, all_min_error);
+    mpi::all_gather(*_mpi_comm, inds_min, all_inds_min);
+
+    int argmin = std::min_element(all_min_error.begin(), all_min_error.end()) - all_min_error.begin();
 
     std::vector<node_ptr> min_nodes(n_dim);
     for(int ii = 0; ii < n_dim; ++ii)
-        min_nodes[ii] = _feat_space.phi_selected()[inds_min[ii]];
+        min_nodes[ii] = _feat_space->phi_selected()[all_inds_min[argmin][ii]];
 
     _models.push_back(Model(_prop, min_nodes));
 }
diff --git a/src/descriptor_identifier/SISSORegressor.hpp b/src/descriptor_identifier/SISSORegressor.hpp
index e1ed044c..43bc2eb6 100644
--- a/src/descriptor_identifier/SISSORegressor.hpp
+++ b/src/descriptor_identifier/SISSORegressor.hpp
@@ -7,6 +7,9 @@
 class SISSORegressor
 {
 protected:
+    std::shared_ptr<FeatureSpace> _feat_space;
+    std::shared_ptr<MPI_Interface> _mpi_comm;
+
     int _n_samp;
     int _n_dim;
     int _lwork;
@@ -21,10 +24,9 @@ protected:
 
     std::vector<Model> _models;
     std::vector<double> _prop;
-    FeatureSpace _feat_space;
 
 public:
-    SISSORegressor(FeatureSpace feat_space, std::vector<double> prop, int n_dim);
+    SISSORegressor(std::shared_ptr<FeatureSpace> feat_space, std::vector<double> prop, int n_dim);
 
     int get_opt_lwork(int n_dim);
     void least_squares(std::vector<int>& inds, double* coeffs);
@@ -35,7 +37,7 @@ public:
     void fit();
     void l0_norm(std::vector<double>& prop, int n_dim);
 
-    inline FeatureSpace& feat_space(){return _feat_space;}
+    inline std::shared_ptr<FeatureSpace> feat_space(){return _feat_space;}
     inline std::vector<double>& prop(){return _prop;}
     inline std::vector<Model>& models(){return _models;}
     inline int n_samp(){return _n_samp;}
diff --git a/src/feature_creation/feature_space/FeatureSpace.cpp b/src/feature_creation/feature_space/FeatureSpace.cpp
index 5abcf78d..b6e315f1 100644
--- a/src/feature_creation/feature_space/FeatureSpace.cpp
+++ b/src/feature_creation/feature_space/FeatureSpace.cpp
@@ -1,122 +1,128 @@
 #include <feature_creation/feature_space/FeatureSpace.hpp>
+BOOST_CLASS_EXPORT(FeatureNode)
+BOOST_CLASS_EXPORT(AddNode)
+BOOST_CLASS_EXPORT(SubNode)
+BOOST_CLASS_EXPORT(AbsDiffNode)
+BOOST_CLASS_EXPORT(MultNode)
+BOOST_CLASS_EXPORT(DivNode)
+BOOST_CLASS_EXPORT(SqNode)
+BOOST_CLASS_EXPORT(SqrtNode)
+BOOST_CLASS_EXPORT(CbNode)
+BOOST_CLASS_EXPORT(CbrtNode)
+BOOST_CLASS_EXPORT(SixPowNode)
+BOOST_CLASS_EXPORT(ExpNode)
+BOOST_CLASS_EXPORT(NegExpNode)
+BOOST_CLASS_EXPORT(LogNode)
+BOOST_CLASS_EXPORT(AbsNode)
+BOOST_CLASS_EXPORT(InvNode)
+BOOST_CLASS_EXPORT(SinNode)
+BOOST_CLASS_EXPORT(CosNode)
 
 FeatureSpace::FeatureSpace(
+    std::shared_ptr<MPI_Interface> mpi_comm,
     std::vector<node_ptr> phi_0,
     std::vector<std::string> allowed_ops,
     int max_phi,
     int n_sis_select,
     double max_abs_feat_val
 ):
-    max_phi_(max_phi),
-    n_sis_select_(n_sis_select),
-    n_samp_(phi_0[0]->n_samp()),
-    n_feat_(phi_0.size()),
-    max_abs_feat_val_(max_abs_feat_val),
-    start_gen_(1, 0),
-    scores_(phi_0.size(), 0.0),
-    D_(0, 0.0),
-    phi_0_(phi_0),
-    phi_(phi_0)
+    _mpi_comm(mpi_comm),
+    _max_phi(max_phi),
+    _n_sis_select(n_sis_select),
+    _n_samp(phi_0[0]->n_samp()),
+    _n_feat(phi_0.size()),
+    _max_abs_feat_val(max_abs_feat_val),
+    _start_gen(1, 0),
+    _scores(phi_0.size(), 0.0),
+    _D(0, 0.0),
+    _phi_0(phi_0),
+    _phi(phi_0)
 {
     for(auto & op : allowed_ops)
     {
         if((op.compare("add") == 0) || (op.compare("sub") == 0) || (op.compare("mult") == 0) || (op.compare("abs_diff") == 0))
-            com_bin_operators_.push_back(allowed_op_maps::binary_operator_map[op]);
+            _com_bin_operators.push_back(allowed_op_maps::binary_operator_map[op]);
         else if((op.compare("div") == 0))
-            bin_operators_.push_back(allowed_op_maps::binary_operator_map[op]);
+            _bin_operators.push_back(allowed_op_maps::binary_operator_map[op]);
         else
-            un_operators_.push_back(allowed_op_maps::unary_operator_map[op]);
+            _un_operators.push_back(allowed_op_maps::unary_operator_map[op]);
     }
     generate_feature_space();
-    scores_.resize(phi_.size());
-    scores_.reserve(phi_.size());
+    _scores.resize(_phi.size());
+    _scores.reserve(_phi.size());
 }
 
 FeatureSpace::FeatureSpace(
+    std::shared_ptr<MPI_Interface> mpi_comm,
     std::vector<node_ptr> phi_0,
     int max_phi,
     int n_sis_select,
     double max_abs_feat_val
 ):
-    max_phi_(max_phi),
-    n_sis_select_(n_sis_select),
-    n_samp_(phi_0[0]->n_samp()),
-    n_feat_(phi_0.size()),
-    max_abs_feat_val_(max_abs_feat_val),
-    start_gen_(1, 0),
-    scores_(phi_0.size(), 0.0),
-    D_(0, 0.0),
-    phi_0_(phi_0),
-    phi_(phi_0)
+    _mpi_comm(mpi_comm),
+    _max_phi(max_phi),
+    _n_sis_select(n_sis_select),
+    _n_samp(phi_0[0]->n_samp()),
+    _n_feat(phi_0.size()),
+    _max_abs_feat_val(max_abs_feat_val),
+    _start_gen(1, 0),
+    _scores(phi_0.size(), 0.0),
+    _D(0, 0.0),
+    _phi_0(phi_0),
+    _phi(phi_0)
 {
     for(auto& op : allowed_op_maps::unary_operator_map)
     {
-        un_operators_.push_back(op.second);
+        _un_operators.push_back(op.second);
     }
     for(auto& op : allowed_op_maps::binary_operator_map)
     {
         if((op.first.compare("div") == 0))
-            bin_operators_.push_back(op.second);
+            _bin_operators.push_back(op.second);
         else
-            com_bin_operators_.push_back(op.second);
+            _com_bin_operators.push_back(op.second);
     }
     generate_feature_space();
-    scores_.resize(phi_.size());
-    scores_.reserve(phi_.size());
+    _scores.resize(_phi.size());
+    _scores.reserve(_phi.size());
 }
 
 FeatureSpace::FeatureSpace(FeatureSpace &o) :
-    max_phi_(o.max_phi_),
-    n_sis_select_(o.n_sis_select_),
-    n_samp_(o.n_samp_),
-    n_feat_(o.n_feat_),
-    max_abs_feat_val_(o.max_abs_feat_val_),
-    start_gen_(o.start_gen_),
-    start_D_ind_(o.start_D_ind_),
-    prop_(o.prop_),
-    scores_(o.scores_),
-    D_(o.D_),
-    un_operators_(o.un_operators_),
-    bin_operators_(o.bin_operators_),
-    com_bin_operators_(o.com_bin_operators_),
-    phi_selected_(o.phi_selected_),
-    phi_(o.phi_),
-    phi_0_(o.phi_0_)
-{}
-
-FeatureSpace::FeatureSpace(FeatureSpace &&o) :
-    max_phi_(o.max_phi_),
-    n_sis_select_(o.n_sis_select_),
-    n_samp_(o.n_samp_),
-    n_feat_(o.n_feat_),
-    max_abs_feat_val_(o.max_abs_feat_val_),
-    start_gen_(o.start_gen_),
-    start_D_ind_(o.start_D_ind_),
-    prop_(o.prop_),
-    scores_(o.scores_),
-    D_(o.D_),
-    un_operators_(o.un_operators_),
-    bin_operators_(o.bin_operators_),
-    com_bin_operators_(o.com_bin_operators_),
-    phi_selected_(o.phi_selected_),
-    phi_(o.phi_),
-    phi_0_(o.phi_0_)
+    _max_phi(o._max_phi),
+    _n_sis_select(o._n_sis_select),
+    _n_samp(o._n_samp),
+    _n_feat(o._n_feat),
+    _max_abs_feat_val(o._max_abs_feat_val),
+    _start_gen(o._start_gen),
+    _start_ind(o._start_ind),
+    _prop(o._prop),
+    _scores(o._scores),
+    _D(o._D),
+    _un_operators(o._un_operators),
+    _bin_operators(o._bin_operators),
+    _com_bin_operators(o._com_bin_operators),
+    _phi_selected(o._phi_selected),
+    _phi(o._phi),
+    _phi_0(o._phi_0)
 {}
 
 void FeatureSpace::generate_feature_space()
 {
-    for(int nn = 1; nn <= max_phi_; ++nn)
+    for(int nn = 1; nn <= _max_phi; ++nn)
     {
         std::vector<node_ptr> next_phi;
-        n_feat_ = phi_.size();
-        int feat_ind = n_feat_;
-        for(int f1 = start_gen_[start_gen_.size()-1]; f1 < phi_.size(); ++f1)
+
+        _n_feat = _phi.size();
+        int feat_ind = _n_feat;
+
+        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]);
+        for(int f1 = start_end[0]; f1 < start_end[1]; ++f1)
         {
-            for(auto& op : un_operators_)
+            for(auto& op : _un_operators)
             {
                 try
                 {
-                    next_phi.push_back(op(phi_[f1], feat_ind));
+                    next_phi.push_back(op(_phi[f1], feat_ind));
                     ++feat_ind;
                 }
                 catch(const InvalidFeatureException& e)
@@ -125,13 +131,13 @@ void FeatureSpace::generate_feature_space()
                 }
             }
 
-            for(auto& op : com_bin_operators_)
+            for(auto& op : _com_bin_operators)
             {
                 for(int f2 = 0; f2 < f1; ++f2)
                 {
                     try
                     {
-                        next_phi.push_back(op(phi_[f1], phi_[f2], feat_ind));
+                        next_phi.push_back(op(_phi[f1], _phi[f2], feat_ind));
                         ++feat_ind;
                     }
                     catch(const InvalidFeatureException& e)
@@ -141,13 +147,13 @@ void FeatureSpace::generate_feature_space()
                 }
             }
 
-            for(auto& op : bin_operators_)
+            for(auto& op : _bin_operators)
             {
                 for(int f2 = 0; f2 < f1; ++f2)
                 {
                     try
                     {
-                        next_phi.push_back(op(phi_[f1], phi_[f2], feat_ind));
+                        next_phi.push_back(op(_phi[f1], _phi[f2], feat_ind));
                         ++feat_ind;
                     }
                     catch(const InvalidFeatureException& e)
@@ -156,7 +162,7 @@ void FeatureSpace::generate_feature_space()
                     }
                     try
                     {
-                        next_phi.push_back(op(phi_[f2], phi_[f1], feat_ind));
+                        next_phi.push_back(op(_phi[f2], _phi[f1], feat_ind));
                         ++feat_ind;
                     }
                     catch(const InvalidFeatureException& e)
@@ -166,40 +172,56 @@ void FeatureSpace::generate_feature_space()
                 }
             }
         }
+        _mpi_comm->barrier();
+        _start_gen.push_back(_phi.size());
 
-        start_gen_.push_back(phi_.size());
-
-        for(auto& feat : next_phi)
-            if(!feat->is_nan() && !feat->is_const())
-                phi_.push_back(feat);
+        std::vector<std::vector<node_ptr>> next_phi_gathered;
+        mpi::all_gather(*_mpi_comm, next_phi, next_phi_gathered);
 
+        for(auto& next_phi_vec : next_phi_gathered)
+            for(auto& feat : next_phi_vec)
+                _phi.push_back(feat);
     }
 }
 
 void FeatureSpace::project_r(std::vector<double>& prop)
 {
-    for(int ff = 0; ff < phi_.size(); ++ff)
-        scores_[ff] = std::abs(util_funcs::r(prop.data(), phi_[ff]->value_ptr(), prop.size()));
+    std::array<int, 2> start_end = _mpi_comm->get_start_end_for_iterator(_phi.size(), 0);
+
+    std::vector<double> scores(start_end[1] - start_end[0], 0.0);
+    for(int ff = start_end[0]; ff < start_end[1]; ++ff)
+        scores[ff - start_end[0]] = std::abs(util_funcs::r(prop.data(), _phi[ff]->value_ptr(), prop.size()));
+
+    std::vector<std::vector<double>> all_scores;
+    mpi::all_gather(*_mpi_comm, scores, all_scores);
+
+    int iter_start = 0;
+    for(auto& score_vec : all_scores)
+    {
+        for(int ii = 0; ii < score_vec.size(); ++ii)
+            _scores[ii + iter_start] = score_vec[ii];
+        iter_start += score_vec.size();
+    }
 }
 
 void FeatureSpace::sis(std::vector<double>& prop)
 {
-    int cur_feat = D_.size() / prop.size();
-    D_.resize(D_.size() + n_sis_select_ * prop.size());
-    D_.reserve(D_.size());
+    int cur_feat = _D.size() / prop.size();
+    _D.resize(_D.size() + _n_sis_select * prop.size());
+    _D.reserve(_D.size());
 
-    phi_selected_.reserve(phi_selected_.size());
+    _phi_selected.reserve(_phi_selected.size());
 
     project_r(prop);
-    std::vector<int> inds = util_funcs::argsort(scores_);
+    std::vector<int> inds = util_funcs::argsort(_scores);
 
     int ii = inds.size() - 1;
-    while((cur_feat != D_.size() / prop.size()) && (ii >= 0))
+    while((cur_feat != _D.size() / prop.size()) && (ii >= 0))
     {
         bool is_valid = true;
         for(int dd = 0; dd < cur_feat; ++dd)
         {
-            if(1.0 - std::abs(util_funcs::r(&D_[dd*prop.size()], phi_[inds[ii]]->value_ptr(), prop.size())) < 1e-12)
+            if(1.0 - std::abs(util_funcs::r(&_D[dd*prop.size()], _phi[inds[ii]]->value_ptr(), prop.size())) < 1e-12)
             {
                 is_valid = false;
                 break;
@@ -207,10 +229,10 @@ void FeatureSpace::sis(std::vector<double>& prop)
         }
         if(is_valid)
         {
-            std::copy_n(phi_[inds[ii]]->value_ptr(), prop.size(), &D_[cur_feat * prop.size()]);
+            std::copy_n(_phi[inds[ii]]->value_ptr(), prop.size(), &_D[cur_feat * prop.size()]);
             ++cur_feat;
 
-            phi_selected_.push_back(phi_[inds[ii]]);
+            _phi_selected.push_back(_phi[inds[ii]]);
         }
         --ii;
     }
diff --git a/src/feature_creation/feature_space/FeatureSpace.hpp b/src/feature_creation/feature_space/FeatureSpace.hpp
index 04957dae..c353f586 100644
--- a/src/feature_creation/feature_space/FeatureSpace.hpp
+++ b/src/feature_creation/feature_space/FeatureSpace.hpp
@@ -1,37 +1,44 @@
 #ifndef FEATURE_SPACE
 #define FEATURE_SPACE
 
+#include <mpi_interface/MPI_Interface.hpp>
 #include <feature_creation/node/FeatureNode.hpp>
 #include <feature_creation/node/operator_nodes/allowed_ops.hpp>
 
+#include <boost/serialization/shared_ptr.hpp>
+
 #include <iostream>
 
+// namespace mpi = boost::mpi;
+
 class FeatureSpace
 {
-    int max_phi_;
-    int n_sis_select_;
-    int n_samp_;
-    int n_feat_;
+    std::shared_ptr<MPI_Interface> _mpi_comm;
+    int _max_phi;
+    int _n_sis_select;
+    int _n_samp;
+    int _n_feat;
 
-    double max_abs_feat_val_;
+    double _max_abs_feat_val;
 
-    std::vector<int> start_gen_;
-    std::vector<int> start_D_ind_;
+    std::vector<int> _start_gen;
+    std::vector<int> _start_ind;
 
-    std::vector<double> prop_;
-    std::vector<double> scores_;
-    std::vector<double> D_;
+    std::vector<double> _prop;
+    std::vector<double> _scores;
+    std::vector<double> _D;
 
-    std::vector<un_op_node_gen> un_operators_;
-    std::vector<bin_op_node_gen> bin_operators_;
-    std::vector<bin_op_node_gen> com_bin_operators_;
+    std::vector<un_op_node_gen> _un_operators;
+    std::vector<bin_op_node_gen> _bin_operators;
+    std::vector<bin_op_node_gen> _com_bin_operators;
 
-    std::vector<node_ptr> phi_selected_;
-    std::vector<node_ptr> phi_;
-    std::vector<node_ptr> phi_0_;
+    std::vector<node_ptr> _phi_selected;
+    std::vector<node_ptr> _phi;
+    std::vector<node_ptr> _phi_0;
 
 public:
     FeatureSpace(
+        std::shared_ptr<MPI_Interface> mpi_comm,
         std::vector<node_ptr> phi_0,
         std::vector<std::string> allowed_ops,
         int max_phi=1,
@@ -40,6 +47,7 @@ public:
     );
 
     FeatureSpace(
+        std::shared_ptr<MPI_Interface> mpi_comm,
         std::vector<node_ptr> phi_0,
         int max_phi=1,
         int n_sis_select=1,
@@ -47,15 +55,15 @@ public:
     );
 
     FeatureSpace(FeatureSpace &o);
-    FeatureSpace(FeatureSpace &&o);
 
     void generate_feature_space();
 
-    inline std::vector<node_ptr> phi_selected(){return phi_selected_;};
-    inline std::vector<node_ptr> phi(){return phi_;};
-    inline std::vector<node_ptr> phi0(){return phi_0_;};
-    inline std::vector<double>& scores(){return scores_;};
-    inline double* D(int ind){return &D_[ind * n_samp_];}
+    inline std::vector<node_ptr> phi_selected(){return _phi_selected;};
+    inline std::vector<node_ptr> phi(){return _phi;};
+    inline std::vector<node_ptr> phi0(){return _phi_0;};
+    inline std::vector<double>& scores(){return _scores;};
+    inline std::shared_ptr<MPI_Interface> mpi_comm(){return _mpi_comm;}
+    inline double* D(int ind){return &_D[ind * _n_samp];}
 
     void project_r(std::vector<double>& prop);
 
diff --git a/src/feature_creation/node/FeatureNode.cpp b/src/feature_creation/node/FeatureNode.cpp
index cd173cfc..b301e74e 100644
--- a/src/feature_creation/node/FeatureNode.cpp
+++ b/src/feature_creation/node/FeatureNode.cpp
@@ -1,5 +1,8 @@
 #include <feature_creation/node/FeatureNode.hpp>
 
+FeatureNode::FeatureNode()
+{}
+
 FeatureNode::FeatureNode(int feat_ind, std::string expr, std::vector<double> value, Unit unit) :
     Node(feat_ind, expr, value, unit)
 {}
@@ -7,5 +10,7 @@ FeatureNode::FeatureNode(int feat_ind, std::string expr, std::vector<double> val
 FeatureNode::FeatureNode(const FeatureNode &o) :
     Node(o)
 {
-    std::copy_n(o._value_ptr.get(), o._n_samp, _value_ptr.get());
+    std::copy_n(o._value_ptr.data(), o._n_samp, _value_ptr.data());
 }
+
+// BOOST_CLASS_EXPORT(FeatureNode)
diff --git a/src/feature_creation/node/FeatureNode.hpp b/src/feature_creation/node/FeatureNode.hpp
index 88b96645..85540ae1 100644
--- a/src/feature_creation/node/FeatureNode.hpp
+++ b/src/feature_creation/node/FeatureNode.hpp
@@ -5,6 +5,7 @@
 #include <utils/enum.hpp>
 #include <memory>
 
+#include <boost/serialization/export.hpp>
 #include <boost/serialization/base_object.hpp>
 
 #include <feature_creation/node/Node.hpp>
@@ -14,7 +15,9 @@ typedef std::function<double(double, double)> binary_op_func;
 
 class FeatureNode: public Node
 {
+    friend class boost::serialization::access;
 public:
+    FeatureNode();
     FeatureNode(int feat_ind, std::string expr, std::vector<double> value, Unit unit);
     FeatureNode(const FeatureNode &o);
 
@@ -22,6 +25,7 @@ public:
 
     inline Unit unit(){return _unit;}
 
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
diff --git a/src/feature_creation/node/Node.cpp b/src/feature_creation/node/Node.cpp
index 1abe2fea..c70d44ff 100644
--- a/src/feature_creation/node/Node.cpp
+++ b/src/feature_creation/node/Node.cpp
@@ -1,5 +1,8 @@
 #include <feature_creation/node/Node.hpp>
 
+Node::Node()
+{}
+
 Node::Node(int feat_ind, std::string expr, std::vector<double> value, Unit unit, NODE_TYPE tag) :
     _n_samp(value.size()),
     _tag(tag),
@@ -7,7 +10,7 @@ Node::Node(int feat_ind, std::string expr, std::vector<double> value, Unit unit,
     _expr(expr),
     _unit(unit),
     // _value_ptr(nullptr),
-    _value_ptr(new double[_n_samp])
+    _value_ptr(_n_samp)
 {
     // _value_ptr = node_value_arrs::access_value_arr(feat_ind);
     // _value_ptr = _value;
@@ -20,7 +23,7 @@ Node::Node(int feat_ind, std::string expr, int n_samp, NODE_TYPE tag) :
     _feat_ind(feat_ind),
     _expr(expr),
     _unit(),
-    _value_ptr(new double[_n_samp])
+    _value_ptr(_n_samp)
 {}
 
 Node::Node(const Node &o) :
@@ -29,8 +32,9 @@ Node::Node(const Node &o) :
     _feat_ind(o._feat_ind),
     _expr(o._expr),
     _unit(o._unit),
-    _value_ptr(new double[o._n_samp])
+    _value_ptr(o._n_samp)
 {
-    std::copy_n(o._value_ptr.get(), o._n_samp, value_ptr());
+    std::copy_n(o._value_ptr.data(), o._n_samp, value_ptr());
 }
+BOOST_SERIALIZATION_ASSUME_ABSTRACT(Node)
 
diff --git a/src/feature_creation/node/Node.hpp b/src/feature_creation/node/Node.hpp
index 41e0811c..338ba4d2 100644
--- a/src/feature_creation/node/Node.hpp
+++ b/src/feature_creation/node/Node.hpp
@@ -9,9 +9,10 @@
 #include <cmath>
 #include <memory>
 
+#include <boost/serialization/assume_abstract.hpp>
 #include <boost/serialization/serialization.hpp>
-#include <boost/serialization/unique_ptr.hpp>
 #include <boost/serialization/string.hpp>
+#include <boost/serialization/unique_ptr.hpp>
 
 typedef std::function<double(double)> unary_op_func;
 typedef std::function<double(double, double)> binary_op_func;
@@ -30,9 +31,10 @@ protected:
     Unit _unit;
 
     // double* _value_ptr;
-    std::unique_ptr<double[]> _value_ptr;
+    std::vector<double> _value_ptr;
     // std::array<double, _n_samp> _value;
 public:
+    Node();
     Node(int feat_ind, std::string expr, std::vector<double> value, Unit unit, NODE_TYPE tag = NODE_TYPE::FEAT);
     Node(int feat_ind, std::string expr, int n_samp, NODE_TYPE tag = NODE_TYPE::FEAT);
     Node(const Node &o);
@@ -54,7 +56,8 @@ public:
     virtual std::string expr() = 0;
     virtual Unit unit() = 0;
 
-    inline double* value_ptr(){return _value_ptr.get();}
+    inline double* value_ptr(){return _value_ptr.data();}
+
 
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
diff --git a/src/feature_creation/node/operator_nodes/OperatorNode.cpp b/src/feature_creation/node/operator_nodes/OperatorNode.cpp
index 0cfdf73a..c895fb58 100644
--- a/src/feature_creation/node/operator_nodes/OperatorNode.cpp
+++ b/src/feature_creation/node/operator_nodes/OperatorNode.cpp
@@ -1,28 +1,17 @@
 #include <feature_creation/node/operator_nodes/OperatorNode.hpp>
 
-OperatorNode::OperatorNode(op_func func_in, std::vector<node_ptr> feats, int feat_ind, std::string expr, NODE_TYPE tag) :
-    Node(feat_ind, expr, feats[0]->n_samp(), tag),
-    _func(func_in),
-    _feats(feats),
-    _feat_values(feats.size())
-{
-    _feat_values.reserve(_feats.size());
-    for(int ff = 0; ff < _feats.size(); ++ff)
-        _feat_values[ff] = _feats[ff]->value_ptr();
+OperatorNode::OperatorNode()
+{}
 
-    set_value();
-}
+OperatorNode::OperatorNode(std::vector<node_ptr> feats, int feat_ind, std::string expr, NODE_TYPE tag) :
+    Node(feat_ind, expr, feats[0]->n_samp(), tag),
+    _feats(feats)
+{}
 
 OperatorNode::OperatorNode(const OperatorNode &o) :
     Node(o),
-    _func(o._func),
-    _feats(o._feats),
-    _feat_values(o._feats.size())
-{
-    _feat_values.reserve(_feats.size());
-    for(int ff = 0; ff < _feats.size(); ++ff)
-        _feat_values[ff] = _feats[ff]->value_ptr();
-}
+    _feats(o._feats)
+{}
 
 std::string OperatorNode::expr()
 {
@@ -33,3 +22,6 @@ std::string OperatorNode::expr()
     return expression.substr(0, expression.size() - 2) + ")";
 }
 
+BOOST_SERIALIZATION_ASSUME_ABSTRACT(OperatorNode)
+// BOOST_CLASS_EXPORTOperatorNode)
+// BOOST_CLASS_EXPORT_GUID(OperatorNode, "OperatorNode")
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/OperatorNode.hpp b/src/feature_creation/node/operator_nodes/OperatorNode.hpp
index 3accfab3..3ef42a74 100644
--- a/src/feature_creation/node/operator_nodes/OperatorNode.hpp
+++ b/src/feature_creation/node/operator_nodes/OperatorNode.hpp
@@ -5,33 +5,35 @@
 #include <feature_creation/node/operator_nodes/functions.hpp>
 
 #include <boost/serialization/base_object.hpp>
+#include <boost/serialization/export.hpp>
+#include <boost/serialization/vector.hpp>
+#include <boost/serialization/shared_ptr.hpp>
 
 class OperatorNode: public Node
 {
+    friend class boost::serialization::access;
 protected:
-    op_func _func;
     std::vector<node_ptr> _feats;
-    std::vector<double*> _feat_values;
 
 public:
-    OperatorNode(op_func func, std::vector<node_ptr> feats, int feat_ind, std::string expr, NODE_TYPE tag = NODE_TYPE::FXN);
+    OperatorNode();
+    OperatorNode(std::vector<node_ptr> feats, int feat_ind, std::string expr, NODE_TYPE tag = NODE_TYPE::FXN);
     OperatorNode(const OperatorNode &o);
 
-    inline op_func func(){return _func;};
 
     virtual std::string expr() = 0;
 
     virtual Unit unit() = 0;
 
-    inline void set_value() {_func(_n_samp, _feat_values, value_ptr());}
+    virtual void set_value() =0;
+
 
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<Node>(*this);
-        ar & _func;
         ar & _feats;
-        ar & _feat_values;
     }
 };
 
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp
index fe319d94..507620d1 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp
@@ -1,15 +1,31 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp>
 
+// BOOST_CLASS_EXPORT(AbsDiffNode)
+
+AbsDiffNode::AbsDiffNode()
+{}
+
 AbsDiffNode::AbsDiffNode(std::vector<node_ptr> feats, int feat_ind) :
-    OperatorNode(allowed_op_funcs::abs_diff, feats, feat_ind, "abs_diff", NODE_TYPE::ABS_DIFF)
+    OperatorNode(feats, feat_ind, "abs_diff", NODE_TYPE::ABS_DIFF)
 {
     if(feats[0]->unit() != feats[1]->unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 AbsDiffNode::AbsDiffNode(node_ptr feat_1, node_ptr feat_2, int feat_ind) :
-    OperatorNode(allowed_op_funcs::abs_diff, {feat_1, feat_2}, feat_ind, "abs_diff", NODE_TYPE::ABS_DIFF)
+    OperatorNode({feat_1, feat_2}, feat_ind, "abs_diff", NODE_TYPE::ABS_DIFF)
 {
     if(feat_1->unit() != feat_2->unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORTAbsDiffNode)
+// BOOST_CLASS_EXPORT_GUID(AbsDiffNode, "AbsDiffNode")
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp
index cdd7eddc..899132a7 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp
@@ -7,6 +7,8 @@ class AbsDiffNode: public OperatorNode
 {
 
 public:
+    AbsDiffNode();
+
     AbsDiffNode(std::vector<node_ptr> feats, int feat_ind);
 
     AbsDiffNode(node_ptr feat_1, node_ptr feat_2, int feat_ind);
@@ -15,6 +17,11 @@ public:
 
     inline std::string expr(){return "|" + _feats[0]->expr() + " - (" + _feats[1]->expr() + ")|";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::abs_diff(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp
index 5876de6f..e2bfc7c5 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp
@@ -1,9 +1,25 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp>
 
-AbsNode::AbsNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::abs, feats, feat_ind, "abs", NODE_TYPE::ABS)
+// BOOST_CLASS_EXPORT(AbsNode)
+
+AbsNode::AbsNode()
 {}
 
+AbsNode::AbsNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "abs", NODE_TYPE::ABS)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 AbsNode::AbsNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::abs, {feat}, feat_ind, "abs", NODE_TYPE::ABS)
-{}
+    OperatorNode({feat}, feat_ind, "abs", NODE_TYPE::ABS)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORTAbsNode)
+// BOOST_CLASS_EXPORT_GUID(AbsNode, "AbsNode")
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp
index 43de4aca..926eb826 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp
@@ -7,6 +7,8 @@ class AbsNode: public OperatorNode
 {
 
 public:
+    AbsNode();
+
     AbsNode(std::vector<node_ptr> feats, int feat_ind);
 
     AbsNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,16 @@ public:
 
     inline std::string expr(){return "|" + _feats[0]->expr() + "|";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::abs(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
         ar & boost::serialization::base_object<OperatorNode>(*this);
-    }};
+    }
+};
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp
index 89a744b8..607d3cc1 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp
@@ -1,15 +1,30 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp>
 
+AddNode::AddNode()
+{}
+
 AddNode::AddNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::add, feats, feat_ind, "add", NODE_TYPE::ADD)
+    OperatorNode(feats, feat_ind, "add", NODE_TYPE::ADD)
 {
     if(feats[0]->unit() != feats[1]->unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+
 }
 
 AddNode::AddNode(node_ptr feat_1, node_ptr feat_2, int feat_ind):
-    OperatorNode(allowed_op_funcs::add, {feat_1, feat_2}, feat_ind, "add", NODE_TYPE::ADD)
+    OperatorNode({feat_1, feat_2}, feat_ind, "add", NODE_TYPE::ADD)
 {
     if(feat_1->unit() != feat_2->unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+
 }
+
+// BOOST_CLASS_EXPORT(AddNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp
index 4ba33743..23e1c513 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp
@@ -7,6 +7,8 @@ class AddNode: public OperatorNode
 {
 
 public:
+    AddNode();
+
     AddNode(std::vector<node_ptr> feats, int feat_ind);
 
     AddNode(node_ptr feat_1, node_ptr feat_2, int feat_ind);
@@ -15,9 +17,15 @@ public:
 
     inline std::string expr(){return "(" + _feats[0]->expr() + " + " + _feats[1]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::add(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp
index 9c213f9e..38a7308d 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp
@@ -1,15 +1,28 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp>
 
+CosNode::CosNode()
+{}
+
 CosNode::CosNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::cos, feats, feat_ind, "cos", NODE_TYPE::COS)
+    OperatorNode(feats, feat_ind, "cos", NODE_TYPE::COS)
 {
     if(feats[0]->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 CosNode::CosNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::cos, {feat}, feat_ind, "cos", NODE_TYPE::COS)
+    OperatorNode({feat}, feat_ind, "cos", NODE_TYPE::COS)
 {
     if(feat->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORT(CosNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp
index d1081b7a..27b215fb 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp
@@ -7,6 +7,8 @@ class CosNode: public OperatorNode
 {
 
 public:
+    CosNode();
+
     CosNode(std::vector<node_ptr> feats, int feat_ind);
 
     CosNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "cos(" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::cos(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp
index 62414d44..1d1d6e38 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp>
 
-CbNode::CbNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::cb, feats, feat_ind, "cb", NODE_TYPE::CB)
+CbNode::CbNode()
 {}
 
+CbNode::CbNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "cb", NODE_TYPE::CB)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 CbNode::CbNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::cb, {feat}, feat_ind, "cb", NODE_TYPE::CB)
-{}
+    OperatorNode({feat}, feat_ind, "cb", NODE_TYPE::CB)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(CbNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp
index 61def2f4..33a0d048 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp
@@ -7,6 +7,8 @@ class CbNode: public OperatorNode
 {
 
 public:
+    CbNode();
+
     CbNode(std::vector<node_ptr> feats, int feat_ind);
 
     CbNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "(" + _feats[0]->expr() + ")^3";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::cb(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp
index 4bc306aa..63191aba 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp>
 
-CbrtNode::CbrtNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::cbrt, feats, feat_ind, "cbrt", NODE_TYPE::CBRT)
+CbrtNode::CbrtNode()
 {}
 
+CbrtNode::CbrtNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "cbrt", NODE_TYPE::CBRT)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 CbrtNode::CbrtNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::cbrt, {feat}, feat_ind, "cbrt", NODE_TYPE::CBRT)
-{}
+    OperatorNode({feat}, feat_ind, "cbrt", NODE_TYPE::CBRT)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(CbrtNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp
index 93e353b0..6a729054 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp
@@ -7,6 +7,8 @@ class CbrtNode: public OperatorNode
 {
 
 public:
+    CbrtNode();
+
     CbrtNode(std::vector<node_ptr> feats, int feat_ind);
 
     CbrtNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "cbrt(" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::cbrt(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp
index 8f0dbd5f..a2e3b698 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp>
 
-DivNode::DivNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::div, feats, feat_ind, "div", NODE_TYPE::DIV)
+DivNode::DivNode()
 {}
 
+DivNode::DivNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "div", NODE_TYPE::DIV)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 DivNode::DivNode(node_ptr feat_1, node_ptr feat_2, int feat_ind):
-    OperatorNode(allowed_op_funcs::div, {feat_1, feat_2}, feat_ind, "div", NODE_TYPE::DIV)
-{}
+    OperatorNode({feat_1, feat_2}, feat_ind, "div", NODE_TYPE::DIV)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(DivNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp
index f4e370cf..bb2cfa8b 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp
@@ -7,6 +7,8 @@ class DivNode: public OperatorNode
 {
 
 public:
+    DivNode();
+
     DivNode(std::vector<node_ptr> feats, int feat_ind);
 
     DivNode(node_ptr feat_1, node_ptr feat_2, int feat_ind);
@@ -15,9 +17,15 @@ public:
 
     inline std::string expr(){return "[(" + _feats[0]->expr() + ") / (" + _feats[1]->expr() + ")]";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::div(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp
index ed1e05d4..55fa4acc 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp
@@ -1,15 +1,28 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp>
 
+ExpNode::ExpNode()
+{}
+
 ExpNode::ExpNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::exp, feats, feat_ind, "exp", NODE_TYPE::EXP)
+    OperatorNode(feats, feat_ind, "exp", NODE_TYPE::EXP)
 {
     if(feats[0]->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 ExpNode::ExpNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::exp, {feat}, feat_ind, "exp", NODE_TYPE::EXP)
+    OperatorNode({feat}, feat_ind, "exp", NODE_TYPE::EXP)
 {
     if(feat->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORT(ExpNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp
index ebe5df5e..fbf6d3e3 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp
@@ -7,6 +7,8 @@ class ExpNode: public OperatorNode
 {
 
 public:
+    ExpNode();
+
     ExpNode(std::vector<node_ptr> feats, int feat_ind);
 
     ExpNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "exp(" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::exp(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp
index fa893191..60869692 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp>
 
-InvNode::InvNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::inv, feats, feat_ind, "inv", NODE_TYPE::INV)
+InvNode::InvNode()
 {}
 
+InvNode::InvNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "inv", NODE_TYPE::INV)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 InvNode::InvNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::inv, {feat}, feat_ind, "inv", NODE_TYPE::INV)
-{}
+    OperatorNode({feat}, feat_ind, "inv", NODE_TYPE::INV)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(InvNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp
index 37c7e54c..cdcbfc3d 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp
@@ -7,6 +7,8 @@ class InvNode: public OperatorNode
 {
 
 public:
+    InvNode();
+
     InvNode(std::vector<node_ptr> feats, int feat_ind);
 
     InvNode(node_ptr feat, int feat_ind);
@@ -15,9 +17,15 @@ public:
 
     inline std::string expr(){return "1.0 / (" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::inv(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp
index ebd60237..05ac20ac 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp
@@ -1,15 +1,28 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp>
 
+LogNode::LogNode()
+{}
+
 LogNode::LogNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::log, feats, feat_ind, "log", NODE_TYPE::LOG)
+    OperatorNode(feats, feat_ind, "log", NODE_TYPE::LOG)
 {
     if(feats[0]->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 LogNode::LogNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::log, {feat}, feat_ind, "log", NODE_TYPE::LOG)
+    OperatorNode({feat}, feat_ind, "log", NODE_TYPE::LOG)
 {
     if(feat->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORT(LogNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp
index ebf95434..a922505e 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp
@@ -7,6 +7,8 @@ class LogNode: public OperatorNode
 {
 
 public:
+    LogNode();
+
     LogNode(std::vector<node_ptr> feats, int feat_ind);
 
     LogNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "log(" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::log(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp
index 134ff85f..7fa29180 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp>
 
-MultNode::MultNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::mult, feats, feat_ind, "mult", NODE_TYPE::MULT)
+MultNode::MultNode()
 {}
 
+MultNode::MultNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "mult", NODE_TYPE::MULT)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 MultNode::MultNode(node_ptr feat_1, node_ptr feat_2, int feat_ind):
-    OperatorNode(allowed_op_funcs::mult, {feat_1, feat_2}, feat_ind, "mult", NODE_TYPE::MULT)
-{}
+    OperatorNode({feat_1, feat_2}, feat_ind, "mult", NODE_TYPE::MULT)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(MultNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp
index 0972c758..35bf6dd7 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp
@@ -7,6 +7,8 @@ class MultNode: public OperatorNode
 {
 
 public:
+    MultNode();
+
     MultNode(std::vector<node_ptr> feats, int feat_ind);
 
     MultNode(node_ptr feat_1, node_ptr feat_2, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "[(" + _feats[0]->expr() + ") * (" + _feats[1]->expr() + ")]";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::mult(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp
index 20460446..c7010fc2 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp
@@ -1,15 +1,28 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp>
 
+NegExpNode::NegExpNode()
+{}
+
 NegExpNode::NegExpNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::neg_exp, feats, feat_ind, "neg_exp", NODE_TYPE::NEG_EXP)
+    OperatorNode(feats, feat_ind, "neg_exp", NODE_TYPE::NEG_EXP)
 {
     if(feats[0]->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 NegExpNode::NegExpNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::neg_exp, {feat}, feat_ind, "neg_exp", NODE_TYPE::NEG_EXP)
+    OperatorNode({feat}, feat_ind, "neg_exp", NODE_TYPE::NEG_EXP)
 {
     if(feat->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORT(NegExpNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp
index ee764879..633d4c25 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp
@@ -7,6 +7,8 @@ class NegExpNode: public OperatorNode
 {
 
 public:
+    NegExpNode();
+
     NegExpNode(std::vector<node_ptr> feats, int feat_ind);
 
     NegExpNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "exp[-1.0*(" + _feats[0]->expr() + ")]";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::neg_exp(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp
index 1e167ef1..acf10d12 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp
@@ -1,15 +1,28 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp>
 
+SinNode::SinNode()
+{}
+
 SinNode::SinNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::sin, feats, feat_ind, "sin", NODE_TYPE::SIN)
+    OperatorNode(feats, feat_ind, "sin", NODE_TYPE::SIN)
 {
     if(feats[0]->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 SinNode::SinNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::sin, {feat}, feat_ind, "sin", NODE_TYPE::SIN)
+    OperatorNode({feat}, feat_ind, "sin", NODE_TYPE::SIN)
 {
     if(feat->unit() != Unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORT(SinNode)
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp
index 8d678304..3ee7effc 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp
@@ -7,6 +7,8 @@ class SinNode: public OperatorNode
 {
 
 public:
+    SinNode();
+
     SinNode(std::vector<node_ptr> feats, int feat_ind);
 
     SinNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "sin(" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::sin(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp
index 7069acea..c8f10003 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp>
 
-SixPowNode::SixPowNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::sixth_pow, feats, feat_ind, "sixth_pow", NODE_TYPE::SIX_POW)
+SixPowNode::SixPowNode()
 {}
 
+SixPowNode::SixPowNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "sixth_pow", NODE_TYPE::SIX_POW)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 SixPowNode::SixPowNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::sixth_pow, {feat}, feat_ind, "sixth_pow", NODE_TYPE::SIX_POW)
-{}
+    OperatorNode({feat}, feat_ind, "sixth_pow", NODE_TYPE::SIX_POW)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(SixPowNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp
index 0ec966c9..3d3917ed 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp
@@ -7,6 +7,8 @@ class SixPowNode: public OperatorNode
 {
 
 public:
+    SixPowNode();
+
     SixPowNode(std::vector<node_ptr> feats, int feat_ind);
 
     SixPowNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "(" + _feats[0]->expr() + ")^6";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::sixth_pow(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp
index 47332683..d3d014e4 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp>
 
-SqNode::SqNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::sq, feats, feat_ind, "sq", NODE_TYPE::SQ)
+SqNode::SqNode()
 {}
 
+SqNode::SqNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "sq", NODE_TYPE::SQ)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 SqNode::SqNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::sq, {feat}, feat_ind, "sq", NODE_TYPE::SQ)
-{}
+    OperatorNode({feat}, feat_ind, "sq", NODE_TYPE::SQ)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(SqNode)
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp
index 78f3676b..3599f849 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp
@@ -7,6 +7,8 @@ class SqNode: public OperatorNode
 {
 
 public:
+    SqNode();
+
     SqNode(std::vector<node_ptr> feats, int feat_ind);
 
     SqNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "(" + _feats[0]->expr() + ")^2";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::sq(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp
index cbcf0c71..76ba3242 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp
@@ -1,9 +1,22 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp>
 
-SqrtNode::SqrtNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::sqrt, feats, feat_ind, "sqrt", NODE_TYPE::SQRT)
+SqrtNode::SqrtNode()
 {}
 
+SqrtNode::SqrtNode(std::vector<node_ptr> feats, int feat_ind):
+    OperatorNode(feats, feat_ind, "sqrt", NODE_TYPE::SQRT)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
 SqrtNode::SqrtNode(node_ptr feat, int feat_ind):
-    OperatorNode(allowed_op_funcs::sqrt, {feat}, feat_ind, "sqrt", NODE_TYPE::SQRT)
-{}
+    OperatorNode({feat}, feat_ind, "sqrt", NODE_TYPE::SQRT)
+{
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
+}
+
+// BOOST_CLASS_EXPORT(SqrtNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp
index fe8d6ae2..27010af5 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp
@@ -7,6 +7,8 @@ class SqrtNode: public OperatorNode
 {
 
 public:
+    SqrtNode();
+
     SqrtNode(std::vector<node_ptr> feats, int feat_ind);
 
     SqrtNode(node_ptr feat, int feat_ind);
@@ -15,10 +17,17 @@ public:
 
     inline std::string expr(){return "sqrt(" + _feats[0]->expr() + ")";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::sqrt(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
+
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp
index de2f5e8d..3246f0ef 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp
@@ -1,15 +1,28 @@
 #include <feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp>
 
+SubNode::SubNode()
+{}
+
 SubNode::SubNode(std::vector<node_ptr> feats, int feat_ind):
-    OperatorNode(allowed_op_funcs::sub, feats, feat_ind, "sub", NODE_TYPE::SUB)
+    OperatorNode(feats, feat_ind, "sub", NODE_TYPE::SUB)
 {
     if(feats[0]->unit() != feats[1]->unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
 
 SubNode::SubNode(node_ptr feat_1, node_ptr feat_2, int feat_ind):
-    OperatorNode(allowed_op_funcs::sub, {feat_1, feat_2}, feat_ind, "sub", NODE_TYPE::SUB)
+    OperatorNode({feat_1, feat_2}, feat_ind, "sub", NODE_TYPE::SUB)
 {
     if(feat_1->unit() != feat_2->unit())
         throw InvalidFeatureException();
+
+    set_value();
+    if(is_nan() || is_const())
+        throw InvalidFeatureException();
 }
+
+// BOOST_CLASS_EXPORT(SubNode)
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp
index 54da5509..eb621a5d 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp
@@ -7,6 +7,8 @@ class SubNode: public OperatorNode
 {
 
 public:
+    SubNode();
+
     SubNode(std::vector<node_ptr> feats, int feat_ind);
 
     SubNode(node_ptr feat_1, node_ptr feat_2, int feat_ind);
@@ -15,9 +17,15 @@ public:
 
     inline std::string expr(){return "[(" + _feats[0]->expr() + ") - (" + _feats[1]->expr() + ")]";}
 
+    inline void set_value()
+    {
+        allowed_op_funcs::sub(_n_samp, _feats, _value_ptr.data());
+    }
+
     template <typename Archive>
     void serialize(Archive& ar, const unsigned int version)
     {
+        // ar.template register_type<OperatorNode>();
         ar & boost::serialization::base_object<OperatorNode>(*this);
     }
 };
diff --git a/src/feature_creation/node/operator_nodes/functions.hpp b/src/feature_creation/node/operator_nodes/functions.hpp
index df037093..6ddd9070 100644
--- a/src/feature_creation/node/operator_nodes/functions.hpp
+++ b/src/feature_creation/node/operator_nodes/functions.hpp
@@ -4,93 +4,93 @@
 #include <algorithm>
 #include <functional>
 
-typedef std::function<void(int, std::vector<double*>&, double*)> op_func;
+typedef std::function<void(int, std::vector<node_ptr>&, double*)> op_func;
 
 namespace allowed_op_funcs
 {
-    inline void add(int size, std::vector<double*>& inputs, double* out)
+    inline void add(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, inputs[1], out, std::plus<double>());
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, inputs[1]->value_ptr(), out, std::plus<double>());
     }
 
-    inline void sub(int size, std::vector<double*>& inputs, double* out)
+    inline void sub(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, inputs[1], out, std::minus<double>());
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, inputs[1]->value_ptr(), out, std::minus<double>());
     }
 
-    inline void abs_diff(int size, std::vector<double*>& inputs, double* out)
+    inline void abs_diff(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, inputs[1], out, [](double in_0, double in_1){return std::abs(in_0 - in_1);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, inputs[1]->value_ptr(), out, [](double in_0, double in_1){return std::abs(in_0 - in_1);});
     }
 
-    inline void mult(int size, std::vector<double*>& inputs, double* out)
+    inline void mult(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, inputs[1], out, std::multiplies<double>());
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, inputs[1]->value_ptr(), out, std::multiplies<double>());
     }
 
-    inline void div(int size, std::vector<double*>& inputs, double* out)
+    inline void div(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, inputs[1], out, std::divides<double>());
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, inputs[1]->value_ptr(), out, std::divides<double>());
     }
 
-    inline void exp(int size, std::vector<double*>& inputs, double* out)
+    inline void exp(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::exp(in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::exp(in_0);});
     }
 
-    inline void neg_exp(int size, std::vector<double*>& inputs, double* out)
+    inline void neg_exp(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::exp(-1.0*in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::exp(-1.0*in_0);});
     }
 
-    inline void sq(int size, std::vector<double*>& inputs, double* out)
+    inline void sq(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::pow(in_0, 2.0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::pow(in_0, 2.0);});
     }
 
-    inline void cb(int size, std::vector<double*>& inputs, double* out)
+    inline void cb(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::pow(in_0, 3.0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::pow(in_0, 3.0);});
     }
 
-    inline void sixth_pow(int size, std::vector<double*>& inputs, double* out)
+    inline void sixth_pow(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::pow(in_0, 6.0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::pow(in_0, 6.0);});
     }
 
-    inline void cbrt(int size, std::vector<double*>& inputs, double* out)
+    inline void cbrt(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::pow(in_0, 1.0/3.0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::pow(in_0, 1.0/3.0);});
     }
 
-    inline void sqrt(int size, std::vector<double*>& inputs, double* out)
+    inline void sqrt(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::sqrt(in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::sqrt(in_0);});
     }
 
-    inline void inv(int size, std::vector<double*>& inputs, double* out)
+    inline void inv(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return 1.0 / in_0;});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return 1.0 / in_0;});
     }
 
-    inline void log(int size, std::vector<double*>& inputs, double* out)
+    inline void log(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::log(in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::log(in_0);});
     }
 
-    inline void sin(int size, std::vector<double*>& inputs, double* out)
+    inline void sin(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::sin(in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::sin(in_0);});
     }
 
-    inline void cos(int size, std::vector<double*>& inputs, double* out)
+    inline void cos(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::cos(in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::cos(in_0);});
     }
 
-    inline void abs(int size, std::vector<double*>& inputs, double* out)
+    inline void abs(int size, std::vector<node_ptr>& inputs, double* out)
     {
-        std::transform(inputs[0], inputs[0] + size, out, [](double in_0){return std::abs(in_0);});
+        std::transform(inputs[0]->value_ptr(), inputs[0]->value_ptr() + size, out, [](double in_0){return std::abs(in_0);});
     }
 };
 
diff --git a/src/inputs/InputParser.cpp b/src/inputs/InputParser.cpp
index d649389f..ee850ae5 100644
--- a/src/inputs/InputParser.cpp
+++ b/src/inputs/InputParser.cpp
@@ -17,7 +17,7 @@ InputParser::InputParser(boost::property_tree::ptree IP, std::string fn) :
         ++_n_samp;
 }
 
-FeatureSpace InputParser::get_feature_space()
+std::shared_ptr<FeatureSpace> InputParser::get_feature_space(std::shared_ptr<MPI_Interface> comm)
 {
     std::string line;
     std::ifstream data_stream;
@@ -85,27 +85,27 @@ FeatureSpace InputParser::get_feature_space()
         throw std::logic_error(msg.substr(0, msg.size() - 2));
     }
 
-    int prop_ind = 0;
-    while((headers[prop_ind] != _prop_key) && (prop_ind < headers.size()))
+    int _propind = 0;
+    while((headers[_propind] != _prop_key) && (_propind < headers.size()))
     {
-        std::cout << _prop_key << '\t' << headers[prop_ind] << std::endl;
-        ++prop_ind;
+        std::cout << _prop_key << '\t' << headers[_propind] << std::endl;
+        ++_propind;
     }
 
-    if(prop_ind >= headers.size())
-        throw std::logic_error("prop_key not found in data_file");
+    if(_propind >= headers.size())
+        throw std::logic_error("_propkey not found in data_file");
     else
-        _prop = data[prop_ind];
+        _prop = data[_propind];
 
     std::vector<node_ptr> phi_0;
     for(int ff = 0; ff < headers.size(); ++ff)
-        if(ff != prop_ind)
+        if(ff != _propind)
             phi_0.push_back(std::make_shared<FeatureNode>(ff, headers[ff], data[ff], units[ff]));
 
     if(_opset.size() == 0)
-        return FeatureSpace(phi_0, _max_rung, _n_sis_select);
+        return std::make_shared<FeatureSpace>(comm, phi_0, _max_rung, _n_sis_select);
     else
-        return FeatureSpace(phi_0, _opset, _max_rung, _n_sis_select);
+        return std::make_shared<FeatureSpace>(comm, phi_0, _opset, _max_rung, _n_sis_select);
 }
 
 void stripComments(std::string& filename)
diff --git a/src/inputs/InputParser.hpp b/src/inputs/InputParser.hpp
index 09f37900..03241afe 100644
--- a/src/inputs/InputParser.hpp
+++ b/src/inputs/InputParser.hpp
@@ -35,7 +35,7 @@ public:
 
     InputParser(boost::property_tree::ptree IP, std::string fn);
 
-    FeatureSpace get_feature_space();
+    std::shared_ptr<FeatureSpace> get_feature_space(std::shared_ptr<MPI_Interface> comm);
 };
 /**
  * @brief      strips comments from the input file
diff --git a/src/main.cpp b/src/main.cpp
index aec6ad35..8a8f12c2 100644
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -6,19 +6,21 @@ int main(int argc, char const *argv[])
 {
     allowed_op_maps::set_node_maps();
 
+    mpi::environment env;
+    std::shared_ptr<MPI_Interface> comm = std::make_shared<MPI_Interface>();
+
     std::string filename;
     double duration= 0.0;
     if (argc < 2)
-    {
-        std::cout << "Provide an input json file" << std::endl;
-        exit(1);
-    }
+        filename = "sisso.json";
     else
-    {
-        // Take in the file name and strip out all comments for the parser
         filename = argv[1];
+
+    if(comm->rank() == 0)
         stripComments(filename);
-    }
+    else
+        filename = "stripped_" + filename;
+    comm->barrier();
 
     //construct the parser and pass it to the inputs
     boost::property_tree::ptree propTree;
@@ -27,14 +29,19 @@ int main(int argc, char const *argv[])
     InputParser IP(propTree, filename);
     boost::filesystem::remove(filename);
 
-    FeatureSpace feat_space = IP.get_feature_space();
+    if(comm->rank() == 0)
+        boost::filesystem::remove(filename);
+    std::array<int, 2> start_end = comm->get_start_end_for_iterator(67, 0);
+    std::shared_ptr<FeatureSpace> feat_space = IP.get_feature_space(comm);
+    feat_space->project_r(IP._prop);
     SISSORegressor sisso(feat_space, IP._prop, IP._n_dim);
     sisso.fit();
 
-    for(auto& model : sisso.models())
-    {
-        std::cout << model.rmse() << std::endl;
-        std::cout << model << '\n' << std::endl;
-    }
+    if(comm->rank() == 0)
+        for(auto& model : sisso.models())
+        {
+            std::cout << model.rmse() << std::endl;
+            std::cout << model << '\n' << std::endl;
+        }
     return 0;
 }
diff --git a/src/mpi_interface/MPI_Interface.cpp b/src/mpi_interface/MPI_Interface.cpp
index d8741ab1..4c18b6f0 100644
--- a/src/mpi_interface/MPI_Interface.cpp
+++ b/src/mpi_interface/MPI_Interface.cpp
@@ -1,10 +1,12 @@
+#include <mpi_interface/MPI_Interface.hpp>
+
 MPI_Interface::MPI_Interface() : boost::mpi::communicator()
 {}
 
-std::array<int, 2> MPI_Interface::get_start_end_from_list(int size, int start)
+std::array<int, 2> MPI_Interface::get_start_end_for_iterator(int length, int start)
 {
-    int els_per_rank = size / size();
-    int remaineder = size % size();
+    int els_per_rank = length / size();
+    int remaineder = length % size();
 
     std::array<int, 2> start_end;
     start_end[0] = start + els_per_rank * rank() + std::min(rank(), remaineder);
diff --git a/src/mpi_interface/MPI_Interface.hpp b/src/mpi_interface/MPI_Interface.hpp
index 7654105f..0024b28d 100644
--- a/src/mpi_interface/MPI_Interface.hpp
+++ b/src/mpi_interface/MPI_Interface.hpp
@@ -21,7 +21,7 @@ public:
      */
     MPI_Interface();
 
-    std::array<int, 2> get_start_end_from_list(int size, int start);
+    std::array<int, 2> get_start_end_for_iterator(int length, int start);
 
     /**
      * @brief      Unique int tag generator
diff --git a/src/utils/math_funcs.hpp b/src/utils/math_funcs.hpp
index 59463a07..a17f5699 100644
--- a/src/utils/math_funcs.hpp
+++ b/src/utils/math_funcs.hpp
@@ -67,6 +67,11 @@ namespace util_funcs
         std::sort(begin, end, [&vec](int i1, int i2){return vec[i1] < vec[i2];});
     }
 
+    inline int factorial(int n)
+    {
+        return (n == 1 || n == 0) ? 1 : factorial(n - 1) * n;
+    }
+
 }
 
 #endif
\ No newline at end of file
-- 
GitLab