diff --git a/Documentation/ProjectFile/prj/time_loop/global_process_coupling/max_iteration/i_max_iteration.md b/Documentation/ProjectFile/prj/time_loop/global_process_coupling/max_iteration/i_max_iteration.md
index 728cc1c46958a097a2a1a7a4bf27c8473a811abf..363af52cb64f2648fbdd326e9a3be63606ec851b 100644
--- a/Documentation/ProjectFile/prj/time_loop/global_process_coupling/max_iteration/i_max_iteration.md
+++ b/Documentation/ProjectFile/prj/time_loop/global_process_coupling/max_iteration/i_max_iteration.md
@@ -1 +1 @@
-Defines the maximum number of the iterations of the un-coupling loop.
\ No newline at end of file
+Defines the maximum iteration number of the coupling loop of the staggered scheme.
diff --git a/ProcessLib/UncoupledProcessesTimeLoop.cpp b/ProcessLib/UncoupledProcessesTimeLoop.cpp
index d0e04085ec06cd0b8ffb68a1d2fe207df01cf4da..da205e99c97d06559d8fe4ee3367fe87bfc4eadc 100644
--- a/ProcessLib/UncoupledProcessesTimeLoop.cpp
+++ b/ProcessLib/UncoupledProcessesTimeLoop.cpp
@@ -612,136 +612,185 @@ bool UncoupledProcessesTimeLoop::loop()
         INFO("=== timestep #%u (t=%gs, dt=%gs) ==============================",
              timestep, t, delta_t);
 
-        bool coupling_iteration_converged = true;
-        // Coupling iteration
-        for (unsigned i = 0; i < _global_coupling_max_iterations; i++)
+        if (is_staggered_coupling)
+            nonlinear_solver_succeeded =
+                solveCoupledEquationSystemsByStaggeredScheme(t, delta_t,
+                                                             timestep);
+        else
+            nonlinear_solver_succeeded =
+                solveUncoupledEquationSystems(t, delta_t, timestep);
+
+        INFO("[time] Timestep #%u took %g s.", timestep,
+             time_timestep.elapsed());
+
+        if (!nonlinear_solver_succeeded)
+            break;
+    }
+
+    // output last timestep
+    if (nonlinear_solver_succeeded)
+    {
+        unsigned pcs_idx = 0;
+        for (auto& spd : _per_process_data)
         {
-            // TODO use process name
-            unsigned pcs_idx = 0;
-            for (auto& spd : _per_process_data)
-            {
-                auto& pcs = spd->process;
-                BaseLib::RunTime time_timestep_process;
-                time_timestep_process.start();
-
-                auto& x = *_process_solutions[pcs_idx];
-                if (i == 0)
-                    pcs.preTimestep(x, t, delta_t);
-
-                if (is_staggered_coupling)
-                {
-                    // Create a coupling term
-                    if (i == 0)
-                        _global_coupling_conv_crit->preFirstIteration();
-                    else
-                        _global_coupling_conv_crit->setNoFirstIteration();
-                    StaggeredCouplingTerm coupled_term(
-                        spd->coupled_processes,
-                        _solutions_of_coupled_processes[pcs_idx], delta_t);
-
-                    nonlinear_solver_succeeded = solveOneTimeStepOneProcess(
-                        x, timestep, t, delta_t, *spd, coupled_term, *_output);
-                }
-                else
-                {
-                    nonlinear_solver_succeeded = solveOneTimeStepOneProcess(
-                        x, timestep, t, delta_t, *spd,
-                        ProcessLib::createVoidStaggeredCouplingTerm(),
-                        *_output);
-                    pcs.postTimestep(x);
-                }
+            auto& pcs = spd->process;
+            auto const& x = *_process_solutions[pcs_idx];
+            _output->doOutputLastTimestep(pcs, spd->process_output, timestep, t,
+                                          x);
+
+            ++pcs_idx;
+        }
+    }
 
-                pcs.computeSecondaryVariable(t, x);
+    return nonlinear_solver_succeeded;
+}
 
