diff --git a/opm/models/discretization/common/fvbasediscretization.hh b/opm/models/discretization/common/fvbasediscretization.hh
index eae22639c53..53eea027230 100644
--- a/opm/models/discretization/common/fvbasediscretization.hh
+++ b/opm/models/discretization/common/fvbasediscretization.hh
@@ -67,8 +67,10 @@
 #include <algorithm>
 #include <array>
 #include <cstddef>
+#include <exception>
 #include <list>
 #include <memory>
+#include <mutex>
 #include <stdexcept>
 #include <sstream>
 #include <string>
@@ -726,49 +728,86 @@ public:
     {
         invalidateIntensiveQuantitiesCache(timeIdx);
 
+        // exceptions must not escape the parallel block below (that calls
+        // std::terminate()); tuck any exception away and rethrow it after the
+        // block, so that e.g. a failed flash in the property evaluation leads
+        // to a time step chop instead of an abort
+        std::mutex exceptionLock;
+        std::exception_ptr exceptionPtr = nullptr;
+
         // loop over all elements...
         ThreadedEntityIterator<GridView, /*codim=*/0> threadedElemIt(gridView_);
 #ifdef _OPENMP
 #pragma omp parallel
 #endif
         {
-            ElementContext elemCtx(simulator_);
-            ElementIterator elemIt = threadedElemIt.beginParallel();
-            for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
-                const Element& elem = *elemIt;
-                elemCtx.updatePrimaryStencil(elem);
-                elemCtx.updatePrimaryIntensiveQuantities(timeIdx);
+            try {
+                ElementContext elemCtx(simulator_);
+                for (ElementIterator elemIt = threadedElemIt.beginParallel();
+                     !threadedElemIt.isFinished(elemIt);
+                     elemIt = threadedElemIt.increment())
+                {
+                    const Element& elem = *elemIt;
+                    elemCtx.updatePrimaryStencil(elem);
+                    elemCtx.updatePrimaryIntensiveQuantities(timeIdx);
+                }
+            }
+            catch (...) {
+                std::lock_guard<std::mutex> take(exceptionLock);
+                exceptionPtr = std::current_exception();
+                threadedElemIt.setFinished();
             }
         }
+
+        if (exceptionPtr) {
+            std::rethrow_exception(exceptionPtr);
+        }
     }
 
     template <class GridViewType>
     void invalidateAndUpdateIntensiveQuantities(unsigned timeIdx, const GridViewType& gridView) const
     {
+        // see the overload above for why exceptions are bridged out of the
+        // parallel block like this
+        std::mutex exceptionLock;
+        std::exception_ptr exceptionPtr = nullptr;
+
         // loop over all elements...
         ThreadedEntityIterator<GridViewType, /*codim=*/0> threadedElemIt(gridView);
 #ifdef _OPENMP
 #pragma omp parallel
 #endif
         {
-
-            ElementContext elemCtx(simulator_);
-            auto elemIt = threadedElemIt.beginParallel();
-            for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
-                if (elemIt->partitionType() != Dune::InteriorEntity) {
-                    continue;
-                }
-                const Element& elem = *elemIt;
-                elemCtx.updatePrimaryStencil(elem);
-                // Mark cache for this element as invalid.
-                const std::size_t numPrimaryDof = elemCtx.numPrimaryDof(timeIdx);
-                for (unsigned dofIdx = 0; dofIdx < numPrimaryDof; ++dofIdx) {
-                    const unsigned globalIndex = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
-                    setIntensiveQuantitiesCacheEntryValidity(globalIndex, timeIdx, false);
+            try {
+                ElementContext elemCtx(simulator_);
+                for (auto elemIt = threadedElemIt.beginParallel();
+                     !threadedElemIt.isFinished(elemIt);
+                     elemIt = threadedElemIt.increment())
+                {
+                    if (elemIt->partitionType() != Dune::InteriorEntity) {
+                        continue;
+                    }
+                    const Element& elem = *elemIt;
+                    elemCtx.updatePrimaryStencil(elem);
+                    // Mark cache for this element as invalid.
+                    const std::size_t numPrimaryDof = elemCtx.numPrimaryDof(timeIdx);
+                    for (unsigned dofIdx = 0; dofIdx < numPrimaryDof; ++dofIdx) {
+                        const unsigned globalIndex = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
+                        setIntensiveQuantitiesCacheEntryValidity(globalIndex, timeIdx, false);
+                    }
+                    // Update for this element.
+                    elemCtx.updatePrimaryIntensiveQuantities(timeIdx);
                 }
-                // Update for this element.
-                elemCtx.updatePrimaryIntensiveQuantities(timeIdx);
             }
+            catch (...) {
+                std::lock_guard<std::mutex> take(exceptionLock);
+                exceptionPtr = std::current_exception();
+                threadedElemIt.setFinished();
+            }
+        }
+
+        if (exceptionPtr) {
+            std::rethrow_exception(exceptionPtr);
         }
     }
 
diff --git a/opm/simulators/flow/Banners.cpp b/opm/simulators/flow/Banners.cpp
index bc1cc80e5e3..d552ad0ae79 100644
--- a/opm/simulators/flow/Banners.cpp
+++ b/opm/simulators/flow/Banners.cpp
@@ -121,7 +121,8 @@ void printFlowTrailer(int nprocs,
                       const double total_setup_time,
                       const double deck_read_time,
                       const SimulatorReport& report,
-                      const std::string_view extra_summary)
+                      const std::string_view extra_summary,
+                      const bool performance_details)
 {
     std::ostringstream ss;
     ss << "\n\n================    End of simulation     ===============\n\n";
@@ -129,7 +130,7 @@ void printFlowTrailer(int nprocs,
     ss << fmt::format("Threads per MPI process: {:9}\n", nthreads);
     ss << fmt::format("Setup time:                 {:9.2f} s\n", total_setup_time);
     ss << fmt::format("  Deck input:               {:9.2f} s\n", deck_read_time);
-    report.reportFullyImplicit(ss);
+    report.reportFullyImplicit(ss, performance_details);
     if (!extra_summary.empty()) {
         if (extra_summary.front() != '\n') {
             ss << '\n';
diff --git a/opm/simulators/flow/Banners.hpp b/opm/simulators/flow/Banners.hpp
index ef342b6d8fc..bb3fce6e776 100644
--- a/opm/simulators/flow/Banners.hpp
+++ b/opm/simulators/flow/Banners.hpp
@@ -46,7 +46,8 @@ void printFlowTrailer(int nprocs,
                       const double total_setup_time,
                       const double deck_read_time,
                       const SimulatorReport& report,
-                      std::string_view extra_summary);
+                      std::string_view extra_summary,
+                      const bool performance_details = false);
 
 } // namespace Opm
 
diff --git a/opm/simulators/flow/ConvergenceOutputConfiguration.cpp b/opm/simulators/flow/ConvergenceOutputConfiguration.cpp
index 74a246ff526..ff0a4e9bf86 100644
--- a/opm/simulators/flow/ConvergenceOutputConfiguration.cpp
+++ b/opm/simulators/flow/ConvergenceOutputConfiguration.cpp
@@ -57,7 +57,7 @@ namespace {
                 fmt::format("Unsupported convergence output "
                             "option value{}: {}\n"
                             "Supported values are \"none\", "
-                            "\"steps\", and \"iterations\"",
+                            "\"steps\", \"iterations\", and \"performance\"",
                             pl, fmt::join(unsupp.begin(), u, ", "))
             };
         }
@@ -65,7 +65,7 @@ namespace {
         throw std::invalid_argument {
             fmt::format("Option {}:\n - Unsupported value{}: {}\n"
                         " - Supported values are \"none\", "
-                        "\"steps\", and \"iterations\"",
+                        "\"steps\", \"iterations\", and \"performance\"",
                         optionName, pl,
                         fmt::join(unsupp.begin(), u, ", "))
         };
@@ -79,11 +79,13 @@ namespace {
         auto opt = std::vector<Option>{};
 
         const auto values = std::unordered_map<std::string, Option> {
-            { "none"      , Option::None       },
-            { "step"      , Option::Steps      }, // Alias for 'steps' (plural)
-            { "steps"     , Option::Steps      },
-            { "iteration" , Option::Iterations }, // Alias for 'iterations' (plural)
-            { "iterations", Option::Iterations },
+            { "none"       , Option::None        },
+            { "step"       , Option::Steps       }, // Alias for 'steps' (plural)
+            { "steps"      , Option::Steps       },
+            { "iteration"  , Option::Iterations  }, // Alias for 'iterations' (plural)
+            { "iterations" , Option::Iterations  },
+            { "perf"       , Option::Performance }, // Alias for 'performance'
+            { "performance", Option::Performance },
         };
 
         auto unsupp = std::vector<std::string>{};
diff --git a/opm/simulators/flow/ConvergenceOutputConfiguration.hpp b/opm/simulators/flow/ConvergenceOutputConfiguration.hpp
index 78f30bd7a14..9103894e63e 100644
--- a/opm/simulators/flow/ConvergenceOutputConfiguration.hpp
+++ b/opm/simulators/flow/ConvergenceOutputConfiguration.hpp
@@ -39,6 +39,9 @@ namespace Opm {
 ///   * "iterations" -- Want additional convergence output pertaining to each
 ///                     non-linar ieration in each timestep.
 ///
+///   * "performance" -- Want additional performance details in the summary
+///                      at the end of the simulation.
+///
 /// Option value "none" overrides all other options.  In other words, if the
 /// user requests "none", then there will be no additional convergence
 /// output, even if there are other options in the option string.
@@ -54,6 +57,7 @@ class ConvergenceOutputConfiguration
         None = 0,
         Steps = 1 << 1,
         Iterations = 1 << 2,
+        Performance = 1 << 3,
     };
 
     /// Constructor
diff --git a/opm/simulators/flow/FlowGenericProblem_impl.hpp b/opm/simulators/flow/FlowGenericProblem_impl.hpp
index e15d5e351c6..bebcfa1963b 100644
--- a/opm/simulators/flow/FlowGenericProblem_impl.hpp
+++ b/opm/simulators/flow/FlowGenericProblem_impl.hpp
@@ -81,9 +81,25 @@ FlowGenericProblem(const EclipseState& eclState,
     // 2. EQLDIMS item 2.  Default value from
     //    opm-common/opm/input/eclipse/share/keywords/000_Eclipse100/E/EQLDIMS
 
-    numPressurePointsEquil_ = Parameters::IsSet<Parameters::NumPressurePointsEquil>()
-        ? Parameters::Get<Parameters::NumPressurePointsEquil>()
-        : eclState.getTableManager().getEqldims().getNumDepthNodesP();
+    if (Parameters::IsSet<Parameters::NumPressurePointsEquil>()) {
+        numPressurePointsEquil_ = Parameters::Get<Parameters::NumPressurePointsEquil>();
+        if (numPressurePointsEquil_ < 1) {
+            throw std::invalid_argument {
+                fmt::format("--num-pressure-points-equil must be at least 1, "
+                            "but {} was given.", numPressurePointsEquil_)
+            };
+        }
+    }
+    else {
+        numPressurePointsEquil_ = eclState.getTableManager().getEqldims().getNumDepthNodesP();
+        if (numPressurePointsEquil_ < 1) {
+            throw std::invalid_argument {
+                fmt::format("EQLDIMS item 2, the number of depth nodes in the "
+                            "equilibration pressure tables, must be at least 1, "
+                            "but {} was given.", numPressurePointsEquil_)
+            };
+        }
+    }
 
     explicitRockCompaction_ = Parameters::Get<Parameters::ExplicitRockCompaction>();
 }
diff --git a/opm/simulators/flow/FlowMain.hpp b/opm/simulators/flow/FlowMain.hpp
index 8ec955096ca..67527dd6a8e 100644
--- a/opm/simulators/flow/FlowMain.hpp
+++ b/opm/simulators/flow/FlowMain.hpp
@@ -29,6 +29,7 @@
 #include <opm/models/utils/start.hh>
 
 #include <opm/simulators/flow/Banners.hpp>
+#include <opm/simulators/flow/ConvergenceOutputConfiguration.hpp>
 #include <opm/simulators/flow/FlowUtils.hpp>
 #include <opm/simulators/flow/NlddReporting.hpp>
 #include <opm/simulators/flow/SimulatorFullyImplicit.hpp>
@@ -435,8 +436,14 @@ namespace Opm {
                 = omp_get_max_threads();
 #endif
 
+            const auto extraConvOutput = ConvergenceOutputConfiguration {
+                Parameters::Get<Parameters::OutputExtraConvergenceInfo>(),
+                R"(OutputExtraConvergenceInfo (--output-extra-convergence-info))"
+            };
+
             printFlowTrailer(mpi_size_, threads, total_setup_time_, deck_read_time_, report,
-                             simulator_->model().simulator().problem().extraTrailerSummary());
+                             simulator_->model().simulator().problem().extraTrailerSummary(),
+                             extraConvOutput.want(ConvergenceOutputConfiguration::Option::Performance));
 
             detail::handleExtraConvergenceOutput(report,
                                                  Parameters::Get<Parameters::OutputExtraConvergenceInfo>(),
diff --git a/opm/simulators/flow/SimulatorFullyImplicit.cpp b/opm/simulators/flow/SimulatorFullyImplicit.cpp
index f540693e6ec..5efb770aa38 100644
--- a/opm/simulators/flow/SimulatorFullyImplicit.cpp
+++ b/opm/simulators/flow/SimulatorFullyImplicit.cpp
@@ -40,7 +40,9 @@ void registerSimulatorParameters()
          "\"none\" gives no extra output and "
          "overrides all other options, "
          "\"steps\" generates an INFOSTEP file, "
-         "\"iterations\" generates an INFOITER file. "
+         "\"iterations\" generates an INFOITER file, "
+         "\"performance\" adds detailed timings to the "
+         "performance summary at the end of the simulation. "
          "Combine options with commas, e.g., "
          "\"steps,iterations\" for multiple outputs.");
     Parameters::Register<Parameters::SaveStep>
diff --git a/opm/simulators/flow/Transmissibility.hpp b/opm/simulators/flow/Transmissibility.hpp
index 9bc3d6c62bc..40e67d491d1 100644
--- a/opm/simulators/flow/Transmissibility.hpp
+++ b/opm/simulators/flow/Transmissibility.hpp
@@ -100,6 +100,24 @@ class Transmissibility {
      */
     Scalar thermalHalfTrans(unsigned insideElemIdx, unsigned outsideElemIdx) const;
 
+    /*!
+     * \brief One-sided (half) transmissibility of the face between two
+     *        cells, seen from the inside cell (NTG applied, before face
+     *        multipliers). Only available when setStoreHalfTrans(true)
+     *        was called before update(); used by the adjoint module for
+     *        the permeability chain rule
+     *        dT/dK_inside = (h_out/(h_in+h_out))^2 dh_in/dK_inside.
+     */
+    Scalar halfTransmissibility(unsigned insideElemIdx, unsigned outsideElemIdx) const;
+
+    //! \brief Enable storage of the one-sided half transmissibilities
+    //!        (call before update()).
+    void setStoreHalfTrans(bool yesno)
+    { storeHalfTrans_ = yesno; }
+
+    bool storeHalfTrans() const
+    { return storeHalfTrans_; }
+
     Scalar thermalHalfTransBoundary(unsigned insideElemIdx, unsigned boundaryFaceIdx) const;
 
     const std::map<std::pair<unsigned, unsigned>, Scalar>& getThermalHalfTransBoundary() const;
@@ -296,6 +314,8 @@ class Transmissibility {
     bool enableDispersivity_;
     bool warnEditNNC_ = true;
     std::unordered_map<std::uint64_t, Scalar> thermalHalfTrans_; //NB this is based on direction map size is ca 2*trans_ (diffusivity_)
+    std::unordered_map<std::uint64_t, Scalar> halfTrans_; // directional, only filled when storeHalfTrans_
+    bool storeHalfTrans_ = false;
     std::unordered_map<std::uint64_t, Scalar> diffusivity_;
     std::unordered_map<std::uint64_t, Scalar> dispersivity_;
 
diff --git a/opm/simulators/flow/Transmissibility_impl.hpp b/opm/simulators/flow/Transmissibility_impl.hpp
index c268d2eed91..0ededd20b43 100644
--- a/opm/simulators/flow/Transmissibility_impl.hpp
+++ b/opm/simulators/flow/Transmissibility_impl.hpp
@@ -138,6 +138,24 @@ thermalHalfTrans(unsigned insideElemIdx, unsigned outsideElemIdx) const
     return thermalHalfTrans_.at(details::directionalIsId(insideElemIdx, outsideElemIdx));
 }
 
+template<class Grid, class GridView, class ElementMapper, class CartesianIndexMapper, class Scalar>
+Scalar Transmissibility<Grid,GridView,ElementMapper,CartesianIndexMapper,Scalar>::
+halfTransmissibility(unsigned insideElemIdx, unsigned outsideElemIdx) const
+{
+    if (!storeHalfTrans_) {
+        // Without this the caller gets a bare std::out_of_range from the empty
+        // map.  One branch on a bool, and this is not a hot path: only the
+        // adjoint code asks for these today.
+        OPM_THROW(std::logic_error,
+                  "One-sided half transmissibilities were not stored. "
+                  "Call setStoreHalfTrans(true) before update() to have them "
+                  "computed; they are off by default because only the adjoint "
+                  "code needs them.");
+    }
+
+    return halfTrans_.at(details::directionalIsId(insideElemIdx, outsideElemIdx));
+}
+
 template<class Grid, class GridView, class ElementMapper, class CartesianIndexMapper, class Scalar>
 Scalar Transmissibility<Grid,GridView,ElementMapper,CartesianIndexMapper,Scalar>::
 thermalHalfTransBoundary(unsigned insideElemIdx, unsigned boundaryFaceIdx) const
@@ -212,6 +230,13 @@ update(bool global, const TransUpdateQuantities update_quantities,
 
     transBoundary_.clear();
 
+    if (storeHalfTrans_) {
+        halfTrans_.clear();
+        if (num_threads == 1) {
+            halfTrans_.reserve(numElements*6*1.05);
+        }
+    }
+
     // if energy is enabled, let's do the same for the "thermal half transmissibilities"
     if ( enableEnergy_ && !onlyTrans) {
         thermalHalfTrans_.clear();
@@ -312,6 +337,8 @@ update(bool global, const TransUpdateQuantities update_quantities,
                                                   MapBuilderInsertionMode::Insert_Or_Assign);
     ThreadSafeMapBuilder thermalHalfTrans(thermalHalfTrans_, num_threads,
                                           MapBuilderInsertionMode::Insert_Or_Assign);
+    ThreadSafeMapBuilder halfTransMap(halfTrans_, num_threads,
+                                      MapBuilderInsertionMode::Insert_Or_Assign);
     ThreadSafeMapBuilder diffusivity(diffusivity_, num_threads,
                                      MapBuilderInsertionMode::Insert_Or_Assign);
     ThreadSafeMapBuilder dispersivity(dispersivity_, num_threads,
@@ -457,6 +484,22 @@ update(bool global, const TransUpdateQuantities update_quantities,
 
                 Scalar trans = computeHalfMean(computeHalfTrans_, permeability_);
 
+                if (storeHalfTrans_) {
+                    // one-sided half transmissibilities (NTG applied),
+                    // for the adjoint permeability chain rule
+                    auto onesided = computeHalf(computeHalfTrans_,
+                                                permeability_[inside.elemIdx],
+                                                permeability_[outside.elemIdx]);
+                    applyNtg_(onesided[0], inside, ntg);
+                    applyNtg_(onesided[1], outside, ntg);
+                    halfTransMap.insert_or_assign(
+                        details::directionalIsId(inside.elemIdx, outside.elemIdx),
+                        onesided[0]);
+                    halfTransMap.insert_or_assign(
+                        details::directionalIsId(outside.elemIdx, inside.elemIdx),
+                        onesided[1]);
+                }
+
                 // apply the full face transmissibility multipliers
                 // for the inside ...
                 if (!pinchActive) {
@@ -559,6 +602,10 @@ update(bool global, const TransUpdateQuantities update_quantities,
 #pragma omp section
 #endif
         dispersivity.finalize();
+#ifdef _OPENMP
+#pragma omp section
+#endif
+        halfTransMap.finalize();
     }
 
     // Potentially overwrite and/or modify transmissibilities based on input from deck
diff --git a/opm/simulators/flow/equil/InitStateEquil_impl.hpp b/opm/simulators/flow/equil/InitStateEquil_impl.hpp
index 0ee10d8b119..db6ccd49292 100644
--- a/opm/simulators/flow/equil/InitStateEquil_impl.hpp
+++ b/opm/simulators/flow/equil/InitStateEquil_impl.hpp
@@ -323,6 +323,11 @@ RK4IVP<Scalar,RHS>::RK4IVP(const RHS& f,
     : N_(N)
     , span_(span)
 {
+    // stepsize() divides by N_ and operator() evaluates interval N_ - 1.  A
+    // non-positive sample count is rejected when the value is loaded in
+    // FlowGenericProblem, where its source is known.
+    assert(N >= 1);
+
     const Scalar h = stepsize();
     const Scalar h2 = h / 2;
     const Scalar h6 = h / 6;
diff --git a/opm/simulators/timestepping/SimulatorReport.cpp b/opm/simulators/timestepping/SimulatorReport.cpp
index 223fe0b3ba7..df602850d1f 100644
--- a/opm/simulators/timestepping/SimulatorReport.cpp
+++ b/opm/simulators/timestepping/SimulatorReport.cpp
@@ -36,7 +36,9 @@ namespace Opm
                                      7.0, 8.0, 9.0, 10.0, 11.0, 12.0,
                                      13, 14, 15, 16, 17, 18, 19,
                                      true, false, false, 20, 21.0, 22.0,
-                                     23, 24, 25, 26, 27, 28, 29, 30};
+                                     23, 24, 25, 26, 27, 28, 29, 30,
+                                     31.0, 32.0, 33.0, 34.0, 35.0, 36.0, 37.0, 38.0, 39.0, 40.0, 41.0, 42.0, 43.0, 44.0, 45.0, 46.0, 47.0, 48.0, 49.0,
+                                     50, 51};
     }
 
     bool SimulatorReportSingle::operator==(const SimulatorReportSingle& rhs) const
@@ -73,7 +75,28 @@ namespace Opm
                this->converged_domains == rhs.converged_domains &&
                this->unconverged_domains == rhs.unconverged_domains &&
                this->accepted_unconverged_domains == rhs.accepted_unconverged_domains &&
-               this->skipped_domains == rhs.skipped_domains;
+               this->skipped_domains == rhs.skipped_domains &&
+               this->props_time == rhs.props_time &&
+               this->convergence_check_time == rhs.convergence_check_time &&
+               this->output_eval_time == rhs.output_eval_time &&
+               this->tracer_solve_time == rhs.tracer_solve_time &&
+               this->temperature_solve_time == rhs.temperature_solve_time &&
+               this->output_disk_write_time == rhs.output_disk_write_time &&
+               this->linear_solve_apply_time == rhs.linear_solve_apply_time &&
+               this->precond_setup_time == rhs.precond_setup_time &&
+               this->precond_apply_time == rhs.precond_apply_time &&
+               this->well_solve_time == rhs.well_solve_time &&
+               this->well_potential_solve_time == rhs.well_potential_solve_time &&
+               this->well_solve_assemble_time == rhs.well_solve_assemble_time &&
+               this->well_solve_linear_solve_time == rhs.well_solve_linear_solve_time &&
+               this->well_control_network_time == rhs.well_control_network_time &&
+               this->gaslift_time == rhs.gaslift_time &&
+               this->well_facility_time == rhs.well_facility_time &&
+               this->group_control_time == rhs.group_control_time &&
+               this->network_balance_time == rhs.network_balance_time &&
+               this->control_well_solve_time == rhs.control_well_solve_time &&
+               this->total_well_potential_iterations == rhs.total_well_potential_iterations &&
+               this->total_network_iterations == rhs.total_network_iterations;
     }
 
     void SimulatorReportSingle::operator+=(const SimulatorReportSingle& sr)
@@ -88,9 +111,30 @@ namespace Opm
         pre_post_time += sr.pre_post_time;
         assemble_time_well += sr.assemble_time_well;
         update_time += sr.update_time;
+        props_time += sr.props_time;
+        convergence_check_time += sr.convergence_check_time;
+        output_eval_time += sr.output_eval_time;
+        tracer_solve_time += sr.tracer_solve_time;
+        temperature_solve_time += sr.temperature_solve_time;
+        output_disk_write_time += sr.output_disk_write_time;
+        linear_solve_apply_time += sr.linear_solve_apply_time;
+        precond_setup_time += sr.precond_setup_time;
+        precond_apply_time += sr.precond_apply_time;
+        well_solve_time += sr.well_solve_time;
+        well_potential_solve_time += sr.well_potential_solve_time;
+        well_solve_assemble_time += sr.well_solve_assemble_time;
+        well_solve_linear_solve_time += sr.well_solve_linear_solve_time;
+        well_control_network_time += sr.well_control_network_time;
+        gaslift_time += sr.gaslift_time;
+        well_facility_time += sr.well_facility_time;
+        group_control_time += sr.group_control_time;
+        network_balance_time += sr.network_balance_time;
+        control_well_solve_time += sr.control_well_solve_time;
         output_write_time += sr.output_write_time;
         total_time += sr.total_time;
         total_well_iterations += sr.total_well_iterations;
+        total_well_potential_iterations += sr.total_well_potential_iterations;
+        total_network_iterations += sr.total_network_iterations;
         total_linearizations += sr.total_linearizations;
         total_newton_iterations += sr.total_newton_iterations;
         relaxed_cnv_acceptances += sr.relaxed_cnv_acceptances;
@@ -132,7 +176,9 @@ namespace Opm
         return val;
     };
 
-    void SimulatorReportSingle::reportFullyImplicit(std::ostream& os, const SimulatorReportSingle* failureReport) const
+    void SimulatorReportSingle::reportFullyImplicit(std::ostream& os,
+                                                    const SimulatorReportSingle* failureReport,
+                                                    bool performance_details) const
     {
         // Disabling this output as it is redundant with the solver_time, now
         // output as "Simulation time". Usually very small difference between them.
@@ -162,6 +208,22 @@ namespace Opm
             }
             os << std::endl;
 
+            // part of the well assembly used for the facility calculations
+            // (further details in the Facility calculations section below)
+            t = well_facility_time + (failureReport ? failureReport->well_facility_time : 0.0);
+            if (performance_details && t > 0.0) {
+                os << fmt::format("      Facility calculations:  {:7.2f} s", t);
+                os << std::endl;
+
+                t = well_control_network_time + (failureReport ? failureReport->well_control_network_time : 0.0);
+                os << fmt::format("        Control/network:      {:7.2f} s", t);
+                os << std::endl;
+
+                t = well_solve_time + (failureReport ? failureReport->well_solve_time : 0.0);
+                os << fmt::format("        Well solves:          {:7.2f} s", t);
+                os << std::endl;
+            }
+
             t = linear_solve_time + (failureReport ? failureReport->linear_solve_time : 0.0);
             os << fmt::format("  Linear solve time:        {:9.2f} s", t);
             if (failureReport) {
@@ -180,6 +242,29 @@ namespace Opm
             }
             os << std::endl;
 
+            t = precond_setup_time + (failureReport ? failureReport->precond_setup_time : 0.0);
+            if (performance_details && t > 0.0) {
+                os << fmt::format("      Precond setup:          {:7.2f} s", t);
+                os << std::endl;
+            }
+
+            if (performance_details) {
+                t = linear_solve_apply_time + (failureReport ? failureReport->linear_solve_apply_time : 0.0);
+                os << fmt::format("    Linear apply:             {:7.2f} s", t);
+                if (failureReport) {
+                  os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                    failureReport->linear_solve_apply_time,
+                                    100*failureReport->linear_solve_apply_time/noZero(t));
+                }
+                os << std::endl;
+            }
+
+            t = precond_apply_time + (failureReport ? failureReport->precond_apply_time : 0.0);
+            if (performance_details && t > 0.0) {
+                os << fmt::format("      Precond apply:          {:7.2f} s", t);
+                os << std::endl;
+            }
+
             if (local_solve_time > 0.0) {
                 t = local_solve_time + (failureReport ? failureReport->local_solve_time : 0.0);
                 os << fmt::format("  Local solve time:           {:7.2f} s", t);
@@ -199,6 +284,29 @@ namespace Opm
                                 100*failureReport->update_time/noZero(t));
             }
             os << std::endl;
+
+            if (performance_details) {
+                t = props_time + (failureReport ? failureReport->props_time : 0.0);
+                os << fmt::format("    Properties evaluation:    {:7.2f} s", t);
+                if (failureReport) {
+                  os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                    failureReport->props_time,
+                                    100*failureReport->props_time/noZero(t));
+                }
+                os << std::endl;
+            }
+
+            if (performance_details) {
+                t = convergence_check_time + (failureReport ? failureReport->convergence_check_time : 0.0);
+                os << fmt::format("    Convergence checks:       {:7.2f} s", t);
+                if (failureReport) {
+                  os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                    failureReport->convergence_check_time,
+                                    100*failureReport->convergence_check_time/noZero(t));
+                }
+                os << std::endl;
+            }
+
             t = pre_post_time + (failureReport ? failureReport->pre_post_time : 0.0);
             os << fmt::format("  Pre/post step:              {:7.2f} s", t);
             if (failureReport) {
@@ -208,9 +316,110 @@ namespace Opm
             }
             os << std::endl;
 
+            if (performance_details) {
+                t = output_eval_time + (failureReport ? failureReport->output_eval_time : 0.0);
+                os << fmt::format("    Output evaluation:        {:7.2f} s", t);
+                if (failureReport) {
+                  os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                    failureReport->output_eval_time,
+                                    100*failureReport->output_eval_time/noZero(t));
+                }
+                os << std::endl;
+            }
+
+            t = tracer_solve_time + (failureReport ? failureReport->tracer_solve_time : 0.0);
+            if (performance_details && t > 0.0) {
+                os << fmt::format("    Tracer solve:             {:7.2f} s", t);
+                if (failureReport) {
+                  os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                    failureReport->tracer_solve_time,
+                                    100*failureReport->tracer_solve_time/noZero(t));
+                }
+                os << std::endl;
+            }
+
+            t = temperature_solve_time + (failureReport ? failureReport->temperature_solve_time : 0.0);
+            if (performance_details && t > 0.0) {
+                os << fmt::format("    Temperature solve:        {:7.2f} s", t);
+                if (failureReport) {
+                  os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                    failureReport->temperature_solve_time,
+                                    100*failureReport->temperature_solve_time/noZero(t));
+                }
+                os << std::endl;
+            }
+
             os << fmt::format("  Output write time:          {:7.2f} s",
                               output_write_time + (failureReport ? failureReport->output_write_time : 0.0));
             os << std::endl;
