diff --git a/src/algebra/CRSMatrixDescriptor.hpp b/src/algebra/CRSMatrixDescriptor.hpp
index 8f123896d486b7d6b67be1350a96183423b4235e..e367db354e356035d18a8742b2410a8781dde09a 100644
--- a/src/algebra/CRSMatrixDescriptor.hpp
+++ b/src/algebra/CRSMatrixDescriptor.hpp
@@ -35,7 +35,7 @@ class CRSMatrixDescriptor
     Array<IndexType> row_map(m_nb_rows + 1);
     row_map[0] = 0;
     for (IndexType i = 0; i < m_nb_rows; ++i) {
-      Assert(non_zeros[i] >= 0);
+      Assert((non_zeros[i] >= 0) and (non_zeros[i] <= m_nb_columns));
       row_map[i + 1] = row_map[i] + non_zeros[i];
     }
 
diff --git a/src/algebra/TinyMatrix.hpp b/src/algebra/TinyMatrix.hpp
index 6c23b068bf1a8c8744e24ffbf60ec88b42b1317b..6bf810cde020681da88afedb6a815e54a563504d 100644
--- a/src/algebra/TinyMatrix.hpp
+++ b/src/algebra/TinyMatrix.hpp
@@ -1,13 +1,12 @@
 #ifndef TINY_MATRIX_HPP
 #define TINY_MATRIX_HPP
 
+#include <algebra/TinyVector.hpp>
+#include <utils/InvalidData.hpp>
 #include <utils/PugsAssert.hpp>
 #include <utils/PugsMacros.hpp>
-
 #include <utils/Types.hpp>
 
-#include <algebra/TinyVector.hpp>
-
 #include <iostream>
 
 template <size_t N, typename T = double>
@@ -278,7 +277,14 @@ class [[nodiscard]] TinyMatrix
   // One does not use the '=default' constructor to avoid (unexpected)
   // performances issues
   PUGS_INLINE
-  constexpr TinyMatrix() noexcept {}
+  constexpr TinyMatrix() noexcept
+  {
+#ifndef NDEBUG
+    for (size_t i = 0; i < N * N; ++i) {
+      m_values[i] = invalid_data_v<T>;
+    }
+#endif   // NDEBUG
+  }
 
   PUGS_INLINE
   constexpr TinyMatrix(ZeroType) noexcept
diff --git a/src/algebra/TinyVector.hpp b/src/algebra/TinyVector.hpp
index f24503aade42275e8b96f6958535fc8439da356c..69ecb14e39d626b1ed0e2709b6fd018af92714da 100644
--- a/src/algebra/TinyVector.hpp
+++ b/src/algebra/TinyVector.hpp
@@ -1,14 +1,13 @@
 #ifndef TINY_VECTOR_HPP
 #define TINY_VECTOR_HPP
 
-#include <iostream>
-
+#include <utils/InvalidData.hpp>
 #include <utils/PugsAssert.hpp>
 #include <utils/PugsMacros.hpp>
-
 #include <utils/Types.hpp>
 
 #include <cmath>
+#include <iostream>
 
 template <size_t N, typename T = double>
 class [[nodiscard]] TinyVector
@@ -211,7 +210,14 @@ class [[nodiscard]] TinyVector
   // One does not use the '=default' constructor to avoid
   // (zero-initialization) performances issues
   PUGS_INLINE
-  constexpr TinyVector() noexcept {}
+  constexpr TinyVector() noexcept
+  {
+#ifndef NDEBUG
+    for (size_t i = 0; i < N; ++i) {
+      m_values[i] = invalid_data_v<T>;
+    }
+#endif   // NDEBUG
+  }
 
   PUGS_INLINE
   constexpr TinyVector(ZeroType) noexcept
diff --git a/src/utils/Array.hpp b/src/utils/Array.hpp
index 9fa02e95a41a5200260fab84cda211412c4747e2..e47e56575fe725f67f6c6815c34c42eac4ce8d44 100644
--- a/src/utils/Array.hpp
+++ b/src/utils/Array.hpp
@@ -1,11 +1,11 @@
 #ifndef ARRAY_HPP
 #define ARRAY_HPP
 
+#include <utils/InvalidData.hpp>
+#include <utils/PugsAssert.hpp>
 #include <utils/PugsMacros.hpp>
 #include <utils/PugsUtils.hpp>
 
-#include <utils/PugsAssert.hpp>
-
 #include <Kokkos_CopyViews.hpp>
 #include <algorithm>
 
@@ -94,6 +94,19 @@ class [[nodiscard]] Array
     } else {
       m_values = Kokkos::View<DataType*>{"anonymous", size};
     }
