diff --git a/Applications/FileIO/AsciiRasterInterface.cpp b/Applications/FileIO/AsciiRasterInterface.cpp
index c0d6495caa2549e7a9f2123fb96861a6543d32d4..04c95d466e8a22ad1812256428d11d6e9cb9f5b8 100644
--- a/Applications/FileIO/AsciiRasterInterface.cpp
+++ b/Applications/FileIO/AsciiRasterInterface.cpp
@@ -13,6 +13,7 @@
 
 #include "AsciiRasterInterface.h"
 
+#include <fstream>
 #include <tuple>
 
 #include "BaseLib/FileTools.h"
diff --git a/Applications/FileIO/AsciiRasterInterface.h b/Applications/FileIO/AsciiRasterInterface.h
index f3b7d09469e6a71bd9c919daf66ae1ad38636ad2..fdcc7610cbda657cfbc92bd8df1e216a6faf07e2 100644
--- a/Applications/FileIO/AsciiRasterInterface.h
+++ b/Applications/FileIO/AsciiRasterInterface.h
@@ -13,7 +13,6 @@
 
 #pragma once
 
-#include <fstream>
 #include <optional>
 #include <string>
 #include <vector>
diff --git a/Applications/FileIO/CsvInterface.cpp b/Applications/FileIO/CsvInterface.cpp
index aa61387717c68ad10b5bfc605b2a7ab8755e2e56..9c067d01635f486417614aee71070aebc319117f 100644
--- a/Applications/FileIO/CsvInterface.cpp
+++ b/Applications/FileIO/CsvInterface.cpp
@@ -14,7 +14,6 @@
 #include "CsvInterface.h"
 
 #include <algorithm>
-#include <fstream>
 #include <numeric>
 
 #include "GeoLib/Point.h"
diff --git a/Applications/FileIO/FEFLOW/FEFLOWGeoInterface.cpp b/Applications/FileIO/FEFLOW/FEFLOWGeoInterface.cpp
index 328cdd4cb4ea0944e378f8e883df78eb8724693c..6f954e7400092578ed845d5c7ff8ce41cd1b9cc6 100644
--- a/Applications/FileIO/FEFLOW/FEFLOWGeoInterface.cpp
+++ b/Applications/FileIO/FEFLOW/FEFLOWGeoInterface.cpp
@@ -13,6 +13,7 @@
 #include <QtXml/QDomDocument>
 #include <boost/algorithm/string/trim.hpp>
 #include <cctype>
+#include <fstream>
 #include <memory>
 
 #include "BaseLib/FileTools.h"
diff --git a/Applications/FileIO/FEFLOW/FEFLOWMeshInterface.cpp b/Applications/FileIO/FEFLOW/FEFLOWMeshInterface.cpp
index 9831ecbec45b73f3f264cc74ae8447e502916a4c..358ba14108378c3a62b3e08b26162d2ca86e647a 100644
--- a/Applications/FileIO/FEFLOW/FEFLOWMeshInterface.cpp
+++ b/Applications/FileIO/FEFLOW/FEFLOWMeshInterface.cpp
@@ -10,6 +10,7 @@
 
 #include <boost/algorithm/string/trim.hpp>
 #include <cctype>
+#include <fstream>
 #include <memory>
 
 #include "Applications/FileIO/FEFLOW/FEFLOWGeoInterface.h"
diff --git a/Applications/FileIO/GocadIO/CoordinateSystem.h b/Applications/FileIO/GocadIO/CoordinateSystem.h
index c0a5aea65b689125cd9aef29cc0fff823a4f5fe7..e29cc65a8bfe22557c6b592b7915e1c5c425e3a6 100644
--- a/Applications/FileIO/GocadIO/CoordinateSystem.h
+++ b/Applications/FileIO/GocadIO/CoordinateSystem.h
@@ -9,7 +9,8 @@
 
 #pragma once
 