+
+            t = output_disk_write_time + (failureReport ? failureReport->output_disk_write_time : 0.0);
+            if (performance_details && t > 0.0) {
+                os << fmt::format("    Actual disk write:        {:7.2f} s", t);
+                os << std::endl;
+            }
+
+        }
+
+        // Summary of all time used in the wells/facility. This cuts across
+        // the simulation phases above: the control/network updates and the
+        // normal well solves happen during the assembly, while the well
+        // potential solves are part of the pre/post step (output
+        // evaluation). The indented items are parts of their parent, but do
+        // not necessarily add up to it.
+        const double total_facility_time =
+            well_facility_time + well_potential_solve_time +
+            (failureReport ? failureReport->well_facility_time +
+                             failureReport->well_potential_solve_time : 0.0);
+        if (performance_details && total_facility_time > 0.0) {
+            os << fmt::format("Facility calculations:      {:9.2f} s", total_facility_time);
+            if (failureReport) {
+              os << fmt::format(" (Wasted: {:2.1f} s; {:2.1f}%)",
+                                failureReport->well_facility_time + failureReport->well_potential_solve_time,
+                                100*(failureReport->well_facility_time + failureReport->well_potential_solve_time)
+                                   /noZero(total_facility_time));
+            }
+            os << std::endl;
+
+            double t = well_control_network_time + (failureReport ? failureReport->well_control_network_time : 0.0);
+            os << fmt::format("  Control/network:            {:7.2f} s", t);
+            os << std::endl;
+
+            t = group_control_time + (failureReport ? failureReport->group_control_time : 0.0);
+            os << fmt::format("    Group/control updates:    {:7.2f} s", t);
+            os << std::endl;
+
+            t = network_balance_time + (failureReport ? failureReport->network_balance_time : 0.0);
+            if (t > 0.0) {
+                os << fmt::format("    Network balance:          {:7.2f} s", t);
+                os << std::endl;
+            }
+
+            t = gaslift_time + (failureReport ? failureReport->gaslift_time : 0.0);
+            if (t > 0.0) {
+                os << fmt::format("    Gas lift optimize:        {:7.2f} s", t);
+                os << std::endl;
+            }
+
+            t = control_well_solve_time + (failureReport ? failureReport->control_well_solve_time : 0.0);
+            os << fmt::format("    Well solves (control):    {:7.2f} s", t);
+            os << std::endl;
+
+            t = well_solve_time + (failureReport ? failureReport->well_solve_time : 0.0);
+            os << fmt::format("  Well solves:                {:7.2f} s", t);
+            os << std::endl;
+
+            t = well_solve_assemble_time + (failureReport ? failureReport->well_solve_assemble_time : 0.0);
+            os << fmt::format("    Assembly:                 {:7.2f} s", t);
+            os << std::endl;
+
+            t = well_solve_linear_solve_time + (failureReport ? failureReport->well_solve_linear_solve_time : 0.0);
+            os << fmt::format("    Linear solve:             {:7.2f} s", t);
+            os << std::endl;
+
+            t = well_potential_solve_time + (failureReport ? failureReport->well_potential_solve_time : 0.0);
+            os << fmt::format("  Potential solves:           {:7.2f} s", t);
+            os << std::endl;
         }
 
         int n = total_linearizations + (failureReport ? failureReport->total_linearizations : 0);