-                INFO(
-                    "[time] Solving process #%u took %g s in timestep #%u "
-                    "uncouling iteration #%u",
-                    pcs_idx, time_timestep.elapsed(), timestep, i);
+bool UncoupledProcessesTimeLoop::solveUncoupledEquationSystems(
+    const double t, const double dt, const std::size_t timestep_id)
+{
+    // TODO use process name
+    unsigned pcs_idx = 0;
+    for (auto& spd : _per_process_data)
+    {
+        auto& pcs = spd->process;
+        BaseLib::RunTime time_timestep_process;
+        time_timestep_process.start();
 
-                if (!nonlinear_solver_succeeded)
-                {
-                    ERR("The nonlinear solver failed in timestep #%u at t = %g "
-                        "s"
-                        " for process #%u.",
-                        timestep, t, pcs_idx);
+        auto& x = *_process_solutions[pcs_idx];
+        pcs.preTimestep(x, t, dt);
 
-                    // save unsuccessful solution
-                    _output->doOutputAlways(pcs, spd->process_output, timestep,
-                                            t, x);
+        const auto nonlinear_solver_succeeded = solveOneTimeStepOneProcess(
+            x, timestep_id, t, dt, *spd,
+            ProcessLib::createVoidStaggeredCouplingTerm(), *_output);
+        pcs.postTimestep(x);
+        pcs.computeSecondaryVariable(t, x);
 
-                    break;
-                }
-                else
-                {
-                    if (!is_staggered_coupling)
-                        _output->doOutput(pcs, spd->process_output, timestep, t,
-                                          x);
-                }
+        INFO("[time] Solving process #%u took %g s in time step #%u ", pcs_idx,
+             time_timestep_process.elapsed(), timestep_id);
 
-                // Check the convergence of the coupling iteration
-                if (is_staggered_coupling)
-                {
-                    auto& x_old = *_solutions_of_last_cpl_iteration[pcs_idx];
-                    MathLib::LinAlg::axpy(x_old, -1.0, x);
-                    _global_coupling_conv_crit->checkResidual(x_old);
-                    coupling_iteration_converged =
-                        coupling_iteration_converged &&
-                        _global_coupling_conv_crit->isSatisfied();
-
-                    if (coupling_iteration_converged)
-                        break;
-                    MathLib::LinAlg::copy(x, x_old);
-                }
+        if (!nonlinear_solver_succeeded)
+        {
+            ERR("The nonlinear solver failed in time step #%u at t = %g "
+                "s"
+                " for process #%u.",
+                timestep_id, t, pcs_idx);
+
+            // save unsuccessful solution
+            _output->doOutputAlways(pcs, spd->process_output, timestep_id, t,
+                                    x);
+
+            return false;
+        }
+        else
+        {
+            _output->doOutput(pcs, spd->process_output, timestep_id, t, x);
+        }
+
+        ++pcs_idx;
+    }
+
+    return true;
+}
 
