diff --git a/Applications/ApplicationsLib/LinearSolverLibrarySetup.h b/Applications/ApplicationsLib/LinearSolverLibrarySetup.h
index 59c94cd317dd620fda5c1493ca00edba7bce144b..e6b48ce2f6ad7506099b3811ca8fd604fe6b47f3 100644
--- a/Applications/ApplicationsLib/LinearSolverLibrarySetup.h
+++ b/Applications/ApplicationsLib/LinearSolverLibrarySetup.h
@@ -19,6 +19,7 @@
 /// The default implementation is empty providing polymorphic behaviour when
 /// using this class.
 
+#include "BaseLib/MPI.h"
 #include "NumLib/DOF/GlobalMatrixProviders.h"
 
 #if defined(USE_PETSC)
@@ -28,9 +29,8 @@ namespace ApplicationsLib
 {
 struct LinearSolverLibrarySetup final
 {
-    LinearSolverLibrarySetup(int argc, char* argv[])
+    LinearSolverLibrarySetup(int argc, char* argv[]) : mpi_setup(argc, argv)
     {
-        MPI_Init(&argc, &argv);
         char help[] = "ogs6 with PETSc \n";
         PetscInitialize(&argc, &argv, nullptr, help);
         MPI_Comm_set_errhandler(PETSC_COMM_WORLD, MPI_ERRORS_RETURN);
@@ -40,8 +40,9 @@ struct LinearSolverLibrarySetup final
     {
         NumLib::cleanupGlobalMatrixProviders();
         PetscFinalize();
-        MPI_Finalize();
     }
+
+    BaseLib::MPI::Setup mpi_setup;
 };
 }    // ApplicationsLib
 #elif defined(USE_LIS)
diff --git a/Applications/Utils/FileConverter/ConvertSHPToGLI.cpp b/Applications/Utils/FileConverter/ConvertSHPToGLI.cpp
index 2fb39d512d68b563c2384c9ab34ab0f7eac93320..6023763a1b083f5aeeb3ab264eca6c76faab5ed4 100644
--- a/Applications/Utils/FileConverter/ConvertSHPToGLI.cpp
+++ b/Applications/Utils/FileConverter/ConvertSHPToGLI.cpp
@@ -14,10 +14,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 // STL
 #include <fstream>
 #include <vector>
@@ -25,6 +21,7 @@
 // ShapeLib
 #include <shapefil.h>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Qt/XmlGmlInterface.h"
 #include "GeoLib/IO/XmlIO/Qt/XmlStnInterface.h"
@@ -185,9 +182,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string fname(shapefile_arg.getValue());
 
@@ -207,9 +202,6 @@ int main(int argc, char* argv[])
             ERR("Shape file contains {:d} polylines.", number_of_elements);
             ERR("This programme only handles only files containing points.");
             SHPClose(hSHP);
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_SUCCESS;
         }
         SHPClose(hSHP);
@@ -304,8 +296,5 @@ int main(int argc, char* argv[])
         ERR("Could not open the database file.");
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/FEFLOW2OGS.cpp b/Applications/Utils/FileConverter/FEFLOW2OGS.cpp
index acd33c06c8998202652ef1f1fdc85b3fb216e0c2..e299c2b56d2b6cfe1c98f0b0247a0d269f2008c4 100644
--- a/Applications/Utils/FileConverter/FEFLOW2OGS.cpp
+++ b/Applications/Utils/FileConverter/FEFLOW2OGS.cpp
@@ -14,12 +14,8 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
-// BaseLib
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "InfoLib/GitInfo.h"
 #ifndef WIN32
@@ -62,9 +58,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // *** read mesh
     INFO("Reading {:s}.", feflow_mesh_arg.getValue());