-#include <fstream>
+#include <iosfwd>
+#include <string>
 
 namespace FileIO
 {
diff --git a/Applications/FileIO/GocadIO/GocadAsciiReader.cpp b/Applications/FileIO/GocadIO/GocadAsciiReader.cpp
index 7bcd07b56a8a947b5a0c1fa2cc7e173c18726d98..ca1d40d3f9ae6b018c9ec5ee7bd8785d6d099ad1 100644
--- a/Applications/FileIO/GocadIO/GocadAsciiReader.cpp
+++ b/Applications/FileIO/GocadIO/GocadAsciiReader.cpp
@@ -9,7 +9,7 @@
 
 #include "GocadAsciiReader.h"
 
-#include <iosfwd>
+#include <fstream>
 
 #include "Applications/FileIO/GocadIO/CoordinateSystem.h"
 #include "BaseLib/FileTools.h"
diff --git a/Applications/FileIO/Legacy/OGSIOVer4.cpp b/Applications/FileIO/Legacy/OGSIOVer4.cpp
index 763cebfd9f454b60cf42955b6b8c0b53171cded6..ec6808de12b580972cdea5eaa8615f198e531182 100644
--- a/Applications/FileIO/Legacy/OGSIOVer4.cpp
+++ b/Applications/FileIO/Legacy/OGSIOVer4.cpp
@@ -14,6 +14,7 @@
 
 #include "OGSIOVer4.h"
 
+#include <fstream>
 #include <iomanip>
 #include <limits>
 #include <sstream>
diff --git a/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp b/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp
index 7bb834e7b28cb00337d7ec819ff5cb77f6653a68..ef556bcaa99f2a5b2ee70ac41fa4913e5eba34b4 100644
--- a/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp
+++ b/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp
@@ -18,6 +18,7 @@
 #include <mpi.h>
 #endif
 
+#include <fstream>
 #include <memory>
 
 #include "BaseLib/FileTools.h"
diff --git a/Applications/Utils/MeshEdit/MeshMapping.cpp b/Applications/Utils/MeshEdit/MeshMapping.cpp
index 311000c52ccc75e5ea222fe9b8f42fc3fcc378ae..e9ce327c60b8cae914092e164df96c6452635a5b 100644
--- a/Applications/Utils/MeshEdit/MeshMapping.cpp
+++ b/Applications/Utils/MeshEdit/MeshMapping.cpp
@@ -11,6 +11,7 @@
 
 #include <tclap/CmdLine.h>
 
+#include <fstream>
 #include <memory>
 #include <string>
 
diff --git a/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp b/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp
index 355c00e3269973ce4554dadae6a2c24898bd2e25..149aaea0c1b3e919a505bea20bdb0a5c9a14ae5c 100644
--- a/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp
+++ b/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp
@@ -17,6 +17,7 @@
 #endif
 
 #include <algorithm>
+#include <fstream>
 #include <memory>
 #include <string>
 #include <vector>
diff --git a/Applications/Utils/ModelPreparation/PartitionMesh/NodeWiseMeshPartitioner.cpp b/Applications/Utils/ModelPreparation/PartitionMesh/NodeWiseMeshPartitioner.cpp
index d2d627a7141cba19145b6ea0dd1bae3e9e25d72d..1942eb37d0aeae8ec708f7f2a59d38f9bafdb819 100644
--- a/Applications/Utils/ModelPreparation/PartitionMesh/NodeWiseMeshPartitioner.cpp
+++ b/Applications/Utils/ModelPreparation/PartitionMesh/NodeWiseMeshPartitioner.cpp
@@ -21,7 +21,6 @@
 #include "BaseLib/Error.h"
 #include "BaseLib/FileTools.h"
 #include "BaseLib/Logging.h"
-#include "BaseLib/Stream.h"
 #include "IntegrationPointDataTools.h"
 #include "MeshLib/Elements/Elements.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
diff --git a/BaseLib/FileTools.cpp b/BaseLib/FileTools.cpp
index cafbd1c6ba95f2de744b534513bbf0f828d13e28..ffe7b4c32e7f92a4114b24fc62bb34d7a4dc7042 100644
--- a/BaseLib/FileTools.cpp
+++ b/BaseLib/FileTools.cpp
@@ -16,6 +16,7 @@
 
 #include <boost/algorithm/string.hpp>
 #include <filesystem>
+#include <fstream>
 #include <typeindex>
 #include <unordered_map>
 
@@ -259,4 +260,69 @@ void removeFiles(std::vector<std::string> const& files)
         removeFile(file);
     }
 }
