diff --git a/Applications/ApplicationsLib/LinearSolverLibrarySetup.h b/Applications/ApplicationsLib/LinearSolverLibrarySetup.h
index f6d6d760554c5adb61b01cbf5c280e87e975de2a..555d1a2a17547ada4e0bea6b811b124342aa1b2f 100644
--- a/Applications/ApplicationsLib/LinearSolverLibrarySetup.h
+++ b/Applications/ApplicationsLib/LinearSolverLibrarySetup.h
@@ -19,7 +19,7 @@
 /// The default implementation is empty providing polymorphic behaviour when
 /// using this class.
 
-#include "AssemblerLib/GlobalMatrixProviders.h"
+#include "NumLib/DOF/GlobalMatrixProviders.h"
 
 #if defined(USE_PETSC)
 #include <petsc.h>
diff --git a/Applications/CLI/CMakeLists.txt b/Applications/CLI/CMakeLists.txt
index 6329b9e1a711d7eb8fe741dbf5280ff6890ba468..085155e363ad21aab906d779aad0feefc0f4ba10 100644
--- a/Applications/CLI/CMakeLists.txt
+++ b/Applications/CLI/CMakeLists.txt
@@ -5,7 +5,6 @@ add_executable(ogs
 
 target_link_libraries(ogs
     ApplicationsLib
-    AssemblerLib
     MeshGeoToolsLib
     ProcessLib
     NumLib
diff --git a/AssemblerLib/CMakeLists.txt b/AssemblerLib/CMakeLists.txt
deleted file mode 100644
index 8bb4b5be726702cce9fb160e006e9971538cdc75..0000000000000000000000000000000000000000
--- a/AssemblerLib/CMakeLists.txt
+++ /dev/null
@@ -1,19 +0,0 @@
-include(${PROJECT_SOURCE_DIR}/scripts/cmake/OGSEnabledElements.cmake)
-
-#Source files grouped by a directory
-GET_SOURCE_FILES(SOURCES_ASSEMBLERLIB)
-set(SOURCES ${SOURCES_ASSEMBLERLIB})
-
-# Create the library
-add_library(AssemblerLib STATIC ${SOURCES})
-
-target_link_libraries(AssemblerLib
-    MeshLib
-)
-
-if(TARGET Eigen)
-    add_dependencies(AssemblerLib Eigen)
-endif()
-if(TARGET Boost)
-    add_dependencies(AssemblerLib Boost)
-endif()
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 1e7efe9521c49e7337e713fe809691426b254bea..073b3b1fc2018958f387aefdc38db74a2c320dd3 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -188,7 +188,6 @@ include_directories( SYSTEM ${CMAKE_CURRENT_SOURCE_DIR}/ThirdParty/tclap/include
 include_directories(SYSTEM ${CMAKE_CURRENT_SOURCE_DIR}/ThirdParty/vtkGUISupportQt)
 
 add_subdirectory( Applications )
-add_subdirectory( AssemblerLib )
 add_subdirectory( BaseLib )
 # TODO This is a hack but we have to make sure that Boost is built first
 if(TARGET Boost)
diff --git a/MathLib/CMakeLists.txt b/MathLib/CMakeLists.txt
index dd6028af94d6df408ae4e2d1dd041ca4c1252a3d..86095a2f16a2454459faa7b70e4f13b3a4fd7634 100644
--- a/MathLib/CMakeLists.txt
+++ b/MathLib/CMakeLists.txt
@@ -49,10 +49,6 @@ add_library(MathLib STATIC ${SOURCES})
 
 set_target_properties(MathLib PROPERTIES LINKER_LANGUAGE CXX)
 
-target_link_libraries(MathLib
-    AssemblerLib
-)
-
 if (CVODE_FOUND)
     target_link_libraries(MathLib ${CVODE_LIBRARIES})
 endif()
diff --git a/AssemblerLib/SerialExecutor.h b/NumLib/Assembler/SerialExecutor.h
similarity index 95%
rename from AssemblerLib/SerialExecutor.h
rename to NumLib/Assembler/SerialExecutor.h
index 37dd1bd1d7917a1effdc968d8c5bfd3ad9f822a8..205ec3513db735f731d229fc301c5f593df1d6f8 100644
--- a/AssemblerLib/SerialExecutor.h
+++ b/NumLib/Assembler/SerialExecutor.h
@@ -10,10 +10,10 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_SERIALEXECUTOR_H_H
-#define ASSEMBLERLIB_SERIALEXECUTOR_H_H
+#ifndef NUMLIB_SERIALEXECUTOR_H_H
+#define NUMLIB_SERIALEXECUTOR_H_H
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 struct SerialExecutor
@@ -107,6 +107,6 @@ struct SerialExecutor
     }
 };
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
 
-#endif  // ASSEMBLERLIB_SERIALEXECUTOR_H_H
+#endif  // NUMLIB_SERIALEXECUTOR_H_H
diff --git a/AssemblerLib/VectorMatrixAssembler.h b/NumLib/Assembler/VectorMatrixAssembler.h
similarity index 93%
rename from AssemblerLib/VectorMatrixAssembler.h
rename to NumLib/Assembler/VectorMatrixAssembler.h
index bc920630a96f5daf36a04a370486445abbc575d1..1a3a4e063a9e54e64e37de32d39be2f79f776a19 100644
--- a/AssemblerLib/VectorMatrixAssembler.h
+++ b/NumLib/Assembler/VectorMatrixAssembler.h
@@ -7,18 +7,17 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_VECTORMATRIXASSEMBLER_H_
-#define ASSEMBLERLIB_VECTORMATRIXASSEMBLER_H_
-
-#include "LocalToGlobalIndexMap.h"
+#ifndef NUMLIB_VECTORMATRIXASSEMBLER_H_
+#define NUMLIB_VECTORMATRIXASSEMBLER_H_
 
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
 #include "NumLib/ODESolver/Types.h"
 
 namespace
 {
-inline AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices
+inline NumLib::LocalToGlobalIndexMap::RowColumnIndices
 getRowColumnIndices(std::size_t const id,
-                    AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+                    NumLib::LocalToGlobalIndexMap const& dof_table,
                     std::vector<GlobalIndexType>& indices)
 {
     assert(dof_table.size() > id);
@@ -33,13 +32,13 @@ getRowColumnIndices(std::size_t const id,
         indices.insert(indices.end(), idcs.begin(), idcs.end());
     }
 
-    return AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices(indices,
+    return NumLib::LocalToGlobalIndexMap::RowColumnIndices(indices,
                                                                  indices);
 }
 
 template <typename Callback, typename GlobalVector, typename... Args>
 void passLocalVector_(Callback& cb, std::size_t const id,
-                      AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+                      NumLib::LocalToGlobalIndexMap const& dof_table,
                       GlobalVector const& x, Args&&... args)
 {
     std::vector<GlobalIndexType> indices;
@@ -57,7 +56,7 @@ void passLocalVector_(Callback& cb, std::size_t const id,
 }
 }
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 /*! Calls the local assemblers of FEM processes and assembles
@@ -215,6 +214,6 @@ private:
     LocalToGlobalIndexMap const& _data_pos;
 };
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
 
-#endif  // ASSEMBLERLIB_VECTORMATRIXASSEMBLER_H_
+#endif  // NUMLIB_VECTORMATRIXASSEMBLER_H_
diff --git a/AssemblerLib/VectorMatrixBuilder.h b/NumLib/Assembler/VectorMatrixBuilder.h
similarity index 85%
rename from AssemblerLib/VectorMatrixBuilder.h
rename to NumLib/Assembler/VectorMatrixBuilder.h
index a2886f0be900e26c26a8bea43d57c1e4025cbd75..fbba5d69459e353b7db87e12d93097d92a1c7f55 100644
--- a/AssemblerLib/VectorMatrixBuilder.h
+++ b/NumLib/Assembler/VectorMatrixBuilder.h
@@ -10,10 +10,10 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_VECTORMATRIXBUILDER_H_
-#define ASSEMBLERLIB_VECTORMATRIXBUILDER_H_
+#ifndef NUMLIB_VECTORMATRIXBUILDER_H_
+#define NUMLIB_VECTORMATRIXBUILDER_H_
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 template <typename MatrixType_, typename VectorType_>
@@ -44,6 +44,6 @@ public:
 
 };
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
 
-#endif  // ASSEMBLERLIB_VECTORMATRIXBUILDER_H_
+#endif  // NUMLIB_VECTORMATRIXBUILDER_H_
diff --git a/NumLib/CMakeLists.txt b/NumLib/CMakeLists.txt
index 7d5a00bb7cf2e6d9777d5443730ab8ab6e9f808a..b23206f9fd9053f15e8aea5b150bd607265fee34 100644
--- a/NumLib/CMakeLists.txt
+++ b/NumLib/CMakeLists.txt
@@ -2,6 +2,11 @@
 GET_SOURCE_FILES(SOURCES_NUMLIB)
 set(SOURCES ${SOURCES_NUMLIB})
 
+GET_SOURCE_FILES(SOURCES_ASSEMBLER Assembler)
+set(SOURCES ${SOURCES} ${SOURCES_ASSEMBLER})
+GET_SOURCE_FILES(SOURCES_DOF DOF)
+set(SOURCES ${SOURCES} ${SOURCES_DOF})
+
 GET_SOURCE_FILES(SOURCES_FEM Fem)
 set(SOURCES ${SOURCES} ${SOURCES_FEM})
 GET_SOURCE_FILES(SOURCES_FEM_COORDINATESMAPPING Fem/CoordinatesMapping)
diff --git a/AssemblerLib/ComponentGlobalIndexDict.h b/NumLib/DOF/ComponentGlobalIndexDict.h
similarity index 93%
rename from AssemblerLib/ComponentGlobalIndexDict.h
rename to NumLib/DOF/ComponentGlobalIndexDict.h
index 21ff6dd89aaf1f62ec0e465d2184edffd3bc27bd..d39598eebcc5941598e25ce907ad056eb1294767 100644
--- a/AssemblerLib/ComponentGlobalIndexDict.h
+++ b/NumLib/DOF/ComponentGlobalIndexDict.h
@@ -10,8 +10,8 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_COMPONENTGLOBALINDEXDICT_H_
-#define ASSEMBLERLIB_COMPONENTGLOBALINDEXDICT_H_
+#ifndef NUMLIB_COMPONENTGLOBALINDEXDICT_H_
+#define NUMLIB_COMPONENTGLOBALINDEXDICT_H_
 
 #include <limits>
 
@@ -23,7 +23,7 @@
 #include "MeshLib/Location.h"
 #include "ProcessLib/NumericsConfig.h"
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 /// \internal
@@ -118,6 +118,6 @@ typedef boost::multi_index::multi_index_container<
     > ComponentGlobalIndexDict;
 
 }    // namespace detail
-}    // namespace AssemblerLib
+}    // namespace NumLib
 
-#endif  // ASSEMBLERLIB_COMPONENTGLOBALINDEXDICT_H_
+#endif  // NUMLIB_COMPONENTGLOBALINDEXDICT_H_
diff --git a/AssemblerLib/ComputeSparsityPattern.cpp b/NumLib/DOF/ComputeSparsityPattern.cpp
similarity index 98%
rename from AssemblerLib/ComputeSparsityPattern.cpp
rename to NumLib/DOF/ComputeSparsityPattern.cpp
index b357ac594e7ba786f8633effd052df496719b0e7..e9b3b97e7385dc6096b25db914e02c5cee92d8f8 100644
--- a/AssemblerLib/ComputeSparsityPattern.cpp
+++ b/NumLib/DOF/ComputeSparsityPattern.cpp
@@ -12,7 +12,7 @@
 #include "LocalToGlobalIndexMap.h"
 #include "MeshLib/NodeAdjacencyTable.h"
 
-namespace AssemblerLib
+namespace NumLib
 {
 MathLib::SparsityPattern computeSparsityPattern(
     LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh)
diff --git a/AssemblerLib/ComputeSparsityPattern.h b/NumLib/DOF/ComputeSparsityPattern.h
similarity index 83%
rename from AssemblerLib/ComputeSparsityPattern.h
rename to NumLib/DOF/ComputeSparsityPattern.h
index f4364c526c44fc387d69a19ef07f0bee08927702..ed0fc431f3924a79b712d38bb828c029bf13c932 100644
--- a/AssemblerLib/ComputeSparsityPattern.h
+++ b/NumLib/DOF/ComputeSparsityPattern.h
@@ -7,8 +7,8 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_COMPUTESPARSITYPATTERN_H
-#define ASSEMBLERLIB_COMPUTESPARSITYPATTERN_H
+#ifndef NUMLIB_COMPUTESPARSITYPATTERN_H
+#define NUMLIB_COMPUTESPARSITYPATTERN_H
 
 #include <vector>
 
@@ -20,7 +20,7 @@ namespace MeshLib
 class Mesh;
 }
 
-namespace AssemblerLib
+namespace NumLib
 {
 class LocalToGlobalIndexMap;
 
@@ -36,5 +36,5 @@ MathLib::SparsityPattern computeSparsityPattern(
     LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh);
 }
 
-#endif // ASSEMBLERLIB_COMPUTESPARSITYPATTERN_H
+#endif // NUMLIB_COMPUTESPARSITYPATTERN_H
 
diff --git a/AssemblerLib/GlobalMatrixProviders.cpp b/NumLib/DOF/GlobalMatrixProviders.cpp
similarity index 100%
rename from AssemblerLib/GlobalMatrixProviders.cpp
rename to NumLib/DOF/GlobalMatrixProviders.cpp
diff --git a/AssemblerLib/GlobalMatrixProviders.h b/NumLib/DOF/GlobalMatrixProviders.h
similarity index 100%
rename from AssemblerLib/GlobalMatrixProviders.h
rename to NumLib/DOF/GlobalMatrixProviders.h
diff --git a/AssemblerLib/LocalToGlobalIndexMap.cpp b/NumLib/DOF/LocalToGlobalIndexMap.cpp
similarity index 97%
rename from AssemblerLib/LocalToGlobalIndexMap.cpp
rename to NumLib/DOF/LocalToGlobalIndexMap.cpp
index e1c55ac7ec89686083799f397fa1bf6dab33fcf8..27432950335095b3ea46762e9debad64452f96ff 100644
--- a/AssemblerLib/LocalToGlobalIndexMap.cpp
+++ b/NumLib/DOF/LocalToGlobalIndexMap.cpp
@@ -9,7 +9,7 @@
 
 #include "LocalToGlobalIndexMap.h"
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 template <typename ElementIterator>
@@ -46,7 +46,7 @@ LocalToGlobalIndexMap::findGlobalIndices(
 
 LocalToGlobalIndexMap::LocalToGlobalIndexMap(
     std::vector<std::unique_ptr<MeshLib::MeshSubsets>>&& mesh_subsets,
-    AssemblerLib::ComponentOrder const order)
+    NumLib::ComponentOrder const order)
     : _mesh_subsets(std::move(mesh_subsets)),
       _mesh_component_map(_mesh_subsets, order)
 {
@@ -71,7 +71,7 @@ LocalToGlobalIndexMap::LocalToGlobalIndexMap(
     std::vector<std::unique_ptr<MeshLib::MeshSubsets>>&& mesh_subsets,
     std::size_t const component_id,
     std::vector<MeshLib::Element*> const& elements,
-    AssemblerLib::MeshComponentMap&& mesh_component_map)
+    NumLib::MeshComponentMap&& mesh_component_map)
     : _mesh_subsets(std::move(mesh_subsets)),
       _mesh_component_map(std::move(mesh_component_map))
 {
@@ -176,4 +176,4 @@ std::ostream& operator<<(std::ostream& os, LocalToGlobalIndexMap const& map)
 }
 #endif  // NDEBUG
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
diff --git a/AssemblerLib/LocalToGlobalIndexMap.h b/NumLib/DOF/LocalToGlobalIndexMap.h
similarity index 93%
rename from AssemblerLib/LocalToGlobalIndexMap.h
rename to NumLib/DOF/LocalToGlobalIndexMap.h
index 657bf8bcc6a3a3ab5018d9f7fdb2bbeb7d11dd7f..9932f447094452c59577fdbd6acf1df5a505d946 100644
--- a/AssemblerLib/LocalToGlobalIndexMap.h
+++ b/NumLib/DOF/LocalToGlobalIndexMap.h
@@ -7,8 +7,8 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_LOCALTOGLOBALINDEXMAP_H_
-#define ASSEMBLERLIB_LOCALTOGLOBALINDEXMAP_H_
+#ifndef NUMLIB_LOCALTOGLOBALINDEXMAP_H_
+#define NUMLIB_LOCALTOGLOBALINDEXMAP_H_
 
 #ifndef NDEBUG
 #include <iosfwd>
@@ -18,11 +18,12 @@
 
 #include <Eigen/Dense>
 
-#include "AssemblerLib/MeshComponentMap.h"
 #include "MathLib/LinAlg/RowColumnIndices.h"
 #include "MeshLib/MeshSubsets.h"
 
-namespace AssemblerLib
+#include "MeshComponentMap.h"
+
+namespace NumLib
 {
 
 /// Row and column indices in global linear algebra objects for each mesh item.
@@ -45,7 +46,7 @@ public:
     /// each mesh element of the given mesh_subsets.
     explicit LocalToGlobalIndexMap(
         std::vector<std::unique_ptr<MeshLib::MeshSubsets>>&& mesh_subsets,
-        AssemblerLib::ComponentOrder const order);
+        NumLib::ComponentOrder const order);
 
     /// Derive a LocalToGlobalIndexMap constrained to a set of mesh subsets and
     /// elements. A new mesh component map will be constructed using the passed
@@ -124,7 +125,7 @@ private:
         std::vector<std::unique_ptr<MeshLib::MeshSubsets>>&& mesh_subsets,
         std::size_t const component_id,
         std::vector<MeshLib::Element*> const& elements,
-        AssemblerLib::MeshComponentMap&& mesh_component_map);
+        NumLib::MeshComponentMap&& mesh_component_map);
 
     template <typename ElementIterator>
     void
@@ -135,7 +136,7 @@ private:
 private:
     /// A vector of mesh subsets for each process variables' components.
     std::vector<std::unique_ptr<MeshLib::MeshSubsets>> const _mesh_subsets;
-    AssemblerLib::MeshComponentMap _mesh_component_map;
+    NumLib::MeshComponentMap _mesh_component_map;
 
     using Table = Eigen::Matrix<LineIndex, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>;
 
@@ -153,6 +154,6 @@ private:
 
 };
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
 
-#endif  // ASSEMBLERLIB_LOCALTOGLOBALINDEXMAP_H_
+#endif  // NUMLIB_LOCALTOGLOBALINDEXMAP_H_
diff --git a/AssemblerLib/MatrixProviderUser.h b/NumLib/DOF/MatrixProviderUser.h
similarity index 95%
rename from AssemblerLib/MatrixProviderUser.h
rename to NumLib/DOF/MatrixProviderUser.h
index 0f0e4b7f0d5a052a3681c71a20b3ab9bd9d3e339..5c90cd2d061c770654833646e77b641c86801c84 100644
--- a/AssemblerLib/MatrixProviderUser.h
+++ b/NumLib/DOF/MatrixProviderUser.h
@@ -14,7 +14,7 @@
 
 #include "MathLib/LinAlg/SparsityPattern.h"
 
-namespace AssemblerLib { class LocalToGlobalIndexMap; }
+namespace NumLib { class LocalToGlobalIndexMap; }
 namespace MeshLib { class Mesh; }
 
 namespace MathLib
@@ -24,7 +24,7 @@ struct MatrixSpecifications
 {
     MatrixSpecifications(std::size_t const nrows_, std::size_t const ncols_,
                          SparsityPattern const*const sparsity_pattern_,
-                         AssemblerLib::LocalToGlobalIndexMap const*const dof_table_,
+                         NumLib::LocalToGlobalIndexMap const*const dof_table_,
                          MeshLib::Mesh const*const mesh_)
         : nrows(nrows_), ncols(ncols_), sparsity_pattern(sparsity_pattern_)
         , dof_table(dof_table_), mesh(mesh_)
@@ -33,7 +33,7 @@ struct MatrixSpecifications
     std::size_t const nrows;
     std::size_t const ncols;
     SparsityPattern const*const sparsity_pattern;
-    AssemblerLib::LocalToGlobalIndexMap const*const dof_table;
+    NumLib::LocalToGlobalIndexMap const*const dof_table;
     MeshLib::Mesh const*const mesh;
 };
 
diff --git a/AssemblerLib/MatrixVectorTraits.cpp b/NumLib/DOF/MatrixVectorTraits.cpp
similarity index 97%
rename from AssemblerLib/MatrixVectorTraits.cpp
rename to NumLib/DOF/MatrixVectorTraits.cpp
index d112e9791a2966b2b6441c9ed80d9a18468bbf94..2cdd4d74e1c16d122ccc071f44424ba8cc1ac55f 100644
--- a/AssemblerLib/MatrixVectorTraits.cpp
+++ b/NumLib/DOF/MatrixVectorTraits.cpp
@@ -7,7 +7,7 @@
  *
  */
 
-#include "AssemblerLib/LocalToGlobalIndexMap.h"
+#include "LocalToGlobalIndexMap.h"
 #include "MatrixProviderUser.h"
 #include "MatrixVectorTraits.h"
 
@@ -97,7 +97,7 @@ newInstance(MatrixSpecifications const& spec)
         spec.dof_table ? spec.dof_table->dofSizeWithoutGhosts() : spec.nrows;
     auto const ncols = spec.dof_table ? nrows : spec.ncols;
 
-    // TODO I guess it is not hard to make AssemblerLib::computeSparsityPattern()
+    // TODO I guess it is not hard to make NumLib::computeSparsityPattern()
     //      work also for NodePartitionedMesh'es.
     assert(((!spec.sparsity_pattern) || spec.sparsity_pattern->size() == 0u) &&
            "The sparsity pattern is not used with PETSc so I rather crash than"
diff --git a/AssemblerLib/MatrixVectorTraits.h b/NumLib/DOF/MatrixVectorTraits.h
similarity index 100%
rename from AssemblerLib/MatrixVectorTraits.h
rename to NumLib/DOF/MatrixVectorTraits.h
diff --git a/AssemblerLib/MeshComponentMap.cpp b/NumLib/DOF/MeshComponentMap.cpp
similarity index 99%
rename from AssemblerLib/MeshComponentMap.cpp
rename to NumLib/DOF/MeshComponentMap.cpp
index 0c6dccd591315e137115c8f5a63c64767a29d26b..cb6ab849ca672e6183f3bc50186fef5406955785 100644
--- a/AssemblerLib/MeshComponentMap.cpp
+++ b/NumLib/DOF/MeshComponentMap.cpp
@@ -18,7 +18,7 @@
 #include "MeshLib/NodePartitionedMesh.h"
 #endif
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 using namespace detail;
@@ -321,4 +321,4 @@ GlobalIndexType MeshComponentMap::getLocalIndex(
 #endif
 }
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
diff --git a/AssemblerLib/MeshComponentMap.h b/NumLib/DOF/MeshComponentMap.h
similarity index 97%
rename from AssemblerLib/MeshComponentMap.h
rename to NumLib/DOF/MeshComponentMap.h
index 221337754d8d190b635e0004cac8e871c50b68cb..4e6e2e71d6d21911cfc007ba5c1e66a52e8a8fd0 100644
--- a/AssemblerLib/MeshComponentMap.h
+++ b/NumLib/DOF/MeshComponentMap.h
@@ -10,8 +10,8 @@
  *
  */
 
-#ifndef ASSEMBLERLIB_MESHCOMPONENTMAP_H_
-#define ASSEMBLERLIB_MESHCOMPONENTMAP_H_
+#ifndef NUMLIB_MESHCOMPONENTMAP_H_
+#define NUMLIB_MESHCOMPONENTMAP_H_
 
 #include "ComponentGlobalIndexDict.h"
 
@@ -20,7 +20,7 @@ namespace MeshLib
     class MeshSubsets;
 }
 
-namespace AssemblerLib
+namespace NumLib
 {
 
 /// Ordering of components in global matrix/vector.
@@ -190,6 +190,6 @@ private:
     std::vector<GlobalIndexType> _ghosts_indices;
 };
 
-}   // namespace AssemblerLib
+}   // namespace NumLib
 
-#endif  // ASSEMBLERLIB_MESHCOMPONENTMAP_H_
+#endif  // NUMLIB_MESHCOMPONENTMAP_H_
diff --git a/AssemblerLib/SimpleMatrixVectorProvider-impl.h b/NumLib/DOF/SimpleMatrixVectorProvider-impl.h
similarity index 100%
rename from AssemblerLib/SimpleMatrixVectorProvider-impl.h
rename to NumLib/DOF/SimpleMatrixVectorProvider-impl.h
diff --git a/AssemblerLib/SimpleMatrixVectorProvider.h b/NumLib/DOF/SimpleMatrixVectorProvider.h
similarity index 100%
rename from AssemblerLib/SimpleMatrixVectorProvider.h
rename to NumLib/DOF/SimpleMatrixVectorProvider.h
diff --git a/AssemblerLib/UnifiedMatrixSetters.cpp b/NumLib/DOF/UnifiedMatrixSetters.cpp
similarity index 100%
rename from AssemblerLib/UnifiedMatrixSetters.cpp
rename to NumLib/DOF/UnifiedMatrixSetters.cpp
diff --git a/AssemblerLib/UnifiedMatrixSetters.h b/NumLib/DOF/UnifiedMatrixSetters.h
similarity index 100%
rename from AssemblerLib/UnifiedMatrixSetters.h
rename to NumLib/DOF/UnifiedMatrixSetters.h
diff --git a/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator-impl.h b/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator-impl.h
index 340ac3ccc6ed8712f854f810008cb27d9c9efc18..c8c35472679456a133338d3d011ef9e534e102f8 100644
--- a/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator-impl.h
+++ b/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator-impl.h
@@ -12,9 +12,9 @@
 #include <logog/include/logog.hpp>
 #include <Eigen/Core>
 
-#include "AssemblerLib/MatrixVectorTraits.h"
-#include "AssemblerLib/SerialExecutor.h"
 #include "MathLib/LinAlg/BLAS.h"
+#include "NumLib/Assembler/SerialExecutor.h"
+#include "NumLib/DOF/MatrixVectorTraits.h"
 #include "NumLib/Function/Interpolation.h"
 #include "LocalLinearLeastSquaresExtrapolator.h"
 
@@ -37,7 +37,7 @@ extrapolate(LocalAssemblers const& local_assemblers, PropertyTag const property)
     using Self = LocalLinearLeastSquaresExtrapolator<
         GlobalVector, PropertyTag, LocalAssembler>;
 
-    AssemblerLib::SerialExecutor::executeMemberDereferenced(
+    NumLib::SerialExecutor::executeMemberDereferenced(
         *this, &Self::extrapolateElement, local_assemblers, property, *counts);
 
     MathLib::BLAS::componentwiseDivide(_nodal_values, _nodal_values, *counts);
@@ -54,7 +54,7 @@ calculateResiduals(LocalAssemblers const& local_assemblers,
     using Self = LocalLinearLeastSquaresExtrapolator<
         GlobalVector, PropertyTag, LocalAssembler>;
 
-    AssemblerLib::SerialExecutor::executeMemberDereferenced(
+    NumLib::SerialExecutor::executeMemberDereferenced(
         *this, &Self::calculateResiudalElement, local_assemblers, property);
 }
 
diff --git a/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h b/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h
index 22a65316e4be595200e69ec564cb5c4b7716913e..b15c46ab4294d40104c54b5acd4ba5cdf332841d 100644
--- a/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h
+++ b/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h
@@ -10,7 +10,7 @@
 #ifndef NUMLIB_LOCAL_LLSQ_EXTRAPOLATOR_H
 #define NUMLIB_LOCAL_LLSQ_EXTRAPOLATOR_H
 
-#include "AssemblerLib/GlobalMatrixProviders.h"
+#include "NumLib/DOF/GlobalMatrixProviders.h"
 #include "Extrapolator.h"
 
 namespace NumLib
@@ -100,7 +100,7 @@ private:
     GlobalVector  _residuals;    //!< extrapolation residuals
 
     //! DOF table used for writing to global vectors.
-    AssemblerLib::LocalToGlobalIndexMap const& _local_to_global;
+    NumLib::LocalToGlobalIndexMap const& _local_to_global;
 
     //! Avoids frequent reallocations.
     Eigen::MatrixXd _local_matrix_cache;
diff --git a/NumLib/ODESolver/EquationSystem.h b/NumLib/ODESolver/EquationSystem.h
index e92d7dbdea5d516a7d9720bfc1abfe4b42128f53..1856bd17abc08fe704932160b70f71ee001ffbbb 100644
--- a/NumLib/ODESolver/EquationSystem.h
+++ b/NumLib/ODESolver/EquationSystem.h
@@ -10,7 +10,7 @@
 #ifndef NUMLIB_EQUATIONSYSTEM_H
 #define NUMLIB_EQUATIONSYSTEM_H
 
-#include "AssemblerLib/MatrixProviderUser.h"
+#include "NumLib/DOF/MatrixProviderUser.h"
 
 namespace NumLib
 {
diff --git a/NumLib/ODESolver/NonlinearSolver-impl.h b/NumLib/ODESolver/NonlinearSolver-impl.h
index b9bcdfc0feb21d418db56701d2e4cc72ac028321..1d356d06e89c9dd8f43bc1abfed7f8da3e83b8ba 100644
--- a/NumLib/ODESolver/NonlinearSolver-impl.h
+++ b/NumLib/ODESolver/NonlinearSolver-impl.h
@@ -14,11 +14,11 @@
 
 #include "BaseLib/ConfigTree.h"
 #include "MathLib/LinAlg/BLAS.h"
-#include "AssemblerLib/GlobalMatrixProviders.h"
+#include "MathLib/LinAlg/VectorNorms.h"
+#include "NumLib/DOF/GlobalMatrixProviders.h"
 
 #include "NonlinearSolver.h"
 
-#include "MathLib/LinAlg/VectorNorms.h"
 
 namespace NumLib
 {
diff --git a/NumLib/ODESolver/ODESystem.h b/NumLib/ODESolver/ODESystem.h
index ee3abeb66296805c151d3fbd8ea5319e72a371e7..546389a2368ae25cea33217c73f943377cac7fbe 100644
--- a/NumLib/ODESolver/ODESystem.h
+++ b/NumLib/ODESolver/ODESystem.h
@@ -10,7 +10,7 @@
 #ifndef NUMLIB_ODESYSTEM_H
 #define NUMLIB_ODESYSTEM_H
 
-#include "AssemblerLib/MatrixVectorTraits.h"
+#include "NumLib/DOF/MatrixVectorTraits.h"
 
 #include "Types.h"
 #include "EquationSystem.h"
diff --git a/NumLib/ODESolver/TimeDiscretization.h b/NumLib/ODESolver/TimeDiscretization.h
index c7da0353e04c90b6d8edea8ce8c53bedaf1cfe78..e60165cb2604c7312e6ec498756ca015c866961e 100644
--- a/NumLib/ODESolver/TimeDiscretization.h
+++ b/NumLib/ODESolver/TimeDiscretization.h
@@ -13,7 +13,7 @@
 #include <vector>
 
 #include "MathLib/LinAlg/BLAS.h"
-#include "AssemblerLib/GlobalMatrixProviders.h"
+#include "NumLib/DOF/GlobalMatrixProviders.h"
 #include "Types.h"
 
 
diff --git a/NumLib/ODESolver/TimeDiscretizedODESystem.h b/NumLib/ODESolver/TimeDiscretizedODESystem.h
index 384b59d2e4534578d77ad1e87d33165abbbf113d..69db5e47688b86f9525db4f115924ce875d9ff23 100644
--- a/NumLib/ODESolver/TimeDiscretizedODESystem.h
+++ b/NumLib/ODESolver/TimeDiscretizedODESystem.h
@@ -12,7 +12,7 @@
 
 #include <memory>
 
-#include "AssemblerLib/UnifiedMatrixSetters.h"
+#include "NumLib/DOF/UnifiedMatrixSetters.h"
 #include "MathLib/LinAlg/ApplyKnownSolution.h"
 #include "ProcessLib/DirichletBc.h"
 
diff --git a/NumLib/ODESolver/TimeLoopSingleODE.h b/NumLib/ODESolver/TimeLoopSingleODE.h
index 21344c5e68c27e8fdbaa192be441f85de6f6fb34..f7c45ad92ce0e61b2263cebbb11efca5ad1690fb 100644
--- a/NumLib/ODESolver/TimeLoopSingleODE.h
+++ b/NumLib/ODESolver/TimeLoopSingleODE.h
@@ -10,7 +10,7 @@
 #ifndef NUMLIB_TIMELOOP_H
 #define NUMLIB_TIMELOOP_H
 
-#include "AssemblerLib/GlobalMatrixProviders.h"
+#include "NumLib/DOF/GlobalMatrixProviders.h"
 #include "TimeDiscretizedODESystem.h"
 #include "NonlinearSolver.h"
 
diff --git a/ProcessLib/CMakeLists.txt b/ProcessLib/CMakeLists.txt
index f60e5716d6a120f0229e5f1dd282e6eb513f68cf..280aac21594086dac29af104c64def5214621b96 100644
--- a/ProcessLib/CMakeLists.txt
+++ b/ProcessLib/CMakeLists.txt
@@ -12,7 +12,6 @@ APPEND_SOURCE_FILES(SOURCES TES)
 add_library(ProcessLib ${SOURCES})
 
 target_link_libraries(ProcessLib
-    AssemblerLib
     MaterialsLib
     MeshGeoToolsLib
     NumLib # for shape matrices
diff --git a/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h b/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
index 331bad5a960439e03dc550aaaad0f2ccf4e9caf6..96fb30f1a55c646945a3a03ff6c15bc9d6ece4c0 100644
--- a/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
+++ b/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
@@ -105,7 +105,7 @@ public:
         }
     }
 
-    void addToGlobal(AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices const& indices,
+    void addToGlobal(NumLib::LocalToGlobalIndexMap::RowColumnIndices const& indices,
         GlobalMatrix& /*M*/, GlobalMatrix& K, GlobalVector& b)
         const override
     {
diff --git a/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h b/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h
index 603c4978c4db0de2bf8d155ec1964a643a1f041c..8f6bc793428ae11088e5b588b3158156ffdfb54e 100644
--- a/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h
+++ b/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h
@@ -12,7 +12,7 @@
 
 #include <cassert>
 
-#include "AssemblerLib/VectorMatrixAssembler.h"
+#include "NumLib/Assembler/VectorMatrixAssembler.h"
 #include "NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h"
 #include "ProcessLib/Process.h"
 #include "ProcessLib/Utils/CreateLocalAssemblers.h"
@@ -78,7 +78,7 @@ private:
     using LocalAssemblerInterface =
         GroundwaterFlowLocalAssemblerInterface<GlobalMatrix, GlobalVector>;
 
-    using GlobalAssembler = AssemblerLib::VectorMatrixAssembler<
+    using GlobalAssembler = NumLib::VectorMatrixAssembler<
             GlobalMatrix, GlobalVector, LocalAssemblerInterface,
             NumLib::ODESystemTag::FirstOrderImplicitQuasilinear>;
 
@@ -88,7 +88,7 @@ private:
         GlobalVector, IntegrationPointValue, LocalAssemblerInterface>;
 
     void initializeConcreteProcess(
-            AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+            NumLib::LocalToGlobalIndexMap const& dof_table,
             MeshLib::Mesh const& mesh,
             unsigned const integration_order) override
     {
diff --git a/ProcessLib/LocalAssemblerInterface.h b/ProcessLib/LocalAssemblerInterface.h
index 063002f80b34d0b29b17752e4aaa46f1a4957794..788ba9fbab39c4cb48b2e632a659e7c0e448b17a 100644
--- a/ProcessLib/LocalAssemblerInterface.h
+++ b/ProcessLib/LocalAssemblerInterface.h
@@ -27,7 +27,7 @@ public:
     virtual void assemble(double const t, std::vector<double> const& local_x) = 0;
 
     virtual void addToGlobal(
-        AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices const&,
+        NumLib::LocalToGlobalIndexMap::RowColumnIndices const&,
         GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b) const = 0;
 
     virtual void assembleJacobian(double const /*t*/,
@@ -39,7 +39,7 @@ public:
         std::abort();
     }
 
-    virtual void addJacobianToGlobal(AssemblerLib::LocalToGlobalIndexMap::
+    virtual void addJacobianToGlobal(NumLib::LocalToGlobalIndexMap::
                                          RowColumnIndices const& /*indices*/,
                                      GlobalMatrix& /*Jac*/) const
     {
diff --git a/ProcessLib/NeumannBc.h b/ProcessLib/NeumannBc.h
index 4d841002a2647b69b92dcefddcbe474c9aaa88c3..f52ad1fd403a05b32e37b4f7dea8bb0f6198bf12 100644
--- a/ProcessLib/NeumannBc.h
+++ b/ProcessLib/NeumannBc.h
@@ -16,7 +16,7 @@
 #include "NumLib/Fem/FiniteElement/TemplateIsoparametric.h"
 #include "NumLib/Fem/ShapeMatrixPolicy.h"
 
-#include "AssemblerLib/VectorMatrixAssembler.h"
+#include "NumLib/Assembler/VectorMatrixAssembler.h"
 #include "MeshLib/MeshSubset.h"
 #include "MeshLib/MeshSearch/NodeSearch.h"
 
@@ -54,7 +54,7 @@ public:
     NeumannBc(
         NeumannBcConfig const& bc,
         unsigned const integration_order,
-        AssemblerLib::LocalToGlobalIndexMap const& local_to_global_index_map,
+        NumLib::LocalToGlobalIndexMap const& local_to_global_index_map,
         std::size_t const variable_id,
         std::size_t const component_id)
         : _function(*bc.getFunction()),
@@ -136,7 +136,7 @@ private:
 
     /// Local dof table, a subset of the global one restricted to the
     /// participating #_elements of the boundary condition.
-    std::unique_ptr<AssemblerLib::LocalToGlobalIndexMap> _local_to_global_index_map;
+    std::unique_ptr<NumLib::LocalToGlobalIndexMap> _local_to_global_index_map;
 
     /// Integration order for integration over the lower-dimensional elements of
     /// the #_function.
@@ -146,7 +146,7 @@ private:
         GlobalMatrix, GlobalVector>;
 
     using GlobalAssembler =
-        AssemblerLib::VectorMatrixAssembler<
+        NumLib::VectorMatrixAssembler<
             GlobalMatrix, GlobalVector, LocalAssembler,
             NumLib::ODESystemTag::NeumannBC>;
 
diff --git a/ProcessLib/NeumannBcAssembler.h b/ProcessLib/NeumannBcAssembler.h
index fe124ef9996aa3659552fd23af1c3e0eddf672c1..3d77a029e6bf7eb3ab9a4bdde6cd7da8173264af 100644
--- a/ProcessLib/NeumannBcAssembler.h
+++ b/ProcessLib/NeumannBcAssembler.h
@@ -27,7 +27,7 @@ public:
 
     virtual void assemble(const double t) = 0;
 
-    virtual void addToGlobal(AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices const&,
+    virtual void addToGlobal(NumLib::LocalToGlobalIndexMap::RowColumnIndices const&,
                              GlobalVector& b) const = 0;
 };
 
@@ -95,7 +95,7 @@ public:
         }
     }
 
-    void addToGlobal(AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices const& indices,
+    void addToGlobal(NumLib::LocalToGlobalIndexMap::RowColumnIndices const& indices,
                      GlobalVector& b) const override
     {
         b.add(indices.rows, *_localRhs);
diff --git a/ProcessLib/NumericsConfig.h b/ProcessLib/NumericsConfig.h
index dcdba10ec651473f5082158303a86f85461b699a..7bd06d6b53d1a642d11653afa32efcaa77a7ee41 100644
--- a/ProcessLib/NumericsConfig.h
+++ b/ProcessLib/NumericsConfig.h
@@ -82,11 +82,11 @@ namespace detail
 // Global vector/matrix builder.
 //
 
-#include "AssemblerLib/VectorMatrixBuilder.h"
+#include "NumLib/Assembler/VectorMatrixBuilder.h"
 namespace detail
 {
 using GlobalVectorMatrixBuilderType =
-        AssemblerLib::VectorMatrixBuilder<
+        NumLib::VectorMatrixBuilder<
             GlobalMatrixType,
             GlobalVectorType>;
 }
@@ -94,10 +94,10 @@ using GlobalVectorMatrixBuilderType =
 //
 // Global executor
 //
-#include "AssemblerLib/SerialExecutor.h"
+#include "NumLib/Assembler/SerialExecutor.h"
 namespace detail
 {
-using GlobalExecutorType = AssemblerLib::SerialExecutor;
+using GlobalExecutorType = NumLib::SerialExecutor;
 }
 
 ///
diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
index 8e142fe0f85c1f8b7333adaad487b715f05ee918..48bc1cd62abf06431bb7468202ecf3ec08b8e3ec 100644
--- a/ProcessLib/Process.h
+++ b/ProcessLib/Process.h
@@ -10,8 +10,8 @@
 #ifndef PROCESS_LIB_PROCESS_H_
 #define PROCESS_LIB_PROCESS_H_
 
-#include "AssemblerLib/ComputeSparsityPattern.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
+#include "NumLib/DOF/ComputeSparsityPattern.h"
 #include "NumLib/ODESolver/ODESystem.h"
 #include "NumLib/ODESolver/TimeDiscretization.h"
 #include "NumLib/ODESolver/NonlinearSolver.h"
@@ -190,7 +190,7 @@ public:
 private:
     /// Process specific initialization called by initialize().
     virtual void initializeConcreteProcess(
-        AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+        NumLib::LocalToGlobalIndexMap const& dof_table,
         MeshLib::Mesh const& mesh,
         unsigned const integration_order) = 0;
 
@@ -232,9 +232,9 @@ private:
         }
 
         _local_to_global_index_map.reset(
-            new AssemblerLib::LocalToGlobalIndexMap(
+            new NumLib::LocalToGlobalIndexMap(
                 std::move(all_mesh_subsets),
-                AssemblerLib::ComponentOrder::BY_LOCATION));
+                NumLib::ComponentOrder::BY_LOCATION));
     }
 
     /// Sets the initial condition values in the solution vector x for a given
@@ -304,7 +304,7 @@ private:
     /// DOF-table.
     void computeSparsityPattern()
     {
-        _sparsity_pattern = std::move(AssemblerLib::computeSparsityPattern(
+        _sparsity_pattern = std::move(NumLib::computeSparsityPattern(
             *_local_to_global_index_map, _mesh));
     }
 
@@ -312,7 +312,7 @@ protected:
     MeshLib::Mesh& _mesh;
     std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_all_nodes;
 
-    std::unique_ptr<AssemblerLib::LocalToGlobalIndexMap>
+    std::unique_ptr<NumLib::LocalToGlobalIndexMap>
         _local_to_global_index_map;
 
     SecondaryVariableCollection<GlobalVector> _secondary_variables;
diff --git a/ProcessLib/ProcessOutput.h b/ProcessLib/ProcessOutput.h
index 0e29556db7bf46f6eb3a1c63e900f9ea6032e07c..9757f6b6e487675d5c687594a6e8e5a837a3f5b1 100644
--- a/ProcessLib/ProcessOutput.h
+++ b/ProcessLib/ProcessOutput.h
@@ -88,7 +88,7 @@ void doProcessOutput(
         std::string const& file_name,
         GlobalVector const& x,
         MeshLib::Mesh& mesh,
-        AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+        NumLib::LocalToGlobalIndexMap const& dof_table,
         std::vector<std::reference_wrapper<ProcessVariable>> const&
         process_variables,
         SecondaryVariableCollection<GlobalVector> secondary_variables,
diff --git a/ProcessLib/ProcessVariable.h b/ProcessLib/ProcessVariable.h
index 6386cbb396e45e6c4042136803e6d08255b8ca65..3a735d3c780aac8b7e3cdb0332d456381b5779e2 100644
--- a/ProcessLib/ProcessVariable.h
+++ b/ProcessLib/ProcessVariable.h
@@ -64,7 +64,7 @@ public:
     void initializeDirichletBCs(
         OutputIterator output_bcs,
         MeshGeoToolsLib::MeshNodeSearcher& searcher,
-        const AssemblerLib::LocalToGlobalIndexMap& dof_table,
+        const NumLib::LocalToGlobalIndexMap& dof_table,
         const unsigned nodal_dof_idx)
     {
         for (auto& bc_config : _dirichlet_bc_configs)
diff --git a/ProcessLib/SecondaryVariable.h b/ProcessLib/SecondaryVariable.h
index 7e3cbfbef0c2ccd7968baa9a8444e5a44f53a172..9983763728d8d3d491c3278e45b60925f9d36279 100644
--- a/ProcessLib/SecondaryVariable.h
+++ b/ProcessLib/SecondaryVariable.h
@@ -14,7 +14,7 @@
 #include "BaseLib/uniqueInsert.h"
 #include "NumLib/Extrapolation/Extrapolator.h"
 
-namespace AssemblerLib { class LocalToGlobalIndexMap; }
+namespace NumLib { class LocalToGlobalIndexMap; }
 
 namespace ProcessLib
 {
@@ -36,7 +36,7 @@ struct SecondaryVariableFunctions final
      */
     using Function = std::function<GlobalVector const&(
         GlobalVector const& x,
-        AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+        NumLib::LocalToGlobalIndexMap const& dof_table,
         std::unique_ptr<GlobalVector>& result_cache)>;
 
     SecondaryVariableFunctions() = default;
@@ -49,7 +49,7 @@ struct SecondaryVariableFunctions final
         // Used to detect nasty implicit conversions.
         static_assert(std::is_same<GlobalVector const&,
             typename std::result_of<F1(
-                GlobalVector const&, AssemblerLib::LocalToGlobalIndexMap const&,
+                GlobalVector const&, NumLib::LocalToGlobalIndexMap const&,
                 std::unique_ptr<GlobalVector>&
                 )>::type>::value,
             "The function eval_field_ does not return a const reference"
@@ -57,7 +57,7 @@ struct SecondaryVariableFunctions final
 
         static_assert(std::is_same<GlobalVector const&,
             typename std::result_of<F2(
-                GlobalVector const&, AssemblerLib::LocalToGlobalIndexMap const&,
+                GlobalVector const&, NumLib::LocalToGlobalIndexMap const&,
                 std::unique_ptr<GlobalVector>&
             )>::type>::value,
             "The function eval_residuals_ does not return a const reference"
@@ -71,7 +71,7 @@ struct SecondaryVariableFunctions final
         // Used to detect nasty implicit conversions.
         static_assert(std::is_same<GlobalVector const&,
             typename std::result_of<F1(
-                GlobalVector const&, AssemblerLib::LocalToGlobalIndexMap const&,
+                GlobalVector const&, NumLib::LocalToGlobalIndexMap const&,
                 std::unique_ptr<GlobalVector>&
                 )>::type>::value,
             "The function eval_field_ does not return a const reference"
@@ -212,7 +212,7 @@ makeExtrapolator(PropertyEnum const property,
 
     auto const eval_field = [property, &extrapolator, &local_assemblers](
             GlobalVector const& /*x*/,
-            AssemblerLib::LocalToGlobalIndexMap const& /*dof_table*/,
+            NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
             std::unique_ptr<GlobalVector>& /*result_cache*/
             ) -> GlobalVector const&
     {
@@ -222,7 +222,7 @@ makeExtrapolator(PropertyEnum const property,
 
     auto const eval_residuals = [property, &extrapolator, &local_assemblers](
             GlobalVector const& /*x*/,
-            AssemblerLib::LocalToGlobalIndexMap const& /*dof_table*/,
+            NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
             std::unique_ptr<GlobalVector>& /*result_cache*/
             ) -> GlobalVector const&
     {
diff --git a/ProcessLib/UniformDirichletBoundaryCondition.h b/ProcessLib/UniformDirichletBoundaryCondition.h
index 80627072198f72420be56692d51040a9860a8468..66438803a4526747838e812d2d4cca459aa26b8f 100644
--- a/ProcessLib/UniformDirichletBoundaryCondition.h
+++ b/ProcessLib/UniformDirichletBoundaryCondition.h
@@ -18,8 +18,8 @@
 #include "NumericsConfig.h" // for GlobalIndexType
 
 #include "BaseLib/ConfigTree.h"
-#include "AssemblerLib/LocalToGlobalIndexMap.h"
 #include "MeshGeoToolsLib/MeshNodeSearcher.h"
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
 
 #include "DirichletBc.h"
 
@@ -56,7 +56,7 @@ public:
     /// object.
     void initialize(
             MeshGeoToolsLib::MeshNodeSearcher& searcher,
-            AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+            NumLib::LocalToGlobalIndexMap const& dof_table,
             std::size_t component_id,
             DirichletBc<GlobalIndexType>& bc)
     {
diff --git a/ProcessLib/Utils/CreateLocalAssemblers.h b/ProcessLib/Utils/CreateLocalAssemblers.h
index 8d61deb840129a95b5fd2fcc33ff15282f8a1e3f..081478cee2aea9f9a4d802e4585822c79501a66a 100644
--- a/ProcessLib/Utils/CreateLocalAssemblers.h
+++ b/ProcessLib/Utils/CreateLocalAssemblers.h
@@ -12,7 +12,7 @@
 #include <vector>
 #include <logog/include/logog.hpp>
 
-#include "AssemblerLib/LocalToGlobalIndexMap.h"
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
 
 #include "LocalDataInitializer.h"
 
@@ -29,7 +29,7 @@ template<unsigned GlobalDim, typename GlobalSetup,
          typename LocalAssemblerInterface,
          typename... ExtraCtorArgs>
 void createLocalAssemblers(
-        AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+        NumLib::LocalToGlobalIndexMap const& dof_table,
         std::vector<MeshLib::Element*> const& mesh_elements,
         unsigned const integration_order,
         std::vector<std::unique_ptr<LocalAssemblerInterface>>& local_assemblers,
@@ -83,7 +83,7 @@ template<typename GlobalSetup,
 void createLocalAssemblers(
         const unsigned dimension,
         std::vector<MeshLib::Element*> const& mesh_elements,
-        AssemblerLib::LocalToGlobalIndexMap const& dof_table,
+        NumLib::LocalToGlobalIndexMap const& dof_table,
         unsigned const integration_order,
         std::vector<std::unique_ptr<LocalAssemblerInterface>>& local_assemblers,
         ExtraCtorArgs&&... extra_ctor_args
diff --git a/ProcessLib/Utils/LocalDataInitializer.h b/ProcessLib/Utils/LocalDataInitializer.h
index 9b464d1bd687d1b8008d2d83d6b7bf45361d815c..0172936204529718d967f800412afd9bc5e638c0 100644
--- a/ProcessLib/Utils/LocalDataInitializer.h
+++ b/ProcessLib/Utils/LocalDataInitializer.h
@@ -16,8 +16,8 @@
 #include <typeinfo>
 #include <unordered_map>
 
-#include "AssemblerLib/LocalToGlobalIndexMap.h"
 #include "MeshLib/Elements/Elements.h"
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
 #include "NumLib/Fem/Integration/GaussIntegrationPolicy.h"
 
 
@@ -133,7 +133,7 @@ public:
     using LADataIntfPtr = std::unique_ptr<LocalAssemblerInterface>;
 
     explicit LocalDataInitializer(
-        AssemblerLib::LocalToGlobalIndexMap const& dof_table)
+        NumLib::LocalToGlobalIndexMap const& dof_table)
         : _dof_table(dof_table)
     {
         // /// Lines and points ///////////////////////////////////
@@ -308,7 +308,7 @@ private:
     /// Mapping of element types to local assembler constructors.
     std::unordered_map<std::type_index, LADataBuilder> _builder;
 
-    AssemblerLib::LocalToGlobalIndexMap const& _dof_table;
+    NumLib::LocalToGlobalIndexMap const& _dof_table;
 };
 
 }   // namespace ProcessLib
diff --git a/Tests/CMakeLists.txt b/Tests/CMakeLists.txt
index 6c3ca22f10f333e19c4efc3a02ba20393510693a..cbd1c253f08f8ecd1ef0760503c4ab7b2d8cbd41 100644
--- a/Tests/CMakeLists.txt
+++ b/Tests/CMakeLists.txt
@@ -8,7 +8,6 @@ endif()
 
 APPEND_SOURCE_FILES(TEST_SOURCES)
 APPEND_SOURCE_FILES(TEST_SOURCES ApplicationsLib)
-APPEND_SOURCE_FILES(TEST_SOURCES AssemblerLib)
 APPEND_SOURCE_FILES(TEST_SOURCES BaseLib)
 APPEND_SOURCE_FILES(TEST_SOURCES FileIO)
 APPEND_SOURCE_FILES(TEST_SOURCES GeoLib)
@@ -24,7 +23,7 @@ if(QT4_FOUND)
 endif()
 
 if(OGS_USE_PETSC OR OGS_USE_MPI)
-    list(REMOVE_ITEM TEST_SOURCES AssemblerLib/TestSerialLinearSolver.cpp)
+    list(REMOVE_ITEM TEST_SOURCES NumLib/TestSerialLinearSolver.cpp)
 endif()
 
 add_executable(testrunner ${TEST_SOURCES})
@@ -33,7 +32,6 @@ set_target_properties(testrunner PROPERTIES FOLDER Testing)
 target_link_libraries(testrunner
     ApplicationsLib
     ApplicationsFileIO
-    AssemblerLib
     FileIO
     GTest
     InSituLib
diff --git a/Tests/AssemblerLib/LocalToGlobalIndexMap.cpp b/Tests/NumLib/LocalToGlobalIndexMap.cpp
similarity index 73%
rename from Tests/AssemblerLib/LocalToGlobalIndexMap.cpp
rename to Tests/NumLib/LocalToGlobalIndexMap.cpp
index d9f6e7e2b77913e9891b2aea0a274b7979afad92..2105ab2404888dcf2436e25eab7f45b4fb5a702c 100644
--- a/Tests/AssemblerLib/LocalToGlobalIndexMap.cpp
+++ b/Tests/NumLib/LocalToGlobalIndexMap.cpp
@@ -12,7 +12,7 @@
 
 #include <gtest/gtest.h>
 
-#include "AssemblerLib/LocalToGlobalIndexMap.h"
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
 #include "MeshLib/MeshGenerators/MeshGenerator.h"
 
 #include "MeshLib/MeshSearch/NodeSearch.h"
@@ -20,11 +20,11 @@
 #include "MeshLib/Mesh.h"
 
 
-class AssemblerLibLocalToGlobalIndexMapTest : public ::testing::Test
+class NumLibLocalToGlobalIndexMapTest : public ::testing::Test
 {
 
 public:
-    AssemblerLibLocalToGlobalIndexMapTest()
+    NumLibLocalToGlobalIndexMapTest()
     {
         mesh.reset(MeshLib::MeshGenerator::generateLineMesh(1.0, mesh_size));
         nodesSubset.reset(new MeshLib::MeshSubset(*mesh, &mesh->getNodes()));
@@ -44,21 +44,21 @@ protected:
     static std::size_t const comp1_id = 1;
     std::vector<std::unique_ptr<MeshLib::MeshSubsets>> components;
 
-    std::unique_ptr<AssemblerLib::LocalToGlobalIndexMap const> dof_map;
+    std::unique_ptr<NumLib::LocalToGlobalIndexMap const> dof_map;
 };
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibLocalToGlobalIndexMapTest, NumberOfRowsByComponent)
+TEST_F(NumLibLocalToGlobalIndexMapTest, NumberOfRowsByComponent)
 #else
-TEST_F(AssemblerLibLocalToGlobalIndexMapTest, DISABLED_NumberOfRowsByComponent)
+TEST_F(NumLibLocalToGlobalIndexMapTest, DISABLED_NumberOfRowsByComponent)
 #endif
 {
     // need to store the size because the components will be moved into the
     // DOF-table.
     std::size_t components_size = components.size();
 
-    dof_map.reset(new AssemblerLib::LocalToGlobalIndexMap(std::move(components),
-        AssemblerLib::ComponentOrder::BY_COMPONENT));
+    dof_map.reset(new NumLib::LocalToGlobalIndexMap(std::move(components),
+        NumLib::ComponentOrder::BY_COMPONENT));
 
     // There must be as many rows as nodes in the input times the number of
     // components.
@@ -66,17 +66,17 @@ TEST_F(AssemblerLibLocalToGlobalIndexMapTest, DISABLED_NumberOfRowsByComponent)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibLocalToGlobalIndexMapTest, NumberOfRowsByLocation)
+TEST_F(NumLibLocalToGlobalIndexMapTest, NumberOfRowsByLocation)
 #else
-TEST_F(AssemblerLibLocalToGlobalIndexMapTest, DISABLED_NumberOfRowsByLocation)
+TEST_F(NumLibLocalToGlobalIndexMapTest, DISABLED_NumberOfRowsByLocation)
 #endif
 {
     // need to store the size because the components will be moved into the
     // DOF-table.
     std::size_t components_size = components.size();
 
-    dof_map.reset(new AssemblerLib::LocalToGlobalIndexMap(std::move(components),
-        AssemblerLib::ComponentOrder::BY_LOCATION));
+    dof_map.reset(new NumLib::LocalToGlobalIndexMap(std::move(components),
+        NumLib::ComponentOrder::BY_LOCATION));
 
     // There must be as many rows as nodes in the input times the number of
     // components.
@@ -84,13 +84,13 @@ TEST_F(AssemblerLibLocalToGlobalIndexMapTest, DISABLED_NumberOfRowsByLocation)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibLocalToGlobalIndexMapTest, SubsetByComponent)
+TEST_F(NumLibLocalToGlobalIndexMapTest, SubsetByComponent)
 #else
-TEST_F(AssemblerLibLocalToGlobalIndexMapTest, DISABLED_SubsetByComponent)
+TEST_F(NumLibLocalToGlobalIndexMapTest, DISABLED_SubsetByComponent)
 #endif
 {
-    dof_map.reset(new AssemblerLib::LocalToGlobalIndexMap(std::move(components),
-        AssemblerLib::ComponentOrder::BY_COMPONENT));
+    dof_map.reset(new NumLib::LocalToGlobalIndexMap(std::move(components),
+        NumLib::ComponentOrder::BY_COMPONENT));
 
     // Select some elements from the full mesh.
     std::array<std::size_t, 3> const ids = {{ 0, 5, 8 }};
@@ -106,7 +106,7 @@ TEST_F(AssemblerLibLocalToGlobalIndexMapTest, DISABLED_SubsetByComponent)
     auto selected_component = std::unique_ptr<MeshLib::MeshSubsets>{
         new MeshLib::MeshSubsets{selected_subset.get()}};
 
-    auto dof_map_subset = std::unique_ptr<AssemblerLib::LocalToGlobalIndexMap>{
+    auto dof_map_subset = std::unique_ptr<NumLib::LocalToGlobalIndexMap>{
         dof_map->deriveBoundaryConstrainedMap(0,  // variable id
                                               1,  // component id
                                               std::move(selected_component),
diff --git a/Tests/AssemblerLib/LocalToGlobalIndexMapMultiComponent.cpp b/Tests/NumLib/LocalToGlobalIndexMapMultiComponent.cpp
similarity index 91%
rename from Tests/AssemblerLib/LocalToGlobalIndexMapMultiComponent.cpp
rename to Tests/NumLib/LocalToGlobalIndexMapMultiComponent.cpp
index 0184d26fc2efa95fe60bc2818a190b91a1561be0..7be34c60b3b4430e63166f4ff1c1b62141bd8ca4 100644
--- a/Tests/AssemblerLib/LocalToGlobalIndexMapMultiComponent.cpp
+++ b/Tests/NumLib/LocalToGlobalIndexMapMultiComponent.cpp
@@ -10,7 +10,7 @@
 
 #include <memory>
 
-#include "AssemblerLib/VectorMatrixAssembler.h"
+#include "NumLib/Assembler/VectorMatrixAssembler.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/Location.h"
 #include "MeshLib/Mesh.h"
@@ -23,15 +23,15 @@
 #include "MeshGeoToolsLib/MeshNodeSearcher.h"
 #include "MeshGeoToolsLib/BoundaryElementsSearcher.h"
 
-namespace AL = AssemblerLib;
+namespace AL = NumLib;
 namespace MeL = MeshLib;
 namespace MGTL = MeshGeoToolsLib;
 
-class AssemblerLibLocalToGlobalIndexMapMultiDOFTest : public ::testing::Test
+class NumLibLocalToGlobalIndexMapMultiDOFTest : public ::testing::Test
 {
 public:
     static const std::size_t mesh_subdivs = 4;
-    AssemblerLibLocalToGlobalIndexMapMultiDOFTest()
+    NumLibLocalToGlobalIndexMapMultiDOFTest()
     {
         mesh.reset(MeL::MeshGenerator::generateRegularQuadMesh(2.0, mesh_subdivs));
         mesh_items_all_nodes.reset(new MeL::MeshSubset(*mesh, &mesh->getNodes()));
@@ -71,7 +71,7 @@ public:
             mesh_items_all_nodes->getIntersectionByNodes(nodes));
     }
 
-    ~AssemblerLibLocalToGlobalIndexMapMultiDOFTest()
+    ~NumLibLocalToGlobalIndexMapMultiDOFTest()
     {
         for (auto e : boundary_elements)
             delete e;
@@ -143,7 +143,7 @@ struct ComputeGlobalIndexByLocation
 
 
 template <AL::ComponentOrder ComponentOrder>
-void AssemblerLibLocalToGlobalIndexMapMultiDOFTest::test(
+void NumLibLocalToGlobalIndexMapMultiDOFTest::test(
     const unsigned num_components,
     const unsigned selected_component,
     std::function<std::size_t(std::size_t, std::size_t)> const
@@ -218,9 +218,9 @@ void assert_equal(AL::LocalToGlobalIndexMap const& dof1, AL::LocalToGlobalIndexM
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, Test1Comp)
+TEST_F(NumLibLocalToGlobalIndexMapMultiDOFTest, Test1Comp)
 #else
-TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_Test1Comp)
+TEST_F(NumLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_Test1Comp)
 #endif
 {
     unsigned const num_components = 1;
@@ -240,9 +240,9 @@ TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_Test1Comp)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, TestMultiCompByComponent)
+TEST_F(NumLibLocalToGlobalIndexMapMultiDOFTest, TestMultiCompByComponent)
 #else
-TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_TestMultiCompByComponent)
+TEST_F(NumLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_TestMultiCompByComponent)
 #endif
 {
     unsigned const num_components = 5;
@@ -253,9 +253,9 @@ TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_TestMultiCompByCo
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, TestMultiCompByLocation)
+TEST_F(NumLibLocalToGlobalIndexMapMultiDOFTest, TestMultiCompByLocation)
 #else
-TEST_F(AssemblerLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_TestMultiCompByLocation)
+TEST_F(NumLibLocalToGlobalIndexMapMultiDOFTest, DISABLED_TestMultiCompByLocation)
 #endif
 {
     unsigned const num_components = 5;
diff --git a/Tests/NumLib/ODEs.h b/Tests/NumLib/ODEs.h
index 352cb5f342f12801eb84ba1d3f5a46efac136fe7..bb3c642e6f985e0e8b6485c5acd23bdfbd272a29 100644
--- a/Tests/NumLib/ODEs.h
+++ b/Tests/NumLib/ODEs.h
@@ -10,7 +10,7 @@
 #ifndef TESTS_NUMLIB_ODES_H
 #define TESTS_NUMLIB_ODES_H
 
-#include "AssemblerLib/UnifiedMatrixSetters.h"
+#include "NumLib/DOF/UnifiedMatrixSetters.h"
 #include "MathLib/LinAlg/BLAS.h"
 #include "NumLib/ODESolver/ODESystem.h"
 
diff --git a/Tests/AssemblerLib/SteadyDiffusion2DExample1.h b/Tests/NumLib/SteadyDiffusion2DExample1.h
similarity index 98%
rename from Tests/AssemblerLib/SteadyDiffusion2DExample1.h
rename to Tests/NumLib/SteadyDiffusion2DExample1.h
index e5a4737c35e369a5b99e588cd33ecf5bbe54eafa..5e72dce9b35b50b7e7ca6ed638feb91127c9ba19 100644
--- a/Tests/AssemblerLib/SteadyDiffusion2DExample1.h
+++ b/Tests/NumLib/SteadyDiffusion2DExample1.h
@@ -49,7 +49,7 @@ template<typename IndexType>struct SteadyDiffusion2DExample1
             // already stored in the _localA matrix.
         }
 
-        void addToGlobal(AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices const& indices,
+        void addToGlobal(NumLib::LocalToGlobalIndexMap::RowColumnIndices const& indices,
                          GlobalMatrix& /*M*/, GlobalMatrix& K, GlobalVector& b) const
         {
             K.add(indices, *_localA);
diff --git a/Tests/NumLib/TestExtrapolation.cpp b/Tests/NumLib/TestExtrapolation.cpp
index cb2709ab07f633ac1360540665e347bb80f27efd..0367ca36cfb198aade563ae1fbb7222715801b08 100644
--- a/Tests/NumLib/TestExtrapolation.cpp
+++ b/Tests/NumLib/TestExtrapolation.cpp
@@ -10,10 +10,10 @@
 #include <random>
 #include <gtest/gtest.h>
 
-#include "AssemblerLib/MatrixProviderUser.h"
-#include "AssemblerLib/MatrixVectorTraits.h"
-#include "AssemblerLib/UnifiedMatrixSetters.h"
-#include "AssemblerLib/VectorMatrixAssembler.h"
+#include "NumLib/DOF/MatrixProviderUser.h"
+#include "NumLib/DOF/MatrixVectorTraits.h"
+#include "NumLib/DOF/UnifiedMatrixSetters.h"
+#include "NumLib/Assembler/VectorMatrixAssembler.h"
 
 #include "MathLib/LinAlg/BLAS.h"
 
@@ -156,7 +156,7 @@ public:
     using GlobalVector = typename GlobalSetup::VectorType;
 
     using LocalAssembler = LocalAssemblerDataInterface<GlobalMatrix, GlobalVector>;
-    using GlobalAssembler = AssemblerLib::VectorMatrixAssembler<
+    using GlobalAssembler = NumLib::VectorMatrixAssembler<
         GlobalMatrix, GlobalVector, LocalAssembler,
         // The exact tag does not matter here.
         NumLib::ODESystemTag::FirstOrderImplicitQuasilinear>;
@@ -175,9 +175,9 @@ public:
         all_mesh_subsets.emplace_back(
                     new MeshLib::MeshSubsets{&_mesh_subset_all_nodes});
 
-        _dof_table.reset(new AssemblerLib::LocalToGlobalIndexMap(
+        _dof_table.reset(new NumLib::LocalToGlobalIndexMap(
               std::move(all_mesh_subsets),
-              AssemblerLib::ComponentOrder::BY_COMPONENT));
+              NumLib::ComponentOrder::BY_COMPONENT));
 
         // Passing _dof_table works, because this process has only one variable
         // and the variable has exactly one component.
@@ -196,7 +196,7 @@ public:
         {
             auto inner_cb = [&loc_asm, property](
                 std::vector<double> const& local_x,
-                AssemblerLib::LocalToGlobalIndexMap::RowColumnIndices const&
+                NumLib::LocalToGlobalIndexMap::RowColumnIndices const&
             ) {
                 loc_asm.interpolateNodalValuesToIntegrationPoints(local_x, property);
             };
@@ -253,7 +253,7 @@ private:
     unsigned const _integration_order;
 
     MeshLib::MeshSubset _mesh_subset_all_nodes;
-    std::unique_ptr<AssemblerLib::LocalToGlobalIndexMap> _dof_table;
+    std::unique_ptr<NumLib::LocalToGlobalIndexMap> _dof_table;
 
     std::unique_ptr<GlobalAssembler> _global_assembler;
     std::vector<std::unique_ptr<LocalAssembler>> _local_assemblers;
diff --git a/Tests/AssemblerLib/TestMeshComponentMap.cpp b/Tests/NumLib/TestMeshComponentMap.cpp
similarity index 84%
rename from Tests/AssemblerLib/TestMeshComponentMap.cpp
rename to Tests/NumLib/TestMeshComponentMap.cpp
index 326593eccea25ec6cff683e5f3b4be7b6ba9b3be..d4ae505f1111e4fe595878d93ee1b6b195847e8e 100644
--- a/Tests/AssemblerLib/TestMeshComponentMap.cpp
+++ b/Tests/NumLib/TestMeshComponentMap.cpp
@@ -14,20 +14,20 @@
 #include <memory>
 #include <vector>
 
-#include "AssemblerLib/MeshComponentMap.h"
+#include "NumLib/DOF/MeshComponentMap.h"
 
 #include "MeshLib/MeshGenerators/MeshGenerator.h"
 #include "MeshLib/MeshSubsets.h"
 
-class AssemblerLibMeshComponentMapTest : public ::testing::Test
+class NumLibMeshComponentMapTest : public ::testing::Test
 {
     public:
     typedef MeshLib::MeshItemType MeshItemType;
     typedef MeshLib::Location Location;
-    typedef AssemblerLib::MeshComponentMap MeshComponentMap;
+    typedef NumLib::MeshComponentMap MeshComponentMap;
 
     public:
-    AssemblerLibMeshComponentMapTest()
+    NumLibMeshComponentMapTest()
         : mesh(nullptr), nodesSubset(nullptr), cmap(nullptr)
     {
         mesh = MeshLib::MeshGenerator::generateLineMesh(1.0, mesh_size);
@@ -38,7 +38,7 @@ class AssemblerLibMeshComponentMapTest : public ::testing::Test
         components.emplace_back(new MeshLib::MeshSubsets{nodesSubset});
     }
 
-    ~AssemblerLibMeshComponentMapTest()
+    ~NumLibMeshComponentMapTest()
     {
         delete cmap;
         delete nodesSubset;
@@ -68,16 +68,16 @@ class AssemblerLibMeshComponentMapTest : public ::testing::Test
 };
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibMeshComponentMapTest, CheckOrderByComponent)
+TEST_F(NumLibMeshComponentMapTest, CheckOrderByComponent)
 #else
-TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_CheckOrderByComponent)
+TEST_F(NumLibMeshComponentMapTest, DISABLED_CheckOrderByComponent)
 #endif
 {
     // - Entries in the vector are arranged in the order of a component type and then node ID
     // - For example, x=[(node 0, comp 0) (node 1, comp 0) ... (node n, comp0), (node 0, comp1) ... ]
 
     cmap = new MeshComponentMap(components,
-        AssemblerLib::ComponentOrder::BY_COMPONENT);
+        NumLib::ComponentOrder::BY_COMPONENT);
 
     ASSERT_EQ(2 * mesh->getNNodes(), cmap->dofSizeWithGhosts());
     for (std::size_t i = 0; i < mesh_size; i++)
@@ -96,16 +96,16 @@ TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_CheckOrderByComponent)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibMeshComponentMapTest, CheckOrderByLocation)
+TEST_F(NumLibMeshComponentMapTest, CheckOrderByLocation)
 #else
-TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_CheckOrderByLocation)
+TEST_F(NumLibMeshComponentMapTest, DISABLED_CheckOrderByLocation)
 #endif
 {
     // - Entries in the vector are arranged in the order of node ID and then a component type
     // - For example, x=[(node 0, comp 0) (node 0, comp 1) ... (node n, comp0), (node n, comp1) ]
 
     cmap = new MeshComponentMap(components,
-        AssemblerLib::ComponentOrder::BY_LOCATION);
+        NumLib::ComponentOrder::BY_LOCATION);
 
     ASSERT_EQ(2 * mesh->getNNodes(), cmap->dofSizeWithGhosts());
     for (std::size_t i = 0; i < mesh_size; i++)
@@ -124,13 +124,13 @@ TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_CheckOrderByLocation)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibMeshComponentMapTest, OutOfRangeAccess)
+TEST_F(NumLibMeshComponentMapTest, OutOfRangeAccess)
 #else
-TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_OutOfRangeAccess)
+TEST_F(NumLibMeshComponentMapTest, DISABLED_OutOfRangeAccess)
 #endif
 {
     cmap = new MeshComponentMap(components,
-        AssemblerLib::ComponentOrder::BY_COMPONENT);
+        NumLib::ComponentOrder::BY_COMPONENT);
 
     ASSERT_EQ(MeshComponentMap::nop, cmap->getGlobalIndex(
         Location(mesh->getID(), MeshItemType::Node, mesh_size + 1), comp0_id));
@@ -143,13 +143,13 @@ TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_OutOfRangeAccess)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibMeshComponentMapTest, SubsetOfNodesByComponent)
+TEST_F(NumLibMeshComponentMapTest, SubsetOfNodesByComponent)
 #else
-TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_SubsetOfNodesByComponent)
+TEST_F(NumLibMeshComponentMapTest, DISABLED_SubsetOfNodesByComponent)
 #endif
 {
     cmap = new MeshComponentMap(components,
-        AssemblerLib::ComponentOrder::BY_COMPONENT);
+        NumLib::ComponentOrder::BY_COMPONENT);
 
     // Select some nodes from the full mesh.
     std::array<std::size_t, 3> const ids = {{ 0, 5, 9 }};
@@ -180,13 +180,13 @@ TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_SubsetOfNodesByComponent)
 }
 
 #ifndef USE_PETSC
-TEST_F(AssemblerLibMeshComponentMapTest, SubsetOfNodesByLocation)
+TEST_F(NumLibMeshComponentMapTest, SubsetOfNodesByLocation)
 #else
-TEST_F(AssemblerLibMeshComponentMapTest, DISABLED_SubsetOfNodesByLocation)
+TEST_F(NumLibMeshComponentMapTest, DISABLED_SubsetOfNodesByLocation)
 #endif
 {
     cmap = new MeshComponentMap(components,
-        AssemblerLib::ComponentOrder::BY_LOCATION);
+        NumLib::ComponentOrder::BY_LOCATION);
 
     // Select some nodes from the full mesh.
     std::array<std::size_t, 3> const ids = {{ 0, 5, 9 }};
diff --git a/Tests/AssemblerLib/TestSerialExecutor.cpp b/Tests/NumLib/TestSerialExecutor.cpp
similarity index 88%
rename from Tests/AssemblerLib/TestSerialExecutor.cpp
rename to Tests/NumLib/TestSerialExecutor.cpp
index 0f47f417ce9cfa52144550ba6ff6d1fff5ac2e8c..0ae72d9ee14658377bda0ab4c95df029a78f4259 100644
--- a/Tests/AssemblerLib/TestSerialExecutor.cpp
+++ b/Tests/NumLib/TestSerialExecutor.cpp
@@ -13,10 +13,10 @@
 #include <vector>
 #include <functional>
 #include <numeric>
-#include "AssemblerLib/SerialExecutor.h"
+#include "NumLib/Assembler/SerialExecutor.h"
 
 template <typename ContainerElement_>
-class AssemblerLibSerialExecutor : public ::testing::Test
+class NumLibSerialExecutor : public ::testing::Test
 {
 public:
     using ContainerElement = ContainerElement_;
@@ -69,9 +69,9 @@ public:
 
 typedef ::testing::Types<int> TestCases;
 
-TYPED_TEST_CASE(AssemblerLibSerialExecutor, TestCases);
+TYPED_TEST_CASE(NumLibSerialExecutor, TestCases);
 
-TYPED_TEST(AssemblerLibSerialExecutor, ContainerArgument)
+TYPED_TEST(NumLibSerialExecutor, ContainerArgument)
 {
     using Elem         = typename TestFixture::ContainerElement;
     using Container    = typename TestFixture::Container;
@@ -84,14 +84,14 @@ TYPED_TEST(AssemblerLibSerialExecutor, ContainerArgument)
                     TestFixture::subtractFromReferenceStatic(value, index, ref_inner);
                 };
 
-            AssemblerLib::SerialExecutor::executeDereferenced(
+            NumLib::SerialExecutor::executeDereferenced(
                 cb_static, ctnr, ref);
         }
     );
 
     TestFixture::test(
         [this](PtrContainer const& ctnr, Container& ref) {
-            AssemblerLib::SerialExecutor::executeMemberDereferenced(
+            NumLib::SerialExecutor::executeMemberDereferenced(
                 *static_cast<TestFixture*>(this), &TestFixture::subtractFromReference,
                 ctnr, ref
             );
diff --git a/Tests/AssemblerLib/TestSerialLinearSolver.cpp b/Tests/NumLib/TestSerialLinearSolver.cpp
similarity index 95%
rename from Tests/AssemblerLib/TestSerialLinearSolver.cpp
rename to Tests/NumLib/TestSerialLinearSolver.cpp
index 402244719e7c896ad2e9e6c7789756d11f1820ef..48fc76250de34753ed5d7fccc5d19f5738185489 100644
--- a/Tests/AssemblerLib/TestSerialLinearSolver.cpp
+++ b/Tests/NumLib/TestSerialLinearSolver.cpp
@@ -16,7 +16,7 @@
 
 #include <gtest/gtest.h>
 
-#include "AssemblerLib/VectorMatrixAssembler.h"
+#include "NumLib/Assembler/VectorMatrixAssembler.h"
 
 #include "MathLib/LinAlg/ApplyKnownSolution.h"
 #include "MathLib/LinAlg/Solvers/GaussAlgorithm.h"
@@ -35,7 +35,7 @@
 #include "../TestTools.h"
 #include "SteadyDiffusion2DExample1.h"
 
-TEST(AssemblerLibSerialLinearSolver, Steady2DdiffusionQuadElem)
+TEST(NumLibSerialLinearSolver, Steady2DdiffusionQuadElem)
 {
     // example
     using Example = SteadyDiffusion2DExample1<GlobalIndexType>;
@@ -58,8 +58,8 @@ TEST(AssemblerLibSerialLinearSolver, Steady2DdiffusionQuadElem)
     // define a mesh item composition in a vector
     std::vector<std::unique_ptr<MeshLib::MeshSubsets>> vec_comp_dis;
     vec_comp_dis.emplace_back(new MeshLib::MeshSubsets{&mesh_items_all_nodes});
-    AssemblerLib::LocalToGlobalIndexMap local_to_global_index_map(
-        std::move(vec_comp_dis), AssemblerLib::ComponentOrder::BY_COMPONENT);
+    NumLib::LocalToGlobalIndexMap local_to_global_index_map(
+        std::move(vec_comp_dis), NumLib::ComponentOrder::BY_COMPONENT);
 
     //--------------------------------------------------------------------------
     // Construct a linear system
@@ -102,7 +102,7 @@ TEST(AssemblerLibSerialLinearSolver, Steady2DdiffusionQuadElem)
 
     // TODO in the future use simpler NumLib::ODESystemTag
     // Local and global assemblers.
-    typedef AssemblerLib::VectorMatrixAssembler<
+    typedef NumLib::VectorMatrixAssembler<
             GlobalMatrix, GlobalVector, LocalAssembler,
             NumLib::ODESystemTag::FirstOrderImplicitQuasilinear> GlobalAssembler;
 
diff --git a/Tests/AssemblerLib/TestVectorMatrixBuilder.cpp b/Tests/NumLib/TestVectorMatrixBuilder.cpp
similarity index 74%
rename from Tests/AssemblerLib/TestVectorMatrixBuilder.cpp
rename to Tests/NumLib/TestVectorMatrixBuilder.cpp
index a47171805fbfc2cc1cf83e3e1c2891717c83c513..dbc1d7811df8ca65e90d58c774836204e245c39f 100644
--- a/Tests/AssemblerLib/TestVectorMatrixBuilder.cpp
+++ b/Tests/NumLib/TestVectorMatrixBuilder.cpp
@@ -10,26 +10,26 @@
 #include <memory>
 #include <gtest/gtest.h>
 
-#include "AssemblerLib/MeshComponentMap.h"
+#include "NumLib/DOF/MeshComponentMap.h"
 
 #include "MeshLib/MeshGenerators/MeshGenerator.h"
 #include "MeshLib/MeshSubsets.h"
 
-#include "AssemblerLib/VectorMatrixBuilder.h"
+#include "NumLib/Assembler/VectorMatrixBuilder.h"
 
 template <typename Builder>
-class AssemblerLibVectorMatrixBuilder : public ::testing::Test
+class NumLibVectorMatrixBuilder : public ::testing::Test
 {
     public:
     typedef MeshLib::MeshItemType MeshItemType;
     typedef MeshLib::Location Location;
-    typedef AssemblerLib::MeshComponentMap MeshComponentMap;
+    typedef NumLib::MeshComponentMap MeshComponentMap;
 
     typedef typename Builder::VectorType VectorType;
     typedef typename Builder::MatrixType MatrixType;
 
     public:
-    AssemblerLibVectorMatrixBuilder()
+    NumLibVectorMatrixBuilder()
         : mesh(nullptr), nodesSubset(nullptr), cmap(nullptr)
     {
         mesh = MeshLib::MeshGenerator::generateLineMesh(1.0, mesh_size);
@@ -40,10 +40,10 @@ class AssemblerLibVectorMatrixBuilder : public ::testing::Test
         components.emplace_back(new MeshLib::MeshSubsets{nodesSubset});
 
         cmap = new MeshComponentMap(components,
-            AssemblerLib::ComponentOrder::BY_COMPONENT);
+            NumLib::ComponentOrder::BY_COMPONENT);
     }
 
-    ~AssemblerLibVectorMatrixBuilder()
+    ~NumLibVectorMatrixBuilder()
     {
         delete cmap;
         delete nodesSubset;
@@ -58,12 +58,12 @@ class AssemblerLibVectorMatrixBuilder : public ::testing::Test
     MeshComponentMap const* cmap;
 };
 
-TYPED_TEST_CASE_P(AssemblerLibVectorMatrixBuilder);
+TYPED_TEST_CASE_P(NumLibVectorMatrixBuilder);
 
 #ifndef USE_PETSC
-TYPED_TEST_P(AssemblerLibVectorMatrixBuilder, createVector)
+TYPED_TEST_P(NumLibVectorMatrixBuilder, createVector)
 #else
-TYPED_TEST_P(AssemblerLibVectorMatrixBuilder, DISABLED_createVector)
+TYPED_TEST_P(NumLibVectorMatrixBuilder, DISABLED_createVector)
 #endif
 {
     typedef typename TestFixture::VectorType V;
@@ -77,9 +77,9 @@ TYPED_TEST_P(AssemblerLibVectorMatrixBuilder, DISABLED_createVector)
 }
 
 #ifndef USE_PETSC
-TYPED_TEST_P(AssemblerLibVectorMatrixBuilder, createMatrix)
+TYPED_TEST_P(NumLibVectorMatrixBuilder, createMatrix)
 #else
-TYPED_TEST_P(AssemblerLibVectorMatrixBuilder, DISABLED_createMatrix)
+TYPED_TEST_P(NumLibVectorMatrixBuilder, DISABLED_createMatrix)
 #endif
 {
     typedef typename TestFixture::MatrixType M;
@@ -94,10 +94,10 @@ TYPED_TEST_P(AssemblerLibVectorMatrixBuilder, DISABLED_createMatrix)
 }
 
 #ifndef USE_PETSC
-REGISTER_TYPED_TEST_CASE_P(AssemblerLibVectorMatrixBuilder,
+REGISTER_TYPED_TEST_CASE_P(NumLibVectorMatrixBuilder,
     createVector, createMatrix);
 #else
-REGISTER_TYPED_TEST_CASE_P(AssemblerLibVectorMatrixBuilder,
+REGISTER_TYPED_TEST_CASE_P(NumLibVectorMatrixBuilder,
     DISABLED_createVector, DISABLED_createMatrix);
 #endif
 
@@ -118,17 +118,17 @@ REGISTER_TYPED_TEST_CASE_P(AssemblerLibVectorMatrixBuilder,
 
 typedef ::testing::Types
     <
-      AssemblerLib::VectorMatrixBuilder<
+      NumLib::VectorMatrixBuilder<
       MathLib::EigenMatrix, MathLib::EigenVector>
 #ifdef USE_LIS
-    , AssemblerLib::VectorMatrixBuilder<
+    , NumLib::VectorMatrixBuilder<
         MathLib::LisMatrix, MathLib::LisVector>
 #endif  // USE_LIS
 #ifdef USE_PETSC
-    , AssemblerLib::VectorMatrixBuilder<
+    , NumLib::VectorMatrixBuilder<
         MathLib::PETScMatrix, MathLib::PETScVector>
 #endif  // USE_PETSC
     > TestTypes;
 
-INSTANTIATE_TYPED_TEST_CASE_P(templated, AssemblerLibVectorMatrixBuilder,
+INSTANTIATE_TYPED_TEST_CASE_P(templated, NumLibVectorMatrixBuilder,
     TestTypes);
diff --git a/Tests/testrunner.cpp b/Tests/testrunner.cpp
index c67c2bfd77bdb087a05f2a91dbef8cdff40671f8..b14fe810ef19a7ea905c8860e3d731d66ecbc8d1 100644
--- a/Tests/testrunner.cpp
+++ b/Tests/testrunner.cpp
@@ -19,7 +19,7 @@
 
 #include "Applications/ApplicationsLib/LogogSetup.h"
 #include "Applications/ApplicationsLib/LinearSolverLibrarySetup.h"
-#include "AssemblerLib/GlobalMatrixProviders.h"
+#include "NumLib/DOF/GlobalMatrixProviders.h"
 #include "BaseLib/BuildInfo.h"
 #include "BaseLib/TemplateLogogFormatterSuppressedGCC.h"