diff --git a/src/common/aka_array_printer.hh b/src/common/aka_array_printer.hh
index e0949fb8b..f66bb7040 100644
--- a/src/common/aka_array_printer.hh
+++ b/src/common/aka_array_printer.hh
@@ -1,101 +1,101 @@
 /**
  * @file   aka_array_printer.hh
  *
  * @author Nicolas Richart
  *
  * @date creation  mer jun 19 2019
  *
  * @brief A Documented file.
  *
  * @section LICENSE
  *
  * Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as  published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A  PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 /* -------------------------------------------------------------------------- */
 #include "aka_array.hh"
 /* -------------------------------------------------------------------------- */
 
 #ifndef __AKANTU_AKA_ARRAY_PRINTER_HH__
 #define __AKANTU_AKA_ARRAY_PRINTER_HH__
 
 namespace akantu {
 
 /* -------------------------------------------------------------------------- */
 template <class container, bool no_explicit = true> class ArrayPrinter {
 public:
   ArrayPrinter(const container & cont) : cont(cont) {}
 
   void printself(std::ostream & stream, int indent = 0) const {
     std::string space(indent, AKANTU_INDENT);
 
     stream << space << "{";
     for (UInt i = 0; i < this->cont.size(); ++i) {
       stream << this->cont[i];
       if (i != this->cont.size() - 1)
         stream << ", ";
     }
-    stream << "}" << std::endl;
+    stream << "}";
   }
 
 private:
   const container & cont;
 };
 
 /* -------------------------------------------------------------------------- */
 template <class T> class ArrayPrinter<Array<T>> {
 public:
   ArrayPrinter(const Array<T> & cont) : cont(cont) {}
 
   void printself(std::ostream & stream, int indent = 0) const {
     std::string space(indent, AKANTU_INDENT);
 
     stream << space << "{";
     for (UInt i = 0; i < this->cont.size(); ++i) {
       stream << "{";
       for (UInt j = 0; j < this->cont.getNbComponent(); ++j) {
         stream << this->cont(i, j);
         if (j != this->cont.getNbComponent() - 1)
           stream << ", ";
       }
       stream << "}";
       if (i != this->cont.size() - 1)
         stream << ", ";
     }
-    stream << "}" << std::endl;
+    stream << "}";
   }
 
 private:
   const Array<T> & cont;
 };
 
 template <class container>
 decltype(auto) make_printer(const container & array) {
   return ArrayPrinter<container>(array);
 }
 
 /* -------------------------------------------------------------------------- */
 template <class T>
 inline std::ostream & operator<<(std::ostream & stream,
                                  const ArrayPrinter<T> & _this) {
   _this.printself(stream);
   return stream;
 }
 
 } // namespace akantu
 
 #endif /* __AKANTU_AKA_ARRAY_PRINTER_HH__ */
diff --git a/src/common/aka_compatibilty_with_cpp_standard.hh b/src/common/aka_compatibilty_with_cpp_standard.hh
index 3afb2ecad..2bcd70e4c 100644
--- a/src/common/aka_compatibilty_with_cpp_standard.hh
+++ b/src/common/aka_compatibilty_with_cpp_standard.hh
@@ -1,156 +1,173 @@
 /**
  * @file   aka_compatibilty_with_cpp_standard.hh
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Fri Jun 18 2010
  * @date last modification: Wed Jan 10 2018
  *
  * @brief  The content of this file is taken from the possible implementations
  * on
  * http://en.cppreference.com
  *
  * @section LICENSE
  *
  * Copyright (©)  2010-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "aka_static_if.hh"
 /* -------------------------------------------------------------------------- */
 #include <tuple>
 #include <type_traits>
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 #ifndef __AKANTU_AKA_COMPATIBILTY_WITH_CPP_STANDARD_HH__
 #define __AKANTU_AKA_COMPATIBILTY_WITH_CPP_STANDARD_HH__
 
 namespace aka {
 
 // Part taken from C++14
 #if __cplusplus < 201402L
 template <bool B, class T = void>
 using enable_if_t = typename enable_if<B, T>::type;
 #else
 template <bool B, class T = void> using enable_if_t = std::enable_if_t<B, T>;
 #endif
 
 // Part taken from C++17
 #if __cplusplus < 201703L
 // bool_constant
 template <bool B> using bool_constant = std::integral_constant<bool, B>;
 namespace {
   template <bool B> constexpr bool bool_constant_v = bool_constant<B>::value;
 }
 /* -------------------------------------------------------------------------- */
 
 // conjunction
 template <class...> struct conjunction : std::true_type {};
 template <class B1> struct conjunction<B1> : B1 {};
 template <class B1, class... Bn>
 struct conjunction<B1, Bn...>
     : std::conditional_t<bool(B1::value), conjunction<Bn...>, B1> {};
 
 // disjunction
 template <class...> struct disjunction : std::false_type {};
 template <class B1> struct disjunction<B1> : B1 {};
 template <class B1, class... Bn>
 struct disjunction<B1, Bn...>
     : std::conditional_t<bool(B1::value), B1, disjunction<Bn...>> {};
 /* -------------------------------------------------------------------------- */
 
 // negations
 template <class B> struct negation : bool_constant<!bool(B::value)> {};
 
 /* -------------------------------------------------------------------------- */
 
 namespace detail {
   template <class T> struct is_reference_wrapper : std::false_type {};
   template <class U>
   struct is_reference_wrapper<std::reference_wrapper<U>> : std::true_type {};
   // template <class T>
   // constexpr bool is_reference_wrapper_v = is_reference_wrapper<T>::value;
 
   // template <class T, class Type, class T1, class... Args>
   // decltype(auto) INVOKE(Type T::*f, T1 && t1, Args &&... args) {
   //   static_assert(std::is_member_function_pointer<decltype(f)>{} and
   //                     std::is_base_of<T, std::decay_t<T1>>{},
   //                 "Does not know what to do with this types");
   //   return (std::forward<T1>(t1).*f)(std::forward<Args>(args)...);
   // }
 
   template <class T, class Type, class T1, class... Args>
   decltype(auto) INVOKE(Type T::*f, T1 && t1, Args &&... args) {
     static_if(std::is_member_function_pointer<decltype(f)>{})
         .then_([&](auto && f) {
           static_if(std::is_base_of<T, std::decay_t<T1>>{})
               .then_([&](auto && f) {
                 return (std::forward<T1>(t1).*f)(std::forward<Args>(args)...);
               })
               .elseif(is_reference_wrapper<std::decay_t<T1>>{})([&](auto && f) {
                 return (t1.get().*f)(std::forward<Args>(args)...);
               })
               .else_([&](auto && f) {
                 return ((*std::forward<T1>(t1)).*
                         f)(std::forward<Args>(args)...);
               })(std::forward<decltype(f)>(f));
         })
         .else_([&](auto && f) {
           static_assert(std::is_member_object_pointer<decltype(f)>::value,
                         "f is not a member object");
           static_assert(sizeof...(args) == 0, "f takes arguments");
           static_if(std::is_base_of<T, std::decay_t<T1>>{})
               .then_([&](auto && f) { return std::forward<T1>(t1).*f; })
               .elseif(std::is_base_of<T, std::decay_t<T1>>{})(
                   [&](auto && f) { return t1.get().*f; })
               .else_([&](auto && f) { return (*std::forward<T1>(t1)).*f; })(
                   std::forward<decltype(f)>(f));
         })(std::forward<decltype(f)>(f));
   }
 
   template <class F, class... Args>
   decltype(auto) INVOKE(F && f, Args &&... args) {
     return std::forward<F>(f)(std::forward<Args>(args)...);
   }
 } // namespace detail
 
 template <class F, class... Args>
 decltype(auto) invoke(F && f, Args &&... args) {
   return detail::INVOKE(std::forward<F>(f), std::forward<Args>(args)...);
 }
 
 namespace detail {
   template <class F, class Tuple, std::size_t... Is>
   constexpr decltype(auto) apply_impl(F && f, Tuple && t,
                                       std::index_sequence<Is...>) {
     return invoke(std::forward<F>(f), std::get<Is>(std::forward<Tuple>(t))...);
   }
 } // namespace detail
 
+/* -------------------------------------------------------------------------- */
 template <class F, class Tuple>
 constexpr decltype(auto) apply(F && f, Tuple && t) {
   return detail::apply_impl(
       std::forward<F>(f), std::forward<Tuple>(t),
       std::make_index_sequence<std::tuple_size<std::decay_t<Tuple>>::value>{});
 }
+/* -------------------------------------------------------------------------- */
+
+template <class InputIt, class UnaryPredicate>
+typename std::iterator_traits<InputIt>::difference_type
+count_if(InputIt first, InputIt last, UnaryPredicate p) {
+  typename std::iterator_traits<InputIt>::difference_type ret = 0;
+  for (; first != last; ++first) {
+    if (p(*first)) {
+      ret++;
+    }
+  }
+  return ret;
+}
 
+#else
+using apply = std::apply;
+using count_if = std::count_if;
 #endif
 } // namespace aka
 
 #endif /* __AKANTU_AKA_COMPATIBILTY_WITH_CPP_STANDARD_HH__ */