+
+template <typename T>
+T readBinaryValue(std::istream& in)
+{
+    T v;
+    in.read(reinterpret_cast<char*>(&v), sizeof(T));
+    return v;
+}
+
+// explicit template instantiation
+template float readBinaryValue<float>(std::istream&);
+template double readBinaryValue<double>(std::istream&);
+
+template <typename T>
+std::vector<T> readBinaryArray(std::string const& filename, std::size_t const n)
+{
+    std::ifstream in(filename.c_str());
+    if (!in)
+    {
+        ERR("readBinaryArray(): Error while reading from file '{:s}'.",
+            filename);
+        ERR("Could not open file '{:s}' for input.", filename);
+        in.close();
+        return std::vector<T>();
+    }
+
+    std::vector<T> result;
+    result.reserve(n);
+
+    for (std::size_t p = 0; in && !in.eof() && p < n; ++p)
+    {
+        result.push_back(BaseLib::readBinaryValue<T>(in));
+    }
+
+    if (result.size() == n)
+    {
+        return result;
+    }
+
+    ERR("readBinaryArray(): Error while reading from file '{:s}'.", filename);
+    ERR("Read different number of values. Expected {:d}, got {:d}.",
+        n,
+        result.size());
+
+    if (!in.eof())
+    {
+        ERR("EOF reached.\n");
+    }
+
+    return std::vector<T>();
+}
+
+// explicit template instantiation
+template std::vector<float> readBinaryArray<float>(std::string const&,
+                                                   std::size_t const);
+
+template <typename T>
+void writeValueBinary(std::ostream& out, T const& val)
+{
+    out.write(reinterpret_cast<const char*>(&val), sizeof(T));
+}
+
+// explicit template instantiation
+template void writeValueBinary<std::size_t>(std::ostream&, std::size_t const&);
+
 }  // end namespace BaseLib
diff --git a/BaseLib/FileTools.h b/BaseLib/FileTools.h
index 504398e1e1924f8f71576ed61dbaaf2df7c7c29a..6376fc1fa9b322320322a9e35215a79c7e13f12d 100644
--- a/BaseLib/FileTools.h
+++ b/BaseLib/FileTools.h
@@ -14,12 +14,10 @@
 
 #pragma once
 
-#include <fstream>
+#include <iosfwd>
 #include <string>
 #include <vector>
 