+
+#ifndef NDEBUG
+    if constexpr (not std::is_const_v<DataType>) {
+      using T = std::decay_t<DataType>;
+      if constexpr (std::is_arithmetic_v<T>) {
+        this->fill(invalid_data_v<T>);
+      } else if constexpr (is_tiny_vector_v<T>) {
+        this->fill(T{});
+      } else if constexpr (is_tiny_matrix_v<T>) {
+        this->fill(T{});
+      }
+    }
+#endif   // NDEBUG
   }
 
   PUGS_INLINE
diff --git a/src/utils/InvalidData.hpp b/src/utils/InvalidData.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..78e354c25a9efcd46a93d4bc98259cc0ece82bde
--- /dev/null
+++ b/src/utils/InvalidData.hpp
@@ -0,0 +1,30 @@
+#ifndef NDEBUG
+
+#ifndef INVALID_DATA_HPP
+#define INVALID_DATA_HPP
+
+#include <limits>
+#include <type_traits>
+
+template <typename DataType>
+struct invalid_data
+{
+  static constexpr DataType value = [] {
+    using T = std::decay_t<DataType>;
+    static_assert(std::is_arithmetic_v<T>);
+    if constexpr (std::is_floating_point_v<T>) {
+      return std::numeric_limits<T>::signaling_NaN();
+    } else if constexpr (std::is_integral_v<T>) {
+      return std::numeric_limits<T>::max() / 2;
+    } else {
+      return T{};
+    }
+  }();
+};
+
+template <typename DataType>
+inline constexpr DataType invalid_data_v = invalid_data<DataType>::value;
+
+#endif   // INVALID_DATA_HPP
+
+#endif   // NDEBUG
diff --git a/src/utils/Table.hpp b/src/utils/Table.hpp
index 8dae32b97b0d15a99256fd38c4f4f7c5ecf29588..5f621ab2719e93ea63d3be610818aab1cb164665 100644
--- a/src/utils/Table.hpp
+++ b/src/utils/Table.hpp
@@ -1,12 +1,14 @@
 #ifndef TABLE_HPP
 #define TABLE_HPP
 
-#include <Kokkos_CopyViews.hpp>
 #include <utils/Array.hpp>
+#include <utils/InvalidData.hpp>
 #include <utils/PugsAssert.hpp>
 #include <utils/PugsMacros.hpp>
 #include <utils/PugsUtils.hpp>
 
+#include <Kokkos_CopyViews.hpp>
+
 template <typename DataType>
 class [[nodiscard]] Table
 {
@@ -109,6 +111,19 @@ class [[nodiscard]] Table
     } else {
       m_values = Kokkos::View<DataType**>{"anonymous", nb_lines, nb_columns};
     }
+
+#ifndef NDEBUG
+    if constexpr (not std::is_const_v<DataType>) {
+      using T = std::decay_t<DataType>;
+      if constexpr (std::is_arithmetic_v<T>) {
+        this->fill(invalid_data_v<T>);
+      } else if constexpr (is_tiny_vector_v<T>) {
+        this->fill(T{});
+      } else if constexpr (is_tiny_matrix_v<T>) {
+        this->fill(T{});
+      }
+    }
+#endif   // NDEBUG
   }
 
   PUGS_INLINE
diff --git a/tests/test_Array.cpp b/tests/test_Array.cpp
index 23d63a4a2dd769c679811ddc00265ba4aa615dac..dd6f57ea61bd52632477cb16fc7f6556757036f5 100644
--- a/tests/test_Array.cpp
+++ b/tests/test_Array.cpp
@@ -243,5 +243,23 @@ TEST_CASE("Array", "[utils]")
     Array<int> b{2 * a.size()};
     REQUIRE_THROWS_AS(copy_to(a, b), AssertError);
   }
+
+  SECTION("checking for nan initialization")
+  {
+    Array<double> array(10);
+
+    for (size_t i = 0; i < array.size(); ++i) {
+      REQUIRE(std::isnan(array[i]));
+    }
+  }
+
+  SECTION("checking for bad initialization")
+  {
+    Array<int> array(10);
+
+    for (size_t i = 0; i < array.size(); ++i) {
+      REQUIRE(array[i] == std::numeric_limits<int>::max() / 2);
+    }
+  }
 #endif   // NDEBUG
 }