-                ++pcs_idx;
+bool UncoupledProcessesTimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
+    const double t, const double dt, const std::size_t timestep_id)
+{
+    // Coupling iteration
+    bool coupling_iteration_converged = true;
+    for (unsigned i = 0; i < _global_coupling_max_iterations; i++)
+    {
+        // TODO use process name
+        bool nonlinear_solver_succeeded = true;
+        unsigned pcs_idx = 0;
+        for (auto& spd : _per_process_data)
+        {
+            auto& pcs = spd->process;
+            BaseLib::RunTime time_timestep_process;
+            time_timestep_process.start();
+
+            auto& x = *_process_solutions[pcs_idx];
+            if (i == 0)
+            {
+                pcs.preTimestep(x, t, dt);
+
+                // Create a coupling term
+                _global_coupling_conv_crit->preFirstIteration();
+            }
+            else
+            {
+                _global_coupling_conv_crit->setNoFirstIteration();
             }
+            StaggeredCouplingTerm coupled_term(
+                spd->coupled_processes,
+                _solutions_of_coupled_processes[pcs_idx], dt);
 
-            if (is_staggered_coupling && coupling_iteration_converged)
-                break;
+            const auto nonlinear_solver_succeeded = solveOneTimeStepOneProcess(
+                x, timestep_id, t, dt, *spd, coupled_term, *_output);
+
+            INFO(
+                "[time] Solving process #%u took %g s in time step #%u "
+                " coupling iteration #%u",
+                pcs_idx, time_timestep_process.elapsed(), timestep_id, i);
 
             if (!nonlinear_solver_succeeded)
             {
+                ERR("The nonlinear solver failed in time step #%u at t = %g "
+                    "s"
+                    " for process #%u.",
+                    timestep_id, t, pcs_idx);
+
+                // save unsuccessful solution
+                _output->doOutputAlways(pcs, spd->process_output, timestep_id,
+                                        t, x);
+
                 break;
             }
-        }
 
-        if (is_staggered_coupling && coupling_iteration_converged)
-        {
-            unsigned pcs_idx = 0;
-            for (auto& spd : _per_process_data)
-            {
-                auto& pcs = spd->process;
-                auto& x = *_process_solutions[pcs_idx];
-                pcs.postTimestep(x);
+            // Check the convergence of the coupling iteration
+            auto& x_old = *_solutions_of_last_cpl_iteration[pcs_idx];
+            MathLib::LinAlg::axpy(x_old, -1.0, x);
+            _global_coupling_conv_crit->checkResidual(x_old);
+            coupling_iteration_converged =
+                coupling_iteration_converged &&
+                _global_coupling_conv_crit->isSatisfied();
 
-                _output->doOutput(pcs, spd->process_output, timestep, t, x);
-                ++pcs_idx;
-            }
+            if (coupling_iteration_converged)
+                break;
+            MathLib::LinAlg::copy(x, x_old);
+
+            ++pcs_idx;
         }
 
-        INFO("[time] Timestep #%u took %g s.", timestep,
-             time_timestep.elapsed());
+        if (coupling_iteration_converged)
+            break;
 
         if (!nonlinear_solver_succeeded)
-            break;
+        {
+            return false;
+        }
     }
 
-    // output last timestep
-    if (nonlinear_solver_succeeded)
+    if (!coupling_iteration_converged)
     {
-        unsigned pcs_idx = 0;
-        for (auto& spd : _per_process_data)
-        {
-            auto& pcs = spd->process;
-            auto const& x = *_process_solutions[pcs_idx];
-            _output->doOutputLastTimestep(pcs, spd->process_output, timestep, t,
-                                          x);
+        WARN(
+            "The coupling iterations reaches its maximum number in time step"
+            "#%u at t = %g s",
+            timestep_id, t);
+    }
 
-            ++pcs_idx;
-        }
+    unsigned pcs_idx = 0;
+    for (auto& spd : _per_process_data)
+    {
+        auto& pcs = spd->process;
+        auto& x = *_process_solutions[pcs_idx];
+        pcs.postTimestep(x);
+        pcs.computeSecondaryVariable(t, x);
+
+        _output->doOutput(pcs, spd->process_output, timestep_id, t, x);
+        ++pcs_idx;
     }
 
-    return nonlinear_solver_succeeded;
+    return true;
 }
 
 UncoupledProcessesTimeLoop::~UncoupledProcessesTimeLoop()
diff --git a/ProcessLib/UncoupledProcessesTimeLoop.h b/ProcessLib/UncoupledProcessesTimeLoop.h
index b9981da30608c65f60be622741965f2dace92dcb..2ce1510f777af715500ed39c494441f5e7bf5169 100644
--- a/ProcessLib/UncoupledProcessesTimeLoop.h
+++ b/ProcessLib/UncoupledProcessesTimeLoop.h
@@ -66,6 +66,33 @@ private:
     /// Solutions of the previous coupling iteration for the convergence
     /// criteria of the coupling iteration.
     std::vector<GlobalVector*> _solutions_of_last_cpl_iteration;
+
+    /**
+     * \brief Member to solver uncoupled systems of equations, which can be
+     *        a single system of equations, or several systems of equations
+     *        without any dependency among the different systems.
+     *
+     * @param t           Current time
+     * @param dt          Time step size
+     * @param timestep_id Index of the time step
+     * @return            true:  if all nonlinear solvers convergence.
+     *                    false: if any of nonlinear solvers divergences.
+     */
+    bool solveUncoupledEquationSystems(const double t, const double dt,
+                                       const std::size_t timestep_id);
+
+    /**
+     * \brief Member to solver coupled systems of equations by the staggered
+     *        scheme.
+     *
+     * @param t           Current time
+     * @param dt          Time step size
+     * @param timestep_id Index of the time step
+     * @return            true:   if all nonlinear solvers convergence.
+     *                    false:  if any of nonlinear solvers divergences.
+     */
+    bool solveCoupledEquationSystemsByStaggeredScheme(
+        const double t, const double dt, const std::size_t timestep_id);
 };
 
 //! Builds an UncoupledProcessesTimeLoop from the given configuration.