-#include "Logging.h"
-
 namespace BaseLib
 {
 /**
@@ -60,10 +58,7 @@ std::string constructFormattedFileName(std::string const& format_specification,
  * \param val   value
  */
 template <typename T>
-void writeValueBinary(std::ostream& out, T const& val)
-{
-    out.write(reinterpret_cast<const char*>(&val), sizeof(T));
-}
+void writeValueBinary(std::ostream& out, T const& val);
 
 template <typename T>
 T swapEndianness(T const& v)
@@ -86,51 +81,11 @@ T swapEndianness(T const& v)
 double swapEndianness(double const& v);
 
 template <typename T>
-T readBinaryValue(std::istream& in)
-{
-    T v;
-    in.read(reinterpret_cast<char*>(&v), sizeof(T));
-    return v;
-}
+T readBinaryValue(std::istream& in);
 
 template <typename T>
-std::vector<T> readBinaryArray(std::string const& filename, std::size_t const n)
-{
-    std::ifstream in(filename.c_str());
-    if (!in)
-    {
-        ERR("readBinaryArray(): Error while reading from file '{:s}'.",
-            filename);
-        ERR("Could not open file '{:s}' for input.", filename);
-        in.close();
-        return std::vector<T>();
-    }
-
-    std::vector<T> result;
-    result.reserve(n);
-
-    for (std::size_t p = 0; in && !in.eof() && p < n; ++p)
-    {
-        result.push_back(BaseLib::readBinaryValue<T>(in));
-    }
-
-    if (result.size() == n)
-    {
-        return result;
-    }
-
-    ERR("readBinaryArray(): Error while reading from file '{:s}'.", filename);
-    ERR("Read different number of values. Expected {:d}, got {:d}.",
-        n,
-        result.size());
-
-    if (!in.eof())
-    {
-        ERR("EOF reached.\n");
-    }
-
-    return std::vector<T>();
-}
+std::vector<T> readBinaryArray(std::string const& filename,
+                               std::size_t const n);
 
 /**
  * Extracts basename from given pathname with extension.
diff --git a/BaseLib/Histogram.cpp b/BaseLib/Histogram.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..7d812272723d5de9e9095b3be158b3e6cd5ceb43
--- /dev/null
+++ b/BaseLib/Histogram.cpp
@@ -0,0 +1,87 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2022, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "Histogram.h"
+
+#include <fstream>
+
+#include "BaseLib/Logging.h"
+
+namespace BaseLib
+{
+
+template <typename T>
+int Histogram<T>::write(std::string const& file_name,
+                        std::string const& data_set_name,
+                        std::string const& param_name) const
+{
+    if (file_name.empty())
+    {
+        ERR("No file name specified.");
+        return 1;
+    }
+
+    std::ofstream out(file_name);
+    if (!out)
+    {
+        ERR("Error writing histogram: Could not open file.");
+        return 1;
+    }
+
+    out << "# Histogram for parameter " << param_name << " of data set "
+        << data_set_name << "\n";
+    std::size_t const n_bins = this->getNumberOfBins();
+    std::vector<std::size_t> const& bin_cnts(this->getBinCounts());
+    double const min(this->getMinimum());
+    double const bin_width(this->getBinWidth());
+
+    for (std::size_t k(0); k < n_bins; k++)
+    {
+        out << min + k * bin_width << " " << bin_cnts[k] << "\n";
+    }
+    out.close();
+    return 0;
+}
+
+template <typename T>
+void Histogram<T>::prettyPrint(std::ostream& os,
+                               const unsigned int line_width) const
+{
+    const std::size_t count_max =
+        *std::max_element(histogram_.begin(), histogram_.end());
+    for (unsigned int bin = 0; bin < nr_bins_; ++bin)
+    {
+        os << "[" << min_ + bin * bin_width_ << ", "
+           << min_ + (bin + 1) * bin_width_ << ")\t";
+        os << histogram_[bin] << "\t";
+
+        const int n_stars = static_cast<int>(
+            std::ceil(line_width * ((double)histogram_[bin] / count_max)));
+        for (int star = 0; star < n_stars; star++)
+        {
+            os << "*";
+        }
+        os << "\n";
+    }
+}
+
+template <typename T>
+std::ostream& operator<<(std::ostream& os, const Histogram<T>& h)
+{
+    os << h.getNumberOfBins() << " " << h.getMinimum() << " " << h.getMaximum()
+       << " ";
+    std::copy(h.getBinCounts().begin(), h.getBinCounts().end(),
+              std::ostream_iterator<std::size_t>(os, " "));
+    return os << std::endl;
+}
+
+template class Histogram<double>;
+template std::ostream& operator<<(std::ostream& os, const Histogram<double>& h);
+}  // namespace BaseLib
diff --git a/BaseLib/Histogram.h b/BaseLib/Histogram.h
index 031b4389556c03371e5b993a75527f276bb6201f..8428c2547d0518a9a981b106c1a90c153f01bc02 100644
--- a/BaseLib/Histogram.h
+++ b/BaseLib/Histogram.h
@@ -14,14 +14,11 @@
 
 #include <algorithm>
 #include <cmath>
-#include <fstream>
 #include <iosfwd>
-#include <iterator>
+#include <string>
 #include <utility>
 #include <vector>
 
-#include "Logging.h"
-
 namespace BaseLib
 {
 /** Basic Histogram implementation.
@@ -115,56 +112,11 @@ public:
     const T& getMaximum() const { return max_; }
     const T& getBinWidth() const { return bin_width_; }
 
-    void prettyPrint(std::ostream& os, const unsigned int line_width = 16) const
-    {
-        const std::size_t count_max =
-            *std::max_element(histogram_.begin(), histogram_.end());
-        for (unsigned int bin = 0; bin < nr_bins_; ++bin)
-        {
-            os << "[" << min_ + bin * bin_width_ << ", "
-               << min_ + (bin + 1) * bin_width_ << ")\t";
-            os << histogram_[bin] << "\t";
-
-            const int n_stars =
-                std::ceil(line_width * ((double)histogram_[bin] / count_max));
-            for (int star = 0; star < n_stars; star++)
-            {
-                os << "*";
-            }
-            os << "\n";
-        }
-    }
+    void prettyPrint(std::ostream& os,
+                     const unsigned int line_width = 16) const;
 
     int write(std::string const& file_name, std::string const& data_set_name,
-              std::string const& param_name) const
-    {
-        if (file_name.empty())
-        {
-            ERR("No file name specified.");
-            return 1;
-        }
-
-        std::ofstream out(file_name);
-        if (!out)
-        {
-            ERR("Error writing histogram: Could not open file.");
-            return 1;
-        }
-
-        out << "# Histogram for parameter " << param_name << " of data set "
-            << data_set_name << "\n";
-        std::size_t const n_bins = this->getNumberOfBins();
-        std::vector<std::size_t> const& bin_cnts(this->getBinCounts());
-        double const min(this->getMinimum());
-        double const bin_width(this->getBinWidth());
-
-        for (std::size_t k(0); k < n_bins; k++)
-        {
-            out << min + k * bin_width << " " << bin_cnts[k] << "\n";
-        }
-        out.close();
-        return 0;
-    }
+              std::string const& param_name) const;
 
 protected:
     /** Initialize class members after constructor call.
@@ -200,12 +152,9 @@ private:
  * number of bins, minimum, maximum, bin0 count, ..., binN-1 count.
  */
 template <typename T>
-std::ostream& operator<<(std::ostream& os, const Histogram<T>& h)
-{
-    os << h.getNumberOfBins() << " " << h.getMinimum() << " " << h.getMaximum()
-       << " ";
-    std::copy(h.getBinCounts().begin(), h.getBinCounts().end(),
-              std::ostream_iterator<T>(os, " "));
-    return os << std::endl;
-}
+std::ostream& operator<<(std::ostream& os, const Histogram<T>& h);
+
+extern template class Histogram<double>;
+extern template std::ostream& operator<<(std::ostream& os,
+                                         const Histogram<double>& h);
 }  // namespace BaseLib