@@ -246,6 +455,39 @@ namespace Opm
             os << fmt::format("Relaxed CNV acceptances:   {:7}", relaxed_cnv_acceptances)
                << std::endl;
         }
+
+        n = total_well_iterations + (failureReport ? failureReport->total_well_iterations : 0);
+        unsigned int np = total_well_potential_iterations +
+            (failureReport ? failureReport->total_well_potential_iterations : 0);
+        if (performance_details && (n > 0 || np > 0)) {
+            os << fmt::format("Overall Well Iterations:   {:7}", n);
+            if (failureReport) {
+              os << fmt::format("      (Wasted: {:5}; {:2.1f}%)",
+                                failureReport->total_well_iterations,
+                                100.0*failureReport->total_well_iterations/noZero(n));
+            }
+            os << std::endl;
+
+            os << fmt::format("Well Potential Iterations: {:7}", np);
+            if (failureReport) {
+              os << fmt::format("      (Wasted: {:5}; {:2.1f}%)",
+                                failureReport->total_well_potential_iterations,
+                                100.0*failureReport->total_well_potential_iterations/noZero(np));
+            }
+            os << std::endl;
+        }
+
+        unsigned int nn = total_network_iterations +
+            (failureReport ? failureReport->total_network_iterations : 0);
+        if (performance_details && nn > 0) {
+            os << fmt::format("Network Balance Iterations:{:7}", nn);
+            if (failureReport) {
+              os << fmt::format("      (Wasted: {:5}; {:2.1f}%)",
+                                failureReport->total_network_iterations,
+                                100.0*failureReport->total_network_iterations/noZero(nn));
+            }
+            os << std::endl;
+        }
     }
 
 