@@ -81,9 +75,6 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         INFO("Could not read mesh from {:s}.", feflow_mesh_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 #ifndef WIN32
@@ -99,8 +90,5 @@ int main(int argc, char* argv[])
     INFO("Writing {:s}.", ogs_mesh_fname);
     MeshLib::IO::writeMeshToFile(*mesh, ogs_mesh_fname);
     INFO("\tDone.");
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/GMSH2OGS.cpp b/Applications/Utils/FileConverter/GMSH2OGS.cpp
index a528076749b32809089eaf5bf03ad9a226cb3c3f..6674dbc902fdf8828b4626c5c575ca517f782fdd 100644
--- a/Applications/Utils/FileConverter/GMSH2OGS.cpp
+++ b/Applications/Utils/FileConverter/GMSH2OGS.cpp
@@ -19,12 +19,8 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
-// BaseLib
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "InfoLib/GitInfo.h"
 #ifndef WIN32
@@ -176,9 +172,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // *** read mesh
     INFO("Reading {:s}.", gmsh_mesh_arg.getValue());
@@ -194,9 +188,6 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         INFO("Could not read mesh from {:s}.", gmsh_mesh_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 #ifndef WIN32
@@ -279,8 +270,5 @@ int main(int argc, char* argv[])
     MeshLib::IO::writeMeshToFile(*mesh, ogs_mesh_arg.getValue());
 
     delete mesh;
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/GocadSGridReader.cpp b/Applications/Utils/FileConverter/GocadSGridReader.cpp
index 0f3e6784c141519e3c80ef9a3b96d3c416148c28..a28124cc5c313cf0e8e36b78d718179614f285bb 100644
--- a/Applications/Utils/FileConverter/GocadSGridReader.cpp
+++ b/Applications/Utils/FileConverter/GocadSGridReader.cpp
@@ -12,16 +12,13 @@
 #include <spdlog/spdlog.h>
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <fstream>
 #include <sstream>
 #include <string>
 
 #include "Applications/FileIO/GocadIO/GenerateFaceSetMeshes.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -61,9 +58,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read the Gocad SGrid
     INFO("Start reading Gocad SGrid.");
@@ -82,8 +77,5 @@ int main(int argc, char* argv[])
 
     INFO("Writing mesh to '{:s}'.", mesh_output_arg.getValue());
     MeshLib::IO::writeMeshToFile(*mesh, mesh_output_arg.getValue());
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/GocadTSurfaceReader.cpp b/Applications/Utils/FileConverter/GocadTSurfaceReader.cpp
index dd93c5024e1aa173a98d928661f29279c39fd536..8bb1391f3e41f3e5af1d40d23f7cabe88ac010f1 100644
--- a/Applications/Utils/FileConverter/GocadTSurfaceReader.cpp
+++ b/Applications/Utils/FileConverter/GocadTSurfaceReader.cpp
@@ -9,11 +9,8 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/GocadIO/GocadAsciiReader.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
 #include "MeshLib/Mesh.h"
@@ -65,17 +62,12 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (export_lines_arg.isSet() && export_surfaces_arg.isSet())
     {
         ERR("Both the 'lines-only'-flag and 'surfaces-only'-flag are set. Only "
             "one is allowed at a time.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 2;
     }
 
@@ -94,9 +86,6 @@ int main(int argc, char* argv[])
     if (!FileIO::Gocad::GocadAsciiReader::readFile(file_name, meshes, t))
     {
         ERR("Error reading file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 1;
     }
     INFO("{:d} meshes found.", meshes.size());
@@ -115,8 +104,5 @@ int main(int argc, char* argv[])
         MeshLib::IO::VtuInterface vtu(mesh.get(), data_mode, compressed);
         vtu.writeToFile(dir + delim + mesh->getName() + ".vtu");
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return 0;
 }
diff --git a/Applications/Utils/FileConverter/Mesh2Raster.cpp b/Applications/Utils/FileConverter/Mesh2Raster.cpp
index 4e32f111d57e707572a6fb0c5ac85e8f0fb556c4..9b090554776f0416475b75b7fa709e443a1fcd77 100644
--- a/Applications/Utils/FileConverter/Mesh2Raster.cpp
+++ b/Applications/Utils/FileConverter/Mesh2Raster.cpp
@@ -9,15 +9,12 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <filesystem>
 #include <fstream>
 #include <memory>
 #include <string>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/AABB.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -52,9 +49,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     INFO("Rasterising mesh...");
     std::unique_ptr<MeshLib::Mesh> const mesh(
@@ -62,18 +57,12 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Error reading mesh file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 1;
     }
     if (mesh->getDimension() != 2)
     {
         ERR("The programme requires a mesh containing two-dimensional elements "
             "(i.e. triangles or quadrilaterals.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 2;
     }
 
@@ -178,8 +167,5 @@ int main(int argc, char* argv[])
     }
     out.close();
     INFO("Result written to {:s}", output_name);
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return 0;
 }
diff --git a/Applications/Utils/FileConverter/Mesh2Shape.cpp b/Applications/Utils/FileConverter/Mesh2Shape.cpp
index b9ffd3f0493f43684b74c9d7fce4982f2e43840c..44bb45b383b51c4a733d30b8df1883fce431d38a 100644
--- a/Applications/Utils/FileConverter/Mesh2Shape.cpp
+++ b/Applications/Utils/FileConverter/Mesh2Shape.cpp
@@ -9,11 +9,8 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/SHPInterface.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/Mesh.h"
@@ -43,22 +40,14 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string const file_name(input_arg.getValue());
     std::unique_ptr<MeshLib::Mesh> const mesh(
         MeshLib::IO::readMeshFromFile(file_name));
     if (FileIO::SHPInterface::write2dMeshToSHP(output_arg.getValue(), *mesh))
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_SUCCESS;
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_FAILURE;
 }
diff --git a/Applications/Utils/FileConverter/NetCdfConverter.cpp b/Applications/Utils/FileConverter/NetCdfConverter.cpp
index 84e3be1b7ab368f37932339a0ab59c6c7906f93c..8009dd19e5f1553d900ee6c375a3f401f872ceb6 100644
--- a/Applications/Utils/FileConverter/NetCdfConverter.cpp
+++ b/Applications/Utils/FileConverter/NetCdfConverter.cpp
@@ -9,10 +9,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 // STL
 #include <cctype>
 #include <iostream>
@@ -26,6 +22,7 @@
 
 #include "BaseLib/FileTools.h"
 #include "BaseLib/Logging.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
 #include "InfoLib/GitInfo.h"
@@ -740,18 +737,13 @@ int main(int argc, char* argv[])
     cmd.add(arg_input);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     NcFile dataset(arg_input.getValue().c_str(), NcFile::read);
 
     if (dataset.isNull())
     {
         ERR("Error opening file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
@@ -761,9 +753,6 @@ int main(int argc, char* argv[])
     {
         ERR("Only one output format can be specified (single-file, multi-file, "
             "or images)");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -783,9 +772,6 @@ int main(int argc, char* argv[])
     if (var.isNull())
     {
         ERR("Variable \"{:s}\" not found in file.", arg_varname.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -797,9 +783,6 @@ int main(int argc, char* argv[])
         if (!assignDimParams(var, dim_idx_map, arg_dim_time, arg_dim1, arg_dim2,
                              arg_dim3))
         {
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -856,15 +839,9 @@ int main(int argc, char* argv[])
     if (!convert(dataset, var, output_name, dim_idx_map, is_time_dep,
                  time_bounds, output, elem_type))
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
     std::cout << "Conversion finished successfully.\n";
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/OGS2VTK.cpp b/Applications/Utils/FileConverter/OGS2VTK.cpp
index 78f360c28831677e5674e49776b26366f2517e62..46bc603a270d4bc555575f66305ed4d88d67c131 100644
--- a/Applications/Utils/FileConverter/OGS2VTK.cpp
+++ b/Applications/Utils/FileConverter/OGS2VTK.cpp
@@ -13,13 +13,10 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 #include <string>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -52,17 +49,12 @@ int main(int argc, char* argv[])
     cmd.add(use_ascii_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh const> mesh(
         MeshLib::IO::readMeshFromFile(mesh_in.getValue()));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("Mesh read: {:d} nodes, {:d} elements.", mesh->getNumberOfNodes(),
@@ -73,8 +65,5 @@ int main(int argc, char* argv[])
 
     MeshLib::IO::writeVtu(*mesh, mesh_out.getValue(), data_mode);
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/PVD2XDMF.cpp b/Applications/Utils/FileConverter/PVD2XDMF.cpp
index 838aa202559ec5b4e0aefe8870fc7c8d3273bbdb..6c114ef6673ced9e09f7d28a6658af1b2f01f796 100644
--- a/Applications/Utils/FileConverter/PVD2XDMF.cpp
+++ b/Applications/Utils/FileConverter/PVD2XDMF.cpp
@@ -9,10 +9,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <array>
 #include <boost/property_tree/ptree.hpp>
 #include <boost/property_tree/xml_parser.hpp>
@@ -20,6 +16,7 @@
 
 #include "BaseLib/FileTools.h"
 #include "BaseLib/Logging.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/MemWatch.h"
 #include "BaseLib/RunTime.h"
 #include "BaseLib/StringTools.h"
@@ -141,9 +138,7 @@ int main(int argc, char* argv[])
     cmd.add(outdir_arg);
 
     cmd.parse(argc, argv);
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     BaseLib::setConsoleLogLevel(log_level_arg.getValue());
 
     auto const pvd_file_dir = BaseLib::extractPath(pvd_file_arg.getValue());
@@ -237,8 +232,5 @@ int main(int argc, char* argv[])
 
         mesh_xdmf_hdf_writer->writeStep(time);
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/Raster2ASC.cpp b/Applications/Utils/FileConverter/Raster2ASC.cpp
index de81adf7d1929dd4b8a32d491478d7f6b257d942..398d8665c11a6ef2099cf069e6b096019f20175a 100644
--- a/Applications/Utils/FileConverter/Raster2ASC.cpp
+++ b/Applications/Utils/FileConverter/Raster2ASC.cpp
@@ -11,10 +11,7 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
 #include "InfoLib/GitInfo.h"
@@ -41,9 +38,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<GeoLib::Raster> raster(
         FileIO::AsciiRasterInterface::readRaster(input_arg.getValue()));
@@ -61,8 +56,5 @@ int main(int argc, char* argv[])
     }
 
     FileIO::AsciiRasterInterface::writeRasterAsASC(*raster, output_name);
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/TIN2VTK.cpp b/Applications/Utils/FileConverter/TIN2VTK.cpp
index 91cc0a1962bf41ebe3f606f76122d14027b97cac..1d2b7290baabbb87ab31c838a49758dd55131d2c 100644
--- a/Applications/Utils/FileConverter/TIN2VTK.cpp
+++ b/Applications/Utils/FileConverter/TIN2VTK.cpp
@@ -9,17 +9,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 // STL
 #include <memory>
 #include <string>
 #include <vector>
 
-// BaseLib
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 
 // GeoLib
@@ -53,9 +49,7 @@ int main(int argc, char* argv[])
     cmd.add(outArg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     INFO("reading the TIN file...");
     const std::string tinFileName(inArg.getValue());
@@ -66,9 +60,6 @@ int main(int argc, char* argv[])
         GeoLib::IO::TINInterface::readTIN(tinFileName, point_vec));
     if (!sfc)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("TIN read:  {:d} points, {:d} triangles", point_vec.size(),
@@ -85,8 +76,5 @@ int main(int argc, char* argv[])
     MeshLib::IO::VtuInterface writer(mesh.get());
     writer.writeToFile(outArg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/TecPlotTools.cpp b/Applications/Utils/FileConverter/TecPlotTools.cpp
index 0d70176cf1972bd2893645627e388d9a0ff5b391..bf209269396a1bb33f9ef92cd83b66116a705e63 100644
--- a/Applications/Utils/FileConverter/TecPlotTools.cpp
+++ b/Applications/Utils/FileConverter/TecPlotTools.cpp
@@ -9,16 +9,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <iostream>
 #include <memory>
 #include <string>
 #include <string_view>
 #include <vector>
 
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "GeoLib/Point.h"
 #include "InfoLib/GitInfo.h"
@@ -480,25 +477,17 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (!input_arg.isSet())
     {
         ERR("No input file given. Please specify TecPlot (*.plt) file");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
     if (convert_arg.getValue() && !output_arg.isSet())
     {
         ERR("No output file given. Please specify OGS mesh (*.vtu) file");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
@@ -506,18 +495,12 @@ int main(int argc, char* argv[])
     if (!in.is_open())
     {
         ERR("Could not open file {:s}.", input_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -2;
     }
 
     if (!convert_arg.isSet() && !split_arg.isSet())
     {
         INFO("Nothing to do. Use -s to split or -c to convert.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 0;
     }
 
@@ -534,8 +517,5 @@ int main(int argc, char* argv[])
     }
 
     in.close();
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return return_val;
 }
diff --git a/Applications/Utils/FileConverter/VTK2OGS.cpp b/Applications/Utils/FileConverter/VTK2OGS.cpp
index 2ec723dc2e79bf35f90701f6a7e72441c7d50a12..80c424df8a1e37bd1a45d44f717755352e4dcc81 100644
--- a/Applications/Utils/FileConverter/VTK2OGS.cpp
+++ b/Applications/Utils/FileConverter/VTK2OGS.cpp
@@ -14,12 +14,9 @@
 // STL
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <string>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/Legacy/MeshIO.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
@@ -47,9 +44,7 @@ int main(int argc, char* argv[])
     cmd.add(mesh_out);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     MeshLib::Mesh* mesh(
         MeshLib::IO::VtuInterface::readVTUFile(mesh_in.getValue()));
@@ -60,8 +55,5 @@ int main(int argc, char* argv[])
     meshIO.setMesh(mesh);
     BaseLib::IO::writeStringToFile(meshIO.writeToString(), mesh_out.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/VTK2TIN.cpp b/Applications/Utils/FileConverter/VTK2TIN.cpp
index b0401d1b199b66373e9b1a76e1518c38e5f3f640..bde167facc2857c2d42ff72106a214945206a3b0 100644
--- a/Applications/Utils/FileConverter/VTK2TIN.cpp
+++ b/Applications/Utils/FileConverter/VTK2TIN.cpp
@@ -9,17 +9,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 // STL
 #include <fstream>
 #include <memory>
 #include <string>
 
-// BaseLib
 #include "BaseLib/Logging.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 
 // GeoLib
@@ -56,9 +52,7 @@ int main(int argc, char* argv[])
     cmd.add(mesh_out);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     std::unique_ptr<MeshLib::Mesh> mesh(
         MeshLib::IO::VtuInterface::readVTUFile(mesh_in.getValue()));
     INFO("Mesh read: {:d} nodes, {:d} elements.", mesh->getNumberOfNodes(),
@@ -74,8 +68,5 @@ int main(int argc, char* argv[])
             mesh_out.getValue());
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/convertGEO.cpp b/Applications/Utils/FileConverter/convertGEO.cpp
index c7f76b7f641356e020ac363d3d3e341a8c62a6f9..3a7d859d519b1299aeb8c17471a887feae2dc2dd 100644
--- a/Applications/Utils/FileConverter/convertGEO.cpp
+++ b/Applications/Utils/FileConverter/convertGEO.cpp
@@ -10,15 +10,12 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <string>
 #include <vector>
 
 #include "Applications/FileIO/readGeometryFromFile.h"
 #include "Applications/FileIO/writeGeometryToFile.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "InfoLib/GitInfo.h"
 
@@ -51,9 +48,7 @@ int main(int argc, char* argv[])
     cmd.add(gmsh_path_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     GeoLib::GEOObjects geoObjects;
     FileIO::readGeometryFromFile(argInputFileName.getValue(), geoObjects,
@@ -64,8 +59,5 @@ int main(int argc, char* argv[])
     FileIO::writeGeometryToFile(geo_names[0], geoObjects,
                                 argOutputFileName.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp b/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp
index 9291eacb3b7bdbfa7c53b025b73bc72d7bde5d50..88032d6354ea69f369c9deafe262a90275eb4e39 100644
--- a/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp
+++ b/Applications/Utils/FileConverter/generateMatPropsFromMatID.cpp
@@ -14,14 +14,11 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <fstream>
 #include <memory>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -50,9 +47,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read mesh
     std::unique_ptr<MeshLib::Mesh> mesh(
@@ -61,9 +56,6 @@ int main(int argc, char* argv[])
     if (!mesh)
     {
         INFO("Could not read mesh from file '{:s}'.", mesh_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -92,9 +84,6 @@ int main(int argc, char* argv[])
     else
     {
         ERR("Could not create property '{:s}' file.", new_matname);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -106,8 +95,5 @@ int main(int argc, char* argv[])
 
     INFO("New files '{:s}' and '{:s}' written.", new_mshname, new_matname);
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/GeoTools/MoveGeometry.cpp b/Applications/Utils/GeoTools/MoveGeometry.cpp
index 33505f1d59475065fed2c5305b8595b99f50f4bc..dd14e7210d22b9be1c2e9d7a59984b4ab4254a30 100644
--- a/Applications/Utils/GeoTools/MoveGeometry.cpp
+++ b/Applications/Utils/GeoTools/MoveGeometry.cpp
@@ -14,10 +14,7 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
 #include "InfoLib/GitInfo.h"
@@ -49,9 +46,7 @@ int main(int argc, char* argv[])
     cmd.add(geo_input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     GeoLib::GEOObjects geo_objects;
     GeoLib::IO::BoostXmlGmlInterface xml(geo_objects);
@@ -59,9 +54,6 @@ int main(int argc, char* argv[])
     {
         if (!xml.readFile(geo_input_arg.getValue()))
         {
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -69,9 +61,6 @@ int main(int argc, char* argv[])
     {
         ERR("Failed to read file `{:s}'.", geo_input_arg.getValue());
         ERR("{:s}", err.what());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -106,8 +95,5 @@ int main(int argc, char* argv[])
     BaseLib::IO::writeStringToFile(xml.writeToString(),
                                    geo_output_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/GeoTools/addDataToRaster.cpp b/Applications/Utils/GeoTools/addDataToRaster.cpp
index 425bfb5d68995ded4d0f5e2ed67ef65c81a03dcc..c78dc0549cc968b4fa5c77f25eb82545c8765116 100644
--- a/Applications/Utils/GeoTools/addDataToRaster.cpp
+++ b/Applications/Utils/GeoTools/addDataToRaster.cpp
@@ -11,16 +11,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <cmath>
 #include <memory>
 #include <numbers>
 #include <numeric>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/AABB.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Point.h"
@@ -153,9 +150,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::array input_points = {
         GeoLib::Point{{ll_x_arg.getValue(), ll_y_arg.getValue(), 0}},
@@ -210,8 +205,5 @@ int main(int argc, char* argv[])
 
     FileIO::AsciiRasterInterface::writeRasterAsASC(*raster,
                                                    out_raster_arg.getValue());
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/GeoTools/createRaster.cpp b/Applications/Utils/GeoTools/createRaster.cpp
index 8d36fc335aa1e31bc76ffac6b557fa1aec426ac7..b3aa8b4b8f03bd945fde2a6f0e6111cf59dbd73b 100644
--- a/Applications/Utils/GeoTools/createRaster.cpp
+++ b/Applications/Utils/GeoTools/createRaster.cpp
@@ -11,12 +11,9 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/AABB.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Point.h"
@@ -71,9 +68,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     GeoLib::RasterHeader header{
         n_cols.getValue(),
@@ -88,8 +83,5 @@ int main(int argc, char* argv[])
     FileIO::AsciiRasterInterface::writeRasterAsASC(raster,
                                                    output_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/GeoTools/generateGeometry.cpp b/Applications/Utils/GeoTools/generateGeometry.cpp
index b84b167e37f1402c53aa0ac1113f331ad4d480ed..56c054e79e08d3304fdc27a8c4b97157bdf73327 100644
--- a/Applications/Utils/GeoTools/generateGeometry.cpp
+++ b/Applications/Utils/GeoTools/generateGeometry.cpp
@@ -11,12 +11,9 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <numeric>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
 #include "GeoLib/Point.h"
@@ -241,9 +238,7 @@ int main(int argc, char* argv[])
     cmd.add(geo_output_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     auto const p0 = GeoLib::Point{x0.getValue(), y0.getValue(), z0.getValue()};
     auto const p1 = GeoLib::Point{x1.getValue(), y1.getValue(), z1.getValue()};
@@ -286,9 +281,6 @@ int main(int argc, char* argv[])
                 std::move(polyline_name.getValue()), geometry_name.getValue(),
                 geometry) == EXIT_FAILURE)
         {
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -298,8 +290,5 @@ int main(int argc, char* argv[])
     BaseLib::IO::writeStringToFile(xml.writeToString(),
                                    geo_output_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/AddElementQuality.cpp b/Applications/Utils/MeshEdit/AddElementQuality.cpp
index 25753925c288c2d753e9bf3b0b122affd6ca3730..08c7087db3b853aa5f9f0b787df975fc2bcf7c4e 100644
--- a/Applications/Utils/MeshEdit/AddElementQuality.cpp
+++ b/Applications/Utils/MeshEdit/AddElementQuality.cpp
@@ -9,13 +9,10 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <array>
 #include <string>
 
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -53,9 +50,7 @@ int main(int argc, char* argv[])
     cmd.add(mesh_in_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read the mesh file
     BaseLib::RunTime run_time;
@@ -64,9 +59,6 @@ int main(int argc, char* argv[])
         mesh_in_arg.getValue(), true /* compute_element_neighbors */));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("Time for reading: {:g} s", run_time.elapsed());
@@ -82,8 +74,5 @@ int main(int argc, char* argv[])
     INFO("Writing mesh '{:s}' ... ", mesh_out_arg.getValue());
     MeshLib::IO::writeMeshToFile(*mesh, mesh_out_arg.getValue());
     INFO("done.");
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/AddFaultToVoxelGrid.cpp b/Applications/Utils/MeshEdit/AddFaultToVoxelGrid.cpp
index e27be76f84654260c39335f24a183166d856f77e..2d4e47d92a53ddd7d2881760a84f021d4cf2751c 100644
--- a/Applications/Utils/MeshEdit/AddFaultToVoxelGrid.cpp
+++ b/Applications/Utils/MeshEdit/AddFaultToVoxelGrid.cpp
@@ -15,10 +15,7 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -64,9 +61,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string const input_name = input_arg.getValue();
     std::string const fault_name = fault_arg.getValue();
@@ -81,17 +76,11 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Input mesh not found...");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     auto const& mat_ids = MeshLib::materialIDs(*mesh);
     if (!mat_ids)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         ERR("Input mesh has no material IDs");
         return EXIT_FAILURE;
     }
@@ -105,15 +94,9 @@ int main(int argc, char* argv[])
     {
         MeshLib::IO::VtuInterface vtu(mesh.get());
         vtu.writeToFile(output_name);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         INFO("The fault was successfully added.");
         return EXIT_SUCCESS;
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     ERR("No fault could be added.");
     return EXIT_FAILURE;
 }
diff --git a/Applications/Utils/MeshEdit/AddLayer.cpp b/Applications/Utils/MeshEdit/AddLayer.cpp
index 7958e404e2de021985d5a768e0f156c75b8f388b..90e2d51efb788096293673bde20d6417dba8e092 100644
--- a/Applications/Utils/MeshEdit/AddLayer.cpp
+++ b/Applications/Utils/MeshEdit/AddLayer.cpp
@@ -11,13 +11,10 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -78,15 +75,10 @@ int main(int argc, char* argv[])
     {
         ERR("It is not possible to set both options '--copy-material-ids' and "
             "'--set-material-id'.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     INFO("Reading mesh '{:s}' ... ", mesh_arg.getValue());
     auto subsfc_mesh = std::unique_ptr<MeshLib::Mesh>(
@@ -94,9 +86,6 @@ int main(int argc, char* argv[])
     if (!subsfc_mesh)
     {
         ERR("Error reading mesh '{:s}'.", mesh_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("done.");
@@ -114,9 +103,6 @@ int main(int argc, char* argv[])
     if (!result)
     {
         ERR("Failure while adding layer.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -124,8 +110,5 @@ int main(int argc, char* argv[])
     MeshLib::IO::writeMeshToFile(*result, mesh_out_arg.getValue());
     INFO("done.");
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/CreateBoundaryConditionsAlongPolylines.cpp b/Applications/Utils/MeshEdit/CreateBoundaryConditionsAlongPolylines.cpp
index 070bcedc532acb904b30592aae6d6e0d53d2aa0f..1debc20fb54ab5c2b8ae3d198d87fa3fc9c41bb1 100644
--- a/Applications/Utils/MeshEdit/CreateBoundaryConditionsAlongPolylines.cpp
+++ b/Applications/Utils/MeshEdit/CreateBoundaryConditionsAlongPolylines.cpp
@@ -12,10 +12,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <fstream>
 #include <map>
 #include <string>
@@ -24,6 +20,7 @@
 #include "Applications/FileIO/readGeometryFromFile.h"
 #include "Applications/FileIO/writeGeometryToFile.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/Point.h"
 #include "InfoLib/GitInfo.h"
@@ -183,9 +180,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // *** read mesh
     INFO("Reading mesh '{:s}' ... ", mesh_arg.getValue());
@@ -216,9 +211,6 @@ int main(int argc, char* argv[])
     {
         ERR("Could not get vector of polylines out of geometry '{:s}'.",
             geo_name);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -251,9 +243,6 @@ int main(int argc, char* argv[])
     if (geo_names.empty())
     {
         ERR("Did not find mesh nodes along polylines.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -313,8 +302,5 @@ int main(int argc, char* argv[])
         BaseLib::dropFileExtension(output_base_fname.getValue()));
     writeBCsAndGeometry(geometry_sets, surface_name, base_fname,
                         bc_type.getValue(), gml_arg.getValue());
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/ExtractBoundary.cpp b/Applications/Utils/MeshEdit/ExtractBoundary.cpp
index 744929e8c64338a9b513091c2ca478988bc58cd4..f643dbd7d1052d5f249414fea16971560036902c 100644
--- a/Applications/Utils/MeshEdit/ExtractBoundary.cpp
+++ b/Applications/Utils/MeshEdit/ExtractBoundary.cpp
@@ -11,16 +11,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <memory>
 #include <string>
 #include <vector>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
@@ -63,18 +60,13 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh const> mesh(MeshLib::IO::readMeshFromFile(
         mesh_in.getValue(), true /* compute_element_neighbors */));
 
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -102,8 +94,5 @@ int main(int argc, char* argv[])
         use_ascii_arg.getValue() ? vtkXMLWriter::Ascii : vtkXMLWriter::Binary;
     MeshLib::IO::writeVtu(*surface_mesh, out_fname, data_mode);
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/ExtractMaterials.cpp b/Applications/Utils/MeshEdit/ExtractMaterials.cpp
index 9edf690aae1eb59b6c4f177e7bc46fa77442d728..050cf13656a9e731ff256df2939684162887676c 100644
--- a/Applications/Utils/MeshEdit/ExtractMaterials.cpp
+++ b/Applications/Utils/MeshEdit/ExtractMaterials.cpp
@@ -12,11 +12,8 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -73,9 +70,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string const input_name = input_arg.getValue();
     std::string const output_name = output_arg.getValue();
@@ -87,9 +82,6 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Error reading input mesh. Aborting...");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -97,9 +89,6 @@ int main(int argc, char* argv[])
     if (mat_ids == nullptr)
     {
         ERR("No material IDs found in mesh. Aborting...");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -107,9 +96,6 @@ int main(int argc, char* argv[])
     if (id_range.first == id_range.second)
     {
         ERR("Mesh only contains one material, no extraction required.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     int min_id, max_id;
@@ -119,9 +105,6 @@ int main(int argc, char* argv[])
         if (min_id < *id_range.first || min_id > *id_range.second)
         {
             ERR("Specified material ID does not exist.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         max_id = min_id;
@@ -160,8 +143,5 @@ int main(int argc, char* argv[])
     {
         ostream.close();
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/ExtractSurface.cpp b/Applications/Utils/MeshEdit/ExtractSurface.cpp
index da959739497657ca995cd29790fadf317770955e..f8a95e304de731ab126bcf5e20caee46ff0fa74c 100644
--- a/Applications/Utils/MeshEdit/ExtractSurface.cpp
+++ b/Applications/Utils/MeshEdit/ExtractSurface.cpp
@@ -11,16 +11,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <memory>
 #include <string>
 #include <vector>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
@@ -75,9 +72,7 @@ int main(int argc, char* argv[])
     cmd.add(mesh_in);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh const> mesh(MeshLib::IO::readMeshFromFile(
         mesh_in.getValue(), true /* compute_element_neighbors */));
@@ -85,18 +80,12 @@ int main(int argc, char* argv[])
     if (!mesh)
     {
         ERR("Error reading mesh file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
     if (mesh->getDimension() != 3)
     {
         ERR("Surfaces can currently only be extracted from 3D meshes.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -122,8 +111,5 @@ int main(int argc, char* argv[])
         use_ascii_arg.getValue() ? vtkXMLWriter::Ascii : vtkXMLWriter::Binary;
     MeshLib::IO::writeVtu(*surface_mesh, out_fname, data_mode);
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/Layers2Grid.cpp b/Applications/Utils/MeshEdit/Layers2Grid.cpp
index 6108bee99fd7ff6ee37ae6c35819c27e6b01ab51..12525d1326c185e082bc66469e391a0520b6588f 100644
--- a/Applications/Utils/MeshEdit/Layers2Grid.cpp
+++ b/Applications/Utils/MeshEdit/Layers2Grid.cpp
@@ -15,11 +15,8 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/IO/readStringListFromFile.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/AABB.h"
 #include "InfoLib/GitInfo.h"
 #include "MathLib/Point3d.h"
@@ -78,18 +75,13 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if ((y_arg.isSet() && !z_arg.isSet()) ||
         ((!y_arg.isSet() && z_arg.isSet())))
     {
         ERR("For equilateral cubes, only x needs to be set. For unequal "
             "cuboids, all three edge lengths (x/y/z) need to be specified.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -104,9 +96,6 @@ int main(int argc, char* argv[])
     if (layer_names.size() < 2)
     {
         ERR("At least two layers are required to create a 3D Mesh");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -124,9 +113,6 @@ int main(int argc, char* argv[])
         if (mesh == nullptr)
         {
             ERR("Input layer '{:s}' not found. Aborting...", layer);
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         layers.emplace_back(mesh);
@@ -141,15 +127,9 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("The VoxelGrid could not be created.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     MeshLib::IO::VtuInterface vtu(mesh.get());
     vtu.writeToFile(output_name);
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/MapGeometryToMeshSurface.cpp b/Applications/Utils/MeshEdit/MapGeometryToMeshSurface.cpp
index d32118f004c8c06014c5410b6f95b5adcf4b227c..106e8305d44f62f4e78d5f96aab0ad1c81c2a6da 100644
--- a/Applications/Utils/MeshEdit/MapGeometryToMeshSurface.cpp
+++ b/Applications/Utils/MeshEdit/MapGeometryToMeshSurface.cpp
@@ -11,14 +11,11 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <cstdlib>
 #include <vector>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
 #include "InfoLib/GitInfo.h"
@@ -60,9 +57,7 @@ int main(int argc, char* argv[])
     cmd.add(output_geometry_fname);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // *** read geometry
     GeoLib::GEOObjects geometries;
@@ -75,9 +70,6 @@ int main(int argc, char* argv[])
         }
         else
         {
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -105,8 +97,5 @@ int main(int argc, char* argv[])
         BaseLib::IO::writeStringToFile(xml_io.writeToString(),
                                        output_geometry_fname.getValue());
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/MeshMapping.cpp b/Applications/Utils/MeshEdit/MeshMapping.cpp
index 6b939c1821d47d9e96f3a97599b963b5ff1d35cc..3d3f3139fda3eec08eefe48f29609afa385644ba 100644
--- a/Applications/Utils/MeshEdit/MeshMapping.cpp
+++ b/Applications/Utils/MeshEdit/MeshMapping.cpp
@@ -15,11 +15,8 @@
 #include <memory>
 #include <string>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/AABB.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
@@ -54,9 +51,7 @@ double getClosestPointElevation(MeshLib::Node const& p,
 
 int main(int argc, char* argv[])
 {
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     TCLAP::CmdLine cmd(
         "Changes the elevation of 2D mesh nodes based on either raster data or "
         "another 2D mesh. In addition, a low pass filter can be applied to "
@@ -107,9 +102,6 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Error reading mesh file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -118,18 +110,12 @@ int main(int argc, char* argv[])
     {
         ERR("Nothing to do. Please choose mapping based on a raster or mesh "
             "file, or to a static value.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
     if (map_raster_arg.isSet() && map_mesh_arg.isSet())
     {
         ERR("Please select mapping based on *either* a mesh or a raster file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -148,9 +134,6 @@ int main(int argc, char* argv[])
         if (!file_stream.good())
         {
             ERR("Opening raster file {} failed.", raster_path);
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         file_stream.close();
@@ -168,9 +151,6 @@ int main(int argc, char* argv[])
         if (ground_truth == nullptr)
         {
             ERR("Error reading mesh file.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
 
@@ -243,15 +223,9 @@ int main(int argc, char* argv[])
 
     if (MeshLib::IO::writeMeshToFile(*mesh, output_arg.getValue()) != 0)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
     INFO("Result successfully written.");
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/MoveMesh.cpp b/Applications/Utils/MeshEdit/MoveMesh.cpp
index a755a792143525429c83aad7152557609c65b3bf..65cefe2c3bb0ca7d1df876d104fc4af237a8aeff 100644
--- a/Applications/Utils/MeshEdit/MoveMesh.cpp
+++ b/Applications/Utils/MeshEdit/MoveMesh.cpp
@@ -11,11 +11,8 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "GeoLib/AABB.h"
 #include "InfoLib/GitInfo.h"
@@ -64,9 +61,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     std::string fname(mesh_arg.getValue());
 
     std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::IO::readMeshFromFile(fname));
@@ -74,9 +69,6 @@ int main(int argc, char* argv[])
     if (!mesh)
     {
         ERR("Could not read mesh from file '{:s}'.", fname);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -107,8 +99,5 @@ int main(int argc, char* argv[])
 
     MeshLib::IO::writeMeshToFile(*mesh, out_fname);
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/NodeReordering.cpp b/Applications/Utils/MeshEdit/NodeReordering.cpp
index cc83457c3a48ef1528d41a41ff9ef977e780cb73..f98f006397521cd019316c530621a0b842c26afb 100644
--- a/Applications/Utils/MeshEdit/NodeReordering.cpp
+++ b/Applications/Utils/MeshEdit/NodeReordering.cpp
@@ -11,16 +11,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <array>
 #include <memory>
 #include <vector>
 
 #include "BaseLib/Algorithm.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -191,18 +188,13 @@ int main(int argc, char* argv[])
     cmd.add(input_mesh_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh> mesh(
         MeshLib::IO::readMeshFromFile(input_mesh_arg.getValue()));
 
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -228,8 +220,5 @@ int main(int argc, char* argv[])
 
     INFO("VTU file written.");
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/PVTU2VTU/PVTU2VTU.cpp b/Applications/Utils/MeshEdit/PVTU2VTU/PVTU2VTU.cpp
index 0073291797d34ea3ebd0b4ac61fdbe8382d15795..753004bf2bcd49485eeb033a6948002406e485cc 100644
--- a/Applications/Utils/MeshEdit/PVTU2VTU/PVTU2VTU.cpp
+++ b/Applications/Utils/MeshEdit/PVTU2VTU/PVTU2VTU.cpp
@@ -24,11 +24,8 @@
 #include <unordered_set>
 #include <vector>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "GeoLib/AABB.h"
 #include "GeoLib/OctTree.h"
@@ -390,9 +387,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (BaseLib::getFileExtension(input_arg.getValue()) != ".pvtu")
     {
@@ -502,9 +497,6 @@ int main(int argc, char* argv[])
     if (!result)
     {
         ERR("Could not write mesh to '{:s}'.", output_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("writing mesh took {} s", writing_timer.elapsed());
@@ -516,8 +508,5 @@ int main(int argc, char* argv[])
     // cleaned.
     merged_mesh.shallowClean();
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/RemoveGhostData.cpp b/Applications/Utils/MeshEdit/RemoveGhostData.cpp
index 1824b8661c7e0d14d3ad1092bd15d75d210302db..1384ac736b2a08ad1a20bf397c06a00e52b8b978 100644
--- a/Applications/Utils/MeshEdit/RemoveGhostData.cpp
+++ b/Applications/Utils/MeshEdit/RemoveGhostData.cpp
@@ -15,11 +15,8 @@
 #include <vtkXMLPUnstructuredGridReader.h>
 #include <vtkXMLUnstructuredGridWriter.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/Logging.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/writeMeshToFile.h"
 
@@ -51,9 +48,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     vtkSmartPointer<vtkXMLPUnstructuredGridReader> reader =
         vtkSmartPointer<vtkXMLPUnstructuredGridReader>::New();
@@ -73,8 +68,5 @@ int main(int argc, char* argv[])
     writer->SetFileName(output_arg.getValue().c_str());
     writer->Write();
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/ReorderMesh.cpp b/Applications/Utils/MeshEdit/ReorderMesh.cpp
index fd5334d5ee156f17d14a18cf75f4dba6d849c3cc..5de0c9879c880ae1ac40eb9be2f1ce603e7233fe 100644
--- a/Applications/Utils/MeshEdit/ReorderMesh.cpp
+++ b/Applications/Utils/MeshEdit/ReorderMesh.cpp
@@ -12,10 +12,7 @@
 #include <memory>
 #include <string>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -136,9 +133,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     const std::string filename(mesh_in_arg.getValue());
 
@@ -147,9 +142,6 @@ int main(int argc, char* argv[])
         std::unique_ptr<MeshLib::Mesh>(MeshLib::IO::readMeshFromFile(filename));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -247,8 +239,5 @@ int main(int argc, char* argv[])
             "transferred to the reordered mesh.");
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/ResetPropertiesInPolygonalRegion.cpp b/Applications/Utils/MeshEdit/ResetPropertiesInPolygonalRegion.cpp
index 76631d3ea973f19b88d0cf9de818b6bed11a343f..ada2f29b6dd0f78adf46abfac93e1dc0d0fc7ed0 100644
--- a/Applications/Utils/MeshEdit/ResetPropertiesInPolygonalRegion.cpp
+++ b/Applications/Utils/MeshEdit/ResetPropertiesInPolygonalRegion.cpp
@@ -16,11 +16,8 @@
 #include <cstdlib>
 #include <vector>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/readGeometryFromFile.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/Polygon.h"
 #include "InfoLib/GitInfo.h"
@@ -99,9 +96,7 @@ int main(int argc, char* argv[])
     cmd.add(gmsh_path_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // *** read geometry
     GeoLib::GEOObjects geometries;
@@ -117,9 +112,6 @@ int main(int argc, char* argv[])
     {
         ERR("Could not get vector of polylines out of geometry '{:s}'.",
             geo_name);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -129,9 +121,6 @@ int main(int argc, char* argv[])
     if (!ply)
     {
         ERR("Polyline '{:s}' not found.", polygon_name_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -140,9 +129,6 @@ int main(int argc, char* argv[])
     {
         ERR("Polyline '{:s}' is not closed, i.e. does not describe a region.",
             polygon_name_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -154,9 +140,6 @@ int main(int argc, char* argv[])
     if (!mesh)
     {
         // error message written already by readMeshFromFile
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     std::string const& property_name(property_name_arg.getValue());
@@ -186,8 +169,5 @@ int main(int argc, char* argv[])
 
     MeshLib::IO::writeMeshToFile(*mesh, mesh_out.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/Vtu2Grid.cpp b/Applications/Utils/MeshEdit/Vtu2Grid.cpp
index 7429e1d9f72e7f45db1944012b1f0dc217891ec1..55a8bc6f5821b1e41b2c64c20b7deba7a3b17ca3 100644
--- a/Applications/Utils/MeshEdit/Vtu2Grid.cpp
+++ b/Applications/Utils/MeshEdit/Vtu2Grid.cpp
@@ -10,6 +10,7 @@
 #include <tclap/CmdLine.h>
 #include <vtkXMLUnstructuredGridReader.h>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/writeMeshToFile.h"
 #include "MeshLib/Mesh.h"
@@ -18,9 +19,6 @@
 #include "MeshToolsLib/MeshGenerators/MeshGenerator.h"
 #include "MeshToolsLib/MeshGenerators/VoxelGridFromMesh.h"
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
 int main(int argc, char* argv[])
 {
     TCLAP::CmdLine cmd(
@@ -64,18 +62,13 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if ((y_arg.isSet() && !z_arg.isSet()) ||
         ((!y_arg.isSet() && z_arg.isSet())))
     {
         ERR("For equilateral cubes, only x needs to be set. For unequal "
             "cuboids, all three edge lengths (x/y/z) need to be specified.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
     using namespace MeshToolsLib::MeshGenerator;
@@ -88,9 +81,6 @@ int main(int argc, char* argv[])
     {
         ERR("A cellsize ({},{},{}) is not allowed to be <= 0", x_size, y_size,
             z_size);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
@@ -113,9 +103,6 @@ int main(int argc, char* argv[])
     {
         ERR("The range ({},{},{}) is not allowed to be < 0", ranges[0],
             ranges[1], ranges[2]);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
     std::array<std::size_t, 3> const dims =
@@ -134,9 +121,6 @@ int main(int argc, char* argv[])
 
     if (!VoxelGridFromMesh::removeUnusedGridCells(mesh, grid))
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -144,14 +128,8 @@ int main(int argc, char* argv[])
 
     if (MeshLib::IO::writeMeshToFile(*grid, output_arg.getValue()) != 0)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/appendLinesAlongPolyline.cpp b/Applications/Utils/MeshEdit/appendLinesAlongPolyline.cpp
index e7ab944bd296bf74ca5e078b603b0ced23f90f28..69e3baa152283e077226a0e0912bb58b7b7fe714 100644
--- a/Applications/Utils/MeshEdit/appendLinesAlongPolyline.cpp
+++ b/Applications/Utils/MeshEdit/appendLinesAlongPolyline.cpp
@@ -11,12 +11,9 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/readGeometryFromFile.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/PolylineVec.h"
 #include "InfoLib/GitInfo.h"
@@ -58,9 +55,7 @@ int main(int argc, char* argv[])
     // parse arguments
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read GEO objects
     GeoLib::GEOObjects geo_objs;
@@ -71,9 +66,6 @@ int main(int argc, char* argv[])
     if (geo_names.empty())
     {
         ERR("No geometries found.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     const GeoLib::PolylineVec* ply_vec(
@@ -81,9 +73,6 @@ int main(int argc, char* argv[])
     if (!ply_vec)
     {
         ERR("Could not find polylines in geometry '{:s}'.", geo_names.front());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -93,9 +82,6 @@ int main(int argc, char* argv[])
     if (!mesh)
     {
         ERR("Mesh file '{:s}' not found", mesh_in.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("Mesh read: {:d} nodes, {:d} elements.", mesh->getNumberOfNodes(),
@@ -109,8 +95,5 @@ int main(int argc, char* argv[])
 
     MeshLib::IO::writeMeshToFile(*new_mesh, mesh_out.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/checkMesh.cpp b/Applications/Utils/MeshEdit/checkMesh.cpp
index 4795937702cf6b75c7a8bc3db00bce38059f22a6..2d31816fd60fb5270ae2dd56f01827c5cb4d66d1 100644
--- a/Applications/Utils/MeshEdit/checkMesh.cpp
+++ b/Applications/Utils/MeshEdit/checkMesh.cpp
@@ -7,18 +7,14 @@
  *              http://www.opengeosys.org/project/license
  */
 
-#include <tclap/CmdLine.h>
-
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <spdlog/fmt/bundled/ranges.h>
+#include <tclap/CmdLine.h>
 
 #include <array>
 #include <string>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/MemWatch.h"
 #include "BaseLib/RunTime.h"
 #include "BaseLib/StringTools.h"
@@ -52,9 +48,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read the mesh file
     BaseLib::MemWatch mem_watch;
@@ -65,9 +59,6 @@ int main(int argc, char* argv[])
         mesh_arg.getValue(), true /* compute_element_neighbors */));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -106,8 +97,5 @@ int main(int argc, char* argv[])
         // Remark: MeshValidation can modify the original mesh
         MeshToolsLib::MeshInformation::writeMeshValidationResults(*mesh);
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/convertToLinearMesh.cpp b/Applications/Utils/MeshEdit/convertToLinearMesh.cpp
index 2eeea03f024b7bb8b6cde9bb0179ebe0203291c3..61569eaaba66bfe1250318530b3ac76511ed0c0f 100644
--- a/Applications/Utils/MeshEdit/convertToLinearMesh.cpp
+++ b/Applications/Utils/MeshEdit/convertToLinearMesh.cpp
@@ -11,13 +11,10 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 #include <string>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -42,25 +39,17 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh> mesh(
         MeshLib::IO::readMeshFromFile(input_arg.getValue()));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (!mesh->hasNonlinearElement())
     {
         ERR("The input mesh is linear. Exit.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -71,8 +60,5 @@ int main(int argc, char* argv[])
     INFO("Save the new mesh into a file");
     MeshLib::IO::writeMeshToFile(*new_mesh, output_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/createLayeredMeshFromRasters.cpp b/Applications/Utils/MeshEdit/createLayeredMeshFromRasters.cpp
index 91ab223b4bc1145f9cb1d90d8f828b7a48a9d860..56805e54f9470572249e6783ae31f8a7c213d731 100644
--- a/Applications/Utils/MeshEdit/createLayeredMeshFromRasters.cpp
+++ b/Applications/Utils/MeshEdit/createLayeredMeshFromRasters.cpp
@@ -18,12 +18,9 @@
 #include <string>
 #include <vector>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/FileTools.h"
 #include "BaseLib/IO/readStringListFromFile.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
@@ -79,9 +76,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (min_thickness_arg.isSet())
     {
@@ -89,9 +84,6 @@ int main(int argc, char* argv[])
         if (min_thickness < 0)
         {
             ERR("Minimum layer thickness must be non-negative value.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -102,17 +94,11 @@ int main(int argc, char* argv[])
     if (!sfc_mesh)
     {
         ERR("Error reading mesh '{:s}'.", mesh_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (sfc_mesh->getDimension() != 2)
     {
         ERR("Input mesh must be a 2D mesh.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("done.");
@@ -122,9 +108,6 @@ int main(int argc, char* argv[])
     if (raster_paths.size() < 2)
     {
         ERR("At least two raster files needed to create 3D mesh.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     std::reverse(raster_paths.begin(), raster_paths.end());
@@ -134,18 +117,12 @@ int main(int argc, char* argv[])
     {
         if (!mapper.createLayers(*sfc_mesh, *rasters, min_thickness))
         {
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
     else
     {
         ERR("Reading raster files.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -160,9 +137,6 @@ int main(int argc, char* argv[])
     if (result_mesh == nullptr)
     {
         ERR("Mapper returned empty result for 'SubsurfaceMesh'.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -172,8 +146,5 @@ int main(int argc, char* argv[])
     MeshLib::IO::writeVtu(*result_mesh, output_name, data_mode);
     INFO("done.");
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/createQuadraticMesh.cpp b/Applications/Utils/MeshEdit/createQuadraticMesh.cpp
index 735c37537c458dd400be9c7355755e55e5fbf9b5..6758cf37146db4e312c87ef78190a21000d745d2 100644
--- a/Applications/Utils/MeshEdit/createQuadraticMesh.cpp
+++ b/Applications/Utils/MeshEdit/createQuadraticMesh.cpp
@@ -9,13 +9,10 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 #include <string>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -44,17 +41,12 @@ int main(int argc, char* argv[])
     cmd.add(add_centre_node_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh> mesh(
         MeshLib::IO::readMeshFromFile(input_arg.getValue()));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -65,8 +57,5 @@ int main(int argc, char* argv[])
     INFO("Save the new mesh into a file");
     MeshLib::IO::writeMeshToFile(*new_mesh, output_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/createTetgenSmeshFromRasters.cpp b/Applications/Utils/MeshEdit/createTetgenSmeshFromRasters.cpp
index 838eb9bde55d8c5fc218a334ea36f0550136e614..4b2a881ff699154c1c2e7cf62e54c5ae57549c25 100644
--- a/Applications/Utils/MeshEdit/createTetgenSmeshFromRasters.cpp
+++ b/Applications/Utils/MeshEdit/createTetgenSmeshFromRasters.cpp
@@ -18,13 +18,10 @@
 #include <string>
 #include <vector>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/TetGenInterface.h"
 #include "BaseLib/FileTools.h"
 #include "BaseLib/IO/readStringListFromFile.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
@@ -85,9 +82,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (min_thickness_arg.isSet())
     {
@@ -95,9 +90,6 @@ int main(int argc, char* argv[])
         if (min_thickness < 0)
         {
             ERR("Minimum layer thickness must be non-negative value.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -108,17 +100,11 @@ int main(int argc, char* argv[])
     if (!sfc_mesh)
     {
         ERR("Error reading mesh '{:s}'.", mesh_arg.getValue());
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (sfc_mesh->getDimension() != 2)
     {
         ERR("Input mesh must be a 2D mesh.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("done.");
@@ -128,9 +114,6 @@ int main(int argc, char* argv[])
     if (raster_paths.size() < 2)
     {
         ERR("At least two raster files needed to create 3D mesh.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     std::reverse(raster_paths.begin(), raster_paths.end());
@@ -148,18 +131,12 @@ int main(int argc, char* argv[])
         if (!lv.createLayers(*sfc_mesh, *rasters, min_thickness))
         {
             ERR("Creating the layers was erroneous.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
     else
     {
         ERR("The raster files are not accessible.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     std::unique_ptr<MeshLib::Mesh> tg_mesh =
@@ -174,14 +151,8 @@ int main(int argc, char* argv[])
     else
     {
         ERR("The tetgen-smesh could not be created.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/editMaterialID.cpp b/Applications/Utils/MeshEdit/editMaterialID.cpp
index 179997796907f5644b346a5f4ae24b3c69e9ff75..61d0892de4625b0e7a1319510e2716f47d1705e0 100644
--- a/Applications/Utils/MeshEdit/editMaterialID.cpp
+++ b/Applications/Utils/MeshEdit/editMaterialID.cpp
@@ -9,12 +9,9 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -68,24 +65,16 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (!replaceArg.isSet() && !condenseArg.isSet() && !specifyArg.isSet())
     {
         INFO("Please select editing mode: -r or -c or -s");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 0;
     }
     if (replaceArg.isSet() && condenseArg.isSet())
     {
         INFO("Please select only one editing mode: -r or -c or -s");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 0;
     }
     if (replaceArg.isSet())
@@ -95,9 +84,6 @@ int main(int argc, char* argv[])
             INFO(
                 "current and new material IDs must be provided for "
                 "replacement");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return 0;
         }
     }
@@ -108,9 +94,6 @@ int main(int argc, char* argv[])
             INFO(
                 "element type and new material IDs must be provided to specify "
                 "elements");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return 0;
         }
     }
@@ -170,8 +153,5 @@ int main(int argc, char* argv[])
     }
     MeshLib::IO::writeMeshToFile(*mesh, mesh_out.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/ipDataToPointCloud.cpp b/Applications/Utils/MeshEdit/ipDataToPointCloud.cpp
index a8a2600b640b2ee262fef44e516cbf69627b1f4b..d1fda1f1f9d93b4c06492c3ccfef13a36fe45f86 100644
--- a/Applications/Utils/MeshEdit/ipDataToPointCloud.cpp
+++ b/Applications/Utils/MeshEdit/ipDataToPointCloud.cpp
@@ -11,10 +11,7 @@
 
 #include <unordered_map>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -255,9 +252,7 @@ int main(int argc, char** argv)
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh const> mesh_in(
         MeshLib::IO::readMeshFromFile(arg_in_file.getValue()));
@@ -273,8 +268,5 @@ int main(int argc, char** argv)
 
     MeshLib::IO::writeMeshToFile(point_cloud, arg_out_file.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/queryMesh.cpp b/Applications/Utils/MeshEdit/queryMesh.cpp
index fd713e8fe74b52b847d04e6e968859c0cd3adb8d..c4a3cb9d3e04517fe983acbd366d98f68aac0441 100644
--- a/Applications/Utils/MeshEdit/queryMesh.cpp
+++ b/Applications/Utils/MeshEdit/queryMesh.cpp
@@ -9,16 +9,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <array>
 #include <memory>
 #include <sstream>
 #include <string>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
@@ -52,9 +49,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     const std::string filename(mesh_arg.getValue());
 
@@ -64,9 +59,6 @@ int main(int argc, char* argv[])
             filename, true /* compute_element_neighbors */));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -155,8 +147,5 @@ int main(int argc, char* argv[])
         out << std::endl;
         INFO("{:s}", out.str());
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/removeMeshElements.cpp b/Applications/Utils/MeshEdit/removeMeshElements.cpp
index 456fe2635f19fa2855c9a339e242044514c81591..bd20620b226dfff18a94a8fe2bcff5c0f1f92dec 100644
--- a/Applications/Utils/MeshEdit/removeMeshElements.cpp
+++ b/Applications/Utils/MeshEdit/removeMeshElements.cpp
@@ -13,12 +13,9 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -182,17 +179,12 @@ int main(int argc, char* argv[])
     cmd.add(mesh_in);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh const> mesh(
         MeshLib::IO::readMeshFromFile(mesh_in.getValue()));
     if (mesh == nullptr)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -229,9 +221,6 @@ int main(int argc, char* argv[])
             !property_name_arg.isSet())
         {
             ERR("Specify a property name for the value/range selected.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
 
@@ -241,9 +230,6 @@ int main(int argc, char* argv[])
         {
             ERR("Specify a value or range ('-min-value' and '-max_value') for "
                 "the property selected.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
 
@@ -264,9 +250,6 @@ int main(int argc, char* argv[])
             {
                 ERR("Specify if the inside or the outside of the selected "
                     "range should be removed.");
-#ifdef USE_PETSC
-                MPI_Finalize();
-#endif
                 return EXIT_FAILURE;
             }
 
@@ -299,9 +282,6 @@ int main(int argc, char* argv[])
         }
         if (aabb_error)
         {
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
 
@@ -324,17 +304,11 @@ int main(int argc, char* argv[])
 
     if (new_mesh == nullptr)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
     // write into a file
     MeshLib::IO::writeMeshToFile(*new_mesh, mesh_out.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/reviseMesh.cpp b/Applications/Utils/MeshEdit/reviseMesh.cpp
index ff9c0d3509ac84a03820ab764f27a845a6e23f3b..ca06a00b084831892f840d6cf2d5e34a29be6b59 100644
--- a/Applications/Utils/MeshEdit/reviseMesh.cpp
+++ b/Applications/Utils/MeshEdit/reviseMesh.cpp
@@ -9,15 +9,12 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <array>
 #include <memory>
 #include <string>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
@@ -59,18 +56,13 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read a mesh file
     std::unique_ptr<MeshLib::Mesh> org_mesh(
         MeshLib::IO::readMeshFromFile(input_arg.getValue()));
     if (!org_mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     INFO("Mesh read: {:d} nodes, {:d} elements.", org_mesh->getNumberOfNodes(),
@@ -92,8 +84,5 @@ int main(int argc, char* argv[])
         MeshLib::IO::writeMeshToFile(*new_mesh, output_arg.getValue());
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshEdit/swapNodeCoordinateAxes.cpp b/Applications/Utils/MeshEdit/swapNodeCoordinateAxes.cpp
index dcc9558a7f1269bd8740919ad3916a803b837d48..0eaadc7013e403b19020184ca3bb47ed82f53d3e 100644
--- a/Applications/Utils/MeshEdit/swapNodeCoordinateAxes.cpp
+++ b/Applications/Utils/MeshEdit/swapNodeCoordinateAxes.cpp
@@ -9,15 +9,12 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <array>
 #include <memory>
 #include <string>
 
 #include "BaseLib/Logging.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -114,17 +111,12 @@ int main(int argc, char* argv[])
     cmd.add(new_order_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     const std::string str_order = new_order_arg.getValue();
     std::array<int, 3> new_order = {{}};
     if (!parseNewOrder(str_order, new_order))
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -132,9 +124,6 @@ int main(int argc, char* argv[])
         MeshLib::IO::readMeshFromFile(input_arg.getValue()));
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -152,8 +141,5 @@ int main(int argc, char* argv[])
     INFO("Save the new mesh into a file");
     MeshLib::IO::writeMeshToFile(*mesh, output_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/AssignRasterDataToMesh.cpp b/Applications/Utils/MeshGeoTools/AssignRasterDataToMesh.cpp
index ceada5695ec3773702a1fdbe69c68812ce0a7c68..338354d0d6eb7bc62e438e33425fa29c228c6944 100644
--- a/Applications/Utils/MeshGeoTools/AssignRasterDataToMesh.cpp
+++ b/Applications/Utils/MeshGeoTools/AssignRasterDataToMesh.cpp
@@ -13,10 +13,7 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
 #include "InfoLib/GitInfo.h"
@@ -74,9 +71,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     bool const create_cell_array(set_cells_arg.isSet());
     bool const create_node_array =
@@ -91,9 +86,6 @@ int main(int argc, char* argv[])
     if (mesh->getDimension() > 2)
     {
         ERR("Method can not be applied to 3D meshes.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -108,9 +100,6 @@ int main(int argc, char* argv[])
         if (!assigned)
         {
             ERR("Error assigning raster data to scalar node array");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         INFO("Created node array {:s}", array_name_arg.getValue());
@@ -124,9 +113,6 @@ int main(int argc, char* argv[])
         if (!assigned)
         {
             ERR("Error assigning raster data to scalar cell array");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         INFO("Created cell array {:s}", array_name_arg.getValue());
@@ -134,8 +120,5 @@ int main(int argc, char* argv[])
 
     MeshLib::IO::VtuInterface vtu(mesh.get());
     vtu.writeToFile(output_name);
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/IntegrateBoreholesIntoMesh.cpp b/Applications/Utils/MeshGeoTools/IntegrateBoreholesIntoMesh.cpp
index 6f0f5b701cd13be145048b1ac4f9b3b3e8b1ecd1..9c88e9215889506dd96bf3d62253eeac9552064c 100644
--- a/Applications/Utils/MeshGeoTools/IntegrateBoreholesIntoMesh.cpp
+++ b/Applications/Utils/MeshGeoTools/IntegrateBoreholesIntoMesh.cpp
@@ -16,10 +16,7 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
 #include "InfoLib/GitInfo.h"
@@ -122,9 +119,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::pair<int, int> mat_limits(0, std::numeric_limits<int>::max());
     std::pair<double, double> elevation_limits(
@@ -134,9 +129,6 @@ int main(int argc, char* argv[])
     {
         ERR("If minimum MaterialID is set, maximum ID must be set, too (and "
             "vice versa).");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (min_id_arg.isSet() && max_id_arg.isSet())
@@ -151,18 +143,12 @@ int main(int argc, char* argv[])
     if (min_id_arg.isSet() && (mat_limits.first < 0 || mat_limits.second < 0))
     {
         ERR("Specified MaterialIDs must have non-negative values.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (min_elevation_arg.isSet() != max_elevation_arg.isSet())
     {
         ERR("If minimum elevation is set, maximum elevation must be set, too "
             "(and vice versa).");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (min_elevation_arg.isSet() && max_elevation_arg.isSet())
@@ -184,9 +170,6 @@ int main(int argc, char* argv[])
     if (!xml_io.readFile(geo_name))
     {
         ERR("Failed to read geometry file `{:s}'.", geo_name);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     std::vector<GeoLib::Point*> const& points =
@@ -197,17 +180,11 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Failed to read input mesh file `{:s}'.", mesh_name);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     if (mesh->getDimension() != 3)
     {
         ERR("Method can only be applied to 3D meshes.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -216,9 +193,6 @@ int main(int argc, char* argv[])
     if (mat_ids == nullptr)
     {
         ERR("Mesh is required to have MaterialIDs");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -256,8 +230,5 @@ int main(int argc, char* argv[])
                                true /* compute_element_neighbors */, props);
     MeshLib::IO::VtuInterface vtu(&result);
     vtu.writeToFile(output_name);
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/Raster2Mesh.cpp b/Applications/Utils/MeshGeoTools/Raster2Mesh.cpp
index c4905205da3b5f9de2be99d36479b576b7899944..423ecc4d961e5818d6e6126e16023bd8407362e0 100644
--- a/Applications/Utils/MeshGeoTools/Raster2Mesh.cpp
+++ b/Applications/Utils/MeshGeoTools/Raster2Mesh.cpp
@@ -13,10 +13,7 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
 #include "InfoLib/GitInfo.h"
@@ -72,9 +69,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string const input_name = input_arg.getValue().c_str();
     std::string const output_name = output_arg.getValue().c_str();
@@ -107,16 +102,10 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Conversion failed.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
     MeshLib::IO::VtuInterface vtu(mesh.get());
     vtu.writeToFile(output_name);
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/VerticalSliceFromLayers.cpp b/Applications/Utils/MeshGeoTools/VerticalSliceFromLayers.cpp
index 671cdaeeae9148c3db68de145dd1817e58845774..f9193a588cdfb36418a70beaecefd21fa5bb8183 100644
--- a/Applications/Utils/MeshGeoTools/VerticalSliceFromLayers.cpp
+++ b/Applications/Utils/MeshGeoTools/VerticalSliceFromLayers.cpp
@@ -16,16 +16,13 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <QCoreApplication>
 
 #include "Applications/FileIO/Gmsh/GMSHInterface.h"
 #include "Applications/FileIO/Gmsh/GmshReader.h"
 #include "BaseLib/FileTools.h"
 #include "BaseLib/IO/readStringListFromFile.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/AABB.h"
 #include "GeoLib/AnalyticalGeometry.h"
 #include "GeoLib/GEOObjects.h"
@@ -418,9 +415,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string const input_name = input_arg.getValue();
     std::string const output_name = output_arg.getValue();
@@ -439,9 +434,6 @@ int main(int argc, char* argv[])
     if (layer_names.size() < 2)
     {
         ERR("At least two layers are required to extract a slice.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -453,9 +445,6 @@ int main(int argc, char* argv[])
     {
         ERR("Less than two geometries could be created from layers. Aborting "
             "extraction...");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -475,9 +464,6 @@ int main(int argc, char* argv[])
     if (mesh == nullptr)
     {
         ERR("Error generating mesh... (GMSH was unable to output mesh)");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     rotateMesh(*mesh, rot_mat, z_shift);
@@ -485,9 +471,6 @@ int main(int argc, char* argv[])
     if (new_mesh == nullptr)
     {
         ERR("Error generating mesh... (GMSH created line mesh)");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -506,8 +489,5 @@ int main(int argc, char* argv[])
     }
     MeshLib::IO::VtuInterface vtu(revised_mesh.get());
     vtu.writeToFile(output_name + ".vtu");
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp b/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp
index 5329c0b9ad9b2bc14bbca3684241043d1ef4ee36..4ca7478f6d46974686f59d28dddc2a979503c829 100644
--- a/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp
+++ b/Applications/Utils/MeshGeoTools/computeSurfaceNodeIDsInPolygonalRegion.cpp
@@ -12,10 +12,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <fstream>
 #include <memory>
@@ -25,6 +21,7 @@
 #include "Applications/FileIO/readGeometryFromFile.h"
 #include "BaseLib/Error.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/Polygon.h"
 #include "InfoLib/GitInfo.h"
@@ -101,9 +98,7 @@ int main(int argc, char* argv[])
     cmd.add(gmsh_path_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh const> mesh(
         MeshLib::IO::readMeshFromFile(mesh_in.getValue()));
@@ -194,8 +189,5 @@ int main(int argc, char* argv[])
     std::for_each(all_sfc_nodes.begin(), all_sfc_nodes.end(),
                   std::default_delete<MeshLib::Node>());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/constructMeshesFromGeometry.cpp b/Applications/Utils/MeshGeoTools/constructMeshesFromGeometry.cpp
index d57974e0476a652f085eda4f5b27e36141c85477..9d99fc3e19d019c7bc758e9f668b57d0a5b3c33e 100644
--- a/Applications/Utils/MeshGeoTools/constructMeshesFromGeometry.cpp
+++ b/Applications/Utils/MeshGeoTools/constructMeshesFromGeometry.cpp
@@ -10,6 +10,7 @@
 
 #include <tclap/CmdLine.h>
 
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
 #include "InfoLib/GitInfo.h"
@@ -19,9 +20,6 @@
 #include "MeshLib/IO/writeMeshToFile.h"
 #include "MeshLib/Mesh.h"
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
 std::unique_ptr<GeoLib::GEOObjects> readGeometry(std::string const& filename)
 {
     auto geo_objects = std::make_unique<GeoLib::GEOObjects>();
@@ -80,9 +78,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     std::unique_ptr<MeshLib::Mesh> mesh{
         MeshLib::IO::readMeshFromFile(mesh_arg.getValue())};
 
@@ -101,9 +97,6 @@ int main(int argc, char* argv[])
         if (!m_ptr)
         {
             ERR("Could not create a mesh for each given geometry.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         if (m_ptr->getNodes().empty())
@@ -117,8 +110,5 @@ int main(int argc, char* argv[])
         MeshLib::IO::writeMeshToFile(*m_ptr, m_ptr->getName() + ".vtu");
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/createIntermediateRasters.cpp b/Applications/Utils/MeshGeoTools/createIntermediateRasters.cpp
index 2a98a25c57d1f6e2231eeba2d9dacf04b82044eb..75a7572a9641c3e81362b0b49899816277e9c509 100644
--- a/Applications/Utils/MeshGeoTools/createIntermediateRasters.cpp
+++ b/Applications/Utils/MeshGeoTools/createIntermediateRasters.cpp
@@ -9,14 +9,11 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <memory>
 #include <string>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
 #include "InfoLib/GitInfo.h"
@@ -50,9 +47,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<GeoLib::Raster> dem1(
         FileIO::AsciiRasterInterface::readRaster(input1_arg.getValue()));
@@ -61,9 +56,6 @@ int main(int argc, char* argv[])
 
     if (dem1 == nullptr || dem2 == nullptr)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 1;
     }
 
@@ -99,9 +91,6 @@ int main(int argc, char* argv[])
 
     if (errors_found)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 2;
     }
 
@@ -120,9 +109,6 @@ int main(int argc, char* argv[])
         if (it2 == dem2->end())
         {
             ERR("Error: File 2 is shorter than File 1.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return 3;
         }
         if (*it1 == h1.no_data || *it2 == h2.no_data)
@@ -147,9 +133,6 @@ int main(int argc, char* argv[])
     if (it2 != dem2->end())
     {
         ERR("Error: File 1 is shorter than File 2.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return 3;
     }
 
@@ -164,8 +147,5 @@ int main(int argc, char* argv[])
             r, basename + std::to_string(i) + ext);
         INFO("Layer {:d} written.", i + 1);
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/geometryToGmshGeo.cpp b/Applications/Utils/MeshGeoTools/geometryToGmshGeo.cpp
index 87777c243a7af24430cd51a2b8a058458fb37ff0..87d735950a42b17ee9d2f8c377cdb1d586d0c717 100644
--- a/Applications/Utils/MeshGeoTools/geometryToGmshGeo.cpp
+++ b/Applications/Utils/MeshGeoTools/geometryToGmshGeo.cpp
@@ -12,11 +12,8 @@
 // ThirdParty
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/Gmsh/GMSHInterface.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
 #include "InfoLib/GitInfo.h"
@@ -70,9 +67,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     GeoLib::GEOObjects geo_objects;
     for (auto const& geometry_name : geo_input_arg.getValue())
@@ -82,9 +77,6 @@ int main(int argc, char* argv[])
         {
             if (!xml.readFile(geometry_name))
             {
-#ifdef USE_PETSC
-                MPI_Finalize();
-#endif
                 return EXIT_FAILURE;
             }
         }
@@ -92,9 +84,6 @@ int main(int argc, char* argv[])
         {
             ERR("Failed to read file '{:s}'.", geometry_name);
             ERR("{:s}", err.what());
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         INFO("Successfully read file '{:s}'.", geometry_name);
@@ -155,13 +144,7 @@ int main(int argc, char* argv[])
                 "for better analysis of the problem.");
         }
         ERR("An error has occurred - programme will be terminated.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/MeshGeoTools/identifySubdomains.cpp b/Applications/Utils/MeshGeoTools/identifySubdomains.cpp
index 46e44c801d5730e72c6e7bca081dec9d4fcc3254..fcfba24983f80302dcadcd15e2446d8e7fc444ab 100644
--- a/Applications/Utils/MeshGeoTools/identifySubdomains.cpp
+++ b/Applications/Utils/MeshGeoTools/identifySubdomains.cpp
@@ -10,10 +10,7 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshGeoToolsLib/IdentifySubdomainMesh.h"
@@ -98,9 +95,7 @@ int main(int argc, char* argv[])
     cmd.add(subdomain_meshes_filenames_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     //
     // The bulk mesh.
@@ -165,9 +160,6 @@ int main(int argc, char* argv[])
     }
     INFO("writing time: {:g} s", writing_time.elapsed());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     INFO("Entire run time: {:g} s", run_time.elapsed());
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/ModelPreparation/ComputeNodeAreasFromSurfaceMesh.cpp b/Applications/Utils/ModelPreparation/ComputeNodeAreasFromSurfaceMesh.cpp
index c29a8bc4bbf82d1f7c42efb6c041d08f63935835..9bad7b3970af3a0b43213c6f93d088e14e957e85 100644
--- a/Applications/Utils/ModelPreparation/ComputeNodeAreasFromSurfaceMesh.cpp
+++ b/Applications/Utils/ModelPreparation/ComputeNodeAreasFromSurfaceMesh.cpp
@@ -11,10 +11,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <fstream>
 #include <memory>
 #include <numeric>
@@ -23,6 +19,7 @@
 
 #include "BaseLib/Error.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/Mesh.h"
@@ -90,9 +87,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh> surface_mesh(
         MeshLib::IO::readMeshFromFile(mesh_in.getValue()));
@@ -112,9 +107,6 @@ int main(int argc, char* argv[])
         if (!orig_node_ids)
         {
             ERR("Fatal error: could not create property.");
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
         orig_node_ids->resize(surface_mesh->getNumberOfNodes());
@@ -149,8 +141,5 @@ int main(int argc, char* argv[])
     writeToFile(id_and_area_fname, csv_fname, ids_and_areas,
                 surface_mesh->getNodes());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/ModelPreparation/PartitionMesh/BinaryToPVTU.cpp b/Applications/Utils/ModelPreparation/PartitionMesh/BinaryToPVTU.cpp
index fd7b3baca9c08d252051cc6f872c28455a526646..81ce322f9ec6a42d72dc9d717a505db3dbf3981f 100644
--- a/Applications/Utils/ModelPreparation/PartitionMesh/BinaryToPVTU.cpp
+++ b/Applications/Utils/ModelPreparation/PartitionMesh/BinaryToPVTU.cpp
@@ -10,7 +10,6 @@
 
 */
 
-#include <mpi.h>
 #include <spdlog/spdlog.h>
 #include <tclap/CmdLine.h>
 #include <vtkMPIController.h>
@@ -18,6 +17,7 @@
 
 #include "BaseLib/CPUTime.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
@@ -72,9 +72,7 @@ int main(int argc, char* argv[])
             OGS_FATAL("spdlog logger error occurred.");
         });
 
-    // init MPI
-    MPI_Init(&argc, &argv);
-
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     // start the timer
     BaseLib::RunTime run_timer;
     run_timer.start();
@@ -111,8 +109,5 @@ int main(int argc, char* argv[])
 
     INFO("Total runtime: {:g} s.", run_timer.elapsed());
     INFO("Total CPU time: {:g} s.", CPU_timer.elapsed());
-
-    MPI_Finalize();
-
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/ModelPreparation/PartitionMesh/PartitionMesh.cpp b/Applications/Utils/ModelPreparation/PartitionMesh/PartitionMesh.cpp
index daa2aa43cd3ea0fe8e9de0e1c61a12f2119dc7ef..3e5316313857119e63917dd33361787ea94387a8 100644
--- a/Applications/Utils/ModelPreparation/PartitionMesh/PartitionMesh.cpp
+++ b/Applications/Utils/ModelPreparation/PartitionMesh/PartitionMesh.cpp
@@ -15,12 +15,9 @@
 #include <spdlog/spdlog.h>
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "BaseLib/CPUTime.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/RunTime.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -102,9 +99,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     BaseLib::setConsoleLogLevel(log_level_arg.getValue());
     spdlog::set_pattern("%^%l:%$ %v");
@@ -144,9 +139,6 @@ int main(int argc, char* argv[])
         INFO("Total runtime: {:g} s.", run_timer.elapsed());
         INFO("Total CPU time: {:g} s.", CPU_timer.elapsed());
 
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_SUCCESS;
     }
 
@@ -192,9 +184,6 @@ int main(int argc, char* argv[])
         {
             INFO("Failed in system calling.");
             INFO("Return value of system call {:d} ", status);
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -244,8 +233,5 @@ int main(int argc, char* argv[])
     INFO("Total runtime: {:g} s.", run_timer.elapsed());
     INFO("Total CPU time: {:g} s.", CPU_timer.elapsed());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/ModelPreparation/convertVtkDataArrayToVtkDataArray.cpp b/Applications/Utils/ModelPreparation/convertVtkDataArrayToVtkDataArray.cpp
index dbb156149c02168c9920bfd374837e8d0e198f85..0fd23a6fe3c1ced13434a4e8591c03562737bb49 100644
--- a/Applications/Utils/ModelPreparation/convertVtkDataArrayToVtkDataArray.cpp
+++ b/Applications/Utils/ModelPreparation/convertVtkDataArrayToVtkDataArray.cpp
@@ -11,15 +11,12 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <cmath>
 #include <memory>
 #include <numeric>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -105,18 +102,13 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh> mesh(
         MeshLib::IO::readMeshFromFile(mesh_arg.getValue()));
 
     if (!mesh)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
@@ -163,16 +155,10 @@ int main(int argc, char* argv[])
     if (!success)
     {
         ERR("{:s}", err_msg);
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
     MeshLib::IO::writeMeshToFile(*mesh, out_mesh_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/ModelPreparation/createNeumannBc.cpp b/Applications/Utils/ModelPreparation/createNeumannBc.cpp
index 7f108478b1f70e54b7fa99ae67114c4d8dea7e28..b69f69c48613b1cd25ffd64797e59e74c9feb4d0 100644
--- a/Applications/Utils/ModelPreparation/createNeumannBc.cpp
+++ b/Applications/Utils/ModelPreparation/createNeumannBc.cpp
@@ -10,12 +10,9 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <fstream>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/IO/readMeshFromFile.h"
@@ -130,9 +127,7 @@ int main(int argc, char* argv[])
     cmd.add(result_file);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read surface mesh
     std::unique_ptr<MeshLib::Mesh> surface_mesh(
@@ -155,9 +150,6 @@ int main(int argc, char* argv[])
     }();
     if (!node_id_pv)
     {
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
 
@@ -192,8 +184,5 @@ int main(int argc, char* argv[])
         result_out << p.first << " " << p.second << "\n";
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/ModelPreparation/scaleProperty.cpp b/Applications/Utils/ModelPreparation/scaleProperty.cpp
index 03eb59846b759d122c7ccdbef0efcd9c5edd0c7d..562ee1d9212647d2ec86dc2ffadab0f53fe88d51 100644
--- a/Applications/Utils/ModelPreparation/scaleProperty.cpp
+++ b/Applications/Utils/ModelPreparation/scaleProperty.cpp
@@ -11,15 +11,12 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <cmath>
 #include <memory>
 #include <numeric>
 
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -72,9 +69,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::unique_ptr<MeshLib::Mesh> mesh(
         MeshLib::IO::readMeshFromFile(mesh_arg.getValue()));
@@ -97,8 +92,5 @@ int main(int argc, char* argv[])
 
     MeshLib::IO::writeMeshToFile(*mesh, out_mesh_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/PostProcessing/Raster2PointCloud.cpp b/Applications/Utils/PostProcessing/Raster2PointCloud.cpp
index 2b22e4fde07524ff0cad6ca773068e6ea674cc42..44fa7de536fe4567393b802aeb96f573063aac21 100644
--- a/Applications/Utils/PostProcessing/Raster2PointCloud.cpp
+++ b/Applications/Utils/PostProcessing/Raster2PointCloud.cpp
@@ -13,9 +13,6 @@
 
 // ThirdParty
 #include <tclap/CmdLine.h>
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
 #include <vtkPolyDataAlgorithm.h>
 #include <vtkSmartPointer.h>
 #include <vtkXMLPolyDataWriter.h>
@@ -23,6 +20,7 @@
 #include "Applications/DataExplorer/VtkVis/VtkGeoImageSource.h"
 #include "Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "InfoLib/GitInfo.h"
 
@@ -80,9 +78,7 @@ int main(int argc, char* argv[])
     cmd.add(input_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     std::string const input_name = input_arg.getValue().c_str();
     std::string const output_name = output_arg.getValue().c_str();
@@ -161,8 +157,5 @@ int main(int argc, char* argv[])
         }
     }
     std::cout << "done." << std::endl;
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/PostProcessing/postLIE.cpp b/Applications/Utils/PostProcessing/postLIE.cpp
index 86199ac4055ea3f0b34b42d7a1b655d95fef5be4..dac283371f5843970749be49b18adba01ae05e4d 100644
--- a/Applications/Utils/PostProcessing/postLIE.cpp
+++ b/Applications/Utils/PostProcessing/postLIE.cpp
@@ -9,10 +9,6 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <boost/property_tree/ptree.hpp>
 #include <boost/property_tree/xml_parser.hpp>
 #include <map>
@@ -20,6 +16,7 @@
 #include <vector>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/IO/readMeshFromFile.h"
 #include "MeshLib/IO/writeMeshToFile.h"
@@ -151,9 +148,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     auto const in_file_ext = BaseLib::getFileExtension(arg_in_file.getValue());
     if (in_file_ext == ".pvd")
@@ -170,8 +165,5 @@ int main(int argc, char* argv[])
         OGS_FATAL("The given file type ({:s}) is not supported.", in_file_ext);
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/SWMMConverter/SWMMConverter.cpp b/Applications/Utils/SWMMConverter/SWMMConverter.cpp
index 7db44c2b4afa29462b28f4c25f050fc2018a8235..be95e7c61282d27ab77600e89affc6404c69bb4e 100644
--- a/Applications/Utils/SWMMConverter/SWMMConverter.cpp
+++ b/Applications/Utils/SWMMConverter/SWMMConverter.cpp
@@ -9,12 +9,9 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include "Applications/FileIO/SWMM/SWMMInterface.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/StringTools.h"
 #include "GeoLib/GEOObjects.h"
 #include "GeoLib/IO/XmlIO/Boost/BoostXmlGmlInterface.h"
@@ -199,18 +196,13 @@ int main(int argc, char* argv[])
     cmd.add(add_system_arg);
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     if (!(geo_output_arg.isSet() || mesh_output_arg.isSet() ||
           csv_output_arg.isSet()))
     {
         ERR("No output format given. Please specify OGS geometry or mesh "
             "output file.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
@@ -219,9 +211,6 @@ int main(int argc, char* argv[])
     {
         ERR("Please specify csv output file for exporting subcatchment or "
             "system parameters.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return -1;
     }
 
@@ -240,8 +229,5 @@ int main(int argc, char* argv[])
                        add_subcatchments_arg.getValue(),
                        add_system_arg.getValue());
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return 0;
 }
diff --git a/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp b/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp
index 9039a44a1349f43e1b3f1167ac065e3d2c16fe83..00b7dfedf5a1e2a10c2f356a28f14b79655450ff 100644
--- a/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp
+++ b/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp
@@ -12,16 +12,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <cmath>
 #include <memory>
 #include <numeric>
 
 #include "BaseLib/FileTools.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/quicksort.h"
 #include "GeoLib/IO/AsciiRasterInterface.h"
 #include "GeoLib/Raster.h"
@@ -94,9 +91,7 @@ int main(int argc, char* argv[])
 
     cmd.parse(argc, argv);
 
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
     // read mesh
     std::unique_ptr<MeshLib::Mesh> dest_mesh(
@@ -138,8 +133,5 @@ int main(int argc, char* argv[])
         MeshLib::IO::writeMeshToFile(*dest_mesh, out_mesh_arg.getValue());
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }
diff --git a/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp b/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp
index 7e3062698a32deaf00debee231fe8011de5fe172..0392559aae06ae71e09f583e7aedac88ec4d206e 100644
--- a/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp
+++ b/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp
@@ -9,16 +9,13 @@
 
 #include <tclap/CmdLine.h>
 
-#ifdef USE_PETSC
-#include <mpi.h>
-#endif
-
 #include <algorithm>
 #include <memory>
 #include <string>
 #include <vector>
 
 #include "BaseLib/Error.h"
+#include "BaseLib/MPI.h"
 #include "BaseLib/Subdivision.h"
 #include "BaseLib/TCLAPCustomOutput.h"
 #include "InfoLib/GitInfo.h"
@@ -161,9 +158,7 @@ int main(int argc, char* argv[])
 
     // parse arguments
     cmd.parse(argc, argv);
-#ifdef USE_PETSC
-    MPI_Init(&argc, &argv);
-#endif
+    BaseLib::MPI::Setup mpi_setup(argc, argv);
     const std::string eleTypeName(eleTypeArg.getValue());
     const MeshLib::MeshElemType eleType =
         MeshLib::String2MeshElemType(eleTypeName);
@@ -194,9 +189,6 @@ int main(int argc, char* argv[])
     if (!isLengthSet)
     {
         ERR("Missing input: Length information is not provided at all.");
-#ifdef USE_PETSC
-        MPI_Finalize();
-#endif
         return EXIT_FAILURE;
     }
     for (unsigned i = 0; i < 3; i++)
@@ -206,9 +198,6 @@ int main(int argc, char* argv[])
             ERR("Missing input: Length for dimension [{:d}] is required but "
                 "missing.",
                 i);
-#ifdef USE_PETSC
-            MPI_Finalize();
-#endif
             return EXIT_FAILURE;
         }
     }
@@ -304,8 +293,5 @@ int main(int argc, char* argv[])
             *(mesh.get()), std::filesystem::path(mesh_out.getValue()));
     }
 
-#ifdef USE_PETSC
-    MPI_Finalize();
-#endif
     return EXIT_SUCCESS;
 }