diff --git a/src/common/aka_error.hh b/src/common/aka_error.hh
index 9dc9e316f..b889408a5 100644
--- a/src/common/aka_error.hh
+++ b/src/common/aka_error.hh
@@ -1,384 +1,389 @@
 /**
  * @file   aka_error.hh
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Mon Jun 14 2010
  * @date last modification: Tue Feb 20 2018
  *
  * @brief  error management and internal exceptions
  *
  * @section LICENSE
  *
  * Copyright (©)  2010-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include <set>
 #include <sstream>
 #include <typeinfo>
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 /* -------------------------------------------------------------------------- */
 #ifndef __AKANTU_ERROR_HH__
 #define __AKANTU_ERROR_HH__
 
 namespace akantu {
 /* -------------------------------------------------------------------------- */
 enum DebugLevel {
   dbl0 = 0,
   dblError = 0,
   dblAssert = 0,
   dbl1 = 1,
   dblException = 1,
   dblCritical = 1,
   dbl2 = 2,
   dblMajor = 2,
   dbl3 = 3,
   dblCall = 3,
   dblSecondary = 3,
   dblHead = 3,
   dbl4 = 4,
   dblWarning = 4,
   dbl5 = 5,
   dblInfo = 5,
   dbl6 = 6,
   dblIn = 6,
   dblOut = 6,
   dbl7 = 7,
   dbl8 = 8,
   dblTrace = 8,
   dbl9 = 9,
   dblAccessory = 9,
   dbl10 = 10,
   dblDebug = 42,
   dbl100 = 100,
   dblDump = 100,
   dblTest = 1337
 };
 
 /* -------------------------------------------------------------------------- */
 #define AKANTU_LOCATION                                                        \
   "(" << __func__ << "(): " << __FILE__ << ":" << __LINE__ << ")"
 
 /* -------------------------------------------------------------------------- */
 namespace debug {
   void setDebugLevel(const DebugLevel & level);
   const DebugLevel & getDebugLevel();
 
   void initSignalHandler();
   std::string demangle(const char * symbol);
   template <class T> std::string demangle() {
     return demangle(typeid(T).name());
   }
-  std::string exec(const std::string & cmd);
+
+  template <class T> std::string demangle(const T & t) {
+    return demangle(typeid(t).name());
+  }
+
+std::string exec(const std::string & cmd);
   void printBacktrace(int sig);
 
   void exit(int status) __attribute__((noreturn));
   /* ------------------------------------------------------------------------ */
   /// exception class that can be thrown by akantu
   class Exception : public std::exception {
     /* ---------------------------------------------------------------------- */
     /* Constructors/Destructors                                               */
     /* ---------------------------------------------------------------------- */
   protected:
     explicit Exception(std::string info = "")
         : _info(std::move(info)), _file("") {}
 
   public:
     //! full constructor
     Exception(std::string info, std::string file, unsigned int line)
         : _info(std::move(info)), _file(std::move(file)), _line(line) {}
 
     //! destructor
     ~Exception() noexcept override = default;
 
     /* ---------------------------------------------------------------------- */
     /*  Methods */
     /* ---------------------------------------------------------------------- */
   public:
     const char * what() const noexcept override { return _info.c_str(); }
 
     virtual const std::string info() const noexcept {
       std::stringstream stream;
       stream << debug::demangle(typeid(*this).name()) << " : " << _info << " ["
              << _file << ":" << _line << "]";
       return stream.str();
     }
 
   public:
     void setInfo(const std::string & info) { _info = info; }
     void setFile(const std::string & file) { _file = file; }
     void setLine(unsigned int line) { _line = line; }
     void setModule(const std::string & module) { _module = module; }
     /* ---------------------------------------------------------------------- */
     /* Class Members                                                          */
     /* ---------------------------------------------------------------------- */
   protected:
     /// exception description and additionals
     std::string _info;
 
   private:
     /// file it is thrown from
     std::string _file;
 
     /// line it is thrown from
     unsigned int _line{0};
 
     /// module in which exception was raised
     std::string _module{"core"};
   };
 
   class CriticalError : public Exception {};
   class AssertException : public Exception {};
   class NotImplementedException : public Exception {};
 
   /// standard output stream operator
   inline std::ostream & operator<<(std::ostream & stream,
                                    const Exception & _this) {
     stream << _this.what();
     return stream;
   }
 
   /* --------------------------------------------------------------------------
    */
   class Debugger {
   public:
     Debugger();
     virtual ~Debugger();
     Debugger(const Debugger &) = default;
     Debugger & operator=(const Debugger &) = default;
 
     void exit(int status) __attribute__((noreturn));
 
     void throwException(const std::string & info, const std::string & file,
                         unsigned int line, bool, const std::string &,
                         const std::string & module) const noexcept(false)
         __attribute__((noreturn));
 
     /*----------------------------------------------------------------------- */
     template <class Except>
     void throwCustomException(const Except & ex, const std::string & info,
                               const std::string & file, unsigned int line,
                               const std::string & module) const noexcept(false)
         __attribute__((noreturn));
     /*----------------------------------------------------------------------- */
     template <class Except>
     void throwCustomException(const Except & ex, const std::string & file,
                               unsigned int line,
                               const std::string & module) const noexcept(false)
         __attribute__((noreturn));
 
     void printMessage(const std::string & prefix, const DebugLevel & level,
                       const std::string & info,
                       const std::string & module) const;
 
     void setOutStream(std::ostream & out) { cout = &out; }
     std::ostream & getOutStream() { return *cout; }
 
   public:
     void setParallelContext(int rank, int size);
 
     void setDebugLevel(const DebugLevel & level);
     const DebugLevel & getDebugLevel() const;
 
     void setLogFile(const std::string & filename);
     std::ostream & getOutputStream();
 
     inline bool testLevel(const DebugLevel & level,
                           const std::string & module = "core") const {
       auto level_reached = (this->level >= (level));
       auto correct_module =
           (level <= dblCritical) or (modules_to_debug.size() == 0) or
           (modules_to_debug.find(module) != modules_to_debug.end());
       return level_reached and correct_module;
     }
 
     void printBacktrace(bool on_off) { this->print_backtrace = on_off; }
     bool printBacktrace() { return this->print_backtrace; }
 
     void addModuleToDebug(const std::string & id) {
       modules_to_debug.insert(id);
     }
     void removeModuleToDebug(const std::string & id) {
       auto it = modules_to_debug.find(id);
       if (it != modules_to_debug.end())
         modules_to_debug.erase(it);
     }
 
     void listModules() {
       for (auto & module_ : modules_to_debug) {
         (*cout) << module_ << std::endl;
       }
     }
 
   private:
     std::string parallel_context;
     std::ostream * cout;
     bool file_open;
     DebugLevel level;
     bool print_backtrace;
     std::set<std::string> modules_to_debug;
   };
 
   extern Debugger debugger;
 } // namespace debug
 
 /* -------------------------------------------------------------------------- */
 #define AKANTU_STRINGIZE_(str) #str
 #define AKANTU_STRINGIZE(str) AKANTU_STRINGIZE_(str)
 /* -------------------------------------------------------------------------- */
 #define AKANTU_DEBUG_MODULE AKANTU_STRINGIZE(AKANTU_MODULE)
 /* -------------------------------------------------------------------------- */
 #define AKANTU_STRINGSTREAM_IN(_str, _sstr)                                    \
   ;                                                                            \
   do {                                                                         \
     std::stringstream _dbg_s_info;                                             \
     _dbg_s_info << _sstr;                                                      \
     _str = _dbg_s_info.str();                                                  \
   } while (false)
 
 /* -------------------------------------------------------------------------- */
 #define AKANTU_EXCEPTION(info) AKANTU_EXCEPTION_(info, false)
 
 #define AKANTU_SILENT_EXCEPTION(info) AKANTU_EXCEPTION_(info, true)
 
 #define AKANTU_EXCEPTION_(info, silent)                                        \
   do {                                                                         \
     std::stringstream _dbg_str;                                                \
     _dbg_str << info;                                                          \
     std::stringstream _dbg_loc;                                                \
     _dbg_loc << AKANTU_LOCATION;                                               \
     ::akantu::debug::debugger.throwException(_dbg_str.str(), __FILE__,         \
                                              __LINE__, silent, _dbg_loc.str(), \
                                              AKANTU_DEBUG_MODULE);             \
   } while (false)
 
 #define AKANTU_CUSTOM_EXCEPTION_INFO(ex, info)                                 \
   do {                                                                         \
     std::stringstream _dbg_str;                                                \
     _dbg_str << info;                                                          \
     ::akantu::debug::debugger.throwCustomException(                            \
         ex, _dbg_str.str(), __FILE__, __LINE__, AKANTU_DEBUG_MODULE);          \
   } while (false)
 
 #define AKANTU_CUSTOM_EXCEPTION(ex)                                            \
   do {                                                                         \
     ::akantu::debug::debugger.throwCustomException(ex, __FILE__, __LINE__,     \
                                                    AKANTU_DEBUG_MODULE);       \
   } while (false)
 
 /* -------------------------------------------------------------------------- */
 #ifdef AKANTU_NDEBUG
 #define AKANTU_DEBUG_TEST(level) (false)
 #define AKANTU_DEBUG_LEVEL_IS_TEST()                                           \
   (::akantu::debug::debugger.testLevel(dblTest, AKANTU_DEBUG_MODULE))
 #define AKANTU_DEBUG(level, info)
 #define AKANTU_DEBUG_(pref, level, info)
 #define AKANTU_DEBUG_IN()
 #define AKANTU_DEBUG_OUT()
 #define AKANTU_DEBUG_INFO(info)
 #define AKANTU_DEBUG_WARNING(info)
 #define AKANTU_DEBUG_TRACE(info)
 #define AKANTU_DEBUG_ASSERT(test, info)
 #define AKANTU_ERROR(info)                                                     \
   AKANTU_CUSTOM_EXCEPTION_INFO(::akantu::debug::CriticalError(), info)
 /* -------------------------------------------------------------------------- */
 #else
 #define AKANTU_DEBUG(level, info) AKANTU_DEBUG_("   ", level, info)
 
 #define AKANTU_DEBUG_(pref, level, info)                                       \
   do {                                                                         \
     std::string _dbg_str;                                                      \
     AKANTU_STRINGSTREAM_IN(_dbg_str, info << " " << AKANTU_LOCATION);          \
     ::akantu::debug::debugger.printMessage(pref, level, _dbg_str,              \
                                            AKANTU_DEBUG_MODULE);               \
   } while (false)
 
 #define AKANTU_DEBUG_TEST(level)                                               \
   (::akantu::debug::debugger.testLevel(level, AKANTU_DEBUG_MODULE))
 
 #define AKANTU_DEBUG_LEVEL_IS_TEST()                                           \
   (::akantu::debug::debugger.testLevel(dblTest))
 
 #define AKANTU_DEBUG_IN()                                                      \
   AKANTU_DEBUG_("==>", ::akantu::dblIn, __func__ << "()")
 
 #define AKANTU_DEBUG_OUT()                                                     \
   AKANTU_DEBUG_("<==", ::akantu::dblOut, __func__ << "()")
 
 #define AKANTU_DEBUG_INFO(info) AKANTU_DEBUG_("---", ::akantu::dblInfo, info)
 
 #define AKANTU_DEBUG_WARNING(info)                                             \
   AKANTU_DEBUG_("/!\\", ::akantu::dblWarning, info)
 
 #define AKANTU_DEBUG_TRACE(info) AKANTU_DEBUG_(">>>", ::akantu::dblTrace, info)
 
 #define AKANTU_DEBUG_ASSERT(test, info)                                        \
   do {                                                                         \
     if (not(test))                                                             \
       AKANTU_CUSTOM_EXCEPTION_INFO(::akantu::debug::AssertException(),         \
                                    "assert [" << #test << "] " << info);       \
   } while (false)
 
 #define AKANTU_ERROR(info)                                                     \
   do {                                                                         \
     AKANTU_DEBUG_("!!! ", ::akantu::dblError, info);                           \
     AKANTU_CUSTOM_EXCEPTION_INFO(::akantu::debug::CriticalError(), info);      \
   } while (false)
 #endif // AKANTU_NDEBUG
 
 #define AKANTU_TO_IMPLEMENT()                                                  \
   AKANTU_CUSTOM_EXCEPTION_INFO(::akantu::debug::NotImplementedException(),     \
                                __func__ << " : not implemented yet !")
 
 /* -------------------------------------------------------------------------- */
 
 namespace debug {
   /* ------------------------------------------------------------------------ */
   template <class Except>
   void
   Debugger::throwCustomException(const Except & ex, const std::string & info,
                                  const std::string & file, unsigned int line,
                                  const std::string & module) const
       noexcept(false) {
     auto & nc_ex = const_cast<Except &>(ex);
     nc_ex.setInfo(info);
     nc_ex.setFile(file);
     nc_ex.setLine(line);
     nc_ex.setModule(module);
     throw ex;
   }
   /* ------------------------------------------------------------------------ */
   template <class Except>
   void Debugger::throwCustomException(const Except & ex,
                                       const std::string & file,
                                       unsigned int line,
                                       const std::string & module) const
       noexcept(false) {
     auto & nc_ex = const_cast<Except &>(ex);
     nc_ex.setFile(file);
     nc_ex.setLine(line);
     nc_ex.setModule(module);
     throw ex;
   }
 } // namespace debug
 } // namespace akantu
 
 #endif /* __AKANTU_ERROR_HH__ */
diff --git a/src/common/aka_random_generator.hh b/src/common/aka_random_generator.hh
index 3d22a4777..41bdf1afc 100644
--- a/src/common/aka_random_generator.hh
+++ b/src/common/aka_random_generator.hh
@@ -1,271 +1,267 @@
 /**
  * @file   aka_random_generator.hh
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Thu Feb 21 2013
  * @date last modification: Wed Nov 08 2017
  *
  * @brief  generic random generator
  *
  * @section LICENSE
  *
  * Copyright (©) 2014-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "aka_array.hh"
 /* -------------------------------------------------------------------------- */
 #include <random>
 /* -------------------------------------------------------------------------- */
 
 #ifndef __AKANTU_AKA_RANDOM_GENERATOR_HH__
 #define __AKANTU_AKA_RANDOM_GENERATOR_HH__
 
 namespace akantu {
 
 /* -------------------------------------------------------------------------- */
 /* List of available distributions                                            */
 /* -------------------------------------------------------------------------- */
 // clang-format off
 #define AKANTU_RANDOM_DISTRIBUTION_TYPES                \
   ((uniform      , std::uniform_real_distribution ))    \
   ((exponential  , std::exponential_distribution  ))    \
   ((gamma        , std::gamma_distribution        ))    \
   ((weibull      , std::weibull_distribution      ))    \
   ((extreme_value, std::extreme_value_distribution))    \
   ((normal       , std::normal_distribution       ))    \
   ((lognormal    , std::lognormal_distribution    ))    \
   ((chi_squared  , std::chi_squared_distribution  ))    \
   ((cauchy       , std::cauchy_distribution       ))    \
   ((fisher_f     , std::fisher_f_distribution     ))    \
   ((student_t    , std::student_t_distribution    ))
 // clang-format on
 
 #define AKANTU_RANDOM_DISTRIBUTION_TYPES_PREFIX(elem) BOOST_PP_CAT(_rdt_, elem)
 #define AKANTU_RANDOM_DISTRIBUTION_PREFIX(s, data, elem)                       \
   AKANTU_RANDOM_DISTRIBUTION_TYPES_PREFIX(BOOST_PP_TUPLE_ELEM(2, 0, elem))
 
 enum RandomDistributionType {
   BOOST_PP_SEQ_ENUM(BOOST_PP_SEQ_TRANSFORM(AKANTU_RANDOM_DISTRIBUTION_PREFIX, _,
                                            AKANTU_RANDOM_DISTRIBUTION_TYPES)),
   _rdt_not_defined
 };
 
 /* -------------------------------------------------------------------------- */
 /* Generator                                                                  */
 /* -------------------------------------------------------------------------- */
 template <typename T> class RandomGenerator {
   /* ------------------------------------------------------------------------ */
+private:
+  static long int _seed;
+  static std::default_random_engine generator;
+  /* ------------------------------------------------------------------------ */
 public:
   inline T operator()() { return generator(); }
 
   /// function to print the contain of the class
-  virtual void printself(std::ostream & stream, int) const {
+  void printself(std::ostream & stream, int) const {
     stream << "RandGenerator [seed=" << _seed << "]";
   }
 
   /* ------------------------------------------------------------------------ */
 public:
   static void seed(long int s) {
     _seed = s;
     generator.seed(_seed);
   }
   static long int seed() { return _seed; }
 
   static constexpr T min() { return generator.min(); }
   static constexpr T max() { return generator.max(); }
-
-  /* ------------------------------------------------------------------------ */
-private:
-  static long int _seed;
-  static std::default_random_engine generator;
 };
 
-// template <typename T>
-// long int RandomGenerator<T>::_seed = 0;
-
 /* -------------------------------------------------------------------------- */
 /* -------------------------------------------------------------------------- */
 /* -------------------------------------------------------------------------- */
 #undef AKANTU_RANDOM_DISTRIBUTION_PREFIX
 
 #define AKANTU_RANDOM_DISTRIBUTION_TYPE_PRINT_CASE(r, data, elem)              \
   case AKANTU_RANDOM_DISTRIBUTION_TYPES_PREFIX(                                \
       BOOST_PP_TUPLE_ELEM(2, 0, elem)): {                                      \
     stream << BOOST_PP_STRINGIZE(AKANTU_RANDOM_DISTRIBUTION_TYPES_PREFIX(      \
         BOOST_PP_TUPLE_ELEM(2, 0, elem)));                                     \
     break;                                                                     \
   }
 
 inline std::ostream & operator<<(std::ostream & stream,
                                  RandomDistributionType type) {
   switch (type) {
     BOOST_PP_SEQ_FOR_EACH(AKANTU_RANDOM_DISTRIBUTION_TYPE_PRINT_CASE, _,
                           AKANTU_RANDOM_DISTRIBUTION_TYPES)
   default:
     stream << UInt(type) << " not a RandomDistributionType";
     break;
   }
   return stream;
 }
 #undef AKANTU_RANDOM_DISTRIBUTION_TYPE_PRINT_CASE
 
 /* -------------------------------------------------------------------------- */
 /* Some Helper                                                                */
 /* -------------------------------------------------------------------------- */
 template <typename T, class Distribution> class RandomDistributionTypeHelper {
   enum { value = _rdt_not_defined };
 };
 
 /* -------------------------------------------------------------------------- */
 #define AKANTU_RANDOM_DISTRIBUTION_TYPE_GET_TYPE(r, data, elem)                \
   template <typename T>                                                        \
       struct RandomDistributionTypeHelper<T, BOOST_PP_TUPLE_ELEM(2, 1, elem) < \
                                                  T>> {                         \
     enum {                                                                     \
       value = AKANTU_RANDOM_DISTRIBUTION_TYPES_PREFIX(                         \
           BOOST_PP_TUPLE_ELEM(2, 0, elem))                                     \
     };                                                                         \
                                                                                \
     static void printself(std::ostream & stream) {                             \
       stream << BOOST_PP_STRINGIZE(BOOST_PP_TUPLE_ELEM(2, 0, elem));           \
     }                                                                          \
   };
 
 BOOST_PP_SEQ_FOR_EACH(AKANTU_RANDOM_DISTRIBUTION_TYPE_GET_TYPE, _,
                       AKANTU_RANDOM_DISTRIBUTION_TYPES)
 
 #undef AKANTU_RANDOM_DISTRIBUTION_TYPE_GET_TYPE
 
 /* -------------------------------------------------------------------------- */
 template <class T> class RandomDistribution {
 public:
   virtual ~RandomDistribution() = default;
   virtual T operator()(RandomGenerator<UInt> & gen) = 0;
   virtual std::unique_ptr<RandomDistribution<T>> make_unique() const = 0;
   virtual void printself(std::ostream & stream, int = 0) const = 0;
 };
 
 template <class T, class Distribution>
 class RandomDistributionProxy : public RandomDistribution<T> {
 public:
   explicit RandomDistributionProxy(Distribution dist)
       : distribution(std::move(dist)) {}
 
   T operator()(RandomGenerator<UInt> & gen) override {
     return distribution(gen);
   }
 
   std::unique_ptr<RandomDistribution<T>> make_unique() const override {
     return std::make_unique<RandomDistributionProxy<T, Distribution>>(
         distribution);
   }
 
   void printself(std::ostream & stream, int = 0) const override {
     RandomDistributionTypeHelper<T, Distribution>::printself(stream);
     stream << " [ " << distribution << " ]";
   }
 
 private:
   Distribution distribution;
 };
 
 /* -------------------------------------------------------------------------- */
 /* RandomParameter                                                            */
 /* -------------------------------------------------------------------------- */
 template <typename T> class RandomParameter {
 public:
   template <class Distribution>
   explicit RandomParameter(T base_value, Distribution dist)
       : base_value(base_value),
         type(RandomDistributionType(
             RandomDistributionTypeHelper<T, Distribution>::value)),
         distribution_proxy(
             std::make_unique<RandomDistributionProxy<T, Distribution>>(
                 std::move(dist))) {}
 
   explicit RandomParameter(T base_value)
       : base_value(base_value),
         type(RandomDistributionType(
             RandomDistributionTypeHelper<
                 T, std::uniform_real_distribution<T>>::value)),
         distribution_proxy(
             std::make_unique<
                 RandomDistributionProxy<T, std::uniform_real_distribution<T>>>(
                 std::uniform_real_distribution<T>(0., 0.))) {}
 
   RandomParameter(const RandomParameter & other)
       : base_value(other.base_value), type(other.type),
         distribution_proxy(other.distribution_proxy->make_unique()) {}
 
   RandomParameter & operator=(const RandomParameter & other) {
     distribution_proxy = other.distribution_proxy->make_unique();
     base_value = other.base_value;
     type = other.type;
     return *this;
   }
 
   virtual ~RandomParameter() = default;
 
   inline void setBaseValue(const T & value) { this->base_value = value; }
   inline T getBaseValue() const { return this->base_value; }
 
   template <template <typename> class Generator, class iterator>
   void setValues(iterator it, iterator end) {
     RandomGenerator<UInt> gen;
     for (; it != end; ++it)
       *it = this->base_value + (*distribution_proxy)(gen);
   }
 
   virtual void printself(std::ostream & stream,
                          __attribute__((unused)) int indent = 0) const {
     stream << base_value;
     stream << " + " << *distribution_proxy;
   }
 
 private:
   /// Value with no random variations
   T base_value;
 
   /// Random distribution type
   RandomDistributionType type;
 
   /// Proxy to store a std random distribution
   std::unique_ptr<RandomDistribution<T>> distribution_proxy;
 };
 
 /* -------------------------------------------------------------------------- */
 template <typename T>
 inline std::ostream & operator<<(std::ostream & stream,
                                  RandomDistribution<T> & _this) {
   _this.printself(stream);
   return stream;
 }
 
 /* -------------------------------------------------------------------------- */
 template <typename T>
 inline std::ostream & operator<<(std::ostream & stream,
                                  RandomParameter<T> & _this) {
   _this.printself(stream);
   return stream;
 }
 
 } // namespace akantu
 
 #endif /* __AKANTU_AKA_RANDOM_GENERATOR_HH__ */
diff --git a/src/fe_engine/element.hh b/src/fe_engine/element.hh
index 4ae23a2cf..1dd218498 100644
--- a/src/fe_engine/element.hh
+++ b/src/fe_engine/element.hh
@@ -1,118 +1,128 @@
 /**
  * @file   element.hh
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Tue Sep 02 2014
  * @date last modification: Tue Jan 23 2018
  *
  * @brief  Element helper class
  *
  * @section LICENSE
  *
  * Copyright (©) 2014-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "aka_common.hh"
 /* -------------------------------------------------------------------------- */
 
 #ifndef __AKANTU_ELEMENT_HH__
 #define __AKANTU_ELEMENT_HH__
 
 namespace akantu {
 
 /* -------------------------------------------------------------------------- */
 /* Element                                                                    */
 /* -------------------------------------------------------------------------- */
 class Element {
 public:
   ElementType type;
   UInt element;
   GhostType ghost_type;
 
   // ElementKind kind;
   // ElementType type{_not_defined};
   // UInt element{0};
   // GhostType ghost_type{_not_ghost};
   // ElementKind kind{_ek_regular};
 
   /* ------------------------------------------------------------------------ */
   /* Methods                                                                  */
   /* ------------------------------------------------------------------------ */
 public:
   inline ElementKind kind() const;
 
   inline bool operator==(const Element & elem) const {
     return std::tie(type, element, ghost_type) ==
            std::tie(elem.type, elem.element, elem.ghost_type);
   }
 
   inline bool operator!=(const Element & elem) const {
     return std::tie(type, element, ghost_type) !=
            std::tie(elem.type, elem.element, elem.ghost_type);
   }
 
   // inline bool operator==(const Element & elem) const {
   //   return ((element == elem.element) && (type == elem.type) &&
   //           (ghost_type == elem.ghost_type) && (kind == elem.kind));
   // }
 
   // inline bool operator!=(const Element & elem) const {
   //   return ((element != elem.element) || (type != elem.type) ||
   //           (ghost_type != elem.ghost_type) || (kind != elem.kind));
   // }
 
   inline bool operator<(const Element & rhs) const;
 };
 
 namespace {
   const Element ElementNull{_not_defined, UInt(-1), _casper};
   //      Element{_not_defined, 0, _casper, _ek_not_defined};
 } // namespace
 
 /* -------------------------------------------------------------------------- */
 inline bool Element::operator<(const Element & rhs) const {
   // bool res =
   //     (rhs == ElementNull) ||
   //     ((this->kind < rhs.kind) ||
   //      ((this->kind == rhs.kind) &&
   //       ((this->ghost_type < rhs.ghost_type) ||
   //        ((this->ghost_type == rhs.ghost_type) &&
   //         ((this->type < rhs.type) ||
   //          ((this->type == rhs.type) && (this->element < rhs.element)))))));
   return ((rhs == ElementNull) ||
           std::tie(ghost_type, type, element) <
               std::tie(rhs.ghost_type, rhs.type, rhs.element));
 }
 
-/// standard output stream operator
-inline std::ostream & operator<<(std::ostream & stream, const Element & _this) {
-  if (_this == ElementNull) {
-    stream << "ElementNull";
-    return stream;
+} // namespace akantu
+
+namespace std {
+inline string to_string(const akantu::Element & _this) {
+  if (_this == akantu::ElementNull) {
+    return "ElementNull";
   }
 
-  stream << "Element [" << _this.type << ", " << _this.element << ", "
-         << _this.ghost_type << "]";
-  return stream;
+  string str = "Element [" + to_string(_this.type) + ", " +
+               to_string(_this.element) + ", " + to_string(_this.ghost_type) +
+               "]";
+  return str;
 }
+} // namespace std
 
+namespace akantu {
+
+/// standard output stream operator
+inline std::ostream & operator<<(std::ostream & stream, const Element & _this) {
+  stream << std::to_string(_this);
+  return stream;
+}
 } // namespace akantu
 
 #endif /* __AKANTU_ELEMENT_HH__ */
diff --git a/src/fe_engine/fe_engine_inline_impl.cc b/src/fe_engine/fe_engine_inline_impl.cc
index b3fec52b3..0d7ff04e4 100644
--- a/src/fe_engine/fe_engine_inline_impl.cc
+++ b/src/fe_engine/fe_engine_inline_impl.cc
@@ -1,199 +1,193 @@
 /**
  * @file   fe_engine_inline_impl.cc
  *
  * @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Tue Jul 20 2010
  * @date last modification: Sun Aug 13 2017
  *
  * @brief  Implementation of the inline functions of the FEEngine Class
  *
  * @section LICENSE
  *
  * Copyright (©)  2010-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "element_class.hh"
 #include "fe_engine.hh"
 #include "mesh.hh"
 /* -------------------------------------------------------------------------- */
 #include "element_type_conversion.hh"
 /* -------------------------------------------------------------------------- */
 
 #ifndef __AKANTU_FE_ENGINE_INLINE_IMPL_CC__
 #define __AKANTU_FE_ENGINE_INLINE_IMPL_CC__
 
 namespace akantu {
 
 /* -------------------------------------------------------------------------- */
 inline Real FEEngine::getElementInradius(const Matrix<Real> & coord,
                                          const ElementType & type) {
-  AKANTU_DEBUG_IN();
-
   Real inradius = 0;
 
 #define GET_INRADIUS(type) inradius = ElementClass<type>::getInradius(coord);
 
   AKANTU_BOOST_ALL_ELEMENT_SWITCH(GET_INRADIUS);
 #undef GET_INRADIUS
 
-  AKANTU_DEBUG_OUT();
   return inradius;
 }
 
 /* -------------------------------------------------------------------------- */
 inline InterpolationType
 FEEngine::getInterpolationType(const ElementType & type) {
   return convertType<ElementType, InterpolationType>(type);
 }
 
 /* -------------------------------------------------------------------------- */
 /// @todo rewrite this function in order to get the cohesive element
 /// type directly from the facet
 #if defined(AKANTU_COHESIVE_ELEMENT)
 inline ElementType FEEngine::getCohesiveElementType(const ElementType & type) {
-  AKANTU_DEBUG_IN();
-
   ElementType ctype;
 #define GET_COHESIVE_TYPE(type)                                                \
   ctype = CohesiveFacetProperty<type>::cohesive_type;
 
   AKANTU_BOOST_ALL_ELEMENT_SWITCH(GET_COHESIVE_TYPE);
 #undef GET_COHESIVE_TYPE
 
-  AKANTU_DEBUG_OUT();
   return ctype;
 }
 #else
 inline ElementType
 FEEngine::getCohesiveElementType(__attribute__((unused))
                                  const ElementType & type_facet) {
   return _not_defined;
 }
 #endif
 
 /* -------------------------------------------------------------------------- */
 #if defined(AKANTU_IGFEM)
 } // akantu
 #include "igfem_helper.hh"
 namespace akantu {
 
 inline Vector<ElementType>
 FEEngine::getIGFEMElementTypes(const ElementType & type) {
 
 #define GET_IGFEM_ELEMENT_TYPES(type)                                          \
   return IGFEMHelper::getIGFEMElementTypes<type>();
 
   AKANTU_BOOST_REGULAR_ELEMENT_SWITCH(GET_IGFEM_ELEMENT_TYPES);
 
 #undef GET_IGFEM_ELEMENT_TYPES
 }
 #endif
 
 /* -------------------------------------------------------------------------- */
 template <typename T>
 void FEEngine::extractNodalToElementField(const Mesh & mesh,
                                           const Array<T> & nodal_f,
                                           Array<T> & elemental_f,
                                           const ElementType & type,
                                           const GhostType & ghost_type,
                                           const Array<UInt> & filter_elements) {
   AKANTU_DEBUG_IN();
 
   UInt nb_nodes_per_element = Mesh::getNbNodesPerElement(type);
   UInt nb_degree_of_freedom = nodal_f.getNbComponent();
   UInt nb_element = mesh.getNbElement(type, ghost_type);
   UInt * conn_val = mesh.getConnectivity(type, ghost_type).storage();
 
   if (filter_elements != empty_filter) {
     nb_element = filter_elements.size();
   }
 
   elemental_f.resize(nb_element);
 
   T * nodal_f_val = nodal_f.storage();
   T * f_val = elemental_f.storage();
 
   UInt * el_conn;
   for (UInt el = 0; el < nb_element; ++el) {
     if (filter_elements != empty_filter)
       el_conn = conn_val + filter_elements(el) * nb_nodes_per_element;
     else
       el_conn = conn_val + el * nb_nodes_per_element;
 
     for (UInt n = 0; n < nb_nodes_per_element; ++n) {
       UInt node = *(el_conn + n);
       std::copy(nodal_f_val + node * nb_degree_of_freedom,
                 nodal_f_val + (node + 1) * nb_degree_of_freedom, f_val);
       f_val += nb_degree_of_freedom;
     }
   }
 
   AKANTU_DEBUG_OUT();
 }
 
 /* -------------------------------------------------------------------------- */
 template <typename T>
 void FEEngine::filterElementalData(const Mesh & mesh, const Array<T> & elem_f,
                                    Array<T> & filtered_f,
                                    const ElementType & type,
                                    const GhostType & ghost_type,
                                    const Array<UInt> & filter_elements) {
   AKANTU_DEBUG_IN();
 
   UInt nb_element = mesh.getNbElement(type, ghost_type);
   if (nb_element == 0) {
     filtered_f.resize(0);
     return;
   }
 
   UInt nb_degree_of_freedom = elem_f.getNbComponent();
   UInt nb_data_per_element = elem_f.size() / nb_element;
 
   if (filter_elements != empty_filter) {
     nb_element = filter_elements.size();
   }
 
   filtered_f.resize(nb_element * nb_data_per_element);
 
   T * elem_f_val = elem_f.storage();
   T * f_val = filtered_f.storage();
 
   UInt el_offset;
   for (UInt el = 0; el < nb_element; ++el) {
     if (filter_elements != empty_filter)
       el_offset = filter_elements(el);
     else
       el_offset = el;
 
     std::copy(elem_f_val +
                   el_offset * nb_data_per_element * nb_degree_of_freedom,
               elem_f_val +
                   (el_offset + 1) * nb_data_per_element * nb_degree_of_freedom,
               f_val);
     f_val += nb_degree_of_freedom * nb_data_per_element;
   }
 
   AKANTU_DEBUG_OUT();
 }
 
 } // namespace akantu
 
 #endif /* __AKANTU_FE_ENGINE_INLINE_IMPL_CC__ */
diff --git a/src/mesh_utils/global_ids_updater.cc b/src/mesh_utils/global_ids_updater.cc
index 834f272f8..b3f55fa96 100644
--- a/src/mesh_utils/global_ids_updater.cc
+++ b/src/mesh_utils/global_ids_updater.cc
@@ -1,134 +1,134 @@
 /**
  * @file   global_ids_updater.cc
  *
  * @author Marco Vocialta <marco.vocialta@epfl.ch>
  *
  * @date creation: Fri Apr 13 2012
  * @date last modification: Fri Dec 08 2017
  *
  * @brief  Functions of the GlobalIdsUpdater
  *
  * @section LICENSE
  *
  * Copyright (©)  2010-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "global_ids_updater.hh"
 #include "element_synchronizer.hh"
 #include "mesh_accessor.hh"
 #include "mesh_utils.hh"
 /* -------------------------------------------------------------------------- */
 #include <numeric>
 /* -------------------------------------------------------------------------- */
 
 namespace akantu {
 
 UInt GlobalIdsUpdater::updateGlobalIDs(UInt local_nb_new_nodes) {
   if (mesh.getCommunicator().getNbProc() == 1)
     return local_nb_new_nodes;
 
   UInt total_nb_new_nodes = this->updateGlobalIDsLocally(local_nb_new_nodes);
 
   if (mesh.isDistributed()) {
     this->synchronizeGlobalIDs();
   }
   return total_nb_new_nodes;
 }
 
 UInt GlobalIdsUpdater::updateGlobalIDsLocally(UInt local_nb_new_nodes) {
   const auto & comm = mesh.getCommunicator();
   Int nb_proc = comm.getNbProc();
   if (nb_proc == 1)
     return local_nb_new_nodes;
 
   /// resize global ids array
   MeshAccessor mesh_accessor(mesh);
   auto && nodes_global_ids = mesh_accessor.getNodesGlobalIds();
   UInt old_nb_nodes = mesh.getNbNodes() - local_nb_new_nodes;
 
   nodes_global_ids.resize(mesh.getNbNodes(), -1);
 
   /// compute the number of global nodes based on the number of old nodes
   Vector<UInt> local_master_nodes(2, 0);
-  for (UInt n = 0; n < old_nb_nodes; ++n)
-    if (mesh.isLocalOrMasterNode(n))
-      ++local_master_nodes(0);
+  auto range_old = arange(old_nb_nodes);
+  local_master_nodes(0) =
+      aka::count_if(range_old.begin(), range_old.end(),
+                    [&](auto && n) { return mesh.isLocalOrMasterNode(n); });
 
   /// compute amount of local or master doubled nodes
-  for (UInt n = old_nb_nodes; n < mesh.getNbNodes(); ++n)
-    if (mesh.isLocalOrMasterNode(n))
-      ++local_master_nodes(1);
+  auto range_new = arange(old_nb_nodes, mesh.getNbNodes());
+  local_master_nodes(1) =
+      aka::count_if(range_new.begin(), range_new.end(),
+                    [&](auto && n) { return mesh.isLocalOrMasterNode(n); });
 
   auto starting_index = local_master_nodes(1);
 
   comm.allReduce(local_master_nodes);
 
   UInt old_global_nodes = local_master_nodes(0);
   UInt total_nb_new_nodes = local_master_nodes(1);
 
   if (total_nb_new_nodes == 0)
     return 0;
 
   /// set global ids of local and master nodes
   comm.exclusiveScan(starting_index);
   starting_index += old_global_nodes;
 
-  for (UInt n = old_nb_nodes; n < mesh.getNbNodes(); ++n) {
+  for (auto n : range_new) {
     if (mesh.isLocalOrMasterNode(n)) {
       nodes_global_ids(n) = starting_index;
       ++starting_index;
-    } else {
-      nodes_global_ids(n) = -1;
     }
   }
 
   mesh_accessor.setNbGlobalNodes(old_global_nodes + total_nb_new_nodes);
   return total_nb_new_nodes;
 }
 
 void GlobalIdsUpdater::synchronizeGlobalIDs() {
   this->reduce = true;
   this->synchronizer.slaveReductionOnce(*this,
                                         SynchronizationTag::_giu_global_conn);
 
 #ifndef AKANTU_NDEBUG
   for (auto node : nodes_flags) {
     auto node_flag = mesh.getNodeFlag(node.first);
     if (node_flag != NodeFlag::_pure_ghost)
       continue;
     auto n = 0u;
 
     for (auto & pair : node.second) {
       if (std::get<1>(pair) == NodeFlag::_pure_ghost)
         ++n;
     }
 
     if (n == node.second.size()) {
       AKANTU_DEBUG_WARNING(
           "The node " << n << "is ghost on all the neighboring processors");
     }
   }
 #endif
 
   this->reduce = false;
   this->synchronizer.synchronizeOnce(*this,
                                      SynchronizationTag::_giu_global_conn);
 }
 
 } // namespace akantu
diff --git a/src/synchronizer/element_synchronizer.cc b/src/synchronizer/element_synchronizer.cc
index 736037373..2af3e42db 100644
--- a/src/synchronizer/element_synchronizer.cc
+++ b/src/synchronizer/element_synchronizer.cc
@@ -1,277 +1,290 @@
 /**
  * @file   element_synchronizer.cc
  *
  * @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
  * @author Dana Christen <dana.christen@epfl.ch>
  * @author Aurelia Isabel Cuba Ramos <aurelia.cubaramos@epfl.ch>
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  * @author Marco Vocialta <marco.vocialta@epfl.ch>
  *
  * @date creation: Wed Sep 01 2010
  * @date last modification: Tue Feb 20 2018
  *
  * @brief  implementation of a  communicator using a static_communicator for
  * real
  * send/receive
  *
  * @section LICENSE
  *
  * Copyright (©)  2010-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "element_synchronizer.hh"
 #include "aka_common.hh"
 #include "mesh.hh"
 #include "mesh_utils.hh"
 /* -------------------------------------------------------------------------- */
 #include <algorithm>
 #include <iostream>
 #include <map>
 /* -------------------------------------------------------------------------- */
 
 namespace akantu {
 
+#if defined(AKANTU_MODULE)
+#define AKANTU_MODULE_SAVE_ AKANTU_MODULE
+#undef AKANTU_MODULE
+#endif
+
+#define AKANTU_MODULE element_synchronizer
+
 /* -------------------------------------------------------------------------- */
 ElementSynchronizer::ElementSynchronizer(Mesh & mesh, const ID & id,
                                          MemoryID memory_id,
                                          bool register_to_event_manager,
                                          EventHandlerPriority event_priority)
     : SynchronizerImpl<Element>(mesh.getCommunicator(), id, memory_id),
       mesh(mesh), element_to_prank("element_to_prank", id, memory_id) {
   AKANTU_DEBUG_IN();
 
   if (register_to_event_manager)
     this->mesh.registerEventHandler(*this, event_priority);
 
   AKANTU_DEBUG_OUT();
 }
 
 /* -------------------------------------------------------------------------- */
 ElementSynchronizer::ElementSynchronizer(const ElementSynchronizer & other,
                                          const ID & id,
                                          bool register_to_event_manager,
                                          EventHandlerPriority event_priority)
     : SynchronizerImpl<Element>(other, id), mesh(other.mesh),
       element_to_prank("element_to_prank", id, other.memory_id) {
   AKANTU_DEBUG_IN();
 
   element_to_prank.copy(other.element_to_prank);
 
   if (register_to_event_manager)
     this->mesh.registerEventHandler(*this, event_priority);
 
   AKANTU_DEBUG_OUT();
 }
 
 /* -------------------------------------------------------------------------- */
 ElementSynchronizer::~ElementSynchronizer() = default;
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::substituteElements(
     const std::map<Element, Element> & old_to_new_elements) {
   auto found_element_end = old_to_new_elements.end();
 
   // substitute old elements with new ones
   for (auto && sr : iterate_send_recv) {
     for (auto && scheme_pair : communications.iterateSchemes(sr)) {
       auto & list = scheme_pair.second;
       for (auto & el : list) {
         auto found_element_it = old_to_new_elements.find(el);
         if (found_element_it != found_element_end)
           el = found_element_it->second;
       }
     }
   }
 }
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::onElementsChanged(
     const Array<Element> & old_elements_list,
     const Array<Element> & new_elements_list, const ElementTypeMapArray<UInt> &,
     const ChangedElementsEvent &) {
   // create a map to link old elements to new ones
   std::map<Element, Element> old_to_new_elements;
 
   for (UInt el = 0; el < old_elements_list.size(); ++el) {
     AKANTU_DEBUG_ASSERT(old_to_new_elements.find(old_elements_list(el)) ==
                             old_to_new_elements.end(),
                         "The same element cannot appear twice in the list");
 
     old_to_new_elements[old_elements_list(el)] = new_elements_list(el);
   }
 
   substituteElements(old_to_new_elements);
 
   communications.invalidateSizes();
 }
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::onElementsRemoved(
     const Array<Element> & element_to_remove,
     const ElementTypeMapArray<UInt> & new_numbering,
     const RemovedElementsEvent &) {
   AKANTU_DEBUG_IN();
 
   this->filterScheme([&](auto && element) {
     return std::find(element_to_remove.begin(), element_to_remove.end(),
                      element) == element_to_remove.end();
   });
 
   this->renumberElements(new_numbering);
 
   communications.invalidateSizes();
 
   AKANTU_DEBUG_OUT();
 }
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::buildElementToPrank() {
   AKANTU_DEBUG_IN();
 
   UInt spatial_dimension = mesh.getSpatialDimension();
   element_to_prank.initialize(mesh, _spatial_dimension = spatial_dimension,
                               _element_kind = _ek_not_defined,
                               _with_nb_element = true, _default_value = rank);
 
   /// assign prank to all ghost elements
   for (auto && scheme : communications.iterateSchemes(_recv)) {
     auto & recv = scheme.second;
     auto proc = scheme.first;
 
     for (auto & element : recv) {
       element_to_prank(element) = proc;
     }
   }
 
   AKANTU_DEBUG_OUT();
 }
 
 /* -------------------------------------------------------------------------- */
 Int ElementSynchronizer::getRank(const Element & element) const {
   if (not element_to_prank.exists(element.type, element.ghost_type)) {
     // Nicolas: Ok This is nasty I know....
     const_cast<ElementSynchronizer *>(this)->buildElementToPrank();
   }
 
   return element_to_prank(element);
 }
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::renumberElements(
     const ElementTypeMapArray<UInt> & new_numbering) {
   for (auto && sr : iterate_send_recv) {
     for (auto && scheme_pair : communications.iterateSchemes(sr)) {
       auto & list = scheme_pair.second;
       for (auto && el : list) {
         if (new_numbering.exists(el.type, el.ghost_type))
           el.element = new_numbering(el);
       }
     }
   }
 }
 
 /* -------------------------------------------------------------------------- */
 UInt ElementSynchronizer::sanityCheckDataSize(const Array<Element> & elements,
                                               const SynchronizationTag & tag,
                                               bool from_comm_desc) const {
   UInt size = SynchronizerImpl<Element>::sanityCheckDataSize(elements, tag,
                                                              from_comm_desc);
 
   // global connectivities;
   size += mesh.getNbNodesPerElementList(elements) * sizeof(UInt);
 
   // barycenters
   size += (elements.size() * mesh.getSpatialDimension() * sizeof(Real));
   return size;
 }
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::packSanityCheckData(
     CommunicationBuffer & buffer, const Array<Element> & elements,
     const SynchronizationTag & /*tag*/) const {
   for (auto && element : elements) {
     Vector<Real> barycenter(mesh.getSpatialDimension());
     mesh.getBarycenter(element, barycenter);
     buffer << barycenter;
 
     const auto & conns = mesh.getConnectivity(element.type, element.ghost_type);
     for (auto n : arange(conns.getNbComponent())) {
       buffer << mesh.getNodeGlobalId(conns(element.element, n));
     }
   }
 }
 
 /* -------------------------------------------------------------------------- */
 void ElementSynchronizer::unpackSanityCheckData(CommunicationBuffer & buffer,
                                                 const Array<Element> & elements,
                                                 const SynchronizationTag & tag,
                                                 UInt proc, UInt rank) const {
   auto spatial_dimension = mesh.getSpatialDimension();
 
   std::set<SynchronizationTag> skip_conn_tags{
       SynchronizationTag::_smmc_facets_conn,
       SynchronizationTag::_giu_global_conn};
 
   bool is_skip_tag_conn = skip_conn_tags.find(tag) != skip_conn_tags.end();
 
   for (auto && element : elements) {
     Vector<Real> barycenter_loc(spatial_dimension);
     mesh.getBarycenter(element, barycenter_loc);
 
     Vector<Real> barycenter(spatial_dimension);
     buffer >> barycenter;
 
     auto dist = barycenter_loc.distance(barycenter);
     if (not Math::are_float_equal(dist, 0.)) {
       AKANTU_EXCEPTION("Unpacking an unknown value for the element "
                        << element << "(barycenter " << barycenter_loc
                        << " != buffer " << barycenter << ") [" << dist
                        << "] - tag: " << tag << " comm from " << proc << " to "
                        << rank);
     }
 
     const auto & conns = mesh.getConnectivity(element.type, element.ghost_type);
     Vector<UInt> global_conn(conns.getNbComponent());
     Vector<UInt> local_global_conn(conns.getNbComponent());
 
     auto is_same = true;
     for (auto n : arange(global_conn.size())) {
       buffer >> global_conn(n);
 
       auto node = conns(element.element, n);
       local_global_conn(n) = mesh.getNodeGlobalId(node);
 
       is_same &= is_skip_tag_conn or mesh.isPureGhostNode(node) or
                  (local_global_conn(n) == global_conn(n));
     }
 
     if (not is_same) {
       AKANTU_DEBUG_WARNING(
           "The connectivity of the element "
           << element << " " << local_global_conn
           << " does not match the connectivity of the equivalent "
              "element on proc "
           << proc << " " << global_conn << " in communication with tag "
           << tag);
     }
   }
 }
 
 /* -------------------------------------------------------------------------- */
 } // namespace akantu
+
+#if defined(AKANTU_MODULE_SAVE_)
+#undef AKANTU_MODULE
+#define AKANTU_MODULE AKANTU_MODULE_SAVE_
+#undef AKANTU_MODULE_SAVE_
+#endif
diff --git a/src/synchronizer/facet_synchronizer.cc b/src/synchronizer/facet_synchronizer.cc
index 9de07e683..7b98455cf 100644
--- a/src/synchronizer/facet_synchronizer.cc
+++ b/src/synchronizer/facet_synchronizer.cc
@@ -1,206 +1,220 @@
 /**
  * @file   facet_synchronizer.cc
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  * @author Marco Vocialta <marco.vocialta@epfl.ch>
  *
  * @date creation: Wed Nov 05 2014
  * @date last modification: Fri Jan 26 2018
  *
  * @brief  Facet synchronizer for parallel simulations with cohesive elments
  *
  * @section LICENSE
  *
  * Copyright (©) 2015-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * Akantu is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Akantu is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Akantu. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "facet_synchronizer.hh"
 /* -------------------------------------------------------------------------- */
 
+#if defined(AKANTU_MODULE)
+#define AKANTU_MODULE_SAVE_ AKANTU_MODULE
+#undef AKANTU_MODULE
+#endif
+
+#define AKANTU_MODULE facet_synchronizer
+
+
 namespace akantu {
 
 /* -------------------------------------------------------------------------- */
 FacetSynchronizer::FacetSynchronizer(
     Mesh & mesh, const ElementSynchronizer & element_synchronizer,
     const ID & id, MemoryID memory_id)
     : ElementSynchronizer(mesh, id, memory_id) {
 
   auto spatial_dimension = mesh.getSpatialDimension();
 
   element_to_prank.initialize(mesh, _spatial_dimension = spatial_dimension - 1,
                               _ghost_type = _ghost, _with_nb_element = true,
                               _default_value = rank);
 
   // Build element to prank
   for (auto && scheme_pair :
        element_synchronizer.communications.iterateSchemes(_recv)) {
     auto proc = std::get<0>(scheme_pair);
     const auto & scheme = std::get<1>(scheme_pair);
 
     for (auto && elem : scheme) {
       const auto & facet_to_element =
           mesh.getSubelementToElement(elem.type, elem.ghost_type);
       Vector<Element> facets = facet_to_element.begin(
           facet_to_element.getNbComponent())[elem.element];
 
       for (UInt f = 0; f < facets.size(); ++f) {
         const auto & facet = facets(f);
         if (facet == ElementNull)
           continue;
 
         if (facet.ghost_type == _not_ghost)
           continue;
 
         auto & facet_rank = element_to_prank(facet);
         if ((proc < UInt(facet_rank)) || (UInt(facet_rank) == rank))
           facet_rank = proc;
       }
     }
   }
 
   ElementTypeMapArray<UInt> facet_global_connectivities(
       "facet_global_connectivities", id, memory_id);
   facet_global_connectivities.initialize(
       mesh, _spatial_dimension = spatial_dimension - 1, _with_nb_element = true,
       _with_nb_nodes_per_element = true);
   mesh.getGlobalConnectivity(facet_global_connectivities);
 
   // \TODO perhaps a global element numbering might be useful here...
   for (auto type : facet_global_connectivities.elementTypes(_spatial_dimension =
                                                                 _all_dimensions,
                    _element_kind = _ek_not_defined, _ghost_type = _not_ghost)) {
     auto & conn = facet_global_connectivities(type, _not_ghost);
     auto conn_view = make_view(conn, conn.getNbComponent());
     std::for_each(conn_view.begin(), conn_view.end(), [&](auto & conn) {
       std::sort(conn.storage(), conn.storage() + conn.size());
     });
   }
 
   /// init facet check tracking
   ElementTypeMapArray<bool> facet_checked("facet_checked", id, memory_id);
   std::map<UInt, ElementTypeMapArray<UInt>> recv_connectivities;
 
   /// Generate the recv scheme and connnectivities to send to the other
   /// processors
   for (auto && scheme_pair :
        element_synchronizer.communications.iterateSchemes(_recv)) {
     facet_checked.initialize(mesh, _spatial_dimension = spatial_dimension - 1,
                              _ghost_type = _ghost, _with_nb_element = true,
                              _default_value = false);
 
     auto proc = scheme_pair.first;
     const auto & elements = scheme_pair.second;
     auto & facet_scheme = communications.createScheme(proc, _recv);
 
     // this creates empty arrays...
     auto & connectivities_for_proc = recv_connectivities[proc];
 
     connectivities_for_proc.setID(
         id + ":connectivities_for_proc:" + std::to_string(proc));
     connectivities_for_proc.initialize(
         mesh, _spatial_dimension = spatial_dimension - 1,
         _with_nb_nodes_per_element = true, _ghost_type = _ghost);
 
     // for every element in the element synchronizer communication scheme,
     // check the facets to see if they should be communicated and create a
     // connectivity array to match with the one other processors might send
     for (auto && element : elements) {
       const auto & facet_to_element =
           mesh.getSubelementToElement(element.type, element.ghost_type);
       Vector<Element> facets = facet_to_element.begin(
           facet_to_element.getNbComponent())[element.element];
 
       for (UInt f = 0; f < facets.size(); ++f) {
         auto & facet = facets(f);
 
         // exclude no valid facets
         if (facet == ElementNull)
           continue;
 
         // exclude _ghost facet from send scheme and _not_ghost from receive
         if (facet.ghost_type != _ghost)
           continue;
 
         // exclude facet from other processors then the one of current
         // interest in case of receive scheme
         if (UInt(element_to_prank(facet)) != proc)
           continue;
 
         auto & checked = facet_checked(facet);
         // skip already checked facets
         if (checked)
           continue;
 
         checked = true;
 
         facet_scheme.push_back(facet);
 
         auto & global_conn =
             facet_global_connectivities(facet.type, facet.ghost_type);
         Vector<UInt> conn =
             global_conn.begin(global_conn.getNbComponent())[facet.element];
         std::sort(conn.storage(), conn.storage() + conn.size());
 
         connectivities_for_proc(facet.type, facet.ghost_type).push_back(conn);
       }
     }
   }
 
   std::vector<CommunicationRequest> send_requests;
   /// do every communication by element type
   for (auto && type : mesh.elementTypes(spatial_dimension - 1)) {
     for (auto && pair : recv_connectivities) {
       auto proc = std::get<0>(pair);
       const auto & connectivities_for_proc = std::get<1>(pair);
       auto && tag = Tag::genTag(proc, type, 1337);
       send_requests.push_back(
           communicator.asyncSend(connectivities_for_proc(type, _ghost), proc,
                                  tag, CommunicationMode::_synchronous));
     }
 
     auto nb_nodes_per_facet = Mesh::getNbNodesPerElement(type);
 
     communicator.receiveAnyNumber<UInt>(
         send_requests,
         [&](auto && proc, auto && message) {
           auto & local_connectivities =
               facet_global_connectivities(type, _not_ghost);
           auto & send_scheme = communications.createScheme(proc, _send);
 
           auto conn_view = make_view(local_connectivities, nb_nodes_per_facet);
           auto conn_begin = conn_view.begin();
           auto conn_end = conn_view.end();
           for (const auto & c_to_match :
                make_view(message, nb_nodes_per_facet)) {
             auto it = std::find(conn_begin, conn_end, c_to_match);
 
             if (it != conn_end) {
               auto facet = Element{type, UInt(it - conn_begin), _not_ghost};
               send_scheme.push_back(facet);
             } else {
               AKANTU_EXCEPTION("No local facet found to send to proc "
                                << proc << " corresponding to " << c_to_match);
             }
           }
         },
         Tag::genTag(rank, type, 1337));
   }
 }
 
 } // namespace akantu
+
+#if defined(AKANTU_MODULE_SAVE_)
+#undef AKANTU_MODULE
+#define AKANTU_MODULE AKANTU_MODULE_SAVE_
+#undef AKANTU_MODULE_SAVE_
+#endif