@@ -388,10 +630,10 @@ namespace Opm
         stepreports.insert(stepreports.end(), sr.stepreports.begin(), sr.stepreports.end());
     }
 
-    void SimulatorReport::reportFullyImplicit(std::ostream& os) const
+    void SimulatorReport::reportFullyImplicit(std::ostream& os, bool performance_details) const
     {
         os << fmt::format("Number of timesteps:     {:9}\n", stepreports.size());
-        success.reportFullyImplicit(os, &failure);
+        success.reportFullyImplicit(os, &failure, performance_details);
     }
 
     void SimulatorReport::reportNLDD(std::ostream& os) const
diff --git a/opm/simulators/timestepping/SimulatorReport.hpp b/opm/simulators/timestepping/SimulatorReport.hpp
index cd56a1a3f6e..7c48c7326b8 100644
--- a/opm/simulators/timestepping/SimulatorReport.hpp
+++ b/opm/simulators/timestepping/SimulatorReport.hpp
@@ -72,6 +72,59 @@ namespace Opm
         int accepted_unconverged_domains = 0;
         int skipped_domains = 0;
 
+        /// Part of update_time used for evaluating properties
+        /// (intensive quantities) after a solution update.
+        double props_time = 0.0;
+        /// Part of update_time used for the convergence checks
+        /// (reservoir and wells).
+        double convergence_check_time = 0.0;
+        /// Part of pre_post_time used for the extra evaluations related
+        /// to output (endTimeStep: summary evaluation, actions etc.).
+        double output_eval_time = 0.0;
+        /// Part of pre_post_time used for solving the tracer equations.
+        double tracer_solve_time = 0.0;
+        /// Part of pre_post_time used for the sequential implicit
+        /// temperature solve.
+        double temperature_solve_time = 0.0;
+        /// Time used by the actual (possibly asynchronous) output writes.
+        double output_disk_write_time = 0.0;
+        /// Part of linear_solve_time spent in the actual solver apply.
+        double linear_solve_apply_time = 0.0;
+        /// Part of linear_solve_setup_time used for creating and updating
+        /// the preconditioner.
+        double precond_setup_time = 0.0;
+        /// Part of linear_solve_apply_time used for applying the
+        /// preconditioner.
+        double precond_apply_time = 0.0;
+        /// Statistics for the separate (standalone) well solves, i.e. the
+        /// local solves of the well equations done outside the global
+        /// linearizations (prepareTimeStep, well testing) and the well
+        /// potential calculations.
+        double well_solve_time = 0.0;
+        double well_potential_solve_time = 0.0;
+        double well_solve_assemble_time = 0.0;
+        double well_solve_linear_solve_time = 0.0;
+        /// Part of assemble_time_well used for updating well controls,
+        /// group targets and balancing the network.
+        double well_control_network_time = 0.0;
+        /// Part of well_control_network_time used for gas lift optimization.
+        double gaslift_time = 0.0;
+        /// Part of assemble_time_well used for the facility calculations,
+        /// i.e. the time step preparation of the wells (including the
+        /// standalone well solves) and the control/network updates.
+        double well_facility_time = 0.0;
+        /// Part of well_control_network_time used for updating and
+        /// communicating the group data, checking the well controls and
+        /// updating the guide rates.
+        double group_control_time = 0.0;
+        /// Part of well_control_network_time used for balancing the network.
+        double network_balance_time = 0.0;
+        /// Part of well_control_network_time used for the well solves
+        /// triggered by control changes (prepareWellsBeforeAssembling).
+        double control_well_solve_time = 0.0;
+        unsigned int total_well_potential_iterations = 0;
+        unsigned int total_network_iterations = 0;
+
         static SimulatorReportSingle serializationTestObject();
 
         bool operator==(const SimulatorReportSingle&) const;
