diff --git a/ProcessLib/CoupledSolutionsForStaggeredScheme.cpp b/ProcessLib/CoupledSolutionsForStaggeredScheme.cpp
index ad547737958d77525cbb657fc3f5f16be02287bb..71766d222cad031446174bc138ce0e8bc3829ace 100644
--- a/ProcessLib/CoupledSolutionsForStaggeredScheme.cpp
+++ b/ProcessLib/CoupledSolutionsForStaggeredScheme.cpp
@@ -36,10 +36,9 @@ std::vector<std::vector<double>> getPreviousLocalSolutions(
     std::vector<std::vector<double>> local_xs_t0;
     local_xs_t0.reserve(number_of_coupled_solutions);
 
-    for (std::size_t i = 0; i < number_of_coupled_solutions; i++)
+    for (auto const& x_t0 : cpl_xs.coupled_xs_t0)
     {
-        auto const& x_t0 = *cpl_xs.coupled_xs_t0[i];
-        local_xs_t0.emplace_back(x_t0.get(indices));
+        local_xs_t0.emplace_back(x_t0->get(indices));
     }
     return local_xs_t0;
 }
@@ -52,10 +51,9 @@ std::vector<std::vector<double>> getCurrentLocalSolutions(
     std::vector<std::vector<double>> local_xs_t1;
     local_xs_t1.reserve(number_of_coupled_solutions);
 
-    for (std::size_t i = 0; i < number_of_coupled_solutions; i++)
+    for (auto const& x_t1 : cpl_xs.coupled_xs)
     {
-        auto const& x_t1 = cpl_xs.coupled_xs[i].get();
-        local_xs_t1.emplace_back(x_t1.get(indices));
+        local_xs_t1.emplace_back(x_t1.get().get(indices));
     }
     return local_xs_t1;
 }
diff --git a/ProcessLib/HT/CreateHTProcess.cpp b/ProcessLib/HT/CreateHTProcess.cpp
index 9d5c31681741f378001a834ce61444d78b9e5dec..cd48556153cf9af9f78cb67c3ba4b36c2e9038b9 100644
--- a/ProcessLib/HT/CreateHTProcess.cpp
+++ b/ProcessLib/HT/CreateHTProcess.cpp
@@ -62,14 +62,11 @@ std::unique_ptr<Process> createHTProcess(
     }
     else  // staggered scheme.
     {
-        std::array<std::string, 2> variable_names = {
-            {"temperature",
-             "pressure"}};  // double-braces required in C++11 (not in C++14)
-
-        for (int i = 0; i < 2; i++)
+        using namespace std::string_literals;
+        for (auto const& variable_name : {"temperature"s, "pressure"s})
         {
             auto per_process_variables =
-                findProcessVariables(variables, pv_config, {variable_names[i]});
+                findProcessVariables(variables, pv_config, {variable_name});
             process_variables.push_back(std::move(per_process_variables));
         }
     }
