diff --git a/CMakeLists.txt b/CMakeLists.txt
index 62e3363416cab5794d69e25e60d70852649fa2a4..16a82546da8a04a714955ab224227f38262a39c4 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -258,21 +258,11 @@ include_directories(${PUGS_SOURCE_DIR}/packages/CLI11/include)
 include_directories(${PUGS_SOURCE_DIR}/packages/PEGTL/include/tao)
 
 # Pugs src
-add_subdirectory(src)
-
-include_directories(src)
-include_directories(src/algebra)
-include_directories(src/language)
-include_directories(src/mesh)
-include_directories(src/output)
-include_directories(src/scheme)
-include_directories(src/utils)
-
-# Pugs generated sources
-include_directories(${PUGS_BINARY_DIR}/src/utils)
+add_subdirectory(${PUGS_SOURCE_DIR}/src)
+include_directories(${PUGS_SOURCE_DIR}/src)
+include_directories(${PUGS_BINARY_DIR}/src)
 
 # Pugs tests
-
 set(CATCH_MODULE_PATH "${PUGS_SOURCE_DIR}/packages/Catch2")
 include("${CATCH_MODULE_PATH}/contrib/ParseAndAddCatchTests.cmake")
 add_subdirectory("${CATCH_MODULE_PATH}")
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index 7557831db4a04bb627024a023c026bf85812f05e..42d75ad7e7146554812cef242214d290903d8f31 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -5,23 +5,18 @@ include_directories(${CMAKE_CURRENT_BINARY_DIR})
 
 # Pugs utils
 add_subdirectory(utils)
-include_directories(utils)
 
 # Pugs language
 add_subdirectory(language)
-include_directories(language)
 
 # Pugs algebra
 #add_subdirectory(algebra)
-include_directories(algebra)
 
 # Pugs mesh
 add_subdirectory(mesh)
-include_directories(mesh)
 
 # Pugs mesh
 #add_subdirectory(mesh)
-include_directories(scheme)
 
 # Pugs output
-include_directories(output)
+#add_subdirectory(output)
diff --git a/src/algebra/CRSMatrix.hpp b/src/algebra/CRSMatrix.hpp
index 740343cd612f22911e187351b7929eedfbf7dcff..97a0a8a8ed44d660f9524eac269abee0324cd136 100644
--- a/src/algebra/CRSMatrix.hpp
+++ b/src/algebra/CRSMatrix.hpp
@@ -1,13 +1,13 @@
 #ifndef CRS_MATRIX_HPP
 #define CRS_MATRIX_HPP
 
-#include <Array.hpp>
 #include <Kokkos_StaticCrsGraph.hpp>
-#include <PugsAssert.hpp>
 
-#include <SparseMatrixDescriptor.hpp>
+#include <utils/Array.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <Vector.hpp>
+#include <algebra/SparseMatrixDescriptor.hpp>
+#include <algebra/Vector.hpp>
 
 #include <iostream>
 
diff --git a/src/algebra/SparseMatrixDescriptor.hpp b/src/algebra/SparseMatrixDescriptor.hpp
index badfb96305c72dbe5777c689a0861dfed299ff51..3533a4017011b0c5c1f5d7f04ae2b893fa8a7d0b 100644
--- a/src/algebra/SparseMatrixDescriptor.hpp
+++ b/src/algebra/SparseMatrixDescriptor.hpp
@@ -1,7 +1,7 @@
 #ifndef SPARSE_MATRIX_DESCRIPTOR_HPP
 #define SPARSE_MATRIX_DESCRIPTOR_HPP
 
-#include <Array.hpp>
+#include <utils/Array.hpp>
 
 #include <map>
 #include <type_traits>
diff --git a/src/algebra/TinyMatrix.hpp b/src/algebra/TinyMatrix.hpp
index 7a67fe0c1f11cd2a08954462b83fe2339a9b4adb..6cff7777f522c348c31ebf5eed9460f6adbfcd63 100644
--- a/src/algebra/TinyMatrix.hpp
+++ b/src/algebra/TinyMatrix.hpp
@@ -1,11 +1,12 @@
 #ifndef TINY_MATRIX_HPP
 #define TINY_MATRIX_HPP
 
-#include <PugsAssert.hpp>
-#include <PugsMacros.hpp>
+#include <utils/PugsAssert.hpp>
+#include <utils/PugsMacros.hpp>
 
-#include <TinyVector.hpp>
-#include <Types.hpp>
+#include <utils/Types.hpp>
+
+#include <algebra/TinyVector.hpp>
 
 #include <iostream>
 
diff --git a/src/algebra/TinyVector.hpp b/src/algebra/TinyVector.hpp
index 187c405ee0a3f2007aeb9d98e2b403cbcf472ee7..a999ec7bcee342dd252d85fd85be8d1364308216 100644
--- a/src/algebra/TinyVector.hpp
+++ b/src/algebra/TinyVector.hpp
@@ -3,10 +3,10 @@
 
 #include <iostream>
 
-#include <PugsAssert.hpp>
-#include <PugsMacros.hpp>
+#include <utils/PugsAssert.hpp>
+#include <utils/PugsMacros.hpp>
 
-#include <Types.hpp>
+#include <utils/Types.hpp>
 
 #include <cmath>
 
diff --git a/src/algebra/Vector.hpp b/src/algebra/Vector.hpp
index 43865dae1827cdb4875b3963127adafe78a8c6f3..ac48d969fc55432d3cd8fe180997af09ccee84d0 100644
--- a/src/algebra/Vector.hpp
+++ b/src/algebra/Vector.hpp
@@ -1,12 +1,12 @@
 #ifndef VECTOR_HPP
 #define VECTOR_HPP
 
-#include <PugsMacros.hpp>
-#include <PugsUtils.hpp>
+#include <utils/PugsMacros.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <PugsAssert.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <Array.hpp>
+#include <utils/Array.hpp>
 
 template <typename DataType>
 class Vector   // LCOV_EXCL_LINE
diff --git a/src/language/CMakeLists.txt b/src/language/CMakeLists.txt
index 2c7d8dd1fea695edc3d689ec2b6f2dbacce61e44..b01587e455a481b01f44db34cd7dd8d75adcb59e 100644
--- a/src/language/CMakeLists.txt
+++ b/src/language/CMakeLists.txt
@@ -1,14 +1,9 @@
-include_directories(${CMAKE_CURRENT_SOURCE_DIR})
-include_directories(${CMAKE_CURRENT_BINARY_DIR})
-
 # ------------------- Source files --------------------
 
 add_library(
   PugsLanguage
   PugsParser.cpp)
 
-#include_directories(${PUGS_SOURCE_DIR}/utils)
-
 # Additional dependencies
 #add_dependencies(PugsMesh)
 
diff --git a/src/language/PugsParser.cpp b/src/language/PugsParser.cpp
index 6efab91fb6b0f63536908a2849601de232a15b05..a413c3fa7f733c47ff64d2b36be97c534b923d5c 100644
--- a/src/language/PugsParser.cpp
+++ b/src/language/PugsParser.cpp
@@ -1,12 +1,12 @@
-#include <PugsParser.hpp>
-
-#include <iostream>
+#include <language/PugsParser.hpp>
 
 #include <rang.hpp>
 
 #include <pegtl.hpp>
 #include <pegtl/analyze.hpp>
 
+#include <iostream>
+
 using namespace TAO_PEGTL_NAMESPACE;
 
 namespace language
diff --git a/src/main.cpp b/src/main.cpp
index 346d6ed623bf50d04f6b285d2c3994337d3448d1..e48f74965267f17f6b2aa849e7739398c1d3044c 100644
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -1,29 +1,29 @@
-#include <PugsUtils.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <rang.hpp>
+#include <mesh/Connectivity.hpp>
+#include <mesh/Mesh.hpp>
 
-#include <Connectivity.hpp>
+#include <scheme/AcousticSolver.hpp>
+#include <scheme/BoundaryCondition.hpp>
 
-#include <AcousticSolver.hpp>
-#include <BoundaryCondition.hpp>
-#include <Mesh.hpp>
+#include <output/VTKWriter.hpp>
 
-#include <VTKWriter.hpp>
+#include <utils/Exceptions.hpp>
+#include <utils/Timer.hpp>
 