@@ -80,7 +133,10 @@ namespace Opm
         /// Print a report suitable for a single simulation step.
         void reportStep(std::ostream& os) const;
         /// Print a report suitable for the end of a fully implicit case, leaving out the pressure/transport time.
-        void reportFullyImplicit(std::ostream& os, const SimulatorReportSingle* failedReport = nullptr) const;
+        /// With performance_details, additional detailed timings are included.
+        void reportFullyImplicit(std::ostream& os,
+                                 const SimulatorReportSingle* failedReport = nullptr,
+                                 bool performance_details = false) const;
         void reportNLDD(std::ostream& os, const SimulatorReportSingle* failedReport = nullptr) const;
         template<class Serializer>
         void serializeOp(Serializer& serializer)
@@ -118,6 +174,27 @@ namespace Opm
             serializer(unconverged_domains);
             serializer(accepted_unconverged_domains);
             serializer(skipped_domains);
+            serializer(props_time);
+            serializer(convergence_check_time);
+            serializer(output_eval_time);
+            serializer(tracer_solve_time);
+            serializer(temperature_solve_time);
+            serializer(output_disk_write_time);
+            serializer(linear_solve_apply_time);
+            serializer(precond_setup_time);
+            serializer(precond_apply_time);
+            serializer(well_solve_time);
+            serializer(well_potential_solve_time);
+            serializer(well_solve_assemble_time);
+            serializer(well_solve_linear_solve_time);
+            serializer(well_control_network_time);
+            serializer(gaslift_time);
+            serializer(well_facility_time);
+            serializer(group_control_time);
+            serializer(network_balance_time);
+            serializer(control_well_solve_time);
+            serializer(total_well_potential_iterations);
+            serializer(total_network_iterations);
         }
     };
 
@@ -132,7 +209,7 @@ namespace Opm
         bool operator==(const SimulatorReport&) const;
         void operator+=(const SimulatorReportSingle& sr);
         void operator+=(const SimulatorReport& sr);
-        void reportFullyImplicit(std::ostream& os) const;
+        void reportFullyImplicit(std::ostream& os, bool performance_details = false) const;
         void reportNLDD(std::ostream& os) const;
         void fullReports(std::ostream& os) const;
 
diff --git a/opm/simulators/wells/BlackoilWellModel.hpp b/opm/simulators/wells/BlackoilWellModel.hpp
index 021da5f35fb..59c9df8e6c4 100644
--- a/opm/simulators/wells/BlackoilWellModel.hpp
+++ b/opm/simulators/wells/BlackoilWellModel.hpp
@@ -615,6 +615,10 @@ template<class Scalar> class WellContributions;
             // and in the well equations.
             void assemble(const double dt);
 
+            // collect the statistics for the standalone well solves gathered
+            // since the last call into last_report_ and reset them.
+            void collectWellSolveStats();
+
             // well controls and network pressures affect each other and are solved in an iterative manner.
             // the function handles one iteration of updating well controls and network pressures.
             // it is possible to decouple the update of well controls and network pressures further.
