diff --git a/src/language/modules/SchemeModule.cpp b/src/language/modules/SchemeModule.cpp
index a6dcec9c8d510e36a8b70547eafbf11775acef8a..542690a4ffcfb73ba174f0014b358548208d830f 100644
--- a/src/language/modules/SchemeModule.cpp
+++ b/src/language/modules/SchemeModule.cpp
@@ -6,6 +6,7 @@
 #include <mesh/Mesh.hpp>
 #include <scheme/DirichletBoundaryConditionDescriptor.hpp>
 #include <scheme/FourierBoundaryConditionDescriptor.hpp>
+#include <scheme/FreeBoundaryConditionDescriptor.hpp>
 #include <scheme/IBoundaryConditionDescriptor.hpp>
 #include <scheme/IBoundaryDescriptor.hpp>
 #include <scheme/NamedBoundaryDescriptor.hpp>
diff --git a/src/mesh/SubItemValuePerItem.hpp b/src/mesh/SubItemValuePerItem.hpp
index 8189fc13a8b6a6b9948735da21456ff3296602d7..66af694c7255da41ff7570d294dd4777acb0164f 100644
--- a/src/mesh/SubItemValuePerItem.hpp
+++ b/src/mesh/SubItemValuePerItem.hpp
@@ -55,8 +55,7 @@ class SubItemValuePerItem
 
    public:
     template <typename IndexType>
-    PUGS_INLINE const DataType&
-    operator[](IndexType i) const noexcept(NO_ASSERT)
+    PUGS_INLINE const DataType& operator[](IndexType i) const noexcept(NO_ASSERT)
     {
       static_assert(std::is_integral_v<IndexType>, "SubView is indexed by integral values");
       Assert(i < m_size);
@@ -64,8 +63,7 @@ class SubItemValuePerItem
     }
 
     template <typename IndexType>
-    PUGS_FORCEINLINE DataType&
-    operator[](IndexType i) noexcept(NO_ASSERT)
+    PUGS_FORCEINLINE DataType& operator[](IndexType i) noexcept(NO_ASSERT)
     {
       static_assert(std::is_integral_v<IndexType>, "SubView is indexed by integral values");
       Assert(i < m_size);
@@ -104,10 +102,10 @@ class SubItemValuePerItem
   PUGS_FORCEINLINE DataType&
   operator()(IndexType item_id, SubIndexType i) const noexcept(NO_ASSERT)
   {
-    static_assert(std::is_same_v<IndexType, ItemId>, "first index must be an ItemId");
+    static_assert(std::is_same_v<IndexType, ItemId>, "first index must be of the correct ItemId type");
     static_assert(std::is_integral_v<SubIndexType>, "second index must be an integral type");
     Assert(this->isBuilt());
-    Assert(i < m_host_row_map(size_t{item_id} + 1) - m_host_row_map(size_t{item_id}));
+    Assert(i + m_host_row_map(size_t{item_id}) < m_host_row_map(size_t{item_id} + 1));
     return m_values[m_host_row_map(size_t{item_id}) + i];
   }
 
@@ -122,8 +120,7 @@ class SubItemValuePerItem
   // Following Kokkos logic, these classes are view and const view does allow
   // changes in data
   template <typename ArrayIndexType>
-  DataType&
-  operator[](const ArrayIndexType& i) const noexcept(NO_ASSERT)
+  DataType& operator[](const ArrayIndexType& i) const noexcept(NO_ASSERT)
   {
     static_assert(std::is_integral_v<ArrayIndexType>, "index must be an integral type");
     Assert(this->isBuilt());
diff --git a/src/scheme/FreeBoundaryConditionDescriptor.hpp b/src/scheme/FreeBoundaryConditionDescriptor.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..b0eb2ceb218d77cae347a44219f13d278a25aee5
--- /dev/null
+++ b/src/scheme/FreeBoundaryConditionDescriptor.hpp
@@ -0,0 +1,46 @@
+#ifndef FREE_BOUNDARY_CONDITION_DESCRIPTOR_HPP
+#define FREE_BOUNDARY_CONDITION_DESCRIPTOR_HPP
+
+#include <scheme/IBoundaryConditionDescriptor.hpp>
+#include <scheme/IBoundaryDescriptor.hpp>
+
+#include <memory>
+
+class FreeBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
+{
+ private:
+  std::ostream&
+  _write(std::ostream& os) const final
+  {
+    os << "free(" << *m_boundary_descriptor << ")";
+    return os;
+  }
+
+  std::shared_ptr<const IBoundaryDescriptor> m_boundary_descriptor;
+
+ public:
+  const IBoundaryDescriptor&
+  boundaryDescriptor() const
+  {
+    return *m_boundary_descriptor;
+  }
+
+  Type
+  type() const final
+  {
+    return Type::free;
+  }
+
+  FreeBoundaryConditionDescriptor(std::shared_ptr<const IBoundaryDescriptor> boundary_descriptor)
+    : m_boundary_descriptor(boundary_descriptor)
+  {
+    ;
+  }
+
+  FreeBoundaryConditionDescriptor(const FreeBoundaryConditionDescriptor&) = delete;
+  FreeBoundaryConditionDescriptor(FreeBoundaryConditionDescriptor&&)      = delete;
+
+  ~FreeBoundaryConditionDescriptor() = default;
+};
+
+#endif   // FREE_BOUNDARY_CONDITION_DESCRIPTOR_HPP
diff --git a/src/scheme/IBoundaryConditionDescriptor.hpp b/src/scheme/IBoundaryConditionDescriptor.hpp
index be25ba73180f001ee264e9087379dc15ee03f378..30a25db66f1278e72eb2bea0ce7a017d8022da65 100644
--- a/src/scheme/IBoundaryConditionDescriptor.hpp
+++ b/src/scheme/IBoundaryConditionDescriptor.hpp
@@ -10,6 +10,7 @@ class IBoundaryConditionDescriptor
   {
     dirichlet,
     fourier,
+    free,
     neumann,
     symmetry
   };