diff --git a/ProcessLib/HT/HTFEM.h b/ProcessLib/HT/HTFEM.h
index 72888f29e4b22b8132f2c3a860c0d02378423d3a..626e71b4f1fd41a62d8dd64975cbe11f2513c219 100644
--- a/ProcessLib/HT/HTFEM.h
+++ b/ProcessLib/HT/HTFEM.h
@@ -104,11 +104,10 @@ protected:
                                     const double fluid_density,
                                     const double specific_heat_capacity_fluid)
     {
-        auto const& material_properties = this->_material_properties;
         auto const specific_heat_capacity_solid =
-            material_properties.specific_heat_capacity_solid(t, pos)[0];
+            _material_properties.specific_heat_capacity_solid(t, pos)[0];
 
-        auto const solid_density = material_properties.density_solid(t, pos)[0];
+        auto const solid_density = _material_properties.density_solid(t, pos)[0];
 
         return solid_density * specific_heat_capacity_solid * (1 - porosity) +
                fluid_density * specific_heat_capacity_fluid * porosity;
@@ -119,23 +118,21 @@ protected:
         const double fluid_density, const double specific_heat_capacity_fluid,
         const GlobalDimVectorType& velocity, const GlobalDimMatrixType& I)
     {
-        auto const& material_properties = this->_material_properties;
-
         auto const thermal_conductivity_solid =
-            material_properties.thermal_conductivity_solid(t, pos)[0];
+            _material_properties.thermal_conductivity_solid(t, pos)[0];
         auto const thermal_conductivity_fluid =
-            material_properties.thermal_conductivity_fluid(t, pos)[0];
+            _material_properties.thermal_conductivity_fluid(t, pos)[0];
         double const thermal_conductivity =
             thermal_conductivity_solid * (1 - porosity) +
             thermal_conductivity_fluid * porosity;
 
-        if (!material_properties.has_fluid_thermal_dispersivity)
+        if (!_material_properties.has_fluid_thermal_dispersivity)
             return thermal_conductivity * I;
 
         double const thermal_dispersivity_longitudinal =
-            material_properties.thermal_dispersivity_longitudinal(t, pos)[0];
+            _material_properties.thermal_dispersivity_longitudinal(t, pos)[0];
         auto const thermal_dispersivity_transversal =
-            material_properties.thermal_dispersivity_transversal(t, pos)[0];
+            _material_properties.thermal_dispersivity_transversal(t, pos)[0];
 
         double const velocity_magnitude = velocity.norm();
         GlobalDimMatrixType const thermal_dispersivity =
diff --git a/ProcessLib/HT/MonolithicHTFEM.h b/ProcessLib/HT/MonolithicHTFEM.h
index 775b56f1370c1078a4b6c17276b91b7a6b57f157..966d03ae82de31745d31728b8427d252e19c31e5 100644
--- a/ProcessLib/HT/MonolithicHTFEM.h
+++ b/ProcessLib/HT/MonolithicHTFEM.h
@@ -108,7 +108,7 @@ public:
         unsigned const n_integration_points =
             this->_integration_method.getNumberOfPoints();
 
-        for (std::size_t ip(0); ip < n_integration_points; ip++)
+        for (unsigned ip(0); ip < n_integration_points; ip++)
         {
             pos.setIntegrationPoint(ip);
 
diff --git a/ProcessLib/HT/StaggeredHTFEM-impl.h b/ProcessLib/HT/StaggeredHTFEM-impl.h
index a5add61c701ced6d0f97d51b55e9fb727d6d50f8..58313561689f331ad57d1f818fdcfe14fa8da324 100644
--- a/ProcessLib/HT/StaggeredHTFEM-impl.h
+++ b/ProcessLib/HT/StaggeredHTFEM-impl.h
@@ -76,7 +76,7 @@ void StaggeredHTFEM<ShapeFunction, IntegrationMethod, GlobalDim>::
     unsigned const n_integration_points =
         this->_integration_method.getNumberOfPoints();
 
-    for (std::size_t ip(0); ip < n_integration_points; ip++)
+    for (unsigned ip(0); ip < n_integration_points; ip++)
     {
         pos.setIntegrationPoint(ip);
 
@@ -202,7 +202,7 @@ void StaggeredHTFEM<ShapeFunction, IntegrationMethod, GlobalDim>::
     unsigned const n_integration_points =
         this->_integration_method.getNumberOfPoints();
 
-    for (std::size_t ip(0); ip < n_integration_points; ip++)
+    for (unsigned ip(0); ip < n_integration_points; ip++)
     {
         pos.setIntegrationPoint(ip);
 
diff --git a/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp b/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp
index a9beea3609ab1502fa6f485215397476f7d375af..46bffb117a7eb5847078fc370aef1dae63d84fbb 100644
--- a/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp
+++ b/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp
@@ -64,14 +64,11 @@ std::unique_ptr<Process> createHydroMechanicsProcess(
     }
     else  // staggered scheme.
     {
-        std::array<std::string, 2> variable_names = {
-            {"pressure",
-             "displacement"}};  // double-braces required in C++11 (not in C++14)
-
-        for (int i = 0; i < 2; i++)
+        using namespace std::string_literals;
+        for (auto const& variable_name : {"pressure"s, "displacement"s})
         {
             auto per_process_variables =
-                findProcessVariables(variables, pv_config, {variable_names[i]});
+                findProcessVariables(variables, pv_config, {variable_name});
             process_variables.push_back(std::move(per_process_variables));
         }
         variable_p = &process_variables[0][0].get();
diff --git a/ProcessLib/Output.cpp b/ProcessLib/Output.cpp
index 90b437556ecab8d52e247f910a2a29016b1d77e9..3efece6337f94aea583d518295633b8ff3743827 100644
--- a/ProcessLib/Output.cpp
+++ b/ProcessLib/Output.cpp
@@ -123,17 +123,19 @@ newInstance(const BaseLib::ConfigTree &config, std::string const& output_directo
     return out;
 }
 
-void Output::addProcess(ProcessLib::Process const& process, const unsigned pcs_idx)
+void Output::addProcess(ProcessLib::Process const& process,
+                        const int process_id)
 {
     auto const filename =
-        BaseLib::joinPaths(_output_directory, _output_file_prefix + "_pcs_" + std::to_string(pcs_idx) + ".pvd");
+        BaseLib::joinPaths(_output_directory, _output_file_prefix + "_pcs_" +
+                            std::to_string(process_id) + ".pvd");
     _single_process_data.emplace(std::piecewise_construct,
                                  std::forward_as_tuple(&process),
                                  std::forward_as_tuple(filename));
 }
 
 Output::SingleProcessData Output::findSingleProcessData(
-    Process const& process, const unsigned process_index) const
+    Process const& process, const unsigned process_id) const
 {
     if (process.isMonolithicSchemeUsed())
     {
@@ -149,7 +151,7 @@ Output::SingleProcessData Output::findSingleProcessData(
     unsigned counter = 0;
     for (auto spd_it=spd_range.first; spd_it!=spd_range.second; ++spd_it)
     {
-        if(counter == process_index)
+        if(counter == process_id)
         {
             return spd_it->second;
         }
@@ -162,7 +164,7 @@ Output::SingleProcessData Output::findSingleProcessData(
 
 
 void Output::doOutputAlways(Process const& process,
-                            const unsigned process_index,
+                            const int process_id,
                             ProcessOutput const& process_output,
                             unsigned timestep,
                             const double t,
@@ -171,16 +173,16 @@ void Output::doOutputAlways(Process const& process,
     BaseLib::RunTime time_output;
     time_output.start();
 
-    auto spd = findSingleProcessData(process, process_index);
+    auto spd = findSingleProcessData(process, process_id);
 
     std::string const output_file_name =
-            _output_file_prefix + "_pcs_" + std::to_string(process_index)
+            _output_file_prefix + "_pcs_" + std::to_string(process_id)
             + "_ts_" + std::to_string(timestep)
             + "_t_"  + std::to_string(t)
             + ".vtu";
     std::string const output_file_path = BaseLib::joinPaths(_output_directory, output_file_name);
 
-    const bool make_out = !(process_index < _single_process_data.size() - 1 &&
+    const bool make_out = !(process_id < _single_process_data.size() - 1 &&
                             !(process.isMonolithicSchemeUsed()));
 
     if (make_out)
@@ -201,7 +203,7 @@ void Output::doOutputAlways(Process const& process,
 }
 
 void Output::doOutput(Process const& process,
-                      const unsigned process_index,
+                      const int process_id,
                       ProcessOutput const& process_output,
                       unsigned timestep,
                       const double t,
@@ -209,7 +211,7 @@ void Output::doOutput(Process const& process,
 {
     if (shallDoOutput(timestep, _repeats_each_steps))
     {
-        doOutputAlways(process, process_index, process_output, timestep, t, x);
+        doOutputAlways(process, process_id, process_output, timestep, t, x);
     }
 #ifdef USE_INSITU
     // Note: last time step may be output twice: here and in
@@ -219,7 +221,7 @@ void Output::doOutput(Process const& process,
 }
 
 void Output::doOutputLastTimestep(Process const& process,
-                                  const unsigned process_index,
+                                  const int process_id,
                                   ProcessOutput const& process_output,
                                   unsigned timestep,
                                   const double t,
@@ -227,7 +229,7 @@ void Output::doOutputLastTimestep(Process const& process,
 {
     if (!shallDoOutput(timestep, _repeats_each_steps))
     {
-        doOutputAlways(process, process_index, process_output, timestep, t, x);
+        doOutputAlways(process, process_id, process_output, timestep, t, x);
     }
 #ifdef USE_INSITU
     InSituLib::CoProcess(process.getMesh(), t, timestep, true);
@@ -235,7 +237,7 @@ void Output::doOutputLastTimestep(Process const& process,
 }
 
 void Output::doOutputNonlinearIteration(Process const& process,
-                                        const unsigned process_index,
+                                        const int process_id,
                                         ProcessOutput const& process_output,
                                         const unsigned timestep, const double t,
                                         GlobalVector const& x,
@@ -249,17 +251,17 @@ void Output::doOutputNonlinearIteration(Process const& process,
     BaseLib::RunTime time_output;
     time_output.start();
 
-    findSingleProcessData(process, process_index);
+    findSingleProcessData(process, process_id);
 
     std::string const output_file_name =
-            _output_file_prefix + "_pcs_" + std::to_string(process_index)
+            _output_file_prefix + "_pcs_" + std::to_string(process_id)
             + "_ts_" + std::to_string(timestep)
             + "_t_"  + std::to_string(t)
             + "_nliter_" + std::to_string(iteration)
             + ".vtu";
     std::string const output_file_path = BaseLib::joinPaths(_output_directory, output_file_name);
     DBUG("output iteration results to %s", output_file_path.c_str());
-    const bool make_out = !(process_index < _single_process_data.size() - 1 &&
+    const bool make_out = !(process_id < _single_process_data.size() - 1 &&
                             !(process.isMonolithicSchemeUsed()));
     doProcessOutput(output_file_path, make_out, _output_file_compression,
                     _output_file_data_mode, t, x, process.getMesh(),
diff --git a/ProcessLib/Output.h b/ProcessLib/Output.h
index 44d054dba5defbb0750852931c45cf7c116ac79a..71dc38a656f60c3db17a687b806890fac3f7a705 100644
--- a/ProcessLib/Output.h
+++ b/ProcessLib/Output.h
@@ -33,11 +33,11 @@ public:
                 const std::string& output_directory);
 
     //! TODO doc. Opens a PVD file for each process.
-    void addProcess(ProcessLib::Process const& process, const unsigned pcs_idx);
+    void addProcess(ProcessLib::Process const& process, const int process_id);
 
     //! Writes output for the given \c process if it should be written in the
     //! given \c timestep.
-    void doOutput(Process const& process, const unsigned process_index,
+    void doOutput(Process const& process, const int process_id,
                   ProcessOutput const& process_output,
                   unsigned timestep, const double t, GlobalVector const& x);
 
@@ -45,7 +45,7 @@ public:
     //! This method is intended for doing output after the last timestep in order
     //! to make sure that its results are written.
     void doOutputLastTimestep(Process const& process,
-                              const unsigned process_index,
+                              const int process_id,
                               ProcessOutput const& process_output,
                               unsigned timestep, const double t,
                               GlobalVector const& x);
@@ -53,14 +53,14 @@ public:
     //! Writes output for the given \c process.
     //! This method will always write.
     //! It is intended to write output in error handling routines.
-    void doOutputAlways(Process const& process, const unsigned process_index,
+    void doOutputAlways(Process const& process, const int process_id,
                         ProcessOutput const& process_output, unsigned timestep,
                         const double t, GlobalVector const& x);
 
     //! Writes output for the given \c process.
     //! To be used for debug output after an iteration of the nonlinear solver.
     void doOutputNonlinearIteration(Process const& process,
-                                    const unsigned process_index,
+                                    const int process_id,
                                     ProcessOutput const& process_output,
                                     const unsigned timestep, const double t,
                                     GlobalVector const& x,
@@ -108,7 +108,7 @@ private:
     std::multimap<Process const*, SingleProcessData> _single_process_data;
 
     SingleProcessData findSingleProcessData(
-        Process const& process, const unsigned process_index) const;
+        Process const& process, const unsigned process_id) const;
 };
 
 }
diff --git a/ProcessLib/PhaseField/CreatePhaseFieldProcess.cpp b/ProcessLib/PhaseField/CreatePhaseFieldProcess.cpp
index 487b9a1bfd5d2e722f557bd6157b2133389bb262..b537383a5b35d7a4cc38b9dcdbd470e6bdf6a47a 100644
--- a/ProcessLib/PhaseField/CreatePhaseFieldProcess.cpp
+++ b/ProcessLib/PhaseField/CreatePhaseFieldProcess.cpp
@@ -65,14 +65,11 @@ std::unique_ptr<Process> createPhaseFieldProcess(
     }
     else  // staggered scheme.
     {
-        std::array<std::string, 2> variable_names = {
-            {"phasefield",
-             "displacement"}};  // double-braces required in C++11 (not in C++14)
-
-        for (int i = 0; i < 2; i++)
+        using namespace std::string_literals;
+        for (auto const& variable_name : {"phasefield"s, "displacement"s})
         {
             auto per_process_variables =
-                findProcessVariables(variables, pv_config, {variable_names[i]});
+                findProcessVariables(variables, pv_config, {variable_name});
             process_variables.push_back(std::move(per_process_variables));
         }
         variable_ph = &process_variables[0][0].get();
diff --git a/ProcessLib/Process.cpp b/ProcessLib/Process.cpp
index 265602869527332994501177f8cf4ac0e0d431bc..b1e8ed075ca86c66f7bcf8b16d375d0ab27b392e 100644
--- a/ProcessLib/Process.cpp
+++ b/ProcessLib/Process.cpp
@@ -37,7 +37,7 @@ Process::Process(
       _coupled_solutions(nullptr),
       _integration_order(integration_order),
       _process_variables(std::move(process_variables)),
-      _boundary_conditions([&](const size_t number_of_processes)
+      _boundary_conditions([&](const std::size_t number_of_processes)
                                -> std::vector<BoundaryConditionCollection> {
           std::vector<BoundaryConditionCollection> pcs_BCs;
           pcs_BCs.reserve(number_of_processes);
@@ -79,13 +79,10 @@ void Process::initialize()
                                                   _integration_order);
 
         std::vector<std::unique_ptr<NodalSourceTerm>> per_process_source_terms;
-        for (std::size_t variable_id = 0;
-             variable_id < per_process_variables.size();
-             variable_id++)
+        for (auto& pv : per_process_variables)
         {
-            ProcessVariable& pv = per_process_variables[variable_id];
-            auto sts = pv.createSourceTerms(*_local_to_global_index_map, 0,
-                                            _integration_order);
+            auto sts = pv.get().createSourceTerms(*_local_to_global_index_map,
+                                                  0, _integration_order);
 
             std::move(sts.begin(), sts.end(),
                       std::back_inserter(per_process_source_terms));
diff --git a/ProcessLib/RichardsComponentTransport/CreateRichardsComponentTransportProcess.cpp b/ProcessLib/RichardsComponentTransport/CreateRichardsComponentTransportProcess.cpp
index 377fa17e6fd9245c9f8bc1187269147730c3ee84..efb8122063c56db24060d0f19a31e88c98824029 100644
--- a/ProcessLib/RichardsComponentTransport/CreateRichardsComponentTransportProcess.cpp
+++ b/ProcessLib/RichardsComponentTransport/CreateRichardsComponentTransportProcess.cpp
@@ -62,14 +62,11 @@ std::unique_ptr<Process> createRichardsComponentTransportProcess(
     }
     else  // staggered scheme.
     {
-        std::array<std::string, 2> variable_names = {
-            {"concentration",
-             "pressure"}};  // double-braces required in C++11 (not in C++14)
-
-        for (int i = 0; i < 2; i++)
+        using namespace std::string_literals;
+        for (auto const& variable_name : {"concentration"s, "pressure"s})
         {
             auto per_process_variables =
-                findProcessVariables(variables, pv_config, {variable_names[i]});
+                findProcessVariables(variables, pv_config, {variable_name});
             process_variables.push_back(std::move(per_process_variables));
         }
     }
diff --git a/ProcessLib/ThermoMechanics/CreateThermoMechanicsProcess.cpp b/ProcessLib/ThermoMechanics/CreateThermoMechanicsProcess.cpp
index 470d7681702b5020b068d35b99662e679be68e2c..4019fe48473b574b1b40b6809efeddaad41c723e 100644
--- a/ProcessLib/ThermoMechanics/CreateThermoMechanicsProcess.cpp
+++ b/ProcessLib/ThermoMechanics/CreateThermoMechanicsProcess.cpp
@@ -65,14 +65,11 @@ std::unique_ptr<Process> createThermoMechanicsProcess(
     }
     else  // staggered scheme.
     {
-        std::array<std::string, 2> variable_names = {
-            {"temperature",
-             "displacement"}};  // double-braces required in C++11 (not in C++14)
-
-        for (int i = 0; i < 2; i++)
+        using namespace std::string_literals;
+        for (auto const& variable_name : {"temperature"s, "displacement"s})
         {
             auto per_process_variables =
-                findProcessVariables(variables, pv_config, {variable_names[i]});
+                findProcessVariables(variables, pv_config, {variable_name});
             process_variables.push_back(std::move(per_process_variables));
         }
         variable_T = &process_variables[0][0].get();