diff --git a/opm/simulators/wells/BlackoilWellModelRescoup.hpp b/opm/simulators/wells/BlackoilWellModelRescoup.hpp
index 40f3a9809f7..d62b119f775 100644
--- a/opm/simulators/wells/BlackoilWellModelRescoup.hpp
+++ b/opm/simulators/wells/BlackoilWellModelRescoup.hpp
@@ -242,6 +242,17 @@ class BlackoilWellModelRescoup {
     ///   activated slaves and gates the master's own iteration.
     bool masterNetworkHasMasterGroupLeavesForSlave_(std::size_t slave_idx) const;
 
+    /// \brief Slave-side: true iff this slave's own deck put the named group
+    ///   into its own surface network as a fixed-pressure node.
+    ///
+    /// This is the condition under which a master-supplied node pressure may
+    /// be imposed on the group at all: the pressure is the fixed pressure of
+    /// the slave group's node in the *slave's* network, so the deck must have
+    /// declared such a node (GRUPNET item 2, or NODEPROP item 2 for an
+    /// extended network).  A slave that declares no network keeps its wells
+    /// on their WCONPROD THP limits.
+    bool slaveGroupIsFixedPressureNodeInOwnNetwork_(const std::string& group_name) const;
+
     BlackoilWellModel<TypeTag>& well_model_;
     BlackoilWellModelNetwork<TypeTag>& network_;
     Simulator& simulator_;
diff --git a/opm/simulators/wells/BlackoilWellModelRescoup_impl.hpp b/opm/simulators/wells/BlackoilWellModelRescoup_impl.hpp
index 74b75457d55..27379d0a55f 100644
--- a/opm/simulators/wells/BlackoilWellModelRescoup_impl.hpp
+++ b/opm/simulators/wells/BlackoilWellModelRescoup_impl.hpp
@@ -182,38 +182,47 @@ receiveMasterGroupNodePressuresFromMaster()
     if (num_pressures > 0) {
         rescoup_slave.receiveMasterGroupNodePressuresFromMaster(num_pressures);
     }
-    // Apply pressures as dynamic THP limits on every producer whose
-    // group has a master-supplied pressure.  Wells in master groups that
-    // are not network leaves are not touched (no entry in the map).
-    // Mirrors the standard local-network apply pattern in
+    // Apply pressures as dynamic THP limits on every producer whose group
+    // has a master-supplied pressure *and* whose group the slave's own deck
+    // placed in its own surface network as a fixed-pressure node.  Wells in
+    // master groups that are not network leaves are not touched either (no
+    // entry in the map).  Mirrors the standard local-network apply pattern in
     // BlackoilWellModelNetworkGeneric::updatePressures(): when the well is
     // currently THP-controlled, also write the new THP into the WellState
     // directly, because setDynamicThpLimit() alone leaves the active
     // control's THP value stale and the subsequent well-solve would
     // converge against the old THP.
     //
-    // TODO (follow-up PR): this THP-direct path is only correct when the
-    // slave has no extended network between its slave group and the wells
-    // (the case exercised here).  When the slave uses an extended network
-    // with the slave group declared as a fixed-pressure node, the master-
-    // supplied pressure should instead be installed as that node's terminal
-    // pressure and the slave's network solver should propagate it down to
-    // the wells.  The plumbing exists -- the solver already reads
+    // The node-pressure exchange is a boundary condition between two
+    // networks: the master's network ends at the master groups, and the
+    // pressure it computes there is the fixed pressure of the corresponding
+    // slave group's node in the *slave's* network.  A slave that declares no
+    // network has no such node, and its wells keep the THP limits from their
+    // own WCONPROD records -- so the messages are still received (the
+    // protocol is unchanged) but nothing is imposed.
+    //
+    // TODO (follow-up PR): when the slave group is a fixed-pressure node with
+    // branches below it, the pressure belongs at the top of the slave's own
+    // network solve rather than directly on the wells; the slave's branches
+    // and their pressure losses then decide the wells' THP limits.  The
+    // plumbing exists -- the solver already reads
     // Network::Node::terminal_pressure() when walking up branches -- but
     // surfacing the master-sent value into the solver needs a runtime
     // override path (e.g. a fixed-pressure override map on
     // BlackoilWellModelNetwork) so we do not mutate the parsed Schedule
-    // network at runtime.  Until then, decks where the slave group is a
-    // fixed-pressure node in the slave's extended network are not handled
-    // correctly.
+    // network at runtime.  Until then the direct write below is exact only
+    // for a slave group that is itself a well group, where there is no
+    // branch between the node and the wells.
     const auto& pressures = rescoup_slave.masterGroupNodePressures();
     if (pressures.empty()) return;
     const auto& summary_state = this->well_model_.summaryState();
     auto& well_state = this->wellState();
     for (auto& well : this->wellContainer()) {
         if (!well->isProducer() || !well->wellEcl().predictionMode()) continue;
-        const auto it = pressures.find(well->wellEcl().groupName());
+        const auto& group_name = well->wellEcl().groupName();
+        const auto it = pressures.find(group_name);
         if (it == pressures.end()) continue;
+        if (!this->slaveGroupIsFixedPressureNodeInOwnNetwork_(group_name)) continue;
         well->setDynamicThpLimit(it->second);
         auto& ws = well_state[well->indexOfWell()];
         if (ws.production_cmode == Well::ProducerCMode::THP) {
@@ -404,6 +413,18 @@ masterNetworkHasMasterGroupLeavesForSlave_(std::size_t slave_idx) const
     return false;
 }
 
+template<typename TypeTag>
+bool
+BlackoilWellModelRescoup<TypeTag>::
+slaveGroupIsFixedPressureNodeInOwnNetwork_(const std::string& group_name) const
+{
+    const int episodeIdx = this->simulator_.episodeIndex();
+    const auto& network = this->schedule()[episodeIdx].network();
+    if (!network.active()) return false;
+    if (!network.has_node(group_name)) return false;
+    return network.node(group_name).terminal_pressure().has_value();
+}
+
 
 } // namespace Opm
 
diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp
index 6a8946e1ade..2359934cf2f 100644
--- a/opm/simulators/wells/BlackoilWellModel_impl.hpp
+++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp
@@ -1181,6 +1181,8 @@ namespace Opm {
         // Timestep initialization: should run once at the start of each timestep.
         if (iterCtx.needsTimestepInit() && this->wellsActive()) {
             OPM_TIMEBLOCK(firstIterationAssemble);
+            Dune::Timer prepareTimer;
+            prepareTimer.start();
             // try-catch is needed here as updateWellControls
             // contains global communication and has either to
             // be reached by all processes or all need to abort
@@ -1193,12 +1195,20 @@ namespace Opm {
             OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
                                            "assemble() failed during well initialization: ",
                                            this->terminal_output_, grid().comm());
+            last_report_.well_facility_time += prepareTimer.stop();
         }
 
+        Dune::Timer networkTimer;
+        networkTimer.start();
         const bool well_group_control_changed = updateWellControlsAndNetwork(
                             /*mandatory_network_balance=*/false,
                             dt,
                             local_deferredLogger);
+        // the control/network update is part of the facility calculations,
+        // together with the time step preparation of the wells above
+        const double network_time = networkTimer.stop();
+        last_report_.well_control_network_time += network_time;
+        last_report_.well_facility_time += network_time;
 
         // even when there is no wells active, the network nodal pressure still need to be updated through updateWellControlsAndNetwork()
         // but there is no need to assemble the well equations
@@ -1210,6 +1220,11 @@ namespace Opm {
         // Pre-compute cell rates to we don't have to do this for every cell during linearization...
         updateCellRates();
 
+        // collect the statistics for the standalone well solves done since
+        // the last assembly (prepareTimeStep, well testing and well
+        // potential calculations)
+        collectWellSolveStats();
+
         // if group or well control changes we don't consider the
         // case converged
         last_report_.well_group_control_changed = well_group_control_changed;
@@ -1219,6 +1234,26 @@ namespace Opm {
 
 
 
+    template<typename TypeTag>
+    void
+    BlackoilWellModel<TypeTag>::
+    collectWellSolveStats()
+    {
+        for (const auto& well : well_container_) {
+            const auto& stats = well->solveStats();
+            last_report_.well_solve_time += stats.solve_time;
+            last_report_.well_potential_solve_time += stats.potential_solve_time;
+            last_report_.well_solve_assemble_time += stats.assemble_time;
+            last_report_.well_solve_linear_solve_time += stats.linear_solve_time;
+            last_report_.total_well_iterations += stats.iterations;
+            last_report_.total_well_potential_iterations += stats.potential_iterations;
+            well->resetSolveStats();
+        }
+    }
+
+
+
+
     template<typename TypeTag>
     bool
     BlackoilWellModel<TypeTag>::
@@ -1275,6 +1310,10 @@ namespace Opm {
         if (this->isRescoupMasterCoupledNetworkIteration_()) {
             this->sendSlaveNetworkLoopTerminationSignal_();
         }
+        // only count the balancing iterations when a network is active
+        if (this->schedule()[simulator_.episodeIndex()].network().active()) {
+            last_report_.total_network_iterations += network_update_iteration;
+        }
         return well_group_control_changed;
     }
 
@@ -1298,16 +1337,26 @@ namespace Opm {
             this->rescoupHelper_.receiveMasterGroupNodePressuresFromMaster();
         }
 #endif
+        Dune::Timer groupTimer;
+        groupTimer.start();
         this->updateAndCommunicateGroupData(reportStepIdx, /*update_wellgrouptarget*/ true);
         // We need to call updateWellControls before we update the network as
         // network updates are only done on thp controlled wells.
         // Note that well controls are allowed to change during updateNetwork
         // and in prepareWellsBeforeAssembling during well solves.
         bool well_group_control_changed = updateWellControls(local_deferredLogger);
+        last_report_.group_control_time += groupTimer.stop();
+
+        Dune::Timer networkBalanceTimer;
+        networkBalanceTimer.start();
         const auto [more_inner_network_update, network_imbalance] =
                 this->network_.update(mandatory_network_balance,
                                       local_deferredLogger,
                                       relax_network_tolerance);
+        // only account for the time when a network is active
+        if (this->schedule()[reportStepIdx].network().active()) {
+            last_report_.network_balance_time += networkBalanceTimer.stop();
+        }
 #ifdef RESERVOIR_COUPLING_ENABLED
         if (this->isReservoirCouplingMaster()) {
             this->rescoupHelper_.maybeExchangeNetworkOuterIterationWithSlaves(more_inner_network_update);
@@ -1322,6 +1371,8 @@ namespace Opm {
                 // the network balancing
                 const bool updatePotentials = (this->network_.shouldBalance(reportStepIdx) ||
                                                mandatory_network_balance);
+                Dune::Timer gliftTimer;
+                gliftTimer.start();
                 alq_updated = gaslift_.maybeDoGasLiftOptimize(simulator_,
                                                           well_container_,
                                                           this->network_.nodePressures(),
@@ -1329,8 +1380,17 @@ namespace Opm {
                                                           this->wellState(),
                                                           this->groupState(),
                                                           local_deferredLogger);
+                // only account for the time when gas lift optimization is in use
+                if (this->schedule().glo(reportStepIdx).active()) {
+                    last_report_.gaslift_time += gliftTimer.stop();
+                }
             }
+            // the well solves triggered by control changes are part of the
+            // control/network update, keep track of the time separately
+            Dune::Timer prepareWellsTimer;
+            prepareWellsTimer.start();
             prepareWellsBeforeAssembling(dt);
+            last_report_.control_well_solve_time += prepareWellsTimer.stop();
         }
         OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
                                        "updateWellControlsAndNetworkIteration() failed: ",
@@ -1340,12 +1400,15 @@ namespace Opm {
         if (alq_updated || BlackoilWellModelGuideRates(*this).
                               guideRateUpdateIsNeeded(reportStepIdx)) {
             const double simulationTime = simulator_.time();
+            Dune::Timer guideRateTimer;
+            guideRateTimer.start();
             // NOTE: For reservoir coupling: Slave group potentials are only communicated
             //    at the start of the time step, see beginTimeStep(). Here, we assume those
             //    potentials remain unchanged during the time step when updating guide rates below.
             this->guide_rate_handler_.updateGuideRates(
                 reportStepIdx, simulationTime, this->wellState(), this->groupState()
             );
+            last_report_.group_control_time += guideRateTimer.stop();
         }
         // we need to re-iterate the network when the well group controls changed or gaslift/alq is changed or
         // the inner iterations are did not converge
diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp
index 70e020b50fc..b29f1ed5d3b 100644
--- a/opm/simulators/wells/MultisegmentWell_impl.hpp
+++ b/opm/simulators/wells/MultisegmentWell_impl.hpp
@@ -309,6 +309,10 @@ namespace Opm
             return;
         }
 
+        // attribute the well solves done below (also those done by well
+        // copies) to the well potential calculations in the solve statistics
+        const auto potential_scope = this->potentialCalculationScope();
+
         debug_cost_counter_ = 0;
         bool converged_implicit = false;
         if (this->param_.local_well_solver_control_switching_) {
@@ -622,7 +626,11 @@ namespace Opm
         // We assemble the well equations, then we check the convergence,
         // which is why we do not put the assembleWellEq here.
         try{
-            const BVectorWell dx_well = this->linSys_.solve();
+            BVectorWell dx_well;
+            {
+                const auto linear_solve_timer = this->solveLinearSolveTimer();
+                dx_well = this->linSys_.solve();
+            }
             updateWellState(simulator, dx_well, groupStateHelper, well_state);
         }
         catch(const NumericalProblem& exp) {
@@ -1557,6 +1565,7 @@ namespace Opm
         std::vector<std::vector<Scalar> > residual_history;
         std::vector<Scalar> measure_history;
         int it = 0;
+        const auto solve_scope = this->solveScope(it);
         // relaxation factor
         Scalar relaxation_factor = 1.;
         bool converged = false;
@@ -1616,7 +1625,10 @@ namespace Opm
 
             BVectorWell dx_well;
             try{
-                dx_well = this->linSys_.solve();
+                {
+                    const auto linear_solve_timer = this->solveLinearSolveTimer();
+                    dx_well = this->linSys_.solve();
+                }
                 updateWellState(simulator, dx_well, groupStateHelper, well_state, relaxation_factor);
                 if constexpr (has_energy) {
                     // segment fluid state is only consumed by the energy equation
@@ -1695,6 +1707,7 @@ namespace Opm
         std::vector<std::vector<Scalar> > residual_history;
         std::vector<Scalar> measure_history;
         int it = 0;
+        const auto solve_scope = this->solveScope(it);
         // relaxation factor
         Scalar relaxation_factor = 1.;
         bool converged = false;
@@ -1809,7 +1822,11 @@ namespace Opm
                 }
             }
             try{
-                const BVectorWell dx_well = this->linSys_.solve();
+                BVectorWell dx_well;
+                {
+                    const auto linear_solve_timer = this->solveLinearSolveTimer();
+                    dx_well = this->linSys_.solve();
+                }
                 updateWellState(simulator, dx_well, groupStateHelper, well_state, relaxation_factor);
                 if constexpr (has_energy) {
                     // segment fluid state is only consumed by the energy equation
@@ -1870,6 +1887,8 @@ namespace Opm
     {
         if (!this->isOperableAndSolvable() && !this->wellIsStopped()) return;
 
+        const auto assemble_timer = this->solveAssembleTimer();
+
         auto& deferred_logger = groupStateHelper.deferredLogger();
 
         // update the upwinding segments
diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp
index de9084e5a4b..529343ef7c9 100644
--- a/opm/simulators/wells/StandardWell_impl.hpp
+++ b/opm/simulators/wells/StandardWell_impl.hpp
@@ -351,6 +351,8 @@ namespace Opm
         // for example, the matrices B C does not need to update if only_wells
         if (!this->isOperableAndSolvable() && !this->wellIsStopped()) return;
 
+        const auto assemble_timer = this->solveAssembleTimer();
+
         // clear all entries
         this->linSys_.clear();
 
@@ -1429,7 +1431,10 @@ namespace Opm
         // which is why we do not put the assembleWellEq here.
         BVectorWell dx_well(1);
         dx_well[0].resize(this->primary_variables_.numWellEq());
-        this->linSys_.solve( dx_well);
+        {
+            const auto linear_solve_timer = this->solveLinearSolveTimer();
+            this->linSys_.solve( dx_well);
+        }
 
         updateWellState(simulator, dx_well, groupStateHelper, well_state);
     }
@@ -1810,6 +1815,10 @@ namespace Opm
             return;
         }
 
+        // attribute the well solves done below (also those done by well
+        // copies) to the well potential calculations in the solve statistics
+        const auto potential_scope = this->potentialCalculationScope();
+
         bool converged_implicit = false;
         // for newly opened wells we dont compute the potentials implicit
         // group controlled wells with defaulted guiderates will have zero targets as
@@ -2404,6 +2413,7 @@ namespace Opm
 
         const int max_iter = this->param_.max_inner_iter_wells_;
         int it = 0;
+        const auto solve_scope = this->solveScope(it);
         bool converged;
         bool relax_convergence = false;
         this->regularize_ = false;
@@ -2471,6 +2481,7 @@ namespace Opm
 
         const int max_iter = this->param_.max_inner_iter_wells_;
         int it = 0;
+        const auto solve_scope = this->solveScope(it);
         bool converged = false;
         bool relax_convergence = false;
         this->regularize_ = false;
diff --git a/opm/simulators/wells/WellInterfaceGeneric.hpp b/opm/simulators/wells/WellInterfaceGeneric.hpp
index a8123741db1..00d4c2ed2e3 100644
--- a/opm/simulators/wells/WellInterfaceGeneric.hpp
+++ b/opm/simulators/wells/WellInterfaceGeneric.hpp
@@ -24,12 +24,15 @@
 #ifndef OPM_WELLINTERFACE_GENERIC_HEADER_INCLUDED
 #define OPM_WELLINTERFACE_GENERIC_HEADER_INCLUDED
 
+#include <dune/common/timer.hh>
+
 #include <opm/input/eclipse/Schedule/Well/Well.hpp>
 #include <opm/simulators/flow/BlackoilModelParameters.hpp>
 #include <opm/simulators/wells/RuntimePerforation.hpp>
 
 #include <ctime>
 #include <map>
+#include <memory>
 #include <optional>
 #include <string>
 #include <vector>
@@ -243,7 +246,107 @@ class WellInterfaceGeneric {
 
     void addPerforations(const std::vector<RuntimePerforation>& perfs);
 
+    /// Statistics for the separate (standalone) well solves, i.e. the local
+    /// nonlinear solves of the single well equations done outside the global
+    /// linearizations (prepareTimeStep, well testing) and during the well
+    /// potential calculations.
+    struct WellSolveStats
+    {
+        double solve_time{0.0};            //!< time used by normal well solves
+        double potential_solve_time{0.0};  //!< time used by well potential solves
+        double assemble_time{0.0};         //!< assembly time in the well solves
+        double linear_solve_time{0.0};     //!< linear solve time in the well solves
+        unsigned int iterations{0};            //!< well iterations in normal solves
+        unsigned int potential_iterations{0};  //!< well iterations in potential solves
+    };
+
+    const WellSolveStats& solveStats() const { return *solve_stats_; }
+    void resetSolveStats() const { *solve_stats_ = WellSolveStats{}; }
+
 protected:
+    /// RAII helper marking a standalone well solve. Accumulates the elapsed
+    /// time and the number of well iterations into the solve statistics on
+    /// scope exit, attributing them to the well potential calculations if
+    /// such a calculation is in progress.
+    class WellSolveScope
+    {
+    public:
+        WellSolveScope(const WellInterfaceGeneric& well, const int& iterations)
+            : well_(well), iterations_(iterations)
+        {
+            well_.in_well_solve_ = true;
+            timer_.start();
+        }
+        ~WellSolveScope()
+        {
+            well_.in_well_solve_ = false;
+            auto& stats = *well_.solve_stats_;
+            if (well_.computing_well_potentials_) {
+                stats.potential_solve_time += timer_.stop();
+                stats.potential_iterations += iterations_;
+            } else {
+                stats.solve_time += timer_.stop();
+                stats.iterations += iterations_;
+            }
+        }
+    private:
+        const WellInterfaceGeneric& well_;
+        const int& iterations_;
+        Dune::Timer timer_{};
+    };
+
+    /// RAII helper accumulating the elapsed time into the given part of the
+    /// solve statistics, but only when inside a standalone well solve.
+    class SolvePartTimer
+    {
+    public:
+        SolvePartTimer(const WellInterfaceGeneric& well, double& target)
+            : well_(well), target_(target)
+        {
+            timer_.start();
+        }
+        ~SolvePartTimer()
+        {
+            if (well_.in_well_solve_) {
+                target_ += timer_.stop();
+            }
+        }
+    private:
+        const WellInterfaceGeneric& well_;
+        double& target_;
+        Dune::Timer timer_{};
+    };
+
+    /// RAII helper marking that a well potential calculation is in progress.
+    class PotentialCalculationScope
+    {
+    public:
+        explicit PotentialCalculationScope(const WellInterfaceGeneric& well)
+            : well_(well), previous_(well.computing_well_potentials_)
+        {
+            well_.computing_well_potentials_ = true;
+        }
+        ~PotentialCalculationScope()
+        {
+            well_.computing_well_potentials_ = previous_;
+        }
+    private:
+        const WellInterfaceGeneric& well_;
+        bool previous_;
+    };
+
+    WellSolveScope solveScope(const int& iterations) const
+    { return WellSolveScope(*this, iterations); }
+
+    SolvePartTimer solveAssembleTimer() const
+    { return SolvePartTimer(*this, solve_stats_->assemble_time); }
+
+    SolvePartTimer solveLinearSolveTimer() const
+    { return SolvePartTimer(*this, solve_stats_->linear_solve_time); }
+
+    PotentialCalculationScope potentialCalculationScope() const
+    { return PotentialCalculationScope(*this); }
+
     bool getAllowCrossFlow() const;
 
     Scalar wmicrobes_() const;
@@ -441,6 +544,14 @@ class WellInterfaceGeneric {
     std::vector<std::string> well_control_log_;
 
     bool changed_to_open_this_step_ = true;
+
+    // statistics for standalone well solves; shared so that well copies used
+    // during e.g. potential calculations contribute to the same statistics
+    std::shared_ptr<WellSolveStats> solve_stats_ = std::make_shared<WellSolveStats>();
+    // true while a well potential calculation is in progress
+    mutable bool computing_well_potentials_ = false;
+    // true while a standalone well solve is in progress
+    mutable bool in_well_solve_ = false;
 };
 
 }
diff --git a/opm/simulators/wells/rescoup/RescoupConstraintsCalculator.cpp b/opm/simulators/wells/rescoup/RescoupConstraintsCalculator.cpp
index e965916701f..55f99e60fa3 100644
--- a/opm/simulators/wells/rescoup/RescoupConstraintsCalculator.cpp
+++ b/opm/simulators/wells/rescoup/RescoupConstraintsCalculator.cpp
@@ -211,7 +211,6 @@ calculateMasterGroupConstraintsAndSendToSlaves()
         this->well_model_,
         this->group_state_helper_
     };
-    this->restoreMasterGroupControlsFromSchedule_();
     // Reset effective-GCW entries from the previous sync step to GCW=1 before
     // excludeInactiveSlaveMasterGroupsFromDistribution_() (and later the
     // Phase 2 cap) repopulate the 0-entries for this step.
@@ -406,9 +405,7 @@ capAndRedistributeProductionTargets_(
     // Step 4: switch ALL master groups to individual control with their
     // allocated targets. I.e. the master completes its own time step, assuming
     // that the production and injection rates of the slave groups remain constant
-    // over the time step. The controls are reset from schedule at the start of the
-    // next sync step (before Phase 1, see calculateMasterGroupConstraintsAndSendToSlaves()
-    // above) so guide rate distribution works correctly.
+    // over the time step.
     for (std::size_t slave_idx = 0; slave_idx < num_slaves; ++slave_idx) {
         const auto& master_groups = rescoup_master.getMasterGroupNamesForSlave(slave_idx);
         for (auto& pc : all_production_constraints[slave_idx]) {
@@ -499,43 +496,6 @@ potentialForProductionCmode_(
     }
 }
 
-// Restore master groups' control modes to their original GCONPROD values
-// (from the Schedule) so they participate correctly in guide rate
-// distribution during Phase 1 (see calculateMasterGroupConstraintsAndSendToSlaves()).
-// The previous sync step's final step may have switched them to individual control.
-//
-// Reading from the Schedule (rather than caching the pre-sync mode) is
-// intentional:
-//   - Master group control modes are treated as read-only outside the
-//     constraint calculation invoked from sendMasterGroupConstraintsToSlaves()
-//     at BlackoilWellModel_impl.hpp:555 (see also the master-group guard at
-//     the top of BlackoilWellModel::updateGroupControls).
-//   - Modes are only changed inside calculateMasterGroupConstraintsAndSendToSlaves(),
-//     which runs once per sync step.  A sync step may be shorter than a
-//     report step, so multiple sync steps can occur within the same report
-//     step.
-//   - The Schedule's cmode is constant within a report step, so reading it
-//     here is equivalent to restoring to the pre-sync state for any sync
-//     step within the report step.  At a report-step boundary, a new
-//     GCONPROD record is picked up automatically.
-template <class Scalar, class IndexTraits>
-void
-RescoupConstraintsCalculator<Scalar, IndexTraits>::
-restoreMasterGroupControlsFromSchedule_()
-{
-    auto& rescoup_master = this->reservoir_coupling_master_;
-    const auto num_slaves = rescoup_master.numSlaves();
-    for (std::size_t slave_idx = 0; slave_idx < num_slaves; ++slave_idx) {
-        const auto& master_groups = rescoup_master.getMasterGroupNamesForSlave(slave_idx);
-        for (const auto& group_name : master_groups) {
-            const auto& group = this->schedule_.getGroup(group_name, this->report_step_idx_);
-            const auto original_cmode = group.productionProperties().cmode;
-            this->group_state_helper_.groupState().production_control(
-                group_name, original_cmode);
-        }
-    }
-}
-
 template <class Scalar, class IndexTraits>
 void
 RescoupConstraintsCalculator<Scalar, IndexTraits>::
diff --git a/tests/test_convergenceoutputconfiguration.cpp b/tests/test_convergenceoutputconfiguration.cpp
index 5547ac449ba..a80e479b3a2 100644
--- a/tests/test_convergenceoutputconfiguration.cpp
+++ b/tests/test_convergenceoutputconfiguration.cpp
@@ -84,6 +84,44 @@ BOOST_AUTO_TEST_CASE(Iterations_Alias)
                         "option value must activate Steps option");
 }
 
+BOOST_AUTO_TEST_CASE(Performance)
+{
+    const auto config = Opm::ConvergenceOutputConfiguration{"performance"};
+    BOOST_CHECK_MESSAGE(config.any(),
+                        "Configuration object with supported "
+                        "option value must activate option");
+    BOOST_CHECK_MESSAGE(config.want(Opm::ConvergenceOutputConfiguration::Option::Performance),
+                        "Configuration object with \"performance\" "
+                        "option value must activate Performance option");
+}
+
+BOOST_AUTO_TEST_CASE(Performance_Alias)
+{
+    const auto config = Opm::ConvergenceOutputConfiguration{"perf"};
+    BOOST_CHECK_MESSAGE(config.any(),
+                        "Configuration object with supported "
+                        "option value must activate option");
+    BOOST_CHECK_MESSAGE(config.want(Opm::ConvergenceOutputConfiguration::Option::Performance),
+                        "Configuration object with \"perf\" "
+                        "option value must activate Performance option");
+}
+
+BOOST_AUTO_TEST_CASE(Performance_Combination)
+{
+    const auto both = Opm::ConvergenceOutputConfiguration{"iterations,performance"};
+    BOOST_CHECK_MESSAGE(both.want(Opm::ConvergenceOutputConfiguration::Option::Iterations),
+                        "Combined option must activate Iterations");
+    BOOST_CHECK_MESSAGE(both.want(Opm::ConvergenceOutputConfiguration::Option::Performance),
+                        "Combined option must activate Performance");
+}
+
+BOOST_AUTO_TEST_CASE(Performance_None_Overrides)
+{
+    const auto config = Opm::ConvergenceOutputConfiguration{"performance,none"};
+    BOOST_CHECK_MESSAGE(!config.any(),
+                        "\"none\" must override the performance option");
+}
+
 BOOST_AUTO_TEST_CASE(Combinations)
 {
     const auto steps_iter = Opm::ConvergenceOutputConfiguration{"steps,iterations"};