diff --git a/tests/test_CRSMatrixDescriptor.cpp b/tests/test_CRSMatrixDescriptor.cpp
index 2b311d84b06dff839b3376d8c72c4439bc6855e0..28276aaadb747f01933fda8ce302cb853af5cf6d 100644
--- a/tests/test_CRSMatrixDescriptor.cpp
+++ b/tests/test_CRSMatrixDescriptor.cpp
@@ -145,12 +145,14 @@ TEST_CASE("CRSMatrixDescriptor", "[algebra]")
   SECTION("incompatible non zero vector size and number of columns")
   {
     Array<int> non_zeros(3);
+    non_zeros.fill(1);
     REQUIRE_THROWS_AS((CRSMatrixDescriptor<int>{2, 4, non_zeros}), AssertError);
   }
 
   SECTION("bad row number")
   {
     Array<int> non_zeros(2);
+    non_zeros.fill(1);
     CRSMatrixDescriptor<int> S{2, 4, non_zeros};
     REQUIRE_THROWS_AS(S(2, 3) = 1, AssertError);
   }
@@ -158,6 +160,7 @@ TEST_CASE("CRSMatrixDescriptor", "[algebra]")
   SECTION("bad column number")
   {
     Array<int> non_zeros(2);
+    non_zeros.fill(2);
     CRSMatrixDescriptor<int> S{2, 4, non_zeros};
     REQUIRE_THROWS_AS(S(0, 5) = 1, AssertError);
   }
diff --git a/tests/test_Table.cpp b/tests/test_Table.cpp
index 5e738867a592c3f06a9fe8ac7bd958dbcf390e4e..5372c697dc8d8bde1e2039141a8d8e226e7a68d3 100644
--- a/tests/test_Table.cpp
+++ b/tests/test_Table.cpp
@@ -199,5 +199,27 @@ TEST_CASE("Table", "[utils]")
       REQUIRE_THROWS_AS(copy_to(a, c), AssertError);
     }
   }
+
+  SECTION("checking for nan initialization")
+  {
+    Table<double> B(3, 4);
+
+    for (size_t i = 0; i < B.numberOfRows(); ++i) {
+      for (size_t j = 0; j < B.numberOfColumns(); ++j) {
+        REQUIRE(std::isnan(B(i, j)));
+      }
+    }
+  }
+
+  SECTION("checking for bad initialization")
+  {
+    Table<int> B(4, 2);
+
+    for (size_t i = 0; i < B.numberOfRows(); ++i) {
+      for (size_t j = 0; j < B.numberOfColumns(); ++j) {
+        REQUIRE(B(i, j) == std::numeric_limits<int>::max() / 2);
+      }
+    }
+  }
 #endif   // NDEBUG
 }
diff --git a/tests/test_TinyMatrix.cpp b/tests/test_TinyMatrix.cpp
index 464e78448192f3dacc248d6c7fde8c63338449bb..3a03b47f71e9e5491485053870f93ab45349d942 100644
--- a/tests/test_TinyMatrix.cpp
+++ b/tests/test_TinyMatrix.cpp
@@ -238,5 +238,27 @@ TEST_CASE("TinyMatrix", "[algebra]")
     REQUIRE_THROWS_AS(constA(3, 0), AssertError);
     REQUIRE_THROWS_AS(constA(0, 3), AssertError);
   }
+
+  SECTION("checking for nan initialization")
+  {
+    TinyMatrix<3, double> B;
+
+    for (size_t i = 0; i < B.numberOfRows(); ++i) {
+      for (size_t j = 0; j < B.numberOfColumns(); ++j) {
+        REQUIRE(std::isnan(B(i, j)));
+      }
+    }
+  }
+
+  SECTION("checking for bad initialization")
+  {
+    TinyMatrix<3, int> B;
+
+    for (size_t i = 0; i < B.numberOfRows(); ++i) {
+      for (size_t j = 0; j < B.numberOfColumns(); ++j) {
+        REQUIRE(B(i, j) == std::numeric_limits<int>::max() / 2);
+      }
+    }
+  }
 #endif   // NDEBUG
 }
diff --git a/tests/test_TinyVector.cpp b/tests/test_TinyVector.cpp
index 499536f9c702d5c3a1f2fcf8c4da7b1ecb461038..c4fb06998505d0eb5e5697835958d9a970d5e01a 100644
--- a/tests/test_TinyVector.cpp
+++ b/tests/test_TinyVector.cpp
@@ -87,5 +87,24 @@ TEST_CASE("TinyVector", "[algebra]")
     const TinyVector<3, int>& const_x = x;
     REQUIRE_THROWS_AS(const_x[-1], AssertError);
   }
+
+  SECTION("checking for nan initialization")
+  {
+    TinyVector<3, double> y;
+
+    for (size_t i = 0; i < y.dimension(); ++i) {
+      REQUIRE(std::isnan(y[i]));
+    }
+  }
+
+  SECTION("checking for bad initialization")
+  {
+    TinyVector<3, int> y;
+
+    for (size_t i = 0; i < y.dimension(); ++i) {
+      REQUIRE(y[i] == std::numeric_limits<int>::max() / 2);
+    }
+  }
+
 #endif   // NDEBUG
 }