-#include <Exceptions.hpp>
-#include <Timer.hpp>
+#include <algebra/TinyMatrix.hpp>
+#include <algebra/TinyVector.hpp>
 
-#include <TinyMatrix.hpp>
-#include <TinyVector.hpp>
+#include <scheme/BoundaryConditionDescriptor.hpp>
 
-#include <BoundaryConditionDescriptor.hpp>
+#include <mesh/MeshNodeBoundary.hpp>
 
-#include <MeshNodeBoundary.hpp>
+#include <language/PugsParser.hpp>
+#include <mesh/GmshReader.hpp>
 
-#include <GmshReader.hpp>
-#include <PugsParser.hpp>
+#include <mesh/SynchronizerManager.hpp>
 
-#include <SynchronizerManager.hpp>
+#include <rang.hpp>
 
 #include <iostream>
 
diff --git a/src/mesh/CMakeLists.txt b/src/mesh/CMakeLists.txt
index ec94fd4b755cf91b30928d9894e9f476860927c6..73d686fb20c158ad2693daa28124081902310f61 100644
--- a/src/mesh/CMakeLists.txt
+++ b/src/mesh/CMakeLists.txt
@@ -1,6 +1,3 @@
-include_directories(${CMAKE_CURRENT_SOURCE_DIR})
-include_directories(${CMAKE_CURRENT_BINARY_DIR})
-
 # ------------------- Source files --------------------
 
 add_library(
@@ -11,8 +8,6 @@ add_library(
   ConnectivityDispatcher.cpp
   SynchronizerManager.cpp)
 
-include_directories(${PUGS_BINARY_DIR}/src/utils)
-
 # Additional dependencies
 #add_dependencies(PugsMesh)
 
diff --git a/src/mesh/CellType.hpp b/src/mesh/CellType.hpp
index c0949ffbe27d6017b0c46a7e21cea188ea6c4abf..d0f6986d63a97677662ae91d980e9af02b114347 100644
--- a/src/mesh/CellType.hpp
+++ b/src/mesh/CellType.hpp
@@ -1,7 +1,8 @@
 #ifndef CELL_TYPE_HPP
 #define CELL_TYPE_HPP
 
-#include <PugsMacros.hpp>
+#include <utils/PugsMacros.hpp>
+
 #include <string_view>
 
 enum class CellType : unsigned short
diff --git a/src/mesh/Connectivity.cpp b/src/mesh/Connectivity.cpp
index 632726a16c1993d986758562fcdf5f7c44ec832c..dc51fe75a6cb773e99d59e235b137c7dd98c1164 100644
--- a/src/mesh/Connectivity.cpp
+++ b/src/mesh/Connectivity.cpp
@@ -1,9 +1,9 @@
-#include <Connectivity.hpp>
-#include <map>
+#include <mesh/Connectivity.hpp>
 
-#include <Messenger.hpp>
+#include <mesh/ConnectivityDescriptor.hpp>
+#include <utils/Messenger.hpp>
 
-#include <ConnectivityDescriptor.hpp>
+#include <map>
 
 template <size_t Dimension>
 Connectivity<Dimension>::Connectivity()
diff --git a/src/mesh/Connectivity.hpp b/src/mesh/Connectivity.hpp
index 6b7e31a4057de6a854c72237ba8482ba7df0abf8..3c319b053eb40a569c660d4c35a0b1235dfa768b 100644
--- a/src/mesh/Connectivity.hpp
+++ b/src/mesh/Connectivity.hpp
@@ -1,37 +1,37 @@
 #ifndef CONNECTIVITY_HPP
 #define CONNECTIVITY_HPP
 
-#include <PugsAssert.hpp>
-#include <PugsMacros.hpp>
+#include <utils/PugsAssert.hpp>
+#include <utils/PugsMacros.hpp>
 
-#include <PugsUtils.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <PugsTraits.hpp>
+#include <utils/PugsTraits.hpp>
 
-#include <TinyVector.hpp>
+#include <algebra/TinyVector.hpp>
 
-#include <ItemValue.hpp>
+#include <mesh/ItemValue.hpp>
 
-#include <IConnectivity.hpp>
+#include <mesh/IConnectivity.hpp>
 
-#include <ConnectivityMatrix.hpp>
-#include <ItemToItemMatrix.hpp>
+#include <mesh/ConnectivityMatrix.hpp>
+#include <mesh/ItemToItemMatrix.hpp>
 
-#include <ConnectivityComputer.hpp>
-#include <SubItemValuePerItem.hpp>
+#include <mesh/ConnectivityComputer.hpp>
+#include <mesh/SubItemValuePerItem.hpp>
 
 #include <algorithm>
 #include <vector>
 
-#include <CellType.hpp>
+#include <mesh/CellType.hpp>
 
-#include <CSRGraph.hpp>
+#include <utils/CSRGraph.hpp>
 
-#include <ItemType.hpp>
-#include <RefId.hpp>
-#include <RefItemList.hpp>
+#include <mesh/ItemType.hpp>
+#include <mesh/RefId.hpp>
+#include <mesh/RefItemList.hpp>
 
-#include <SynchronizerManager.hpp>
+#include <mesh/SynchronizerManager.hpp>
 
 #include <iostream>
 #include <set>
diff --git a/src/mesh/ConnectivityComputer.cpp b/src/mesh/ConnectivityComputer.cpp
index 8e1d5d52690e39a779dff6b4cac96fd64b068dbc..a47b491160bb3c27a935925083ba0fb230bfc6c2 100644
--- a/src/mesh/ConnectivityComputer.cpp
+++ b/src/mesh/ConnectivityComputer.cpp
@@ -1,8 +1,7 @@
-#include <Connectivity.hpp>
+#include <mesh/Connectivity.hpp>
 
-#include <ConnectivityComputer.hpp>
-
-#include <Exceptions.hpp>
+#include <mesh/ConnectivityComputer.hpp>
+#include <utils/Exceptions.hpp>
 
 #include <iostream>
 #include <map>
diff --git a/src/mesh/ConnectivityComputer.hpp b/src/mesh/ConnectivityComputer.hpp
index b56b1895ef74bec0afc55bbeeea7b6f7b911dfe6..c6e2a569bc00e0c620224ee9ce9037e4085824e7 100644
--- a/src/mesh/ConnectivityComputer.hpp
+++ b/src/mesh/ConnectivityComputer.hpp
@@ -1,8 +1,8 @@
 #ifndef CONNECTIVITY_COMPUTER_HPP
 #define CONNECTIVITY_COMPUTER_HPP
 
-#include <ConnectivityMatrix.hpp>
-#include <SubItemValuePerItem.hpp>
+#include <mesh/ConnectivityMatrix.hpp>
+#include <mesh/SubItemValuePerItem.hpp>
 
 class ConnectivityComputer
 {
diff --git a/src/mesh/ConnectivityDescriptor.hpp b/src/mesh/ConnectivityDescriptor.hpp
index 8b190be0c01abfd266508e64a22d4be1ff629b71..352a664ea061f38e006332af2abbfa052fad9192 100644
--- a/src/mesh/ConnectivityDescriptor.hpp
+++ b/src/mesh/ConnectivityDescriptor.hpp
@@ -1,9 +1,10 @@
 #ifndef CONNECTIVITY_DESCRIPTOR_HPP
 #define CONNECTIVITY_DESCRIPTOR_HPP
 
-#include <ItemOfItemType.hpp>
-#include <PugsTraits.hpp>
-#include <RefItemList.hpp>
+#include <mesh/ItemOfItemType.hpp>
+#include <mesh/RefItemList.hpp>
+
+#include <utils/PugsTraits.hpp>
 
 #include <vector>
 
diff --git a/src/mesh/ConnectivityDispatcher.cpp b/src/mesh/ConnectivityDispatcher.cpp
index 0136c264c66b9cd00a2be037a1019d1c83886383..ebba8dbe16312c61695f378291d79bfaf82d06d3 100644
--- a/src/mesh/ConnectivityDispatcher.cpp
+++ b/src/mesh/ConnectivityDispatcher.cpp
@@ -1,7 +1,7 @@
-#include <ConnectivityDispatcher.hpp>
-#include <Partitioner.hpp>
+#include <mesh/ConnectivityDispatcher.hpp>
+#include <mesh/ItemOfItemType.hpp>
 
-#include <ItemOfItemType.hpp>
+#include <utils/Partitioner.hpp>
 
 #include <iostream>
 #include <unordered_map>
diff --git a/src/mesh/ConnectivityDispatcher.hpp b/src/mesh/ConnectivityDispatcher.hpp
index 0e28f02e87ebef116ba02713450fe9e0a27cdc04..3de9c1c592a0dcc447104ce032c9d6efa9daceca 100644
--- a/src/mesh/ConnectivityDispatcher.hpp
+++ b/src/mesh/ConnectivityDispatcher.hpp
@@ -1,11 +1,12 @@
 #ifndef CONNECTIVITY_DISPATCHER_HPP
 #define CONNECTIVITY_DISPATCHER_HPP
 
-#include <ItemValue.hpp>
-#include <ItemValueUtils.hpp>
-#include <Mesh.hpp>
+#include <mesh/ItemValue.hpp>
+#include <mesh/ItemValueUtils.hpp>
+#include <mesh/Mesh.hpp>
+
+#include <mesh/ConnectivityDescriptor.hpp>
 
-#include <ConnectivityDescriptor.hpp>
 #include <unordered_map>
 
 template <int Dimension>
diff --git a/src/mesh/ConnectivityMatrix.hpp b/src/mesh/ConnectivityMatrix.hpp
index 0ddb21a964e18e9adae051e536c5a3f266c0e8bb..bef24bde58b5300d2285a7929bb06e45df00de1b 100644
--- a/src/mesh/ConnectivityMatrix.hpp
+++ b/src/mesh/ConnectivityMatrix.hpp
@@ -1,9 +1,10 @@
 #ifndef CONNECTIVITY_MATRIX_HPP
 #define CONNECTIVITY_MATRIX_HPP
 
-#include <Array.hpp>
+#include <utils/Array.hpp>
+#include <utils/PugsUtils.hpp>
+
 #include <Kokkos_StaticCrsGraph.hpp>
-#include <PugsUtils.hpp>
 
 class ConnectivityMatrix
 {
diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp
index 25a224caff4bdbfa3d7203a7a6891fadba5a756f..cf1706c9c54c727f884ae151817e486f17404201 100644
--- a/src/mesh/GmshReader.cpp
+++ b/src/mesh/GmshReader.cpp
@@ -1,29 +1,25 @@
-#include <GmshReader.hpp>
-#include <PugsMacros.hpp>
+#include <mesh/GmshReader.hpp>
 
-#include <fstream>
-#include <iostream>
-#include <rang.hpp>
-#include <set>
-
-#include <CellType.hpp>
-#include <Connectivity.hpp>
+#include <utils/PugsMacros.hpp>
 
-#include <Mesh.hpp>
-#include <MeshData.hpp>
+#include <mesh/CellType.hpp>
+#include <mesh/Connectivity.hpp>
+#include <mesh/ConnectivityDispatcher.hpp>
+#include <mesh/ItemValueUtils.hpp>
+#include <mesh/Mesh.hpp>
+#include <mesh/MeshData.hpp>
+#include <mesh/RefItemList.hpp>
 
-#include <Messenger.hpp>
-#include <RefItemList.hpp>
+#include <utils/ArrayUtils.hpp>
 
-#include <ArrayUtils.hpp>
-#include <ItemValueUtils.hpp>
+#include <utils/Exceptions.hpp>
 
-#include <ConnectivityDispatcher.hpp>
-
-#include <Exceptions.hpp>
+#include <rang.hpp>
 
+#include <fstream>
 #include <iomanip>
 #include <iostream>
+#include <set>
 #include <sstream>
 
 #include <map>
diff --git a/src/mesh/GmshReader.hpp b/src/mesh/GmshReader.hpp
index d6d5f3887eff282833ec169cd410bd3026188181..3501cafdf8c735ab73edf808ca7c886ec572ae34 100644
--- a/src/mesh/GmshReader.hpp
+++ b/src/mesh/GmshReader.hpp
@@ -1,12 +1,12 @@
 #ifndef GMSH_READER_HPP
 #define GMSH_READER_HPP
 
-#include <Array.hpp>
-#include <TinyVector.hpp>
+#include <algebra/TinyVector.hpp>
 
-#include <Mesh.hpp>
+#include <utils/Array.hpp>
 
-#include <RefId.hpp>
+#include <mesh/Mesh.hpp>
+#include <mesh/RefId.hpp>
 
 #include <array>
 #include <fstream>
diff --git a/src/mesh/IConnectivity.hpp b/src/mesh/IConnectivity.hpp
index d548be0a4177036ca082e990d50aa222ba44368f..76a85e0e3d8aae38ba97234395452cbba645db60 100644
--- a/src/mesh/IConnectivity.hpp
+++ b/src/mesh/IConnectivity.hpp
@@ -1,9 +1,10 @@
 #ifndef ICONNECTIVITY_HPP
 #define ICONNECTIVITY_HPP
 
-#include <ConnectivityMatrix.hpp>
-#include <ItemOfItemType.hpp>
-#include <ItemType.hpp>
+#include <mesh/ConnectivityMatrix.hpp>
+#include <mesh/ItemOfItemType.hpp>
+#include <mesh/ItemType.hpp>
+
 #include <memory>
 
 class IConnectivity : public std::enable_shared_from_this<IConnectivity>
diff --git a/src/mesh/ItemId.hpp b/src/mesh/ItemId.hpp
index de7790668e17261edf0446cea870ff1c1f9202a9..8863c1a2867773344b8ed6907ebc6786bff0a7f6 100644
--- a/src/mesh/ItemId.hpp
+++ b/src/mesh/ItemId.hpp
@@ -1,8 +1,9 @@
 #ifndef ITEM_ID_HPP
 #define ITEM_ID_HPP
 
-#include <ItemType.hpp>
-#include <PugsMacros.hpp>
+#include <mesh/ItemType.hpp>
+
+#include <utils/PugsMacros.hpp>
 
 template <ItemType item_type>
 class ItemIdT
diff --git a/src/mesh/ItemOfItemType.hpp b/src/mesh/ItemOfItemType.hpp
index 9921eea68ef8393b11bf386a93086bf3c1d31786..8d27a599bf80fd754ffb73b2756b03ce4193c1b0 100644
--- a/src/mesh/ItemOfItemType.hpp
+++ b/src/mesh/ItemOfItemType.hpp
@@ -1,7 +1,7 @@
 #ifndef ITEM_OF_ITEM_TYPE_HPP
 #define ITEM_OF_ITEM_TYPE_HPP
 
-#include <ItemType.hpp>
+#include <mesh/ItemType.hpp>
 
 template <ItemType sub_item_t, ItemType item_t>
 struct ItemOfItemType
diff --git a/src/mesh/ItemToItemMatrix.hpp b/src/mesh/ItemToItemMatrix.hpp
index 04b8a7b7487cb22789f4ff48c09292b403476695..a3fc307e03bf928203510514bb2f33cafd3d453c 100644
--- a/src/mesh/ItemToItemMatrix.hpp
+++ b/src/mesh/ItemToItemMatrix.hpp
@@ -1,11 +1,11 @@
 #ifndef ITEM_TO_ITEM_MATRIX_HPP
 #define ITEM_TO_ITEM_MATRIX_HPP
 
-#include <ItemId.hpp>
-#include <ItemType.hpp>
+#include <mesh/ConnectivityMatrix.hpp>
+#include <mesh/ItemId.hpp>
+#include <mesh/ItemType.hpp>
 
-#include <ConnectivityMatrix.hpp>
-#include <PugsUtils.hpp>
+#include <utils/PugsUtils.hpp>
 
 template <ItemType SourceItemType, ItemType TargetItemType>
 class ItemToItemMatrix
diff --git a/src/mesh/ItemType.hpp b/src/mesh/ItemType.hpp
index 1f4cbbd745514fc589fad56ad342a353dd581a94..21aebfc5d46b3d9c29712dfb54795826f076c9a3 100644
--- a/src/mesh/ItemType.hpp
+++ b/src/mesh/ItemType.hpp
@@ -1,6 +1,8 @@
 #ifndef ITEM_TYPE_HPP
 #define ITEM_TYPE_HPP
 
+#include <utils/PugsMacros.hpp>
+
 #include <limits>
 #include <string>
 #include <utility>
diff --git a/src/mesh/ItemValue.hpp b/src/mesh/ItemValue.hpp
index 81301f7a3d9a33e66906c20a48d984c22139bf8a..eb9f23113c3d2dc6694a06fd1af080eaf4705a49 100644
--- a/src/mesh/ItemValue.hpp
+++ b/src/mesh/ItemValue.hpp
@@ -1,14 +1,13 @@
 #ifndef ITEM_VALUE_HPP
 #define ITEM_VALUE_HPP
 
-#include <PugsAssert.hpp>
+#include <utils/Array.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <Array.hpp>
+#include <mesh/IConnectivity.hpp>
+#include <mesh/ItemId.hpp>
+#include <mesh/ItemType.hpp>
 
-#include <ItemId.hpp>
-#include <ItemType.hpp>
-
-#include <IConnectivity.hpp>
 #include <memory>
 
 template <typename DataType, ItemType item_type, typename ConnectivityPtr = std::shared_ptr<const IConnectivity>>
diff --git a/src/mesh/ItemValueUtils.hpp b/src/mesh/ItemValueUtils.hpp
index 666d4a0e9bf99b69b5d2b92e4588f84da30e8715..b9570fc9c2bb62d26a8ab2ee1a6b0f9709ba60ac 100644
--- a/src/mesh/ItemValueUtils.hpp
+++ b/src/mesh/ItemValueUtils.hpp
@@ -1,13 +1,12 @@
 #ifndef ITEM_VALUE_UTILS_HPP
 #define ITEM_VALUE_UTILS_HPP
 
-#include <ItemValue.hpp>
-#include <Messenger.hpp>
+#include <utils/Messenger.hpp>
 
-#include <Connectivity.hpp>
-
-#include <Synchronizer.hpp>
-#include <SynchronizerManager.hpp>
+#include <mesh/Connectivity.hpp>
+#include <mesh/ItemValue.hpp>
+#include <mesh/Synchronizer.hpp>
+#include <mesh/SynchronizerManager.hpp>
 
 #include <iostream>
 
diff --git a/src/mesh/Mesh.hpp b/src/mesh/Mesh.hpp
index a8e9411007783bb5a6406510c38c2bd40eb03d54..afafc7a2115d201055c93e0952c3a1e5216e05c2 100644
--- a/src/mesh/Mesh.hpp
+++ b/src/mesh/Mesh.hpp
@@ -1,10 +1,10 @@
 #ifndef MESH_HPP
 #define MESH_HPP
 
-#include <ItemValue.hpp>
-#include <TinyVector.hpp>
+#include <algebra/TinyVector.hpp>
 
-#include <CSRGraph.hpp>
+#include <mesh/ItemValue.hpp>
+#include <utils/CSRGraph.hpp>
 
 #include <memory>
 
diff --git a/src/mesh/MeshData.hpp b/src/mesh/MeshData.hpp
index 17a1ac4c38a7eee7dc90a5584f8d91f4a7413d73..2a1406c447f4eeb2d8e2b40b0ae3f229650d107e 100644
--- a/src/mesh/MeshData.hpp
+++ b/src/mesh/MeshData.hpp
@@ -1,11 +1,11 @@
 #ifndef MESH_DATA_HPP
 #define MESH_DATA_HPP
 
-#include <PugsUtils.hpp>
-#include <TinyVector.hpp>
+#include <algebra/TinyVector.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <ItemValue.hpp>
-#include <SubItemValuePerItem.hpp>
+#include <mesh/ItemValue.hpp>
+#include <mesh/SubItemValuePerItem.hpp>
 
 #include <map>
 
diff --git a/src/mesh/MeshNodeBoundary.hpp b/src/mesh/MeshNodeBoundary.hpp
index fdb771eb559a154f2a2ea9d5ceb3c127d9f9c663..8a3d4c913b96c46ca7d7b93573b41eedc0312d24 100644
--- a/src/mesh/MeshNodeBoundary.hpp
+++ b/src/mesh/MeshNodeBoundary.hpp
@@ -1,20 +1,20 @@
 #ifndef MESH_NODE_BOUNDARY_HPP
 #define MESH_NODE_BOUNDARY_HPP
 
-#include <Array.hpp>
-#include <ItemValue.hpp>
+#include <utils/Array.hpp>
 
-#include <Kokkos_Vector.hpp>
-#include <TinyVector.hpp>
+#include <algebra/TinyVector.hpp>
 
-#include <RefItemList.hpp>
+#include <mesh/ItemValue.hpp>
+#include <mesh/RefItemList.hpp>
 
-#include <ConnectivityMatrix.hpp>
-#include <IConnectivity.hpp>
-#include <iostream>
+#include <mesh/ConnectivityMatrix.hpp>
+#include <mesh/IConnectivity.hpp>
 
-#include <Exceptions.hpp>
-#include <Messenger.hpp>
+#include <utils/Exceptions.hpp>
+#include <utils/Messenger.hpp>
+
+#include <Kokkos_Vector.hpp>
 
 #include <iostream>
 
diff --git a/src/mesh/RefItemList.hpp b/src/mesh/RefItemList.hpp
index a3b9bb71ff08eb214e831cda6ce4d8e0b6115ff1..594b5ad25dc4a9769cd1d4645202cceeec6273e6 100644
--- a/src/mesh/RefItemList.hpp
+++ b/src/mesh/RefItemList.hpp
@@ -1,9 +1,10 @@
 #ifndef REF_ITEM_LIST_HPP
 #define REF_ITEM_LIST_HPP
 
-#include <Array.hpp>
-#include <ItemId.hpp>
-#include <RefId.hpp>
+#include <utils/Array.hpp>
+
+#include <mesh/ItemId.hpp>
+#include <mesh/RefId.hpp>
 
 template <ItemType item_type>
 class RefItemList
diff --git a/src/mesh/SubItemValuePerItem.hpp b/src/mesh/SubItemValuePerItem.hpp
index 132a0fb83d3d6c9a04a21f59a99e8beccb765d18..d09af2536a526318eb88541572e48e1ab5ab791a 100644
--- a/src/mesh/SubItemValuePerItem.hpp
+++ b/src/mesh/SubItemValuePerItem.hpp
@@ -3,17 +3,17 @@
 
 #include <Kokkos_StaticCrsGraph.hpp>
 
-#include <PugsAssert.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <Array.hpp>
+#include <utils/Array.hpp>
 
-#include <ItemId.hpp>
+#include <mesh/ItemId.hpp>
 
-#include <ConnectivityMatrix.hpp>
-#include <IConnectivity.hpp>
+#include <mesh/ConnectivityMatrix.hpp>
+#include <mesh/IConnectivity.hpp>
 
-#include <ItemOfItemType.hpp>
-#include <ItemType.hpp>
+#include <mesh/ItemOfItemType.hpp>
+#include <mesh/ItemType.hpp>
 
 #include <memory>
 
diff --git a/src/mesh/Synchronizer.hpp b/src/mesh/Synchronizer.hpp
index 7c2f97a6c08326c556ded90d2492fc12a22d9f8f..69b7c8b15bd530941836e366c6db8768b50d93ee 100644
--- a/src/mesh/Synchronizer.hpp
+++ b/src/mesh/Synchronizer.hpp
@@ -1,8 +1,10 @@
 #ifndef SYNCHRONIZER_HPP
 #define SYNCHRONIZER_HPP
 
-#include <Connectivity.hpp>
-#include <ItemValue.hpp>
+#include <mesh/Connectivity.hpp>
+#include <mesh/ItemValue.hpp>
+
+#include <utils/Messenger.hpp>
 
 #include <iostream>
 #include <map>
diff --git a/src/mesh/SynchronizerManager.cpp b/src/mesh/SynchronizerManager.cpp
index b02ae8300e83703dcd64e09d8cccb1e51b2060f1..8e2346c793053bd4e8799978542eda2d6e554dc7 100644
--- a/src/mesh/SynchronizerManager.cpp
+++ b/src/mesh/SynchronizerManager.cpp
@@ -1,8 +1,7 @@
-#include <PugsAssert.hpp>
-#include <SynchronizerManager.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <Messenger.hpp>
-#include <Synchronizer.hpp>
+#include <mesh/Synchronizer.hpp>
+#include <mesh/SynchronizerManager.hpp>
 
 SynchronizerManager* SynchronizerManager::m_instance{nullptr};
 
diff --git a/src/mesh/SynchronizerManager.hpp b/src/mesh/SynchronizerManager.hpp
index da01e98458a43f619932a3ad5b85ecead22e28f5..e2f89cbdad0a41cc260abcee39e602f11168ba31 100644
--- a/src/mesh/SynchronizerManager.hpp
+++ b/src/mesh/SynchronizerManager.hpp
@@ -1,8 +1,8 @@
 #ifndef SYNCHRONIZER_MANAGER_HPP
 #define SYNCHRONIZER_MANAGER_HPP
 
-#include <PugsAssert.hpp>
-#include <PugsMacros.hpp>
+#include <utils/PugsAssert.hpp>
+#include <utils/PugsMacros.hpp>
 
 #include <map>
 #include <memory>
diff --git a/src/output/OutputNamedItemValueSet.hpp b/src/output/OutputNamedItemValueSet.hpp
index a893ab31cc01c9ec6e44900db229ad08dd478160..7dcc23e3d0aadaa95cb3a4333e38b73d1edad9a5 100644
--- a/src/output/OutputNamedItemValueSet.hpp
+++ b/src/output/OutputNamedItemValueSet.hpp
@@ -1,9 +1,10 @@
 #ifndef OUTPUT_NAMED_ITEM_VALUE_SET_HPP
 #define OUTPUT_NAMED_ITEM_VALUE_SET_HPP
 
-#include <ItemValue.hpp>
-#include <TinyMatrix.hpp>
-#include <TinyVector.hpp>
+#include <mesh/ItemValue.hpp>
+
+#include <algebra/TinyMatrix.hpp>
+#include <algebra/TinyVector.hpp>
 
 #include <map>
 #include <string>
diff --git a/src/output/VTKWriter.hpp b/src/output/VTKWriter.hpp
index c48da737abf5fda84ef575c2f340a3ceb5dc8fa3..17cb61ebc8ab235d6af53fdb1ac2cdbf8804cd15 100644
--- a/src/output/VTKWriter.hpp
+++ b/src/output/VTKWriter.hpp
@@ -1,14 +1,15 @@
 #ifndef VTK_WRITER_HPP
 #define VTK_WRITER_HPP
 
-#include <IConnectivity.hpp>
-#include <TinyVector.hpp>
+#include <mesh/IConnectivity.hpp>
 
-#include <ItemValue.hpp>
-#include <Messenger.hpp>
-#include <OutputNamedItemValueSet.hpp>
+#include <algebra/TinyVector.hpp>
 
-#include <Exceptions.hpp>
+#include <mesh/ItemValue.hpp>
+#include <output/OutputNamedItemValueSet.hpp>
+#include <utils/Messenger.hpp>
+
+#include <utils/Exceptions.hpp>
 
 #include <fstream>
 #include <iomanip>
diff --git a/src/scheme/AcousticSolver.hpp b/src/scheme/AcousticSolver.hpp
index 3fb130013b121dab13ea4d639ab57bffd99e7255..1310ddc4c86389d62d17435820528f8ef22a540f 100644
--- a/src/scheme/AcousticSolver.hpp
+++ b/src/scheme/AcousticSolver.hpp
@@ -3,23 +3,25 @@
 
 #include <rang.hpp>
 
-#include <ArrayUtils.hpp>
+#include <utils/ArrayUtils.hpp>
 
-#include <BlockPerfectGas.hpp>
-#include <PugsAssert.hpp>
+#include <scheme/BlockPerfectGas.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <BoundaryCondition.hpp>
-#include <FiniteVolumesEulerUnknowns.hpp>
-#include <Mesh.hpp>
-#include <MeshData.hpp>
-#include <TinyMatrix.hpp>
-#include <TinyVector.hpp>
+#include <scheme/BoundaryCondition.hpp>
+#include <scheme/FiniteVolumesEulerUnknowns.hpp>
 
-#include <ItemValueUtils.hpp>
-#include <Messenger.hpp>
-#include <SubItemValuePerItem.hpp>
+#include <mesh/Mesh.hpp>
+#include <mesh/MeshData.hpp>
 
-#include <Exceptions.hpp>
+#include <algebra/TinyMatrix.hpp>
+#include <algebra/TinyVector.hpp>
+
+#include <mesh/ItemValueUtils.hpp>
+#include <mesh/SubItemValuePerItem.hpp>
+
+#include <utils/Exceptions.hpp>
+#include <utils/Messenger.hpp>
 
 #include <iostream>
 
diff --git a/src/scheme/BlockPerfectGas.hpp b/src/scheme/BlockPerfectGas.hpp
index 69f1c24e8064c5201e371f271503f39eebea26d2..1a51da28ed344a12c601081ff9bf61721854a631 100644
--- a/src/scheme/BlockPerfectGas.hpp
+++ b/src/scheme/BlockPerfectGas.hpp
@@ -1,7 +1,7 @@
 #ifndef BLOCK_PERFECTGAS_HPP
 #define BLOCK_PERFECTGAS_HPP
 
-#include <ItemValue.hpp>
+#include <mesh/ItemValue.hpp>
 
 struct BlockPerfectGas
 {
diff --git a/src/scheme/BoundaryCondition.hpp b/src/scheme/BoundaryCondition.hpp
index 23aaaa8995d4cda0e118200932f45d44f20065c6..9e1c65cd3b45a71cc996d419acfc17cdc33f82c2 100644
--- a/src/scheme/BoundaryCondition.hpp
+++ b/src/scheme/BoundaryCondition.hpp
@@ -1,13 +1,13 @@
 #ifndef BOUNDARY_CONDITION_HANDLER_HPP
 #define BOUNDARY_CONDITION_HANDLER_HPP
 
-#include <memory>
-#include <vector>
+#include <utils/Array.hpp>
 
-#include <Array.hpp>
+#include <mesh/MeshNodeBoundary.hpp>
+#include <mesh/RefItemList.hpp>
 
-#include <MeshNodeBoundary.hpp>
-#include <RefItemList.hpp>
+#include <memory>
+#include <vector>
 
 class BoundaryCondition
 {
diff --git a/src/scheme/BoundaryConditionDescriptor.hpp b/src/scheme/BoundaryConditionDescriptor.hpp
index 78a2ba947095672f452f4ed72d0498fcc6757321..b73cd6fca1211c2adbb48d7a177e13c28831e0f1 100644
--- a/src/scheme/BoundaryConditionDescriptor.hpp
+++ b/src/scheme/BoundaryConditionDescriptor.hpp
@@ -1,7 +1,8 @@
 #ifndef BOUNDARY_CONDITION_DESCRIPTOR_HPP
 #define BOUNDARY_CONDITION_DESCRIPTOR_HPP
 
-#include <RefId.hpp>
+#include <mesh/RefId.hpp>
+
 #include <memory>
 
 class BoundaryDescriptor
diff --git a/src/scheme/FiniteVolumesEulerUnknowns.hpp b/src/scheme/FiniteVolumesEulerUnknowns.hpp
index d433e74af51fba3998cba8c2d65b7f0ae37eca35..6997b63a05e4c0a26002977d944fff21c72d9bf3 100644
--- a/src/scheme/FiniteVolumesEulerUnknowns.hpp
+++ b/src/scheme/FiniteVolumesEulerUnknowns.hpp
@@ -1,8 +1,8 @@
 #ifndef FINITE_VOLUMES_EULER_UNKNOWNS_HPP
 #define FINITE_VOLUMES_EULER_UNKNOWNS_HPP
 
-#include <ItemValue.hpp>
-#include <TinyVector.hpp>
+#include <algebra/TinyVector.hpp>
+#include <mesh/ItemValue.hpp>
 
 template <typename TMeshData>
 class FiniteVolumesEulerUnknowns
diff --git a/src/utils/Array.hpp b/src/utils/Array.hpp
index eaf368b3e7207e6ce849ff111a0dfff6ab113c6d..2d6668c68db33aa0816f9f83ef0c95d95ede0e64 100644
--- a/src/utils/Array.hpp
+++ b/src/utils/Array.hpp
@@ -1,10 +1,10 @@
 #ifndef ARRAY_HPP
 #define ARRAY_HPP
 
-#include <PugsMacros.hpp>
-#include <PugsUtils.hpp>
+#include <utils/PugsMacros.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <PugsAssert.hpp>
+#include <utils/PugsAssert.hpp>
 
 #include <Kokkos_CopyViews.hpp>
 #include <algorithm>
diff --git a/src/utils/ArrayUtils.hpp b/src/utils/ArrayUtils.hpp
index 54f5497bc7b96b8f8b9f23644d9280e9c7fb7eb9..6b5b09685ce6266cce13fdce10930750111909da 100644
--- a/src/utils/ArrayUtils.hpp
+++ b/src/utils/ArrayUtils.hpp
@@ -1,10 +1,10 @@
 #ifndef ARRAY_UTILS_HPP
 #define ARRAY_UTILS_HPP
 
-#include <PugsMacros.hpp>
-#include <PugsUtils.hpp>
+#include <utils/PugsMacros.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <Types.hpp>
+#include <utils/Types.hpp>
 
 template <typename DataType, template <typename> typename ArrayT>
 std::remove_const_t<DataType>
diff --git a/src/utils/BacktraceManager.cpp b/src/utils/BacktraceManager.cpp
index 80d9473444436e6c20ef9e47fbc81b733806047e..a07f34d572947ef709e565a86c52a07fb6247d94 100644
--- a/src/utils/BacktraceManager.cpp
+++ b/src/utils/BacktraceManager.cpp
@@ -1,4 +1,4 @@
-#include <BacktraceManager.hpp>
+#include <utils/BacktraceManager.hpp>
 
 #include <iomanip>
 #include <rang.hpp>
diff --git a/src/utils/BuildInfo.cpp b/src/utils/BuildInfo.cpp
index 2f46fb37174f83e7cf88d2f7fd299cfc91d3919c..13b8633c2664add44d91b11b0a6d55afe79d2955 100644
--- a/src/utils/BuildInfo.cpp
+++ b/src/utils/BuildInfo.cpp
@@ -1,6 +1,6 @@
-#include <BuildInfo.hpp>
-#include <pugs_build_info.hpp>
-#include <pugs_config.hpp>
+#include <utils/BuildInfo.hpp>
+#include <utils/pugs_build_info.hpp>
+#include <utils/pugs_config.hpp>
 
 #include <sstream>
 
diff --git a/src/utils/CMakeLists.txt b/src/utils/CMakeLists.txt
index 9c4e305844acec64d11c5cbf39b77659e41215c9..21b313acd63d3396793bb67af7406443fe902344 100644
--- a/src/utils/CMakeLists.txt
+++ b/src/utils/CMakeLists.txt
@@ -1,6 +1,3 @@
-include_directories(${CMAKE_CURRENT_SOURCE_DIR})
-include_directories(${CMAKE_CURRENT_BINARY_DIR})
-
 # ------------------- Source files --------------------
 
 add_library(
@@ -73,10 +70,6 @@ list(
   ${PUGS_BINARY_DIR}/pugs_config.hpp
   )
 
-include_directories(${CMAKE_CURRENT_BINARY_DIR})
-include_directories(${PUGS_BINARY_DIR})
-include_directories(${PUGS_SOURCE_DIR}/packages/rang/include)
-
 # Additional dependencies
 add_dependencies(
   PugsUtils
diff --git a/src/utils/CSRGraph.hpp b/src/utils/CSRGraph.hpp
index 52a7989a9922bbccfc4c7e8633934d47581f4d39..f2e8ea2865e5ea69a1deda594caee2f57e207dfe 100644
--- a/src/utils/CSRGraph.hpp
+++ b/src/utils/CSRGraph.hpp
@@ -1,7 +1,7 @@
 #ifndef CSR_GRAPH_HPP
 #define CSR_GRAPH_HPP
 
-#include <Array.hpp>
+#include <utils/Array.hpp>
 
 class CSRGraph
 {
diff --git a/src/utils/CastArray.hpp b/src/utils/CastArray.hpp
index b7587abcb8a8991a88296e5c952f9126a03f5e48..26397a1605f2384efab5f4ff7fcf9bd1c78241e8 100644
--- a/src/utils/CastArray.hpp
+++ b/src/utils/CastArray.hpp
@@ -1,13 +1,12 @@
 #ifndef CAST_ARRAY_HPP
 #define CAST_ARRAY_HPP
 
-#include <Array.hpp>
-#include <PugsTraits.hpp>
+#include <utils/Array.hpp>
+#include <utils/Exceptions.hpp>
+#include <utils/PugsTraits.hpp>
 
 #include <iostream>
 
-#include <Exceptions.hpp>
-
 template <typename DataType, typename CastDataType>
 class CastArray
 {
diff --git a/src/utils/ConsoleManager.cpp b/src/utils/ConsoleManager.cpp
index 83e4c84cb628fb48e15e68c21c76042e62ddc788..aa6d23b042a9f7aa7365d6d38668183a154de800 100644
--- a/src/utils/ConsoleManager.cpp
+++ b/src/utils/ConsoleManager.cpp
@@ -1,4 +1,4 @@
-#include <ConsoleManager.hpp>
+#include <utils/ConsoleManager.hpp>
 
 #include <rang.hpp>
 
diff --git a/src/utils/Exceptions.cpp b/src/utils/Exceptions.cpp
index fa32207161371a2a231d1fb9cd411cb70475882a..6c7603e833e9fb4c6bc3fdf315c1c80dcbc78b7f 100644
--- a/src/utils/Exceptions.cpp
+++ b/src/utils/Exceptions.cpp
@@ -1,8 +1,9 @@
-#include <Exceptions.hpp>
-#include <sstream>
+#include <utils/Exceptions.hpp>
 
 #include <rang.hpp>
 
+#include <sstream>
+
 NormalError::NormalError(std::string_view error_msg)
   : std::runtime_error([&] {
       std::ostringstream os;
diff --git a/src/utils/FPEManager.cpp b/src/utils/FPEManager.cpp
index b0295f77d6b95abf42e91fcf8efc772a9329bb67..2bf268392c5a37dbf4bf058cf3f137ff6c7e13d5 100644
--- a/src/utils/FPEManager.cpp
+++ b/src/utils/FPEManager.cpp
@@ -1,7 +1,8 @@
-#include <FPEManager.hpp>
-#include <PugsMacros.hpp>
+#include <utils/FPEManager.hpp>
+
+#include <utils/PugsMacros.hpp>
+#include <utils/pugs_config.hpp>
 
-#include <pugs_config.hpp>
 #include <rang.hpp>
 
 #ifdef PUGS_HAS_FENV_H
diff --git a/src/utils/Messenger.cpp b/src/utils/Messenger.cpp
index bb79a3719166c517657c539b60ea005e06e89119..e3fff49d02f94649f641f5a273f9f607cf27edd9 100644
--- a/src/utils/Messenger.cpp
+++ b/src/utils/Messenger.cpp
@@ -1,4 +1,4 @@
-#include <Messenger.hpp>
+#include <utils/Messenger.hpp>
 
 #include <iostream>
 
diff --git a/src/utils/Messenger.hpp b/src/utils/Messenger.hpp
index 85a95711ac066fc9f91e10ba11a32ffebcd91f49..b1acd5e0f11c5c6e6640a1e049a51087231be7d1 100644
--- a/src/utils/Messenger.hpp
+++ b/src/utils/Messenger.hpp
@@ -1,24 +1,24 @@
 #ifndef MESSENGER_HPP
 #define MESSENGER_HPP
 
-#include <PugsAssert.hpp>
-#include <PugsMacros.hpp>
+#include <utils/PugsAssert.hpp>
+#include <utils/PugsMacros.hpp>
 
-#include <Array.hpp>
-#include <ArrayUtils.hpp>
-#include <CastArray.hpp>
+#include <utils/Array.hpp>
+#include <utils/ArrayUtils.hpp>
+#include <utils/CastArray.hpp>
 
 #include <type_traits>
 #include <vector>
 
-#include <pugs_config.hpp>
+#include <utils/pugs_config.hpp>
+
 #ifdef PUGS_HAS_MPI
 #include <mpi.h>
 #endif   // PUGS_HAS_MPI
 
-#include <PugsTraits.hpp>
-
-#include <Exceptions.hpp>
+#include <utils/Exceptions.hpp>
+#include <utils/PugsTraits.hpp>
 
 #include <iostream>
 
diff --git a/src/utils/Partitioner.cpp b/src/utils/Partitioner.cpp
index 2e22aa4f27fc2e74fe24538e88a81b9fdf1cdfb1..43ee64851b8d63e80afaf4d6960fe8015c1fa675 100644
--- a/src/utils/Partitioner.cpp
+++ b/src/utils/Partitioner.cpp
@@ -1,6 +1,7 @@
-#include <Messenger.hpp>
-#include <Partitioner.hpp>
-#include <pugs_config.hpp>
+#include <utils/Partitioner.hpp>
+
+#include <utils/Messenger.hpp>
+#include <utils/pugs_config.hpp>
 
 #ifdef PUGS_HAS_MPI
 
@@ -11,7 +12,7 @@
 #include <iostream>
 #include <vector>
 
-#include <Exceptions.hpp>
+#include <utils/Exceptions.hpp>
 
 Array<int>
 Partitioner::partition(const CSRGraph& graph)
diff --git a/src/utils/Partitioner.hpp b/src/utils/Partitioner.hpp
index 33e9f85f61dc2b07b530727c4e9144a5e00d01a3..30f50723acff089cd1a21fa4f289030f5bbdf7a3 100644
--- a/src/utils/Partitioner.hpp
+++ b/src/utils/Partitioner.hpp
@@ -1,7 +1,7 @@
 #ifndef PARTITIONER_HPP
 #define PARTITIONER_HPP
 
-#include <CSRGraph.hpp>
+#include <utils/CSRGraph.hpp>
 
 class Partitioner
 {
diff --git a/src/utils/PugsAssert.hpp b/src/utils/PugsAssert.hpp
index 5c0a2cc62420234b7b68aa804df4d24c80cee622..837d637612b411a5cf8fc951fa5c6ccbdf41d4df 100644
--- a/src/utils/PugsAssert.hpp
+++ b/src/utils/PugsAssert.hpp
@@ -1,7 +1,7 @@
 #ifndef PUGS_ASSERT_HPP
 #define PUGS_ASSERT_HPP
 
-#include <PugsMacros.hpp>
+#include <utils/PugsMacros.hpp>
 
 #include <iostream>
 #include <rang.hpp>
diff --git a/src/utils/PugsUtils.cpp b/src/utils/PugsUtils.cpp
index 9aa20750ee42367f09d131972a7c66a0ac20218d..ec024977f7cd21b7ebaca795f379b0eea6ad06ec 100644
--- a/src/utils/PugsUtils.cpp
+++ b/src/utils/PugsUtils.cpp
@@ -1,17 +1,17 @@
-#include <PugsUtils.hpp>
+#include <utils/PugsUtils.hpp>
 
-#include <Kokkos_Core.hpp>
+#include <utils/BuildInfo.hpp>
+#include <utils/RevisionInfo.hpp>
 
-#include <BuildInfo.hpp>
-#include <RevisionInfo.hpp>
+#include <utils/Messenger.hpp>
 
-#include <Messenger.hpp>
+#include <utils/ConsoleManager.hpp>
+#include <utils/FPEManager.hpp>
+#include <utils/SignalManager.hpp>
 
 #include <rang.hpp>
 
-#include <ConsoleManager.hpp>
-#include <FPEManager.hpp>
-#include <SignalManager.hpp>
+#include <Kokkos_Core.hpp>
 
 #include <CLI/CLI.hpp>
 
diff --git a/src/utils/PugsUtils.hpp b/src/utils/PugsUtils.hpp
index 339f5de531d694ce27f921e38f3336e863569bae..6a64867dbc3dc7c5fe21a9cf254efb704cb1655a 100644
--- a/src/utils/PugsUtils.hpp
+++ b/src/utils/PugsUtils.hpp
@@ -2,7 +2,7 @@
 #define PUGS_UTILS_HPP
 
 #include <Kokkos_Core.hpp>
-#include <PugsMacros.hpp>
+#include <utils/PugsMacros.hpp>
 
 #include <functional>
 #include <string>
diff --git a/src/utils/RevisionInfo.cpp b/src/utils/RevisionInfo.cpp
index 4d3965b33d570a329d7c8ca8e3d1cd284ccb7fcf..f31e67950593b6a8fecfaeb4e1216af0329fb4ad 100644
--- a/src/utils/RevisionInfo.cpp
+++ b/src/utils/RevisionInfo.cpp
@@ -1,6 +1,7 @@
-#include <RevisionInfo.hpp>
-#include <pugs_git_revision.hpp>
-#include <pugs_version.hpp>
+#include <utils/RevisionInfo.hpp>
+
+#include <utils/pugs_git_revision.hpp>
+#include <utils/pugs_version.hpp>
 
 std::string
 RevisionInfo::version()
diff --git a/src/utils/SignalManager.cpp b/src/utils/SignalManager.cpp
index a9981848d01aa0cb18facff197e708842d84ec19..2ff1499b4f520d96a960b1ffbaab4ac518ed942f 100644
--- a/src/utils/SignalManager.cpp
+++ b/src/utils/SignalManager.cpp
@@ -1,19 +1,17 @@
-#include <SignalManager.hpp>
+#include <utils/SignalManager.hpp>
 
-#include <PugsAssert.hpp>
+#include <utils/BacktraceManager.hpp>
+#include <utils/ConsoleManager.hpp>
+#include <utils/Messenger.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <BacktraceManager.hpp>
-#include <ConsoleManager.hpp>
+#include <utils/pugs_config.hpp>
 
-#include <csignal>
-#include <unistd.h>
-
-#include <pugs_config.hpp>
 #include <rang.hpp>
 
-#include <Messenger.hpp>
-
+#include <csignal>
 #include <iostream>
+#include <unistd.h>
 
 bool SignalManager::s_pause_on_error = false;
 
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index 79892a849a652a56020da0d7748e1f473101b187..b1670a12cb90f646e73cceccb215c912dd1b6767 100644
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -1,6 +1,9 @@
-include_directories("${PUGS_BINARY_DIR}/src/utils")
 set(EXECUTABLE_OUTPUT_PATH ${PUGS_BINARY_DIR})
 
+include_directories(${PUGS_SOURCE_DIR}/src)
+include_directories(${PUGS_BINARY_DIR}/src)
+include_directories(${PUGS_SOURCE_DIR}/test)
+
 add_executable (unit_tests
   test_main.cpp
   test_Array.cpp
diff --git a/tests/mpi_test_Messenger.cpp b/tests/mpi_test_Messenger.cpp
index 7dae0088d7453e71f4b91489a5beee1956d5c891..d7a78259521f322dc9e3ab95ff6a014d5c3de100 100644
--- a/tests/mpi_test_Messenger.cpp
+++ b/tests/mpi_test_Messenger.cpp
@@ -1,10 +1,9 @@
 #include <catch2/catch.hpp>
 
-#include <Array.hpp>
-#include <Messenger.hpp>
+#include <utils/Array.hpp>
+#include <utils/Messenger.hpp>
 
-#include <fstream>
-#include <pugs_config.hpp>
+#include <utils/pugs_config.hpp>
 
 #ifdef PUGS_HAS_MPI
 #include <mpi.h>
diff --git a/tests/mpi_test_main.cpp b/tests/mpi_test_main.cpp
index ec01b679be37c65f2bdbbf5e7a99618f814f76fe..3b40c79cf4f37a8424f7c7d4ad9fefd1171fbdc0 100644
--- a/tests/mpi_test_main.cpp
+++ b/tests/mpi_test_main.cpp
@@ -2,7 +2,8 @@
 #include <catch2/catch.hpp>
 
 #include <Kokkos_Core.hpp>
-#include <Messenger.hpp>
+
+#include <utils/Messenger.hpp>
 
 #include <cstdlib>
 
diff --git a/tests/test_Array.cpp b/tests/test_Array.cpp
index e9f7238ef256bbc4f53901153902c9060722cc72..948cfb554216c41ea62a84c3b7632d986d167602 100644
--- a/tests/test_Array.cpp
+++ b/tests/test_Array.cpp
@@ -1,8 +1,8 @@
 #include <catch2/catch.hpp>
 
-#include <Array.hpp>
-#include <PugsAssert.hpp>
-#include <Types.hpp>
+#include <utils/Array.hpp>
+#include <utils/PugsAssert.hpp>
+#include <utils/Types.hpp>
 
 #include <deque>
 #include <list>
diff --git a/tests/test_ArrayUtils.cpp b/tests/test_ArrayUtils.cpp
index 0f54e23a0d5a1de69f22cab9d252cfc915389067..6a7774e1de458a535db66af5000150fad820609e 100644
--- a/tests/test_ArrayUtils.cpp
+++ b/tests/test_ArrayUtils.cpp
@@ -1,11 +1,11 @@
 #include <catch2/catch.hpp>
 
-#include <Array.hpp>
-#include <ArrayUtils.hpp>
-#include <PugsAssert.hpp>
+#include <utils/Array.hpp>
+#include <utils/ArrayUtils.hpp>
+#include <utils/PugsAssert.hpp>
 
-#include <TinyMatrix.hpp>
-#include <TinyVector.hpp>
+#include <algebra/TinyMatrix.hpp>
+#include <algebra/TinyVector.hpp>
 
 // Instantiate to ensure full coverage is performed
 template class Array<int>;
diff --git a/tests/test_BiCGStab.cpp b/tests/test_BiCGStab.cpp
index 8769200f65f7b0f393a1fde3838a68307e798730..a40fa8a4f9f9607ce7eea53b2c4a1c5685a4281d 100644
--- a/tests/test_BiCGStab.cpp
+++ b/tests/test_BiCGStab.cpp
@@ -1,7 +1,7 @@
 #include <catch2/catch.hpp>
 
-#include <BiCGStab.hpp>
-#include <CRSMatrix.hpp>
+#include <algebra/BiCGStab.hpp>
+#include <algebra/CRSMatrix.hpp>
 
 TEST_CASE("BiCGStab", "[algebra]")
 {
diff --git a/tests/test_CRSMatrix.cpp b/tests/test_CRSMatrix.cpp
index 4d555b2810467edc761a1758299887c2d851a1b7..651eb2350aa0fe95a4895a7ee257cd303486325e 100644
--- a/tests/test_CRSMatrix.cpp
+++ b/tests/test_CRSMatrix.cpp
@@ -1,6 +1,6 @@
 #include <catch2/catch.hpp>
 
-#include <CRSMatrix.hpp>
+#include <algebra/CRSMatrix.hpp>
 
 // Instantiate to ensure full coverage is performed
 template class SparseMatrixDescriptor<int, uint8_t>;
diff --git a/tests/test_ItemType.cpp b/tests/test_ItemType.cpp
index 8560f8790c3c380503f18189a609ca219a545438..88e3ad877024f9c4ea4bd3304d8d3414f438abe0 100644
--- a/tests/test_ItemType.cpp
+++ b/tests/test_ItemType.cpp
@@ -1,7 +1,7 @@
-#include <PugsMacros.hpp>
 #include <catch2/catch.hpp>
 
-#include <ItemType.hpp>
+#include <mesh/ItemType.hpp>
+#include <utils/PugsMacros.hpp>
 
 TEST_CASE("ItemType", "[connectivity]")
 {
diff --git a/tests/test_PCG.cpp b/tests/test_PCG.cpp
index 97bbe1ce25076d05f24b3e250d779d799012beab..d0939d327602f440af82d449c64fb329afa8a662 100644
--- a/tests/test_PCG.cpp
+++ b/tests/test_PCG.cpp
@@ -1,7 +1,7 @@
 #include <catch2/catch.hpp>
 
-#include <CRSMatrix.hpp>
-#include <PCG.hpp>
+#include <algebra/CRSMatrix.hpp>
+#include <algebra/PCG.hpp>
 
 TEST_CASE("PCG", "[algebra]")
 {
diff --git a/tests/test_PugsAssert.cpp b/tests/test_PugsAssert.cpp
index ccbc472034ebab1f511235efe3a967b44b18a39a..9f155fc631173a0a4c6b2ca5cfffb2f55df3d52b 100644
--- a/tests/test_PugsAssert.cpp
+++ b/tests/test_PugsAssert.cpp
@@ -1,6 +1,7 @@
 #include <catch2/catch.hpp>
 
-#include <PugsAssert.hpp>
+#include <utils/PugsAssert.hpp>
+
 #include <string>
 
 TEST_CASE("PugsAssert", "[utils]")
diff --git a/tests/test_RevisionInfo.cpp b/tests/test_RevisionInfo.cpp
index 8bf2cf49ff218e6ebd91cd7a39f34522f2fd0108..e3becb6cbd3c7b60b2c92ccc7e24ae82573d217e 100644
--- a/tests/test_RevisionInfo.cpp
+++ b/tests/test_RevisionInfo.cpp
@@ -1,9 +1,9 @@
 #include <catch2/catch.hpp>
 
-#include <RevisionInfo.hpp>
+#include <utils/RevisionInfo.hpp>
 
-#include <pugs_git_revision.hpp>
-#include <pugs_version.hpp>
+#include <utils/pugs_git_revision.hpp>
+#include <utils/pugs_version.hpp>
 
 TEST_CASE("RevisionInfo", "[utils]")
 {
diff --git a/tests/test_SparseMatrixDescriptor.cpp b/tests/test_SparseMatrixDescriptor.cpp
index c34352f61cfdf81b7e9c2b3e3753562b29a94de0..039ba08242013e339afb1a1ba37725a319603cb0 100644
--- a/tests/test_SparseMatrixDescriptor.cpp
+++ b/tests/test_SparseMatrixDescriptor.cpp
@@ -1,7 +1,8 @@
 #include <catch2/catch.hpp>
 
-#include <PugsAssert.hpp>
-#include <SparseMatrixDescriptor.hpp>
+#include <utils/PugsAssert.hpp>
+
+#include <algebra/SparseMatrixDescriptor.hpp>
 
 // Instantiate to ensure full coverage is performed
 template class SparseMatrixDescriptor<int, uint8_t>;
diff --git a/tests/test_TinyMatrix.cpp b/tests/test_TinyMatrix.cpp
index fe45b86e0da57e884c5c0e70d4c6417aecb66edd..d764a88fc83ccbfb683a7037cece1042eea0f265 100644
--- a/tests/test_TinyMatrix.cpp
+++ b/tests/test_TinyMatrix.cpp
@@ -1,9 +1,11 @@
 #include <catch2/catch.hpp>
 
 #include <Kokkos_Core.hpp>
-#include <PugsAssert.hpp>
-#include <TinyMatrix.hpp>
-#include <Types.hpp>
+
+#include <utils/PugsAssert.hpp>
+#include <utils/Types.hpp>
+
+#include <algebra/TinyMatrix.hpp>
 
 // Instantiate to ensure full coverage is performed
 template class TinyMatrix<1, int>;
diff --git a/tests/test_TinyVector.cpp b/tests/test_TinyVector.cpp
index 79ffb373e550921c69a08eef02178018be573a63..f64e2da09b410a53c1549efd09602b959beef2ce 100644
--- a/tests/test_TinyVector.cpp
+++ b/tests/test_TinyVector.cpp
@@ -1,7 +1,8 @@
 #include <catch2/catch.hpp>
 
-#include <PugsAssert.hpp>
-#include <TinyVector.hpp>
+#include <utils/PugsAssert.hpp>
+
+#include <algebra/TinyVector.hpp>
 
 // Instantiate to ensure full coverage is performed
 template class TinyVector<1, int>;
diff --git a/tests/test_Vector.cpp b/tests/test_Vector.cpp
index 15e02d0eb2c7f9c0ad28af2cdb365ed790bd563f..9b6730f0576233cd8c155a7328fbe6172a049e51 100644
--- a/tests/test_Vector.cpp
+++ b/tests/test_Vector.cpp
@@ -1,7 +1,8 @@
 #include <catch2/catch.hpp>
 
-#include <PugsAssert.hpp>
-#include <Vector.hpp>
+#include <utils/PugsAssert.hpp>
+
+#include <algebra/Vector.hpp>
 
 // Instantiate to ensure full coverage is performed
 template class Vector<int>;