added timing blocks to most important parts

ebos/eclgenerictracermodel.hh

@@ -29,7 +29,7 @@
 #define EWOMS_ECL_GENERIC_TRACER_MODEL_HH
 
 #include <opm/grid/common/CartesianIndexMapper.hpp>
-
+#include <opm/common/ErrorMacros.hpp>
 #include <opm/models/blackoil/blackoilmodel.hh>
 
 #include <opm/simulators/linalg/matrixblock.hh>

ebos/eclproblem.hh

@@ -972,6 +972,7 @@ public:
      */
     void beginEpisode()
     {
+        OPM_TIMEBLOCK(beginEpisode);
         // Proceed to the next report step
         auto& simulator = this->simulator();
         int episodeIdx = simulator.episodeIndex();
@@ -1035,6 +1036,7 @@ public:
      */
     void beginTimeStep()
     {
+        OPM_TIMEBLOCK(beginTimeStep);
         int episodeIdx = this->episodeIndex();
 
         this->beginTimeStep_(enableExperiments,
@@ -1056,9 +1058,11 @@ public:
 
         // the derivatives may have change
         bool invalidateIntensiveQuantities = invalidateFromMaxWaterSat || invalidateFromMinPressure || invalidateFromHyst || invalidateFromMaxOilSat;
-        if (invalidateIntensiveQuantities)
+        if (invalidateIntensiveQuantities){
+            OPM_TIMEBLOCK(beginTimeStepInvalidateIntensiveQuantities);
             this->model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
-
+        }
+        
         if constexpr (getPropValue<TypeTag, Properties::EnablePolymer>())
             updateMaxPolymerAdsorption_();
 
@@ -1074,6 +1078,7 @@ public:
      */
     void beginIteration()
     {
+        OPM_TIMEBLOCK(beginIteration);
         wellModel_.beginIteration();
         if (enableAquifers_)
             aquiferModel_.beginIteration();
@@ -1084,6 +1089,7 @@ public:
      */
     void endIteration()
     {
+        OPM_TIMEBLOCK(endIteration);
         wellModel_.endIteration();
         if (enableAquifers_)
             aquiferModel_.endIteration();
@@ -1094,6 +1100,7 @@ public:
      */
     void endTimeStep()
     {
+        OPM_TIMEBLOCK(endTimeStep);
 #ifndef NDEBUG
         if constexpr (getPropValue<TypeTag, Properties::EnableDebuggingChecks>()) {
             // in debug mode, we don't care about performance, so we check if the model does
@@ -1166,6 +1173,7 @@ public:
      */
     void endEpisode()
     {
+        OPM_TIMEBLOCK(endEpisode);
         auto& simulator = this->simulator();
         auto& schedule = simulator.vanguard().schedule();
 
@@ -1190,6 +1198,7 @@ public:
      */
     void writeOutput(bool verbose = true)
     {
+        OPM_TIMEBLOCK(writeOutput);
         // use the generic code to prepare the output fields and to
         // write the desired VTK files.
         ParentType::writeOutput(verbose);
@@ -1200,6 +1209,7 @@ public:
     }
 
     void finalizeOutput() {
+        OPM_TIMEBLOCK(finalizeOutput);
         // this will write all pending output to disk
         // to avoid corruption of output files
         eclWriter_.reset();
@@ -1460,6 +1470,7 @@ public:
         FluidState &fluidState,
         unsigned globalSpaceIdx) const
     {
+        OPM_TIMEBLOCK_LOCAL(updateRelperms);
         {
             // calculate relative permeabilities. note that we store the result into the
             // mobility_ class attribute. the division by the phase viscosity happens later.
@@ -1557,6 +1568,7 @@ public:
                   unsigned spaceIdx,
                   unsigned timeIdx) const
     {
+        OPM_TIMEBLOCK(eclProblemBoundary);
         if (!context.intersection(spaceIdx).boundary())
             return;
 
@@ -1766,6 +1778,7 @@ public:
                 unsigned globalDofIdx,
                 unsigned timeIdx) const
     {
+        OPM_TIMEBLOCK(eclProblemSource);
         rate = 0.0;
 
         wellModel_.computeTotalRatesForDof(rate, globalDofIdx);
@@ -1847,6 +1860,7 @@ public:
 
     const InitialFluidState boundaryFluidState(unsigned globalDofIdx, const int directionId) const
     {
+        OPM_TIMEBLOCK_LOCAL(boundaryFluidState);
         FaceDir::DirEnum dir = FaceDir::FromIntersectionIndex(directionId);
         const auto& dirichlet = dirichlet_(dir)[globalDofIdx];
         if(std::get<0>(dirichlet) == BCComponent::NONE)
@@ -1936,6 +1950,7 @@ public:
      */
     Scalar nextTimeStepSize() const
     {
+        OPM_TIMEBLOCK(nexTimeStepSize);
         // allow external code to do the timestepping
         if (this->nextTimeStepSize_ > 0.0)
             return this->nextTimeStepSize_;
@@ -1962,7 +1977,7 @@ public:
     template <class LhsEval>
     LhsEval rockCompPoroMultiplier(const IntensiveQuantities& intQuants, unsigned elementIdx) const
     {
-
+        OPM_TIMEBLOCK_LOCAL(rockCompPoroMultiplier);
         if (this->rockCompPoroMult_.empty() && this->rockCompPoroMultWc_.empty())
             return 1.0;
 
@@ -2002,6 +2017,7 @@ public:
     template <class LhsEval>
     LhsEval rockCompTransMultiplier(const IntensiveQuantities& intQuants, unsigned elementIdx) const
     {
+        OPM_TIMEBLOCK_LOCAL(rockCompTransMultiplier);
         if (this->rockCompTransMult_.empty() && this->rockCompTransMultWc_.empty())
             return 1.0;
 
@@ -2034,6 +2050,7 @@ public:
 
     std::pair<bool, RateVector> boundaryCondition(const unsigned int globalSpaceIdx, const int directionId)
     {
+        OPM_TIMEBLOCK_LOCAL(boundaryCondition);
         if (!nonTrivialBoundaryConditions_) {
             return { false, RateVector(0.0) };
         }
@@ -2053,6 +2070,7 @@ private:
     void updateProperty_(const std::string& failureMsg,
                          UpdateFunc func)
     {
+        OPM_TIMEBLOCK(updateProperty);
         ElementContext elemCtx(this->simulator());
         const auto& vanguard = this->simulator().vanguard();
         OPM_BEGIN_PARALLEL_TRY_CATCH();
@@ -2070,6 +2088,7 @@ private:
     // update the parameters needed for DRSDT and DRVDT
     void updateCompositionChangeLimits_()
     {
+        OPM_TIMEBLOCK(updateCompositionChangeLimits);
         // update the "last Rs" values for all elements, including the ones in the ghost
         // and overlap regions
         int episodeIdx = this->episodeIndex();
@@ -2140,6 +2159,7 @@ private:
 
     bool updateMaxOilSaturation_()
     {
+        OPM_TIMEBLOCK(updateMaxOilSaturation);
         int episodeIdx = this->episodeIndex();
 
         // we use VAPPARS
@@ -2160,6 +2180,7 @@ private:
 
     bool updateMaxWaterSaturation_()
     {
+        OPM_TIMEBLOCK(updateMaxWaterSaturation);
         // water compaction is activated in ROCKCOMP
         if (this->maxWaterSaturation_.empty())
             return false;
@@ -2178,6 +2199,7 @@ private:
 
     bool updateMinPressure_()
     {
+        OPM_TIMEBLOCK(updateMinPressure);
         // IRREVERS option is used in ROCKCOMP
         if (this->minOilPressure_.empty())
             return false;
@@ -2195,6 +2217,7 @@ private:
 
     void readMaterialParameters_()
     {
+        OPM_TIMEBLOCK(readMaterialParameters);
         const auto& simulator = this->simulator();
         const auto& vanguard = simulator.vanguard();
         const auto& eclState = vanguard.eclState();

opm/simulators/flow/BlackoilModelEbos.hpp

@@ -631,6 +631,7 @@ namespace Opm {
         /// Apply an update to the primary variables.
         void updateSolution(const BVector& dx)
         {
+            OPM_TIMEBLOCK(updateSolution);
             auto& ebosNewtonMethod = ebosSimulator_.model().newtonMethod();
             SolutionVector& solution = ebosSimulator_.model().solution(/*timeIdx=*/0);
 
@@ -642,7 +643,10 @@ namespace Opm {
                                                     // residual
 
             // if the solution is updated, the intensive quantities need to be recalculated
-            ebosSimulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
+            {
+                OPM_TIMEBLOCK(invalidateAndUpdateIntensiveQuantities);
+                ebosSimulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
+            }
         }
 
         /// Return true if output to cout is wanted.
@@ -659,6 +663,7 @@ namespace Opm {
                                                        std::vector< Scalar >& maxCoeff,
                                                        std::vector< Scalar >& B_avg)
         {
+            OPM_TIMEBLOCK(convergenceReduction);
             // Compute total pore volume (use only owned entries)
             double pvSum = pvSumLocal;
             double numAquiferPvSum = numAquiferPvSumLocal;
@@ -718,6 +723,7 @@ namespace Opm {
                                     std::vector<Scalar>& maxCoeff,
                                     std::vector<Scalar>& B_avg)
         {
+            OPM_TIMEBLOCK(localConvergenceData);
             double pvSumLocal = 0.0;
             double numAquiferPvSumLocal = 0.0;
             const auto& ebosModel = ebosSimulator_.model();
@@ -848,6 +854,7 @@ namespace Opm {
         ///        of a numerical aquifer.
         double computeCnvErrorPv(const std::vector<Scalar>& B_avg, double dt)
         {
+            OPM_TIMEBLOCK(computeCnvErrorPv);
             double errorPV{};
             const auto& ebosModel = ebosSimulator_.model();
             const auto& ebosProblem = ebosSimulator_.problem();
@@ -895,6 +902,7 @@ namespace Opm {
                                                   std::vector<Scalar>& B_avg,
                                                   std::vector<Scalar>& residual_norms)
         {
+            OPM_TIMEBLOCK(getReservoirConvergence);
             typedef std::vector< Scalar > Vector;
 
             const int numComp = numEq;
@@ -1004,13 +1012,16 @@ namespace Opm {
                                          const int iteration,
                                          std::vector<double>& residual_norms)
         {
+            OPM_TIMEBLOCK(getConvergence);
             // Get convergence reports for reservoir and wells.
             std::vector<Scalar> B_avg(numEq, 0.0);
             auto report = getReservoirConvergence(timer.simulationTimeElapsed(),
                                                   timer.currentStepLength(),
                                                   iteration, B_avg, residual_norms);
-            report += wellModel().getWellConvergence(B_avg, /*checkWellGroupControls*/report.converged());
+            {
-
+                OPM_TIMEBLOCK(getWellConvergence);
+                report += wellModel().getWellConvergence(B_avg, /*checkWellGroupControls*/report.converged());
+            }
             return report;
         }
 
@@ -1033,6 +1044,7 @@ namespace Opm {
         std::vector<std::vector<double> >
         computeFluidInPlace(const std::vector<int>& /*fipnum*/) const
         {
+            OPM_TIMEBLOCK(computeFluidInPlace);
             //assert(true)
             //return an empty vector
             std::vector<std::vector<double> > regionValues(0, std::vector<double>(0,0.0));

opm/simulators/linalg/ISTLSolverEbos.hpp

@@ -237,6 +237,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
               converged_(false),
               matrix_()
         {
+            OPM_TIMEBLOCK(istlSolverEbos);
             const bool on_io_rank = (simulator.gridView().comm().rank() == 0);
 #if HAVE_MPI
             comm_.reset( new CommunicationType( simulator_.vanguard().grid().comm() ) );
@@ -311,6 +312,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
 
         void prepare(const SparseMatrixAdapter& M, Vector& b)
         {
+            OPM_TIMEBLOCK(istlSolverEbosPrepare);        
             static bool firstcall = true;
 #if HAVE_MPI
             if (firstcall) {
@@ -370,6 +372,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
 
         bool solve(Vector& x)
         {
+            OPM_TIMEBLOCK(istlSolverEbosSolve);        
             calls_ += 1;
             // Write linear system if asked for.
             const int verbosity = prm_.get<int>("verbosity", 0);
@@ -396,6 +399,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
                                   x, result))
 #endif
             {
+                OPM_TIMEBLOCK(flexibleSolverApply);
                 assert(flexibleSolver_.solver_);
                 flexibleSolver_.solver_->apply(x, *rhs_, result);
             }
@@ -448,7 +452,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
 
         void prepareFlexibleSolver()
         {
-
+            OPM_TIMEBLOCK(prepareFlexibleSolver);        
             if (shouldCreateSolver()) {
                 std::function<Vector()> trueFunc =
                     [this]
@@ -460,7 +464,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
                     auto wellOp = std::make_unique<WellModelOperator>(simulator_.problem().wellModel());
                     flexibleSolver_.wellOperator_ = std::move(wellOp);
                 }
-
+                OPM_TIMEBLOCK(flexibleSolverCreate);
                 flexibleSolver_.create(getMatrix(),
                                        isParallel(),
                                        prm_,
@@ -470,6 +474,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
             }
             else
             {
+                OPM_TIMEBLOCK(flexibleSolverUpdate);        
                 flexibleSolver_.pre_->update();
             }
         }
@@ -527,6 +532,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
         // conservation equations, ignoring all other terms.
         Vector getTrueImpesWeights(int pressureVarIndex) const
         {
+            OPM_TIMEBLOCK(getTrueImpesWeights);
             Vector weights(rhs_->size());
             ElementContext elemCtx(simulator_);
             Amg::getTrueImpesWeights(pressureVarIndex, weights,