diff --git a/ProcessLib/GroundwaterFlowProcess.h b/ProcessLib/GroundwaterFlowProcess.h
index 1f6530c42592e9a7e2481125bed8558aff8f875d..df7efde53b018db15c94c9d392d3fcc675372e0e 100644
--- a/ProcessLib/GroundwaterFlowProcess.h
+++ b/ProcessLib/GroundwaterFlowProcess.h
@@ -22,12 +22,10 @@
 #include "AssemblerLib/LocalAssemblerBuilder.h"
 #include "AssemblerLib/VectorMatrixAssembler.h"
 #include "AssemblerLib/LocalDataInitializer.h"
-#include "AssemblerLib/ComputeSparsityPattern.h"
 
 #include "FileIO/VtkIO/VtuInterface.h"
 
 #include "MathLib/LinAlg/ApplyKnownSolution.h"
-#include "MathLib/LinAlg/SetMatrixSparsity.h"
 
 #include "MeshLib/MeshSubset.h"
 #include "MeshLib/MeshSubsets.h"
@@ -218,10 +216,8 @@ public:
             new AssemblerLib::LocalToGlobalIndexMap(_all_mesh_subsets, AssemblerLib::ComponentOrder::BY_COMPONENT));
 
         DBUG("Compute sparsity pattern");
-        _sparsity_pattern = std::move(
-            AssemblerLib::computeSparsityPattern(
-                *_local_to_global_index_map, this->_mesh));
-
+        Process<GlobalSetup>::computeSparsityPattern(
+            *_local_to_global_index_map);
 
         // create global vectors and linear solver
         Process<GlobalSetup>::createLinearSolver(*_local_to_global_index_map,
@@ -258,8 +254,6 @@ public:
     {
         DBUG("Assemble GroundwaterFlowProcess.");
 
-        this->_A->setZero();
-        MathLib::setMatrixSparsity(*this->_A, _sparsity_pattern);
         *this->_rhs = 0;   // This resets the whole vector.
 
         // Call global assembler for each local assembly item.
@@ -372,8 +366,6 @@ private:
     } _dirichlet_bc;
 
     std::vector<NeumannBc<GlobalSetup>*> _neumann_bcs;
-
-    AssemblerLib::SparsityPattern _sparsity_pattern;
 };
 
 }   // namespace ProcessLib
diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
index a54dcf5c45d987d4f9ae61fcc508791ccabb9095..c3e0b1d0c808787c649cabc4834601a3a975cdc6 100644
--- a/ProcessLib/Process.h
+++ b/ProcessLib/Process.h
@@ -12,8 +12,10 @@
 
 #include <string>
 
+#include "AssemblerLib/ComputeSparsityPattern.h"
 #include "AssemblerLib/LocalToGlobalIndexMap.h"
 #include "BaseLib/ConfigTree.h"
+#include "MathLib/LinAlg/SetMatrixSparsity.h"
 
 #ifdef USE_PETSC
 #include "MeshLib/NodePartitionedMesh.h"
@@ -45,6 +47,9 @@ public:
 
 	bool solve(const double delta_t)
 	{
+		_A->setZero();
+		MathLib::setMatrixSparsity(*_A, _sparsity_pattern);
+
 		bool const result = assemble(delta_t);
 
 		_linear_solver->solve(*_rhs, *_x);
@@ -84,6 +89,14 @@ protected:
 		    *_A, solver_name, _linear_solver_options.get()));
 	}
 
+	/// Computes and stores global matrix' sparsity pattern from given
+	/// DOF-table.
+	void computeSparsityPattern(
+	    AssemblerLib::LocalToGlobalIndexMap const& local_to_global_index_map)
+	{
+		_sparsity_pattern = std::move(AssemblerLib::computeSparsityPattern(
+		    local_to_global_index_map, _mesh));
+	}
 
 protected:
 	MeshLib::Mesh& _mesh;
@@ -96,6 +109,8 @@ protected:
 	std::unique_ptr<typename GlobalSetup::MatrixType> _A;
 	std::unique_ptr<typename GlobalSetup::VectorType> _rhs;
 	std::unique_ptr<typename GlobalSetup::VectorType> _x;
+
+	AssemblerLib::SparsityPattern _sparsity_pattern;
 };
 
 }  // namespace ProcessLib