diff --git a/third-party/akantu_iterators/include/aka_compatibilty_with_cpp_standard.hh b/third-party/akantu_iterators/include/aka_compatibilty_with_cpp_standard.hh
index f82a03c0f..ed6f90b0a 100644
--- a/third-party/akantu_iterators/include/aka_compatibilty_with_cpp_standard.hh
+++ b/third-party/akantu_iterators/include/aka_compatibilty_with_cpp_standard.hh
@@ -1,200 +1,211 @@
 /**
  * @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 <iterator>
 #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)> {};
 
 /* -------------------------------------------------------------------------- */
 // invoke
 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, 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)...);
 }
 
 /* -------------------------------------------------------------------------- */
 // apply
 namespace detail {
   template <class F, class Tuple, std::size_t... Is>
   constexpr decltype(auto) apply_impl(F && f, Tuple && t,
                                       std::index_sequence<Is...> /*unused*/) {
     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>{});
 }
 
 /* -------------------------------------------------------------------------- */
 // count_if
 template <class InputIt, class UnaryPredicate>
-typename std::iterator_traits<InputIt>::difference_type
-count_if(InputIt first, InputIt last, UnaryPredicate p) {
+auto count_if(InputIt first, InputIt last, UnaryPredicate p) ->
+    typename std::iterator_traits<InputIt>::difference_type {
   typename std::iterator_traits<InputIt>::difference_type ret = 0;
   for (; first != last; ++first) {
     if (p(*first)) {
       ret++;
     }
   }
   return ret;
 }
 
 namespace detail {
   template <class T, class Tuple, std::size_t... I>
-  constexpr T make_from_tuple_impl(Tuple && t, std::index_sequence<I...>) {
+  constexpr auto make_from_tuple_impl(Tuple && t,
+                                      std::index_sequence<I...> /*unused*/)
+      -> T {
     return T(std::get<I>(std::forward<Tuple>(t))...);
   }
 } // namespace detail
 
-template <class T, class Tuple> constexpr T make_from_tuple(Tuple && t) {
+template <class T, class Tuple>
+constexpr auto make_from_tuple(Tuple && t) -> T {
   return detail::make_from_tuple_impl<T>(
       std::forward<Tuple>(t),
       std::make_index_sequence<
           std::tuple_size<std::remove_reference_t<Tuple>>::value>{});
 }
 
 #else
 template <bool B> using bool_constant = std::bool_constant<B>;
 template <bool B> using bool_constant_v = std::bool_constant_v<B>;
 
 template <class... Args> using conjunction = std::conjunction<Args...>;
 template <class... Args> using disjunction = std::disjunction<Args...>;
 template <class B> using negation = std::negation<B>;
 
 using invoke = std::invoke;
 using apply = std::apply;
 
 using count_if = std::count_if;
 #endif
 
 template <typename cat1, typename cat2>
 using is_iterator_category_at_least =
     std::is_same<std::common_type_t<cat1, cat2>, cat2>;
 
+template<typename T>
+struct size_type {
+  using type = typename std::decay_t<T>::size_type;
+};
+
+template<typename T>
+using size_type_t = typename size_type<T>::type;
+
 } // namespace aka
 
 #endif /* AKANTU_AKA_COMPATIBILTY_WITH_CPP_STANDARD_HH */
diff --git a/third-party/akantu_iterators/include/aka_tuple_tools.hh b/third-party/akantu_iterators/include/aka_tuple_tools.hh
index c80f1cd01..0a33cf133 100644
--- a/third-party/akantu_iterators/include/aka_tuple_tools.hh
+++ b/third-party/akantu_iterators/include/aka_tuple_tools.hh
@@ -1,364 +1,544 @@
 /**
  * @file   aka_tuple_tools.hh
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Fri Aug 11 2017
  * @date last modification: Mon Jan 29 2018
  *
  * @brief  iterator interfaces
  *
  * @section LICENSE
  *
  * Copyright 2019 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_compatibilty_with_cpp_standard.hh"
 #include "aka_str_hash.hh"
 /* -------------------------------------------------------------------------- */
 #include <tuple>
 /* -------------------------------------------------------------------------- */
 
 #ifndef AKANTU_AKA_TUPLE_TOOLS_HH
 #define AKANTU_AKA_TUPLE_TOOLS_HH
 
 #ifndef AKANTU_ITERATORS_NAMESPACE
 #define AKANTU_ITERATORS_NAMESPACE akantu
 #endif
 
 namespace AKANTU_ITERATORS_NAMESPACE {
 
 namespace tuple {
   /* ---------------------------------------------------------------------- */
   template <typename tag, typename type> struct named_tag {
     using _tag = tag;   ///< key
     using _type = type; ///< value type
 
     template <
         typename T,
         std::enable_if_t<not std::is_same<named_tag, T>::value> * = nullptr>
-    explicit named_tag(T && value) // NOLINT
+    named_tag(T && value) // NOLINT
         : _value(std::forward<T>(value)) {}
 
     type _value;
   };
 
   namespace details {
 /* ---------------------------------------------------------------------- */
 #if (defined(__GNUC__) || defined(__GNUG__))
 #pragma GCC diagnostic push
 #pragma GCC diagnostic ignored "-Weffc++"
 #endif
     template <typename tag> struct named_tag_proxy {
       using _tag = tag;
 
       template <typename T> decltype(auto) operator=(T && value) {
         return named_tag<_tag, T>{std::forward<T>(value)};
       }
     };
 #if (defined(__GNUC__) || defined(__GNUG__))
 #pragma GCC diagnostic pop
 #endif
   } // namespace details
 
   /* ---------------------------------------------------------------------- */
   template <typename T> struct is_named_tag : public std::false_type {};
   template <typename tag>
   struct is_named_tag<details::named_tag_proxy<tag>> : public std::true_type {};
   template <typename tag, typename type>
   struct is_named_tag<named_tag<tag, type>> : public std::true_type {};
   /* ---------------------------------------------------------------------- */
 
   template <class... Params>
   struct named_tuple : public std::tuple<typename Params::_type...> {
     using Names_t = std::tuple<typename Params::_tag...>;
     using parent = std::tuple<typename Params::_type...>;
 
     named_tuple(Params &&... params)
         : parent(std::forward<typename Params::_type>(params._value)...) {}
 
     named_tuple(typename Params::_type &&... args)
         : parent(std::forward<typename Params::_type>(args)...) {}
 
   private:
     template <typename tag, std::size_t Idx,
               std::enable_if_t<Idx == sizeof...(Params)> * = nullptr>
-    static constexpr std::size_t get_element_index() noexcept {
-      return -1;
+    static constexpr std::size_t get_element_index_() noexcept {
+      return std::size_t(-1);
     }
 
     template <typename tag, std::size_t Idx,
               std::enable_if_t<(Idx < sizeof...(Params))> * = nullptr>
-    static constexpr std::size_t get_element_index() noexcept {
+    static constexpr std::size_t get_element_index_() noexcept {
       using _tag = std::tuple_element_t<Idx, Names_t>;
       return (std::is_same<_tag, tag>::value)
                  ? Idx
-                 : get_element_index<tag, Idx + 1>();
+                 : get_element_index_<tag, Idx + 1>();
     }
 
   public:
+    template <typename NT>
+    static constexpr std::size_t get_element_index() noexcept {
+      return get_element_index_<typename NT::_tag, 0>();
+    }
+
+    template <typename NT>
+    static constexpr std::size_t get_element_index(NT && /*unused*/) noexcept {
+      return get_element_index_<typename NT::_tag, 0>();
+    }
+
     template <typename NT,
               std::enable_if_t<is_named_tag<NT>::value> * = nullptr>
     constexpr decltype(auto) get(NT && /*unused*/) noexcept {
-      const auto index = get_element_index<typename NT::_tag, 0>();
-      static_assert((index != -1), "wrong named_tag");
+      const auto index = get_element_index<NT>();
+      static_assert((index != std::size_t(-1)), "wrong named_tag");
       return (std::get<index>(*this));
     }
 
     template <typename NT,
               std::enable_if_t<is_named_tag<NT>::value> * = nullptr>
     constexpr decltype(auto) get(NT && /*unused*/) const noexcept {
-      const auto index = get_element_index<typename NT::_tag, 0>();
-      static_assert((index != -1), "wrong named_tag");
-      return std::get<index>(*this);
+      const auto index = get_element_index<NT>();
+      static_assert((index != std::size_t(-1)), "wrong named_tag");
+      return (std::get<index>(*this));
+    }
+
+    template <typename NT,
+              std::enable_if_t<is_named_tag<NT>::value> * = nullptr>
+    constexpr auto has(NT && /*unused*/) const noexcept -> bool {
+      const auto index = get_element_index<NT>();
+      return (index != std::size_t(-1));
     }
   };
 
   /* ---------------------------------------------------------------------- */
   template <typename T> struct is_named_tuple : public std::false_type {};
   template <typename... Params>
   struct is_named_tuple<named_tuple<Params...>> : public std::true_type {};
   /* ---------------------------------------------------------------------- */
 
   template <typename... Params>
   constexpr decltype(auto) make_named_tuple(Params &&... params) noexcept {
     return named_tuple<Params...>(std::forward<Params>(params)...);
   }
 
   template <typename tag> constexpr decltype(auto) make_named_tag() noexcept {
     return details::named_tag_proxy<tag>{};
   }
 
   template <size_t HashCode> constexpr decltype(auto) get() {
     return make_named_tag<std::integral_constant<size_t, HashCode>>();
   }
 
   template <size_t HashCode, class Tuple>
   constexpr decltype(auto) get(Tuple && tuple) {
     return tuple.get(get<HashCode>());
   }
 
+  template <size_t HashCode, class Tuple> constexpr bool has(Tuple && tuple) {
+    return tuple.has(get<HashCode>());
+  }
+
   template <typename Param, typename Tuple,
             std::enable_if_t<is_named_tag<Param>::value> * = nullptr>
   constexpr decltype(auto) get(Tuple && tuple) noexcept {
     return tuple.template get<typename Param::hash>();
   }
 
+  template <size_t I, class Tuple> struct tuple_name_tag {
+    using type = std::tuple_element_t<I, typename Tuple::Names_t>;
+  };
+
+  template <size_t I, class Tuple>
+  using tuple_name_tag_t = typename tuple_name_tag<I, Tuple>::type;
+
+  template <size_t I, class Tuple> struct tuple_element {
+    using type = std::tuple_element_t<I, typename Tuple::parent>;
+  };
+
+  template <size_t I, class Tuple>
+  using tuple_element_t = typename tuple_element<I, Tuple>::type;
+
 #if defined(__INTEL_COMPILER)
 // intel warnings here
 #elif defined(__clang__)
 // clang warnings here
 #pragma clang diagnostic push
 #pragma clang diagnostic ignored "-Wgnu-string-literal-operator-template"
 #elif (defined(__GNUC__) || defined(__GNUG__))
 // gcc warnings here
 #pragma GCC diagnostic push
 #pragma GCC diagnostic ignored "-Wpedantic"
 #endif
   /// this is a GNU exstension
   /// http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2013/n3599.html
   template <class CharT, CharT... chars>
   constexpr decltype(auto) operator"" _n() {
     return make_named_tag<std::integral_constant<
         std::size_t, string_literal<CharT, chars...>::hash>>();
   }
 #if defined(__INTEL_COMPILER)
 #elif defined(__clang__)
 #pragma clang diagnostic pop
 #elif (defined(__GNUC__) || defined(__GNUG__))
 #pragma GCC diagnostic pop
 #endif
 
-  /* ------------------------------------------------------------------------ */
+  /* ------------------------------------------------------------------------
+   */
   namespace details {
     template <std::size_t N> struct Foreach {
       template <class Tuple>
       static inline bool not_equal(Tuple && a, Tuple && b) {
         if (std::get<N - 1>(std::forward<Tuple>(a)) ==
             std::get<N - 1>(std::forward<Tuple>(b))) {
           return false;
         }
         return Foreach<N - 1>::not_equal(std::forward<Tuple>(a),
                                          std::forward<Tuple>(b));
       }
     };
 
-    /* ---------------------------------------------------------------------- */
+    /* ----------------------------------------------------------------------
+     */
     template <> struct Foreach<0> {
       template <class Tuple>
       static inline bool not_equal(Tuple && a, Tuple && b) {
         return std::get<0>(std::forward<Tuple>(a)) !=
                std::get<0>(std::forward<Tuple>(b));
       }
     };
 
     template <typename... Ts>
     decltype(auto) make_tuple_no_decay(Ts &&... args) {
       return std::tuple<Ts...>(std::forward<Ts>(args)...);
     }
 
     template <typename... Names, typename... Ts>
     decltype(auto) make_named_tuple_no_decay(std::tuple<Names...> /*unused*/,
                                              Ts &&... args) {
       return named_tuple<named_tag<Names, Ts>...>(std::forward<Ts>(args)...);
     }
 
     template <class F, class Tuple, std::size_t... Is>
     void foreach_impl(F && func, Tuple && tuple,
                       std::index_sequence<Is...> && /*unused*/) {
       (void)std::initializer_list<int>{
           (std::forward<F>(func)(std::get<Is>(std::forward<Tuple>(tuple))),
            0)...};
     }
 
     template <class F, class Tuple, std::size_t... Is>
     decltype(auto) transform_impl(F && func, Tuple && tuple,
                                   std::index_sequence<Is...> && /*unused*/) {
       return make_tuple_no_decay(
           std::forward<F>(func)(std::get<Is>(std::forward<Tuple>(tuple)))...);
     }
 
     template <class F, class Tuple, std::size_t... Is>
     decltype(auto)
     transform_named_impl(F && func, Tuple && tuple,
                          std::index_sequence<Is...> && /*unused*/) {
       return make_named_tuple_no_decay(
           typename std::decay_t<Tuple>::Names_t{},
           std::forward<F>(func)(std::get<Is>(std::forward<Tuple>(tuple)))...);
     }
 
+    template <class Tuple, std::size_t... Is>
+    decltype(auto) flatten(Tuple && tuples,
+                           std::index_sequence<Is...> /*unused*/) {
+      return std::tuple_cat(std::get<Is>(tuples)...);
+    }
+
+    template <class Tuple, class... Vals, std::size_t... Is>
+    decltype(auto) append_impl(std::index_sequence<Is...> && /*unused*/,
+                               Tuple && tuple, Vals &&... other_vals) {
+      return make_tuple_no_decay(std::get<Is>(std::forward<Tuple>(tuple))...,
+                                 std::forward<Vals>(other_vals)...);
+    }
+
+    template <class Tuple, class... Vals, std::size_t... Is>
+    decltype(auto) append_named_impl(std::index_sequence<Is...> && /*unused*/,
+                                     Tuple && tuple, Vals &&... other_vals) {
+      return make_named_tuple_no_decay(
+          std::tuple<tuple_name_tag<Is, Tuple>..., Vals...>{},
+          std::get<Is>(std::forward<Tuple>(tuple))...,
+          std::forward<Vals>(other_vals)...);
+    }
+
+    template <std::size_t nth, class Tuple, std::size_t... Is_before,
+              std::size_t... Is_after>
+    decltype(auto) remove_impl(std::index_sequence<Is_before...> && /*unused*/,
+                               std::index_sequence<Is_after...> && /*unused*/,
+                               Tuple && tuple) {
+      return make_tuple_no_decay(
+          std::get<Is_before>(std::forward<Tuple>(tuple))...,
+          std::get<Is_after + nth + 1>(std::forward<Tuple>(tuple))...);
+    }
+
+    template <std::size_t nth, class Tuple, std::size_t... Is_before,
+              std::size_t... Is_after>
+    decltype(auto)
+    remove_named_impl(std::index_sequence<Is_before...> && /*unused*/,
+                      std::index_sequence<Is_after...> && /*unused*/,
+                      Tuple && tuple) {
+      return make_named_tuple_no_decay(
+          std::tuple<tuple::tuple_name_tag_t<Is_before, std::decay_t<Tuple>>...,
+                     tuple::tuple_name_tag_t<Is_after + nth + 1,
+                                             std::decay_t<Tuple>>...>{},
+          std::get<Is_before>(std::forward<Tuple>(tuple))...,
+          std::get<Is_after + nth + 1>(std::forward<Tuple>(tuple))...);
+    }
+
+    template <std::size_t nth, class Tuple, class Value,
+              std::size_t... Is_before, std::size_t... Is_after>
+    decltype(auto) replace_impl(std::index_sequence<Is_before...> && /*unused*/,
+                                std::index_sequence<Is_after...> && /*unused*/,
+                                Tuple && tuple, Value && value) {
+      return make_tuple_no_decay(
+          std::get<Is_before>(std::forward<Tuple>(tuple))...,
+          std::forward<Value>(value),
+          std::get<Is_after + nth + 1>(std::forward<Tuple>(tuple))...);
+    }
+
+    template <std::size_t nth, size_t HashCode, class Tuple, class Value,
+              std::size_t... Is_before, std::size_t... Is_after>
+    decltype(auto)
+    replace_named_impl(std::index_sequence<Is_before...> && /*unused*/,
+                       std::index_sequence<Is_after...> && /*unused*/,
+                       Tuple && tuple, Value && value) {
+      using tuple_type = std::decay_t<Tuple>;
+      return make_named_tuple_no_decay(
+          std::tuple<
+              tuple::tuple_name_tag_t<Is_before, tuple_type>...,
+              decltype(get<HashCode>()),
+              tuple::tuple_name_tag_t<Is_after + nth + 1, tuple_type>...>{},
+          std::get<Is_before>(std::forward<Tuple>(tuple))...,
+          std::forward<Value>(value),
+          std::get<Is_after + nth + 1>(std::forward<Tuple>(tuple))...);
+    }
   } // namespace details
 
   /* ------------------------------------------------------------------------ */
   template <class Tuple,
             std::enable_if_t<not is_named_tuple<std::decay_t<Tuple>>::value> * =
                 nullptr>
   bool are_not_equal(Tuple && a, Tuple && b) {
     return details::Foreach<std::tuple_size<std::decay_t<Tuple>>::value>::
         not_equal(std::forward<Tuple>(a), std::forward<Tuple>(b));
   }
 
   template <
       class Tuple,
       std::enable_if_t<is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
   bool are_not_equal(Tuple && a, Tuple && b) {
     return details::Foreach<
         std::tuple_size<typename std::decay_t<Tuple>::parent>::value>::
         not_equal(std::forward<Tuple>(a), std::forward<Tuple>(b));
   }
 
+  /* ------------------------------------------------------------------------ */
   template <class F, class Tuple,
             std::enable_if_t<not is_named_tuple<std::decay_t<Tuple>>::value> * =
                 nullptr>
   void foreach (F && func, Tuple && tuple) {
     return details::foreach_impl(
         std::forward<F>(func), std::forward<Tuple>(tuple),
         std::make_index_sequence<
             std::tuple_size<std::decay_t<Tuple>>::value>{});
   }
 
   template <
       class F, class Tuple,
       std::enable_if_t<is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
   void foreach (F && func, Tuple && tuple) {
     return details::foreach_impl(
         std::forward<F>(func), std::forward<Tuple>(tuple),
         std::make_index_sequence<
             std::tuple_size<typename std::decay_t<Tuple>::parent>::value>{});
   }
 
+  /* ------------------------------------------------------------------------ */
   template <class F, class Tuple,
             std::enable_if_t<not is_named_tuple<std::decay_t<Tuple>>::value> * =
                 nullptr>
   decltype(auto) transform(F && func, Tuple && tuple) {
     return details::transform_impl(
         std::forward<F>(func), std::forward<Tuple>(tuple),
         std::make_index_sequence<
             std::tuple_size<std::decay_t<Tuple>>::value>{});
   }
 
   template <
       class F, class Tuple,
       std::enable_if_t<is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
   decltype(auto) transform(F && func, Tuple && tuple) {
     return details::transform_named_impl(
         std::forward<F>(func), std::forward<Tuple>(tuple),
         std::make_index_sequence<
             std::tuple_size<typename std::decay_t<Tuple>::parent>::value>{});
   }
 
-  namespace details {
-    template <class Tuple, std::size_t... Is>
-    decltype(auto) flatten(Tuple && tuples,
-                           std::index_sequence<Is...> /*unused*/) {
-      return std::tuple_cat(std::get<Is>(tuples)...);
-    }
-  } // namespace details
-
+  /* ------------------------------------------------------------------------ */
   template <class Tuple> decltype(auto) flatten(Tuple && tuples) {
     return details::flatten(std::forward<Tuple>(tuples),
                             std::make_index_sequence<
                                 std::tuple_size<std::decay_t<Tuple>>::value>());
   }
 
+  /* ------------------------------------------------------------------------ */
   template <typename... Ts> struct cat {
     using type = decltype(std::tuple_cat(std::declval<Ts>()...));
   };
 
   template <typename... T> using cat_t = typename cat<T...>::type;
 
-  /* --------------------------------------------------------------------------
-   */
+  /* ------------------------------------------------------------------------ */
   template <template <typename> class Pred, typename... Ts> struct filter {};
 
   template <template <typename> class Pred, typename T>
   struct filter<Pred, std::tuple<T>> {
     using type =
         std::conditional_t<Pred<T>::value, std::tuple<T>, std::tuple<>>;
   };
 
   template <template <typename> class Pred, typename T, typename... Ts>
   struct filter<Pred, std::tuple<T, Ts...>> {
     using type = cat_t<typename filter<Pred, std::tuple<T>>::type,
                        typename filter<Pred, std::tuple<Ts...>>::type>;
   };
 
   template <template <typename> class Pred, typename... Ts>
   using filter_t = typename filter<Pred, Ts...>::type;
 
+  /* ------------------------------------------------------------------------ */
+  template <class Tuple, class... Vals,
+            std::enable_if_t<not is_named_tuple<std::decay_t<Tuple>>::value> * =
+                nullptr>
+  decltype(auto) append(Tuple && tuple, Vals &&... vals) {
+    return details::append_impl(
+        std::make_index_sequence<std::tuple_size<Tuple>::value>{},
+        std::forward<decltype(tuple)>(tuple), std::forward<Vals>(vals)...);
+  }
+
+  template <
+      class Tuple, class... Vals,
+      std::enable_if_t<is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
+  decltype(auto) append(Tuple && tuple, Vals &&... vals) {
+    return details::append_named_impl(
+        std::make_index_sequence<
+            std::tuple_size<typename std::decay_t<Tuple>::parent>::value>{},
+        std::forward<decltype(tuple)>(tuple), std::forward<Vals>(vals)...);
+  }
+
+  /* ------------------------------------------------------------------------ */
+  template <size_t nth, class Tuple,
+            std::enable_if_t<not is_named_tuple<std::decay_t<Tuple>>::value> * =
+                nullptr>
+  decltype(auto) remove(Tuple && tuple) {
+    return details::remove_impl<nth>(
+        std::make_index_sequence<nth>{},
+        std::make_index_sequence<std::tuple_size<Tuple>::value - nth - 1>{},
+        std::forward<Tuple>(tuple));
+  }
+
+  template <
+      size_t HashCode, class Tuple,
+      std::enable_if_t<is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
+  decltype(auto) remove(Tuple && tuple) {
+    constexpr auto nth =
+        tuple.template get_element_index(tuple::get<HashCode>());
+    return details::remove_named_impl<nth>(
+        std::make_index_sequence<nth>{},
+        std::make_index_sequence<
+            std::tuple_size<typename std::decay_t<Tuple>::parent>::value - nth -
+            1>{},
+        std::forward<Tuple>(tuple));
+  }
+
+  template <size_t nth, class Tuple, class Value,
+            std::enable_if_t<not is_named_tuple<std::decay_t<Tuple>>::value> * =
+                nullptr>
+  decltype(auto) replace(Tuple && tuple, Value && value) {
+    return details::replace_impl<nth>(
+        std::make_index_sequence<nth>{},
+        std::make_index_sequence<std::tuple_size<Tuple>::value - nth - 1>{},
+        std::forward<Tuple>(tuple), std::forward<Value>(value));
+  }
+
+  template <
+      size_t HashCode, class Tuple, class Value,
+      std::enable_if_t<is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
+  decltype(auto) replace(Tuple && tuple, Value && value) {
+    constexpr auto nth =
+        tuple.template get_element_index(tuple::get<HashCode>());
+    return details::replace_named_impl<nth, HashCode>(
+        std::make_index_sequence<nth>{},
+        std::make_index_sequence<
+            std::tuple_size<typename std::decay_t<Tuple>::parent>::value - nth -
+            1>{},
+        std::forward<Tuple>(tuple), std::forward<Value>(value));
+  }
 } // namespace tuple
 
 } // namespace AKANTU_ITERATORS_NAMESPACE
 
+namespace aka {
+template <typename tag, typename type_>
+struct size_type<AKANTU_ITERATORS_NAMESPACE::tuple::named_tag<tag, type_>> {
+  using type = typename std::decay_t<type_>::size_type;
+};
+} // namespace aka
+
 /* -------------------------------------------------------------------------- */
 #include <iterator>
 /* -------------------------------------------------------------------------- */
 
 namespace std {
 template <typename tag, typename type>
 struct iterator_traits<
     ::AKANTU_ITERATORS_NAMESPACE::tuple::named_tag<tag, type>> {
   using iterator_category = typename type::iterator_category;
   using value_type = typename type::value_type;
   using difference_type = typename type::difference_type;
   using pointer = typename type::pointer;
   using reference = typename type::reference;
 };
 } // namespace std
 #endif /* AKANTU_AKA_TUPLE_TOOLS_HH */
diff --git a/third-party/akantu_iterators/include/iterators/aka_arange_iterator.hh b/third-party/akantu_iterators/include/iterators/aka_arange_iterator.hh
index fc02ecbbd..351a208b4 100644
--- a/third-party/akantu_iterators/include/iterators/aka_arange_iterator.hh
+++ b/third-party/akantu_iterators/include/iterators/aka_arange_iterator.hh
@@ -1,177 +1,179 @@
 /**
  * @file   aka_arange_iterator.hh
  *
  * @author Nicolas Richart
  *
  * @date creation  jeu déc 12 2019
  *
  * @brief implementation of arange
  *
  * @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 <cstddef>
 #include <iterator>
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 #ifndef AKA_ARANGE_ITERATOR_HH
 #define AKA_ARANGE_ITERATOR_HH
 
 #ifndef AKANTU_ITERATORS_NAMESPACE
 #define AKANTU_ITERATORS_NAMESPACE akantu
 #endif
 
 namespace AKANTU_ITERATORS_NAMESPACE {
 
 namespace containers {
 
   template <class Iterator> class Range {
   public:
+    using size_type = std::size_t;//typename std::iterator_traits<Iterator>::difference_type;
     using iterator = Iterator;
     // ugly trick
     using const_iterator = Iterator;
 
     explicit Range(Iterator && it1, Iterator && it2)
         : iterators(std::forward<Iterator>(it1), std::forward<Iterator>(it2)) {}
 
     decltype(auto) begin() const { return std::get<0>(iterators); }
     decltype(auto) begin() { return std::get<0>(iterators); }
 
     decltype(auto) end() const { return std::get<1>(iterators); }
     decltype(auto) end() { return std::get<1>(iterators); }
 
   private:
     std::tuple<Iterator, Iterator> iterators;
   };
 } // namespace containers
 
 template <class Iterator>
 decltype(auto) range(Iterator && it1, Iterator && it2) {
   return containers::Range<Iterator>(std::forward<Iterator>(it1),
                                      std::forward<Iterator>(it2));
 }
 /* -------------------------------------------------------------------------- */
 /* Arange                                                                     */
 /* -------------------------------------------------------------------------- */
 namespace iterators {
   template <class T> class ArangeIterator {
   public:
     using value_type = T;
     using pointer = T *;
     using reference = T &;
     using difference_type = size_t;
     using iterator_category = std::forward_iterator_tag;
 
     constexpr ArangeIterator(T value, T step) : value(value), step(step) {}
     constexpr ArangeIterator(const ArangeIterator &) = default;
 
-    constexpr ArangeIterator & operator++() {
+    constexpr auto operator++() -> ArangeIterator & {
       value += step;
       return *this;
     }
 
-    constexpr T operator*() const { return value; }
+    constexpr auto operator*() const -> T { return value; }
 
-    constexpr bool operator==(const ArangeIterator & other) const {
+    constexpr auto operator==(const ArangeIterator & other) const -> bool {
       return (value == other.value) and (step == other.step);
     }
 
-    constexpr bool operator!=(const ArangeIterator & other) const {
+    constexpr auto operator!=(const ArangeIterator & other) const -> bool {
       return not operator==(other);
     }
 
   private:
     T value{0};
     const T step{1};
   };
 } // namespace iterators
 
 namespace containers {
   template <class T> class ArangeContainer {
   public:
     using iterator = iterators::ArangeIterator<T>;
     using const_iterator = iterators::ArangeIterator<T>;
+    using size_type = T;
 
     constexpr ArangeContainer(T start, T stop, T step = 1)
         : start(start), stop((stop - start) % step == 0
                                  ? stop
                                  : start + (1 + (stop - start) / step) * step),
           step(step) {}
     explicit constexpr ArangeContainer(T stop) : ArangeContainer(0, stop, 1) {}
 
-    constexpr T operator[](std::size_t i) {
+    constexpr auto operator[](std::size_t i) -> T {
       T val = start + i * step;
       assert(val < stop && "i is out of range");
       return val;
     }
 
-    constexpr T size() { return (stop - start) / step; }
+    constexpr auto size() -> T { return (stop - start) / step; }
 
-    constexpr iterator begin() { return iterator(start, step); }
-    constexpr iterator end() { return iterator(stop, step); }
+    constexpr auto begin() -> iterator { return iterator(start, step); }
+    constexpr auto end() -> iterator { return iterator(stop, step); }
 
   private:
     const T start{0}, stop{0}, step{1};
   };
 } // namespace containers
 
 template <class T,
           typename = std::enable_if_t<std::is_integral<std::decay_t<T>>::value>>
 inline decltype(auto) arange(const T & stop) {
   return containers::ArangeContainer<T>(stop);
 }
 
 template <class T1, class T2,
           typename = std::enable_if_t<
               std::is_integral<std::common_type_t<T1, T2>>::value>>
 inline constexpr decltype(auto) arange(const T1 & start, const T2 & stop) {
   return containers::ArangeContainer<std::common_type_t<T1, T2>>(start, stop);
 }
 
 template <class T1, class T2, class T3,
           typename = std::enable_if_t<
               std::is_integral<std::common_type_t<T1, T2, T3>>::value>>
 inline constexpr decltype(auto) arange(const T1 & start, const T2 & stop,
                                        const T3 & step) {
   return containers::ArangeContainer<std::common_type_t<T1, T2, T3>>(
       start, stop, step);
 }
 } // namespace AKANTU_ITERATORS_NAMESPACE
 
 namespace std {
 template <class T>
 struct iterator_traits<
     ::AKANTU_ITERATORS_NAMESPACE::iterators::ArangeIterator<T>> {
 private:
   using iterator_type =
       typename ::AKANTU_ITERATORS_NAMESPACE::iterators::ArangeIterator<T>;
 
 public:
   using iterator_category = typename iterator_type::iterator_category;
   using value_type = typename iterator_type::value_type;
   using difference_type = typename iterator_type::difference_type;
   using pointer = typename iterator_type::pointer;
   using reference = typename iterator_type::reference;
 };
 
 } // namespace std
 
 #endif /* AKA_ARANGE_ITERATOR_HH */
diff --git a/third-party/akantu_iterators/include/iterators/aka_enumerate_iterator.hh b/third-party/akantu_iterators/include/iterators/aka_enumerate_iterator.hh
index d4ca5b94b..7630b8136 100644
--- a/third-party/akantu_iterators/include/iterators/aka_enumerate_iterator.hh
+++ b/third-party/akantu_iterators/include/iterators/aka_enumerate_iterator.hh
@@ -1,155 +1,162 @@
 /**
  * @file   aka_enumerate_iterator.hh
  *
  * @author Nicolas Richart
  *
  * @date creation  jeu déc 12 2019
  *
  * @brief implementation of enumerate.
  *
  * @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 "iterators/aka_zip_iterator.hh"
 /* -------------------------------------------------------------------------- */
 #include <iterator>
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 #ifndef AKA_ENUMERATE_ITERATOR_HH
 #define AKA_ENUMERATE_ITERATOR_HH
 
 #ifndef AKANTU_ITERATORS_NAMESPACE
 #define AKANTU_ITERATORS_NAMESPACE akantu
 #endif
 
 namespace AKANTU_ITERATORS_NAMESPACE {
 
 /* -------------------------------------------------------------------------- */
 namespace iterators {
-  template <class Iterator> class EnumerateIterator {
+  template <class Iterator, class size_type_> class EnumerateIterator {
   public:
     using value_type =
         std::tuple<std::size_t,
                    typename std::iterator_traits<Iterator>::value_type>;
-    using difference_type = std::size_t;
+    using difference_type = size_type_;
     using pointer =
-        std::tuple<std::size_t,
+        std::tuple<size_type_,
                    typename std::iterator_traits<Iterator>::pointer>;
     using reference =
-        std::tuple<std::size_t,
+        std::tuple<size_type_,
                    typename std::iterator_traits<Iterator>::reference>;
     using iterator_category = std::input_iterator_tag;
 
   public:
     explicit EnumerateIterator(Iterator && iterator) : iterator(iterator) {}
 
     // input iterator ++it
     EnumerateIterator & operator++() {
       ++iterator;
       ++index;
       return *this;
     }
 
     // input iterator it++
     EnumerateIterator operator++(int) {
       auto cpy = *this;
       this->operator++();
       return cpy;
     }
 
     // input iterator it != other_it
     bool operator!=(const EnumerateIterator & other) const {
       return iterator != other.iterator;
     }
 
     // input iterator dereference *it
     decltype(auto) operator*() {
       return std::tuple_cat(std::make_tuple(index), *iterator);
     }
 
     bool operator==(const EnumerateIterator & other) const {
       return not this->operator!=(other);
     }
 
   private:
     Iterator iterator;
-    size_t index{0};
+    size_type_ index{0};
   };
 
-  template <class Iterator>
-  inline constexpr decltype(auto) enumerate(Iterator && iterator) {
-    return EnumerateIterator<Iterator>(std::forward<Iterator>(iterator));
+  template <class Iterator, class size_type_ = std::size_t>
+  inline constexpr decltype(auto) enumerate(Iterator && iterator,
+                                            size_type_ /*size*/) {
+    return EnumerateIterator<Iterator, size_type_>(
+        std::forward<Iterator>(iterator));
   }
 
 } // namespace iterators
 
 namespace containers {
   template <class... Containers> class EnumerateContainer {
+    using ZipContainer_t = ZipContainer<Containers...>;
+    using size_type = typename ZipContainer_t::size_type;
+
   public:
     explicit EnumerateContainer(Containers &&... containers)
         : zip_container(std::forward<Containers>(containers)...) {}
 
     decltype(auto) begin() {
-      return iterators::enumerate(zip_container.begin());
+      return iterators::enumerate(zip_container.begin(), size_type{});
     }
 
     decltype(auto) begin() const {
-      return iterators::enumerate(zip_container.begin());
+      return iterators::enumerate(zip_container.begin(), size_type{});
     }
 
-    decltype(auto) end() { return iterators::enumerate(zip_container.end()); }
+    decltype(auto) end() {
+      return iterators::enumerate(zip_container.end(), size_type{});
+    }
 
     decltype(auto) end() const {
-      return iterators::enumerate(zip_container.end());
+      return iterators::enumerate(zip_container.end(), size_type{});
     }
 
   private:
-    ZipContainer<Containers...> zip_container;
+    ZipContainer_t zip_container;
   };
 } // namespace containers
 
 template <class... Container>
 inline constexpr decltype(auto) enumerate(Container &&... container) {
   return containers::EnumerateContainer<Container...>(
       std::forward<Container>(container)...);
 }
 
 } // namespace AKANTU_ITERATORS_NAMESPACE
 
 namespace std {
-template <class Iterator>
-struct iterator_traits<
-    ::AKANTU_ITERATORS_NAMESPACE::iterators::EnumerateIterator<Iterator>> {
+template <class Iterator, class size_type>
+struct iterator_traits<::AKANTU_ITERATORS_NAMESPACE::iterators::
+                           EnumerateIterator<Iterator, size_type>> {
 private:
   using iterator_type =
       typename ::AKANTU_ITERATORS_NAMESPACE::iterators::EnumerateIterator<
-          Iterator>;
+          Iterator, size_type>;
 
 public:
   using iterator_category = typename iterator_type::iterator_category;
   using value_type = typename iterator_type::value_type;
   using difference_type = typename iterator_type::difference_type;
   using pointer = typename iterator_type::pointer;
   using reference = typename iterator_type::reference;
 };
 } // namespace std
 
 #endif /* AKA_ENUMERATE_ITERATOR_HH */
diff --git a/third-party/akantu_iterators/include/iterators/aka_filter_iterator.hh b/third-party/akantu_iterators/include/iterators/aka_filter_iterator.hh
index 5e7ef559e..5c72c177e 100644
--- a/third-party/akantu_iterators/include/iterators/aka_filter_iterator.hh
+++ b/third-party/akantu_iterators/include/iterators/aka_filter_iterator.hh
@@ -1,160 +1,162 @@
 /**
  * @file   aka_filter_iterator.hh
  *
  * @author Nicolas Richart
  *
  * @date creation  jeu déc 12 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 <iterator>
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 #ifndef AKA_FILTER_ITERATOR_H
 #define AKA_FILTER_ITERATOR_H
 
 #ifndef AKANTU_ITERATORS_NAMESPACE
 #define AKANTU_ITERATORS_NAMESPACE akantu
 #endif
 
 namespace AKANTU_ITERATORS_NAMESPACE {
 
 /* -------------------------------------------------------------------------- */
 namespace iterators {
   template <class filter_iterator_t, class container_iterator_t>
   class FilterIterator {
   public:
     using value_type =
         decltype(std::declval<container_iterator_t>().operator[](0));
     using difference_type = typename std::decay_t<filter_iterator_t>::difference_type;
     using pointer = std::decay_t<value_type> *;
     using reference = value_type &;
     using iterator_category = typename std::decay_t<filter_iterator_t>::iterator_category;
 
     FilterIterator(filter_iterator_t filter_it,
                    container_iterator_t && container_begin)
         : filter_it(std::move(filter_it)),
           container_begin(std::move(container_begin)) {}
 
     FilterIterator(const FilterIterator &) = default;
 
     FilterIterator & operator++() {
       ++filter_it;
       return *this;
     }
 
     decltype(auto) operator*() {
       auto container_it = this->container_begin + *this->filter_it;
       return *container_it;
     }
 
     decltype(auto) operator*() const {
       auto container_it = this->container_begin + *this->filter_it;
       return *container_it;
     }
 
     bool operator==(const FilterIterator & other) const {
       return (filter_it == other.filter_it) and
              (container_begin == other.container_begin);
     }
 
     bool operator!=(const FilterIterator & other) const {
       return filter_it != other.filter_it;
     }
 
   private:
     filter_iterator_t filter_it;
     container_iterator_t container_begin;
   };
 
   template <class filter_iterator_t, class container_iterator_t>
   decltype(auto) make_filter_iterator(filter_iterator_t && filter_it,
                                       container_iterator_t && container_begin) {
     return FilterIterator<filter_iterator_t, container_iterator_t>(
         std::forward<filter_iterator_t>(filter_it),
         std::forward<container_iterator_t>(container_begin));
   }
 } // namespace iterators
 
 namespace containers {
   template <class filter_t, class Container> class FilterAdaptor {
   public:
+    using size_type = typename std::decay_t<Container>::size_type;
+
     FilterAdaptor(filter_t && filter, Container && container)
         : filter(std::forward<filter_t>(filter)),
           container(std::forward<Container>(container)) {
       static_assert(
           std::is_same<typename decltype(container.begin())::iterator_category,
                        std::random_access_iterator_tag>::value,
           "Containers must all have random iterators");
     }
 
     decltype(auto) begin() const {
       return iterators::make_filter_iterator(filter.begin(), container.begin());
     }
 
     decltype(auto) begin() {
       return iterators::make_filter_iterator(filter.begin(), container.begin());
     }
 
     decltype(auto) end() const {
       return iterators::make_filter_iterator(filter.end(), container.begin());
     }
 
     decltype(auto) end() {
       return iterators::make_filter_iterator(filter.end(), container.begin());
     }
 
   private:
     filter_t filter;
     Container container;
   };
 } // namespace containers
 
 template <class filter_t, class Container>
 decltype(auto) filter(filter_t && filter, Container && container) {
   return containers::FilterAdaptor<filter_t, Container>(
       std::forward<filter_t>(filter), std::forward<Container>(container));
 }
 
 } // namespace AKANTU_ITERATORS_NAMESPACE
 
 namespace std {
 template <class filter_iterator_t, class container_iterator_t>
 struct iterator_traits<::AKANTU_ITERATORS_NAMESPACE::iterators::FilterIterator<
     filter_iterator_t, container_iterator_t>> {
 private:
   using iterator_type =
       typename ::AKANTU_ITERATORS_NAMESPACE::iterators::FilterIterator<
           filter_iterator_t, container_iterator_t>;
 
 public:
   using iterator_category = typename iterator_type::iterator_category;
   using value_type = typename iterator_type::value_type;
   using difference_type = typename iterator_type::difference_type;
   using pointer = typename iterator_type::pointer;
   using reference = typename iterator_type::reference;
 };
 } // namespace std
 
 #endif // AKA_FILTER_ITERATOR_H
diff --git a/third-party/akantu_iterators/include/iterators/aka_transform_iterator.hh b/third-party/akantu_iterators/include/iterators/aka_transform_iterator.hh
index 7830f74ef..753f65f5e 100644
--- a/third-party/akantu_iterators/include/iterators/aka_transform_iterator.hh
+++ b/third-party/akantu_iterators/include/iterators/aka_transform_iterator.hh
@@ -1,156 +1,163 @@
 /**
  * @file   aka_transform_iterator.hh
  *
  * @author Nicolas Richart
  *
  * @date creation  jeu déc 12 2019
  *
  * @brief transform adaptors
  *
  * @section LICENSE
  *
  * Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * akantu-iterators 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-iterators 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-iterators. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 /* -------------------------------------------------------------------------- */
 #include "aka_compatibilty_with_cpp_standard.hh"
 /* -------------------------------------------------------------------------- */
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 #ifndef AKA_TRANSFORM_ITERATOR_H
 #define AKA_TRANSFORM_ITERATOR_H
 
 #ifndef AKANTU_ITERATORS_NAMESPACE
 #define AKANTU_ITERATORS_NAMESPACE akantu
 #endif
 
 namespace AKANTU_ITERATORS_NAMESPACE {
 
 namespace iterators {
   template <class iterator_t, class operator_t>
   class transform_adaptor_iterator {
   public:
     using value_type = decltype(std::declval<operator_t>()(
         std::declval<typename iterator_t::value_type>()));
     using difference_type = typename iterator_t::difference_type;
     using pointer = std::decay_t<value_type> *;
     using reference = value_type &;
     using iterator_category = typename iterator_t::iterator_category;
 
     transform_adaptor_iterator(iterator_t it, operator_t op)
         : it(std::move(it)), op(op) {}
     transform_adaptor_iterator(const transform_adaptor_iterator &) = default;
 
     transform_adaptor_iterator & operator++() {
       ++it;
       return *this;
     }
 
     decltype(auto) operator*() { return op(std::forward<decltype(*it)>(*it)); }
 
     bool operator==(const transform_adaptor_iterator & other) const {
       return (it == other.it);
     }
 
     bool operator!=(const transform_adaptor_iterator & other) const {
       return not operator==(other);
     }
 
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::random_access_iterator_tag>::value> * = nullptr>
     difference_type operator-(const transform_adaptor_iterator & other) {
       return other - *this;
     }
 
   private:
     iterator_t it;
     operator_t op;
   };
 
   template <class iterator_t, class operator_t>
   decltype(auto) make_transform_adaptor_iterator(iterator_t it, operator_t op) {
     return transform_adaptor_iterator<iterator_t, operator_t>(
         it, std::forward<operator_t>(op));
   }
 
 } // namespace iterators
 
 namespace containers {
   template <class container_t, class operator_t>
   class TransformIteratorAdaptor {
   public:
     // using const_iterator = typename
     // std::decay_t<container_t>::const_iterator; using iterator = typename
     // std::decay_t<container_t>::iterator;
+    using size_type = typename std::decay_t<container_t>::size_type;
 
     TransformIteratorAdaptor(container_t && cont, operator_t op)
         : cont(std::forward<container_t>(cont)),
           op(std::forward<operator_t>(op)) {}
 
     decltype(auto) begin() const {
       return iterators::make_transform_adaptor_iterator(cont.begin(), op);
     }
     decltype(auto) begin() {
       return iterators::make_transform_adaptor_iterator(cont.begin(), op);
     }
 
     decltype(auto) end() const {
       return iterators::make_transform_adaptor_iterator(cont.end(), op);
     }
     decltype(auto) end() {
       return iterators::make_transform_adaptor_iterator(cont.end(), op);
     }
 
   private:
     container_t cont;
     operator_t op;
   };
 } // namespace containers
 
 template <class container_t, class operator_t>
 decltype(auto) make_transform_adaptor(container_t && cont, operator_t && op) {
   return containers::TransformIteratorAdaptor<container_t, operator_t>(
       std::forward<container_t>(cont), std::forward<operator_t>(op));
 }
 
 template <class container_t>
 decltype(auto) make_keys_adaptor(container_t && cont) {
   return make_transform_adaptor(
       std::forward<container_t>(cont),
       [](auto && pair) -> const auto & { return pair.first; });
 }
 
 template <class container_t>
 decltype(auto) make_values_adaptor(container_t && cont) {
   return make_transform_adaptor(
       std::forward<container_t>(cont),
       [](auto && pair) -> auto & { return pair.second; });
 }
 
 template <class container_t>
 decltype(auto) make_dereference_adaptor(container_t && cont) {
   return make_transform_adaptor(
       std::forward<container_t>(cont),
       [](auto && value) -> decltype(*value) { return *value; });
 }
 
+template <typename Type, typename Size>
+decltype(auto) broadcast(Type && data, Size size) {
+  return make_transform_adaptor(
+      arange(size), [data](auto && value) -> decltype(auto) { return data; });
+}
+
 } // namespace AKANTU_ITERATORS_NAMESPACE
 
 #endif // AKA_TRANSFORM_ITERATOR_H
diff --git a/third-party/akantu_iterators/include/iterators/aka_zip_iterator.hh b/third-party/akantu_iterators/include/iterators/aka_zip_iterator.hh
index 9b8e7b589..c10517bae 100644
--- a/third-party/akantu_iterators/include/iterators/aka_zip_iterator.hh
+++ b/third-party/akantu_iterators/include/iterators/aka_zip_iterator.hh
@@ -1,304 +1,445 @@
 /**
  * @file   aka_zip_iterator.hh
  *
  * @author Nicolas Richart
  *
  * @date creation  jeu déc 12 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_compatibilty_with_cpp_standard.hh"
 #include "aka_tuple_tools.hh"
 /* -------------------------------------------------------------------------- */
 #include <iterator>
 #include <tuple>
 #include <utility>
 /* -------------------------------------------------------------------------- */
 
 #ifndef AKA_ZIP_ITERATOR_HH
 #define AKA_ZIP_ITERATOR_HH
 
 #ifndef AKANTU_ITERATORS_NAMESPACE
 #define AKANTU_ITERATORS_NAMESPACE akantu
 #endif
 
 namespace AKANTU_ITERATORS_NAMESPACE {
 
 /* -------------------------------------------------------------------------- */
 namespace iterators {
   namespace details {
     template <bool enable> struct CopyAssignmentEnabler {};
     template <> struct CopyAssignmentEnabler<false> {
       CopyAssignmentEnabler() = default;
       CopyAssignmentEnabler(const CopyAssignmentEnabler &) = default;
       CopyAssignmentEnabler(CopyAssignmentEnabler &&) = default;
       CopyAssignmentEnabler & operator=(const CopyAssignmentEnabler &) = delete;
       CopyAssignmentEnabler & operator=(CopyAssignmentEnabler &&) = default;
     };
 
     template <bool enable> struct MoveAssignmentEnabler {};
     template <> struct MoveAssignmentEnabler<false> {
       MoveAssignmentEnabler() = default;
       MoveAssignmentEnabler(const MoveAssignmentEnabler &) = default;
       MoveAssignmentEnabler(MoveAssignmentEnabler &&) = default;
       MoveAssignmentEnabler & operator=(const MoveAssignmentEnabler &) = delete;
       MoveAssignmentEnabler & operator=(MoveAssignmentEnabler &&) = default;
     };
 
   } // namespace details
 
   /* ------------------------------------------------------------------------ */
   template <template <class...> class Tuple, class... Iterators>
   class ZipIterator_
       : public details::CopyAssignmentEnabler<
             aka::conjunction<std::is_copy_assignable<Iterators>...,
                              std::is_copy_constructible<Iterators>...>::value>,
         public details::MoveAssignmentEnabler<
             aka::conjunction<std::is_move_assignable<Iterators>...,
                              std::is_move_constructible<Iterators>...>::value> {
   private:
     using tuple_t = Tuple<Iterators...>;
 
   public:
     using value_type =
         Tuple<typename std::iterator_traits<Iterators>::value_type...>;
     using difference_type = std::common_type_t<
         typename std::iterator_traits<Iterators>::difference_type...>;
-    using pointer =
-        Tuple<typename std::iterator_traits<Iterators>::pointer...>;
+    using pointer = Tuple<typename std::iterator_traits<Iterators>::pointer...>;
     using reference =
         Tuple<typename std::iterator_traits<Iterators>::reference...>;
     using iterator_category = // std::input_iterator_tag;
         std::common_type_t<
             typename std::iterator_traits<Iterators>::iterator_category...>;
 
-//    using nb_iterators = sizeof...(Iterators);
+    //    using nb_iterators = sizeof...(Iterators);
 
   public:
     explicit ZipIterator_(tuple_t iterators)
         : iterators(std::move(iterators)) {}
 
     /* ---------------------------------------------------------------------- */
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::bidirectional_iterator_tag>::value> * = nullptr>
     ZipIterator_ & operator--() {
       tuple::foreach ([](auto && it) { --it; }, iterators);
       return *this;
     }
 
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::bidirectional_iterator_tag>::value> * = nullptr>
     ZipIterator_ operator--(int) {
       auto cpy = *this;
       this->operator--();
       return cpy;
     }
 
     // input iterator ++it
     ZipIterator_ & operator++() {
       tuple::foreach ([](auto && it) { ++it; }, iterators);
       return *this;
     }
 
     // input iterator it++
     ZipIterator_ operator++(int) {
       auto cpy = *this;
       this->operator++();
       return cpy;
     }
 
     // input iterator it != other_it
     bool operator!=(const ZipIterator_ & other) const {
       // return tuple::are_not_equal(iterators, other.iterators);
       return std::get<0>(iterators) !=
              std::get<0>(other.iterators); // helps the compiler to optimize
     }
 
     // input iterator dereference *it
     decltype(auto) operator*() {
       return tuple::transform([](auto && it) -> decltype(auto) { return *it; },
                               iterators);
     }
 
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::random_access_iterator_tag>::value> * = nullptr>
     difference_type operator-(const ZipIterator_ & other) {
       return std::get<0>(this->iterators) - std::get<0>(other.iterators);
     }
 
     // random iterator it[idx]
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::random_access_iterator_tag>::value> * = nullptr>
     decltype(auto) operator[](std::size_t idx) {
       return tuple::transform(
           [idx](auto && it) -> decltype(auto) { return it[idx]; }, iterators);
     }
 
     // random iterator it + n
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::random_access_iterator_tag>::value> * = nullptr>
     decltype(auto) operator+(std::size_t n) {
       return ZipIterator_(std::forward<tuple_t>(tuple::transform(
           [n](auto && it) -> decltype(auto) { return it + n; }, iterators)));
     }
 
     // random iterator it - n
     template <class iterator_category_ = iterator_category,
               std::enable_if_t<aka::is_iterator_category_at_least<
                   iterator_category_,
                   std::random_access_iterator_tag>::value> * = nullptr>
     decltype(auto) operator-(std::size_t n) {
       return ZipIterator_(std::forward<tuple_t>(tuple::transform(
           [n](auto && it) -> decltype(auto) { return it - n; }, iterators)));
     }
 
     template <
         class iterator_category_ = iterator_category,
         std::enable_if_t<aka::is_iterator_category_at_least<
             iterator_category_, std::forward_iterator_tag>::value> * = nullptr>
     bool operator==(const ZipIterator_ & other) const {
       return not tuple::are_not_equal(iterators, other.iterators);
     }
 
   private:
     tuple_t iterators;
   };
 
   template <class... Iterators>
   using ZipIterator = ZipIterator_<std::tuple, Iterators...>;
 
   template <class... Iterators>
   using NamedZipIterator = ZipIterator_<tuple::named_tuple, Iterators...>;
 } // namespace iterators
 
 /* -------------------------------------------------------------------------- */
 template <class... Iterators>
 decltype(auto) zip_iterator(std::tuple<Iterators...> && iterators_tuple) {
   auto zip = iterators::ZipIterator<Iterators...>(
       std::forward<decltype(iterators_tuple)>(iterators_tuple));
   return zip;
 }
 
 template <class... Iterators>
 decltype(auto)
 zip_iterator(tuple::named_tuple<Iterators...> && iterators_tuple) {
   auto zip = iterators::NamedZipIterator<Iterators...>(
       std::forward<decltype(iterators_tuple)>(iterators_tuple));
   return zip;
 }
 
 /* -------------------------------------------------------------------------- */
 namespace containers {
   template <template <class...> class Tuple, class... Containers>
   class ZipContainer_ {
+  public:
     using containers_t = Tuple<Containers...>;
+    using size_type = std::common_type_t<aka::size_type_t<Containers>...>;
 
-  public:
     explicit ZipContainer_(Containers &&... containers)
         : containers(std::forward<Containers>(containers)...) {}
 
     decltype(auto) begin() const {
       return zip_iterator(
           tuple::transform([](auto && c) { return c.begin(); },
                            std::forward<containers_t>(containers)));
     }
 
     decltype(auto) end() const {
       return zip_iterator(
           tuple::transform([](auto && c) { return c.end(); },
                            std::forward<containers_t>(containers)));
     }
 
     decltype(auto) begin() {
       return zip_iterator(
           tuple::transform([](auto && c) { return c.begin(); },
                            std::forward<containers_t>(containers)));
     }
 
     decltype(auto) end() {
       return zip_iterator(
           tuple::transform([](auto && c) { return c.end(); },
                            std::forward<containers_t>(containers)));
     }
 
+  protected:
+    template <class... OtherContainers, std::size_t... Is>
+    decltype(auto) append(std::index_sequence<Is...> && /*unused*/,
+                          OtherContainers &&... other_containers) {
+      return ZipContainer_<Tuple, Containers..., OtherContainers...>(
+          std::get<Is>(std::forward<containers_t>(containers))...,
+          std::forward<OtherContainers>(
+              std::forward<OtherContainers>(other_containers))...);
+    }
+
+    template <class OtherContainers, std::size_t... Is>
+    decltype(auto) extend(std::index_sequence<Is...> && /*unused*/,
+                          OtherContainers && other_containers) {
+      return append(
+          std::make_index_sequence<sizeof...(Containers)>{},
+          std::get<Is>(std::forward<containers_t>(other_containers))...);
+    }
+
+    // template <std::size_t nth, std::size_t... Is_before,
+    //           std::size_t... Is_after>
+    // decltype(auto) remove(std::index_sequence<Is_before...> && /*unused*/,
+    //                       std::index_sequence<Is_after...> && /*unused*/) {
+    //   using tuple::tuple_element_t;
+    //   return ZipContainer_<
+    //       Tuple, tuple_element_t<Is_before, containers_t>...,
+    //       tuple_element_t<Is_after + nth + 1, containers_t>...>(
+    //       std::forward<tuple_element_t<Is_before, containers_t>>(
+    //           std::get<Is_before>(containers))...,
+    //       std::forward<tuple_element_t<Is_after + nth + 1, containers_t>>(
+    //           std::get<Is_after + nth + 1>(containers))...);
+    // }
+
+  public:
+    template <class... OtherContainers>
+    decltype(auto) append(OtherContainers &&... other_containers) {
+      static_assert(
+          tuple::is_named_tuple<containers_t>::value ==
+              aka::conjunction<
+                  tuple::is_named_tag<std::decay_t<OtherContainers>>...>::value,
+          "Cannot append named tag in non named zip or vice "
+          "versa");
+      return append(std::make_index_sequence<sizeof...(Containers)>{},
+                    std::forward<OtherContainers>(other_containers)...);
+    }
+
+    template <class... OtherContainers>
+    decltype(auto)
+    extend(const ZipContainer_<Tuple, OtherContainers...> & other) {
+      return extend(std::make_index_sequence<sizeof...(OtherContainers)>{},
+                    other.containers);
+    }
+
+    // template <size_t Tag, class Containers_t = containers_t,
+    //           std::enable_if_t<tuple::is_named_tuple<Containers_t>::value> *
+    //           =
+    //               nullptr>
+    // decltype(auto) remove() {
+    //   constexpr auto nth =
+    //       containers.template get_element_index(tuple::get<Tag>());
+    //   return remove<nth>(
+    //       std::make_index_sequence<nth>{},
+    //       std::make_index_sequence<sizeof...(Containers) - nth - 1>{});
+    // }
+
+    template <size_t Tag> decltype(auto) remove();
+
+    template <size_t Tag, class Container>
+    decltype(auto) replace(Container && cont);
+
   private:
+    template <template <class...> class OtherTuple, class... OtherContainers>
+    friend class ZipContainer_;
+
     containers_t containers;
   };
 
   template <class... Containers>
   using ZipContainer = ZipContainer_<std::tuple, Containers...>;
 
   template <class... Containers>
   using NamedZipContainer = ZipContainer_<tuple::named_tuple, Containers...>;
 } // namespace containers
 
 /* -------------------------------------------------------------------------- */
-template <class... Containers> decltype(auto) zip(Containers &&... conts) {
+template <class... Containers,
+          std::enable_if_t<not aka::conjunction<tuple::is_named_tag<
+              std::decay_t<Containers>>...>::value> * = nullptr>
+decltype(auto) zip(Containers &&... conts) {
   return containers::ZipContainer<Containers...>(
       std::forward<Containers>(conts)...);
 }
 
-template <class... NamedContainers>
-decltype(auto) named_zip(NamedContainers &&... conts) {
-  return containers::NamedZipContainer<NamedContainers...>(
-      std::forward<NamedContainers>(conts)...);
+template <class... Containers,
+          std::enable_if_t<aka::conjunction<tuple::is_named_tag<
+              std::decay_t<Containers>>...>::value> * = nullptr>
+decltype(auto) zip(Containers &&... conts) {
+  return containers::NamedZipContainer<Containers...>(
+      std::forward<Containers>(conts)...);
 }
 
 /* -------------------------------------------------------------------------- */
-template <class... zip_container_t>
-decltype(auto) make_zip_cat(zip_container_t &&... cont) {
-  return make_transform_adaptor(
-      zip(std::forward<zip_container_t>(cont)...),
-      [](auto && value) { return tuple::flatten(value); });
+template <class zip_container_1_t, class zip_container_2_t>
+decltype(auto) zip_cat(zip_container_1_t && cont1, zip_container_2_t && cont2) {
+  return std::forward<zip_container_1_t>(cont1).extend(
+      std::forward<zip_container_2_t>(cont2));
+}
+
+template <class zip_container_t, class... container_t>
+decltype(auto) zip_append(zip_container_t && zip_container,
+                          container_t &&... cont) {
+  return std::forward<zip_container_t>(zip_container)
+      .append(std::forward<container_t>(cont)...);
+}
+
+template <size_t Tag, class zip_container_t, class container_t>
+decltype(auto) zip_replace(zip_container_t && zip_container,
+                           container_t && cont) {
+  return std::forward<zip_container_t>(zip_container)
+      .template replace<Tag>(std::forward<container_t>(cont));
+}
+
+template <size_t Tag, class zip_container_t>
+decltype(auto) zip_remove(zip_container_t && zip_container) {
+  return std::forward<zip_container_t>(zip_container).template remove<Tag>();
 }
 
+namespace details {
+  template <class Tuple, std::size_t... Is,
+            std::enable_if_t<
+                tuple::is_named_tuple<std::decay_t<Tuple>>::value> * = nullptr>
+  decltype(auto)
+  make_zip_from_tuple_impl(std::index_sequence<Is...> && /*unused*/,
+                           Tuple && tuple) {
+    return zip(
+        std::get<tuple::tuple_name_tag<Is, std::decay_t<Tuple>>>() =
+            std::forward<tuple::tuple_element_t<Is, std::decay_t<Tuple>>>(
+                std::get<Is>(tuple))...);
+  }
+} // namespace details
+
+template <class Tuple, std::enable_if_t<tuple::is_named_tuple<
+                           std::decay_t<Tuple>>::value> * = nullptr>
+decltype(auto) make_zip_from_tuple(Tuple && tuple) {
+  return details::make_zip_from_tuple_impl(
+      std::make_index_sequence<
+          std::tuple_size<typename std::decay_t<Tuple>::parent>::value>{},
+      std::forward<Tuple>(tuple));
+}
+
+/* -------------------------------------------------------------------------- */
+/* -------------------------------------------------------------------------- */
+namespace containers {
+  template <template <class...> class Tuple, class... Containers>
+  template <std::size_t nth>
+  decltype(auto) ZipContainer_<Tuple, Containers...>::remove() {
+    return make_zip_from_tuple(
+        tuple::remove<nth>(std::forward<containers_t>(containers)));
+  }
+
+  template <template <class...> class Tuple, class... Containers>
+  template <std::size_t nth, class Container>
+  decltype(auto)
+  ZipContainer_<Tuple, Containers...>::replace(Container && cont) {
+    return make_zip_from_tuple(tuple::replace<nth>(
+        std::forward<containers_t>(containers), std::forward<Container>(cont)));
+  }
+
+} // namespace containers
+/* -------------------------------------------------------------------------- */
+/* -------------------------------------------------------------------------- */
+
 } // namespace AKANTU_ITERATORS_NAMESPACE
 
 namespace std {
-template <template <class ...> class Tuple, typename... Its>
+template <template <class...> class Tuple, typename... Its>
 struct iterator_traits<
     ::AKANTU_ITERATORS_NAMESPACE::iterators::ZipIterator_<Tuple, Its...>> {
 private:
   using iterator_type =
-      typename ::AKANTU_ITERATORS_NAMESPACE::iterators::ZipIterator_<Tuple, Its...>;
+      typename ::AKANTU_ITERATORS_NAMESPACE::iterators::ZipIterator_<Tuple,
+                                                                     Its...>;
 
 public:
   using iterator_category = typename iterator_type::iterator_category;
   using value_type = typename iterator_type::value_type;
   using difference_type = typename iterator_type::difference_type;
   using pointer = typename iterator_type::pointer;
   using reference = typename iterator_type::reference;
 };
 } // namespace std
 
 #endif /* AKA_ZIP_ITERATOR_HH */
diff --git a/third-party/akantu_iterators/test/test_iterators.cc b/third-party/akantu_iterators/test/test_iterators.cc
index fd42557f0..1486778ea 100644
--- a/third-party/akantu_iterators/test/test_iterators.cc
+++ b/third-party/akantu_iterators/test/test_iterators.cc
@@ -1,395 +1,486 @@
 /**
  * @file   test_zip_iterator.cc
  *
  * @author Nicolas Richart <nicolas.richart@epfl.ch>
  *
  * @date creation: Fri Jul 21 2017
  * @date last modification: Fri Dec 08 2017
  *
  * @brief  test the zip container and iterator
  *
  * @section LICENSE
  *
  * Copyright (©) 2016-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
  * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
  *
  * akantu-iterators 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-iterators 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-iterators. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 
 /* -------------------------------------------------------------------------- */
 #include "aka_iterators.hh"
 /* -------------------------------------------------------------------------- */
 #include <gtest/gtest.h>
 #include <vector>
 /* -------------------------------------------------------------------------- */
 
 using namespace aka;
 
 /* -------------------------------------------------------------------------- */
 
 // Non Trivial class that counts moves and copies
 template <class T> class A {
 public:
   A() = default;
   A(T a) : a(a){};
   A(const A & other)
       : a(other.a), copy_counter(other.copy_counter + 1),
         move_counter(other.move_counter) {}
   A & operator=(const A & other) {
     if (this != &other) {
       a = other.a;
       copy_counter = other.copy_counter + 1;
     }
     return *this;
   }
 
   A(A && other)
       : a(std::move(other.a)), copy_counter(std::move(other.copy_counter)),
         move_counter(std::move(other.move_counter) + 1) {}
 
   A & operator=(A && other) {
     if (this != &other) {
       a = std::move(other.a);
       copy_counter = std::move(other.copy_counter);
       move_counter = std::move(other.move_counter) + 1;
     }
     return *this;
   }
 
   A & operator*=(const T & b) {
     a *= b;
     return *this;
   }
 
   T a{};
   size_t copy_counter{0};
   size_t move_counter{0};
 };
 
 template <typename T> struct C {
   struct iterator {
     using reference = A<T>;
     using difference_type = void;
     using iterator_category = std::input_iterator_tag;
     using value_type = A<T>;
     using pointer = A<T> *;
 
     iterator(T pos) : pos(std::move(pos)) {}
 
     A<T> operator*() { return A<int>(pos); }
     bool operator!=(const iterator & other) const { return pos != other.pos; }
     bool operator==(const iterator & other) const { return pos == other.pos; }
     iterator & operator++() {
       ++pos;
       return *this;
     }
     T pos;
   };
 
   C(T begin_, T end_) : begin_(std::move(begin_)), end_(std::move(end_)) {}
 
   iterator begin() { return iterator(begin_); }
   iterator end() { return iterator(end_); }
 
   T begin_, end_;
 };
 
 class TestZipFixutre : public ::testing::Test {
 protected:
   void SetUp() override {
     a.reserve(size);
     b.reserve(size);
 
     for (size_t i = 0; i < size; ++i) {
       a.emplace_back(i);
       b.emplace_back(i + size);
     }
   }
 
   template <typename A, typename B>
   void check(A && a, B && b, size_t pos, size_t nb_copy, size_t nb_move) {
     EXPECT_EQ(pos, a.a);
     EXPECT_EQ(nb_copy, a.copy_counter);
     EXPECT_EQ(nb_move, a.move_counter);
 
     EXPECT_FLOAT_EQ(pos + this->size, b.a);
     EXPECT_EQ(nb_copy, b.copy_counter);
     EXPECT_EQ(nb_move, b.move_counter);
   }
 
 protected:
   size_t size{20};
   std::vector<A<int>> a{};
   std::vector<A<float>> b{};
 };
 
 TEST_F(TestZipFixutre, SimpleTest) {
   size_t i = 0;
   for (auto && pair : zip(this->a, this->b)) {
     this->check(std::get<0>(pair), std::get<1>(pair), i, 0, 0);
     ++i;
   }
 }
 
 TEST_F(TestZipFixutre, ConstTest) {
   size_t i = 0;
   const auto & ca = this->a;
   const auto & cb = this->b;
   for (auto && pair : zip(ca, cb)) {
     this->check(std::get<0>(pair), std::get<1>(pair), i, 0, 0);
     EXPECT_EQ(true,
               std::is_const<
                   std::remove_reference_t<decltype(std::get<0>(pair))>>::value);
     EXPECT_EQ(true,
               std::is_const<
                   std::remove_reference_t<decltype(std::get<1>(pair))>>::value);
     ++i;
   }
 }
 
 TEST_F(TestZipFixutre, MixteTest) {
   size_t i = 0;
   const auto & cb = this->b;
   for (auto && pair : zip(a, cb)) {
     this->check(std::get<0>(pair), std::get<1>(pair), i, 0, 0);
     EXPECT_EQ(false,
               std::is_const<
                   std::remove_reference_t<decltype(std::get<0>(pair))>>::value);
     EXPECT_EQ(true,
               std::is_const<
                   std::remove_reference_t<decltype(std::get<1>(pair))>>::value);
     ++i;
   }
 }
 
 TEST_F(TestZipFixutre, MoveTest) {
   size_t i = 0;
   for (auto && pair :
        zip(C<int>(0, this->size), C<int>(this->size, 2 * this->size))) {
     this->check(std::get<0>(pair), std::get<1>(pair), i, 0, 1);
     ++i;
   }
 }
 
 TEST_F(TestZipFixutre, Bidirectional) {
   auto _zip = zip(a, b);
   auto begin = _zip.begin();
 
   auto it = begin;
   ++it;
   EXPECT_EQ(begin, --it);
 
   it = begin;
   EXPECT_EQ(begin, it++);
   EXPECT_EQ(begin, --it);
 
   auto it2 = it = begin;
   ++it;
   ++it2;
   EXPECT_EQ(it2, it--);
   EXPECT_EQ(begin, it);
 }
 
 TEST_F(TestZipFixutre, RandomAccess) {
   auto _zip = zip(a, b);
   auto begin = _zip.begin();
   auto end = _zip.end();
 
   auto && val5 = begin[5];
   this->check(std::get<0>(val5), std::get<1>(val5), 5, 0, 0);
 
   auto && val13 = begin[13];
   this->check(std::get<0>(val13), std::get<1>(val13), 13, 0, 0);
 
   EXPECT_EQ(end - begin, a.size());
   auto it = ++begin;
   EXPECT_EQ(begin + 1, ++it);
   EXPECT_EQ(begin, it - 1);
 }
 
 TEST_F(TestZipFixutre, Cat) {
   size_t i = 0;
   for (auto && data : make_zip_cat(zip(a, b), zip(a, b))) {
     this->check(std::get<0>(data), std::get<1>(data), i, 0, 0);
     this->check(std::get<2>(data), std::get<3>(data), i, 0, 0);
     ++i;
   }
 }
 
+/* -------------------------------------------------------------------------- */
+TEST_F(TestZipFixutre, SimpleNamedTest) {
+  size_t i = 0;
+  for (auto && pair :
+       named_zip(tuple::get<"a"_h>() = this->a, tuple::get<"b"_h> = this->b)) {
+    this->check(tuple::get<"a"_h>(pair), tuple::get<"b"_h>(pair), i, 0, 0);
+    ++i;
+  }
+}
+
+TEST_F(TestZipFixutre, ConstTest) {
+  size_t i = 0;
+  const auto & ca = this->a;
+  const auto & cb = this->b;
+  for (auto && pair : named_zip(tuple::get<"a"_h>() = ca, tuple::get<"b"_h> = cb) {
+    this->check(tuple::get<"a"_h>(pair), tuple::get<"b">(pair), i, 0, 0);
+    EXPECT_EQ(
+        true,
+        std::is_const<
+            std::remove_reference_t<decltype(tuple::get<"a"_h>(pair))>>::value);
+    EXPECT_EQ(
+        true,
+        std::is_const<
+            std::remove_reference_t<decltype(tuple::get<"b"_h>(pair))>>::value);
+    ++i;
+  }
+}
+
+TEST_F(TestZipFixutre, MixteTest) {
+  size_t i = 0;
+  const auto & cb = this->b;
+  for (auto && pair : zip(tuple::get<"a"_h>() = a, tuple::get<"b"_h>() = cb)) {
+    this->check(std::get<0>(pair), std::get<1>(pair), i, 0, 0);
+    EXPECT_EQ(
+        false,
+        std::is_const<
+            std::remove_reference_t<decltype(tuple::get<"a"_h>(pair))>>::value);
+    EXPECT_EQ(
+        true,
+        std::is_const<
+            std::remove_reference_t<decltype(tuple::get<"b"_h>(pair))>>::value);
+    ++i;
+  }
+}
+
+TEST_F(TestZipFixutre, MoveNamedTest) {
+  size_t i = 0;
+  for (auto && pair :
+       named_zip(tuple::get<"a"_h>() = C<int>(0, this->size),
+       tuple::get<"b"_h>() = C<int>(this->size, 2 * this->size))) {
+    this->check(tuple::get<"a"_h>(pair), tuple::get<"b"_h>(pair), i, 0, 1);
+    ++i;
+  }
+}
+
+TEST_F(TestZipFixutre, BidirectionalNamed) {
+  auto _zip = named_zip(tuple::get<"a"_h>() = a, tuple::get<"b"_h>() = b);
+  auto begin = _zip.begin();
+
+  auto it = begin;
+  ++it;
+  EXPECT_EQ(begin, --it);
+
+  it = begin;
+  EXPECT_EQ(begin, it++);
+  EXPECT_EQ(begin, --it);
+
+  auto it2 = it = begin;
+  ++it;
+  ++it2;
+  EXPECT_EQ(it2, it--);
+  EXPECT_EQ(begin, it);
+}
+
+TEST_F(TestZipFixutre, RandomAccessNamed) {
+  auto _zip = zip(tuple::get<"a"_h>() = a, tuple::get<"b"_h>() = b);
+  auto begin = _zip.begin();
+  auto end = _zip.end();
+
+  auto && val5 = begin[5];
+  this->check(tuple::get<"a"_h>(val5), tuple::get<"b"_h>(val5), 5, 0, 0);
+
+  auto && val13 = begin[13];
+  this->check(tuple::get<"a"_h>(val13), tuple::get<"b"_h>(val13), 13, 0, 0);
+
+  EXPECT_EQ(end - begin, a.size());
+  auto it = ++begin;
+  EXPECT_EQ(begin + 1, ++it);
+  EXPECT_EQ(begin, it - 1);
+}
+
 TEST(TestNamedZipFixutre, Simple) {
   std::vector<int> a{0, 10, 20, 30, 40};
   std::vector<int> b{0, 1, 2, 3, 4};
 
   using namespace tuple;
   for (auto && data : named_zip(get<"a"_h>() = a, get<"b"_h>() = b)) {
     auto & a = tuple::get<"a"_h>(data);
     auto & b = tuple::get<"b"_h>(data);
     b *= 10;
     EXPECT_EQ(b, a);
   }
 
-   for (auto && data : named_zip(get<"a"_h>() = a, get<"b"_h>() = b)) {
+  for (auto && data : named_zip(get<"a"_h>() = a, get<"b"_h>() = b)) {
     auto & a = tuple::get<"a"_h>(data);
     auto & b = tuple::get<"b"_h>(data);
     EXPECT_EQ(b, a);
   }
 }
 
 /* -------------------------------------------------------------------------- */
 TEST(TestArangeIterator, Stop) {
   size_t ref_i = 0;
   for (auto i : arange(10)) {
     EXPECT_EQ(ref_i, i);
     ++ref_i;
   }
 }
 
 TEST(TestArangeIterator, StartStop) {
   size_t ref_i = 1;
   for (auto i : arange(1, 10)) {
     EXPECT_EQ(ref_i, i);
     ++ref_i;
   }
 }
 
 TEST(TestArangeIterator, StartStopStep) {
   size_t ref_i = 1;
   for (auto i : arange(1, 22, 2)) {
     EXPECT_EQ(ref_i, i);
     ref_i += 2;
   }
 }
 
 TEST(TestArangeIterator, StartStopStepZipped) {
   int ref_i1 = -1, ref_i2 = 1;
   for (auto && i : zip(arange(-1, -10, -1), arange(1, 18, 2))) {
     EXPECT_EQ(ref_i1, std::get<0>(i));
     EXPECT_EQ(ref_i2, std::get<1>(i));
     ref_i1 += -1;
     ref_i2 += 2;
   }
 }
 
 /* -------------------------------------------------------------------------- */
 TEST(TestEnumerateIterator, SimpleTest) {
   std::vector<int> a{0, 10, 20, 30, 40};
   std::vector<int> b{0, 2, 4, 6, 8};
   for (auto && data : enumerate(a, b)) {
     EXPECT_EQ(std::get<0>(data) * 10, std::get<1>(data));
     EXPECT_EQ(std::get<0>(data) * 2, std::get<2>(data));
   }
 }
 
 /* -------------------------------------------------------------------------- */
 TEST(TestTransformAdaptor, Keys) {
   std::map<std::string, int> map{
       {"1", 1}, {"2", 2}, {"3", 3}, {"3", 3}, {"4", 4}};
 
   char counter = '1';
   for (auto && key : make_keys_adaptor(map)) {
     EXPECT_EQ(counter, key[0]);
     ++counter;
   }
 }
 
 TEST(TestTransformAdaptor, Values) {
   std::map<std::string, int> map{
       {"1", 1}, {"2", 2}, {"3", 3}, {"3", 3}, {"4", 4}};
 
   int counter = 1;
   for (auto && value : make_values_adaptor(map)) {
     EXPECT_EQ(counter, value);
     ++counter;
   }
 }
 
 static int plus1(int value) { return value + 1; }
 
 struct Plus {
   Plus(int a) : a(a) {}
   int operator()(int b) { return a + b; }
 
 private:
   int a{0};
 };
 
 TEST(TestTransformAdaptor, Lambda) {
   auto && container = arange(10);
 
   for (auto && data :
        zip(container, make_transform_adaptor(container, [](auto && value) {
              return value + 1;
            }))) {
     EXPECT_EQ(std::get<0>(data) + 1, std::get<1>(data));
   }
 }
 
 TEST(TestTransformAdaptor, LambdaLambda) {
   std::map<std::string, int> map{
       {"1", 1}, {"2", 2}, {"3", 3}, {"3", 3}, {"4", 4}};
 
   int counter = 1;
   for (auto && data : make_transform_adaptor(
            make_values_adaptor(map), [](auto && value) { return value + 1; })) {
     EXPECT_EQ(counter + 1, data);
     ++counter;
   }
 
   auto && container = arange(10);
 
   for (auto && data :
        zip(container, make_transform_adaptor(container, [](auto && value) {
              return value + 1;
            }))) {
     EXPECT_EQ(std::get<0>(data) + 1, std::get<1>(data));
   }
 }
 
 TEST(TestTransformAdaptor, Function) {
   auto && container = arange(10);
 
   for (auto && data :
        zip(container, make_transform_adaptor(container, plus1))) {
     EXPECT_EQ(std::get<0>(data) + 1, std::get<1>(data));
   }
 }
 
 TEST(TestTransformAdaptor, Functor) {
   auto && container = arange(10);
 
   for (auto && data :
        zip(container, make_transform_adaptor(container, Plus(1)))) {
     EXPECT_EQ(std::get<0>(data) + 1, std::get<1>(data));
   }
 }
 
 /* -------------------------------------------------------------------------- */
 TEST(TestFilteredIterator, Simple) {
   std::vector<int> values{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
   std::vector<int> filter_{1, 3, 4, 10, 8, 6};
   for (auto && data : zip(filter_, filter(filter_, values))) {
     EXPECT_EQ(std::get<0>(data), std::get<1>(data));
   }
 }
 
 /* -------------------------------------------------------------------------- */
 TEST(TestFilteredIterator, Temporary) {
   std::vector<int> filter_{1, 3, 4, 10, 8, 6};
   for (auto && data :
        zip(filter_, filter(filter_, std::vector<int>{0, 1, 2, 3, 4, 5, 6, 7, 8,
                                                      9, 10}))) {
     EXPECT_EQ(std::get<0>(data), std::get<1>(data));
   }
 }