diff --git a/BaseLib/PrjProcessing.cpp b/BaseLib/PrjProcessing.cpp
index 744f730386cc828fab90b9b2a0d24974275a0df1..7e912c6dec4be2b35e6964a89dab0b18705a38aa 100644
--- a/BaseLib/PrjProcessing.cpp
+++ b/BaseLib/PrjProcessing.cpp
@@ -17,6 +17,7 @@
 #include <filesystem>
 #include <fstream>
 #include <regex>
+#include <sstream>
 
 #include "DisableFPE.h"
 #include "Error.h"
diff --git a/BaseLib/PrjProcessing.h b/BaseLib/PrjProcessing.h
index 138ee7d0faa5f7acc6db61c22fec98362731e86d..e08f64f60d6ba553c16d742b99a9374edea42927 100644
--- a/BaseLib/PrjProcessing.h
+++ b/BaseLib/PrjProcessing.h
@@ -10,7 +10,7 @@
 
 #pragma once
 
-#include <sstream>
+#include <iosfwd>
 #include <vector>
 
 namespace BaseLib
diff --git a/BaseLib/Stream.h b/BaseLib/Stream.h
deleted file mode 100644
index a05b6d0fa17f12b45863a7f8c2ec57ea1a75636b..0000000000000000000000000000000000000000
--- a/BaseLib/Stream.h
+++ /dev/null
@@ -1,28 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date   2014-03-07
- * \brief  Helper tools for debugging
- *
- * \file
- * \copyright
- * Copyright (c) 2012-2022, OpenGeoSys Community (http://www.opengeosys.org)
- *            Distributed under a Modified BSD License.
- *              See accompanying file LICENSE.txt or
- *              http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <algorithm>
-#include <iterator>
-#include <ostream>
-#include <vector>
-
-template <typename T>
-std::ostream& operator<<(std::ostream& os, const std::vector<T>& v)
-{
-    std::copy(v.begin(), v.end(), std::ostream_iterator<T>(os, " "));
-    os << "\n";
-    return os;
-}
diff --git a/MaterialLib/MPL/Properties/TemperatureDependentFraction.h b/MaterialLib/MPL/Properties/TemperatureDependentFraction.h
index e7a2f7932d4bced1fe571c1e611f3f677756169b..a68d71a27f7aaa665a470137a5cee3e63598a6da 100644
--- a/MaterialLib/MPL/Properties/TemperatureDependentFraction.h
+++ b/MaterialLib/MPL/Properties/TemperatureDependentFraction.h
@@ -10,9 +10,6 @@
  */
 #pragma once
 
-#include <fstream>
-#include <iostream>
-
 #include "MaterialLib/MPL/Property.h"
 #include "MaterialLib/MPL/Utils/SigmoidFunction.h"
 namespace MaterialPropertyLib
diff --git a/MaterialLib/SolidModels/Ehlers.h b/MaterialLib/SolidModels/Ehlers.h
index d984878cf1397093c1fa8409c7c63f3c80cc43c1..e333206fa6989ab2431bbb1644c892af90c3d37a 100644
--- a/MaterialLib/SolidModels/Ehlers.h
+++ b/MaterialLib/SolidModels/Ehlers.h
@@ -19,7 +19,7 @@
 #pragma once
 
 #ifndef NDEBUG
-#include <ostream>
+#include <iosfwd>
 #endif
 
 #include "BaseLib/Error.h"
diff --git a/MathLib/LinAlg/Eigen/EigenMatrix.cpp b/MathLib/LinAlg/Eigen/EigenMatrix.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..be1aa4c9f66e957c21b4152657115bb0d90e0fed
--- /dev/null
+++ b/MathLib/LinAlg/Eigen/EigenMatrix.cpp
@@ -0,0 +1,40 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2022, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "EigenMatrix.h"
+
+#include <fstream>
+
+namespace MathLib
+{
+
+/// printout this matrix for debugging
+void EigenMatrix::write(const std::string& filename) const
+{
+    std::ofstream of(filename);
+    if (of)
+    {
+        write(of);
+    }
+}
+
+/// printout this matrix for debugging
+void EigenMatrix::write(std::ostream& os) const
+{
+    for (int k = 0; k < mat_.outerSize(); ++k)
+    {
+        for (RawMatrixType::InnerIterator it(mat_, k); it; ++it)
+        {
+            os << it.row() << " " << it.col() << ": " << it.value() << "\n";
+        }
+    }
+    os << std::endl;
+}
+}  // namespace MathLib
diff --git a/MathLib/LinAlg/Eigen/EigenMatrix.h b/MathLib/LinAlg/Eigen/EigenMatrix.h
index 17e53c2f83ec42b30b8c7d27221a0bd13d71a506..b551cbc658404fbacb72a5e865cfcef8d3373cdb 100644
--- a/MathLib/LinAlg/Eigen/EigenMatrix.h
+++ b/MathLib/LinAlg/Eigen/EigenMatrix.h
@@ -11,10 +11,9 @@
 #pragma once
 
 #include <Eigen/Sparse>
-#include <fstream>
+#include <iosfwd>
 #include <string>
 
-#include "EigenVector.h"
 #include "MathLib/LinAlg/RowColumnIndices.h"
 #include "MathLib/LinAlg/SetMatrixSparsity.h"
 
@@ -137,27 +136,10 @@ public:
     static constexpr bool isAssembled() { return true; }
 
     /// printout this matrix for debugging
-    void write(const std::string& filename) const
-    {
-        std::ofstream of(filename);
-        if (of)
-        {
-            write(of);
-        }
-    }
+    void write(const std::string& filename) const;
 
     /// printout this matrix for debugging
-    void write(std::ostream& os) const
-    {
-        for (int k = 0; k < mat_.outerSize(); ++k)
-        {
-            for (RawMatrixType::InnerIterator it(mat_, k); it; ++it)
-            {
-                os << it.row() << " " << it.col() << ": " << it.value() << "\n";
-            }
-        }
-        os << std::endl;
-    }
+    void write(std::ostream& os) const;
 
     RawMatrixType& getRawMatrix() { return mat_; }
     const RawMatrixType& getRawMatrix() const { return mat_; }
diff --git a/MathLib/LinAlg/Eigen/EigenVector.cpp b/MathLib/LinAlg/Eigen/EigenVector.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..07a3b22abc7c2520b33ce0c101b6f7149d2ba465
--- /dev/null
+++ b/MathLib/LinAlg/Eigen/EigenVector.cpp
@@ -0,0 +1,35 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2022, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "EigenVector.h"
+
+#include "EigenMapTools.h"
+
+#ifndef NDEBUG
+#include <fstream>
+#endif
+
+namespace MathLib
+{
+void EigenVector::copyValues(std::vector<double>& u) const
+{
+    u.resize(size());
+    toVector(u) = vec_;
+}
+
+#ifndef NDEBUG
+void EigenVector::write(const std::string& filename) const
+{
+    std::ofstream os(filename);
+    os << vec_;
+}
+#endif
+
+}  // namespace MathLib
diff --git a/MathLib/LinAlg/Eigen/EigenVector.h b/MathLib/LinAlg/Eigen/EigenVector.h
index 7854fd4d3900c88e8d52326011c6e8860e1b8b8f..50828cf5309dacb4bc7905248abd8a8f9514a826 100644
--- a/MathLib/LinAlg/Eigen/EigenVector.h
+++ b/MathLib/LinAlg/Eigen/EigenVector.h
@@ -11,16 +11,13 @@
 #pragma once
 
 #include <vector>
+
 #ifndef NDEBUG
-#include <fstream>
 #include <string>
 #endif
 
-#include <Eigen/Core>
 #include <Eigen/Sparse>
 
-#include "EigenMapTools.h"
-
 namespace MathLib
 {
 /// Global vector based on Eigen vector
@@ -103,19 +100,11 @@ public:
     /// Copy local entries to a vector.
     /// \param u a vector for the values of local entries. It will be resized to
     /// hold the current vector data.
-    void copyValues(std::vector<double>& u) const
-    {
-        u.resize(size());
-        toVector(u) = vec_;
-    }
+    void copyValues(std::vector<double>& u) const;
 
 #ifndef NDEBUG
-    /// printout this equation for debugging
-    void write(const std::string& filename) const
-    {
-        std::ofstream os(filename);
-        os << vec_;
-    }
+    /// write this vector to a file for debugging
+    void write(const std::string& filename) const;
 #endif
 
     /// return a raw Eigen vector object
diff --git a/MeshLib/IO/Legacy/MeshIO.cpp b/MeshLib/IO/Legacy/MeshIO.cpp
index 3eb1c38874cbff667f62a576eca39afc01fbb273..6cc4a378c3ac9ce2648d5d202e650e525dee4887 100644
--- a/MeshLib/IO/Legacy/MeshIO.cpp
+++ b/MeshLib/IO/Legacy/MeshIO.cpp
@@ -18,6 +18,7 @@
 
 #include "MeshIO.h"
 
+#include <fstream>
 #include <iomanip>
 #include <memory>
 #include <sstream>
diff --git a/MeshLib/PropertyVector.h b/MeshLib/PropertyVector.h
index dea884ae9539c00868b71b6c688891fd7c27108b..0eb345ea0dcdb490796e9f7c291fa2c3e4411699 100644
--- a/MeshLib/PropertyVector.h
+++ b/MeshLib/PropertyVector.h
@@ -11,7 +11,6 @@
 
 #include <cassert>
 #include <iterator>
-#include <ostream>
 #include <string>
 #include <utility>
 #include <vector>
@@ -238,23 +237,6 @@ public:
         return p[component];
     }
 
-#ifndef NDEBUG
-    std::ostream& print(std::ostream& os) const
-    {
-        os << "\nPropertyVector<T*> at address: " << this << ":\n";
-        os << "\tmapping (" << size() << "):\n";
-        std::copy(this->cbegin(), this->cend(),
-                  std::ostream_iterator<std::size_t>(os, " "));
-        os << "\n\tvalues (" << _values.size() << "):\n";
-        for (std::size_t k(0); k < _values.size(); k++)
-        {
-            os << "val: " << *(_values[k]) << ", address: " << _values[k]
-               << "\n";
-        }
-        return os;
-    }
-#endif
-
 protected:
     /// @brief The constructor taking meta information for the data.
     /// @param n_prop_groups number of different property values
diff --git a/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h b/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h
index 33193496d7227c59a7899aa07f961e44867b8b86..556eff411197a2864f81f5d1ad769281e4f236d9 100644
--- a/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h
+++ b/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h
@@ -13,7 +13,7 @@
 #pragma once
 
 #include <Eigen/Core>
-#include <ostream>
+#include <iosfwd>
 
 namespace NumLib
 {
diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/Python/BHEInflowPythonBoundaryCondition.h b/ProcessLib/BoundaryConditionAndSourceTerm/Python/BHEInflowPythonBoundaryCondition.h
index b5fbebb6c22701f977dadbf09c771eff78c3e342..2ba00a2bdf466442ccf9a4d1d918396f12210633 100644
--- a/ProcessLib/BoundaryConditionAndSourceTerm/Python/BHEInflowPythonBoundaryCondition.h
+++ b/ProcessLib/BoundaryConditionAndSourceTerm/Python/BHEInflowPythonBoundaryCondition.h
@@ -13,7 +13,6 @@
 #include <pybind11/pybind11.h>
 #endif  // OGS_USE_PYTHON
 #include <algorithm>
-#include <iostream>
 #include <vector>
 
 #include "BHEInflowPythonBoundaryConditionPythonSideInterface.h"
diff --git a/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h b/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h
index b7c3e33d48c42c462e81a4473f6be6096c2dc0ef..7596c26d78e1e4f7de50596b32efdd4e22615140 100644
--- a/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h
+++ b/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h
@@ -11,7 +11,7 @@
 
 #include <Eigen/LU>
 #include <cassert>
-#include <ostream>
+#include <iosfwd>
 
 #include "MaterialLib/MPL/Medium.h"
 #include "MaterialLib/MPL/Utils/FormEigenTensor.h"
diff --git a/Tests/NumLib/ODEs.h b/Tests/NumLib/ODEs.h
index d60572c38de3ea6a2f09924c3bc20afcaac91ca1..c3df0a4e1931bfcb0a281374108bd12ea5105f90 100644
--- a/Tests/NumLib/ODEs.h
+++ b/Tests/NumLib/ODEs.h
@@ -10,13 +10,11 @@
 #pragma once
 
 #include <boost/math/constants/constants.hpp>
+
 #include "MathLib/LinAlg/LinAlg.h"
 #include "MathLib/LinAlg/UnifiedMatrixSetters.h"
 #include "NumLib/ODESolver/ODESystem.h"
 
-// debug
-//#include <iostream>
-
 template <class Ode>
 class ODETraits;