diff --git a/opm/autodiff/BlackoilModelBase.hpp b/opm/autodiff/BlackoilModelBase.hpp index 3f4e01e66..f32aa46fb 100644 --- a/opm/autodiff/BlackoilModelBase.hpp +++ b/opm/autodiff/BlackoilModelBase.hpp @@ -37,6 +37,7 @@ #include #include #include +#include #include @@ -175,7 +176,7 @@ namespace Opm { /// \param[in, out] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables template - IterationReport nonlinearIteration(const int iteration, + SimulatorReport nonlinearIteration(const int iteration, const SimulatorTimerInterface& timer, NonlinearSolverType& nonlinear_solver, ReservoirState& reservoir_state, @@ -194,8 +195,7 @@ namespace Opm { /// \param[in] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables /// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep - /// \return well iterations. - IterationReport + SimulatorReport assemble(const ReservoirState& reservoir_state, WellState& well_state, const bool initial_assembly); @@ -396,7 +396,7 @@ namespace Opm { assembleMassBalanceEq(const SolutionState& state); - IterationReport + SimulatorReport solveWellEq(const std::vector& mob_perfcells, const std::vector& b_perfcells, const ReservoirState& reservoir_state, diff --git a/opm/autodiff/BlackoilModelBase_impl.hpp b/opm/autodiff/BlackoilModelBase_impl.hpp index dd68db30a..42f00cf07 100644 --- a/opm/autodiff/BlackoilModelBase_impl.hpp +++ b/opm/autodiff/BlackoilModelBase_impl.hpp @@ -50,8 +50,10 @@ #include #include #include - #include + +#include + #include #include #include @@ -233,7 +235,7 @@ typedef Eigen::Array template - IterationReport + SimulatorReport BlackoilModelBase:: nonlinearIteration(const int iteration, const SimulatorTimerInterface& timer, @@ -241,6 +243,10 @@ typedef Eigen::Array - IterationReport + SimulatorReport BlackoilModelBase:: assemble(const ReservoirState& reservoir_state, WellState& well_state, @@ -710,6 +746,8 @@ typedef Eigen::Array mob_perfcells; @@ -767,7 +804,7 @@ typedef Eigen::Array cq_s; @@ -783,7 +820,7 @@ typedef Eigen::Array - IterationReport + SimulatorReport BlackoilModelBase:: solveWellEq(const std::vector& mob_perfcells, const std::vector& b_perfcells, @@ -1072,9 +1109,11 @@ typedef Eigen::Array #include +#include #include #include #include @@ -64,6 +65,7 @@ #include #include +#include #include #include @@ -247,12 +249,16 @@ namespace Opm { /// \param[in, out] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables template - IterationReport nonlinearIteration(const int iteration, + SimulatorReport nonlinearIteration(const int iteration, const SimulatorTimerInterface& timer, NonlinearSolverType& nonlinear_solver, ReservoirState& reservoir_state, WellState& well_state) { + SimulatorReport report; + Dune::Timer perfTimer; + + perfTimer.start(); if (iteration == 0) { // For each iteration we store in a vector the norms of the residual of // the mass balance for each active phase, the well flux and the well equations. @@ -264,16 +270,30 @@ namespace Opm { // reset intensive quantities cache useless other options are set // further down invalidateIntensiveQuantitiesCache_ = true; + report.total_linearizations = 1; + + try { + report += assemble(timer, iteration, reservoir_state, well_state); + report.assemble_time += perfTimer.stop(); + } + catch (...) { + report.assemble_time += perfTimer.stop(); + // todo (?): make the report an attribute of the class + throw; // continue throwing the stick + } - IterationReport iter_report = assemble(timer, iteration, reservoir_state, well_state); std::vector residual_norms; - const bool converged = getConvergence(timer, iteration,residual_norms); + perfTimer.reset(); + perfTimer.start(); + report.converged = getConvergence(timer, iteration,residual_norms); + report.update_time += perfTimer.stop(); residual_norms_history_.push_back(residual_norms); - bool must_solve = (iteration < nonlinear_solver.minIter()) || (!converged); - // don't solve if we have reached the maximum number of iteration. - must_solve = must_solve && (iteration < nonlinear_solver.maxIter()); + bool must_solve = iteration < nonlinear_solver.minIter() || !report.converged; if (must_solve) { + perfTimer.reset(); + perfTimer.start(); + report.total_newton_iterations = 1; // enable single precision for solvers when dt is smaller then 20 days //residual_.singlePrecision = (unit::convert::to(dt, unit::day) < 20.) ; @@ -284,7 +304,20 @@ namespace Opm { BVector x(nc); BVector xw(nw); - solveJacobianSystem(x, xw); + try { + solveJacobianSystem(x, xw); + report.linear_solve_time += perfTimer.stop(); + report.total_linear_iterations += linearIterationsLastSolve(); + } + catch (...) { + report.linear_solve_time += perfTimer.stop(); + report.total_linear_iterations += linearIterationsLastSolve(); + // todo (?): make the report an attribute of the class + throw; // re-throw up + } + + perfTimer.reset(); + perfTimer.start(); // Stabilize the nonlinear update. bool isOscillate = false; @@ -301,24 +334,20 @@ namespace Opm { } nonlinear_solver.stabilizeNonlinearUpdate(x, dx_old_, current_relaxation_); - // Apply the update, applying model-dependent - // limitations and chopping of the update. + // Apply the update, with considering model-dependent limitations and + // chopping of the update. updateState(x,reservoir_state); wellModel().updateWellState(xw, well_state); - - // since the solution was changed, the cache for the intensive quantities are invalid - // ebosSimulator_.model().invalidateIntensiveQuantitiesCache(/*timeIdx=*/0); + report.update_time += perfTimer.stop(); } - - if( converged && (iteration >= nonlinear_solver.minIter()) ) - { - // in case of convergence we do not need to reset intensive quantities + else { + // if the solution is not updated, we do not need to recalculate the + // intensive quantities in the next iteration. + assert(report.converged); invalidateIntensiveQuantitiesCache_ = false ; } - const bool failed = false; // Not needed in this model. - const int linear_iters = must_solve ? linearIterationsLastSolve() : 0; - return IterationReport{ failed, converged, linear_iters, iter_report.well_iterations }; + return report; } void printIf(int c, double x, double y, double eps, std::string type) { @@ -346,30 +375,31 @@ namespace Opm { /// \param[in] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables /// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep - IterationReport assemble(const SimulatorTimerInterface& timer, + SimulatorReport assemble(const SimulatorTimerInterface& timer, const int iterationIdx, const ReservoirState& reservoir_state, WellState& well_state) { using namespace Opm::AutoDiffGrid; + SimulatorReport report; + // -------- Mass balance equations -------- assembleMassBalanceEq(timer, iterationIdx, reservoir_state); // -------- Well equations ---------- double dt = timer.currentStepLength(); - IterationReport iter_report; try { - iter_report = wellModel().assemble(ebosSimulator_, iterationIdx, dt, well_state); + report = wellModel().assemble(ebosSimulator_, iterationIdx, dt, well_state); } catch ( const Dune::FMatrixError& e ) { OPM_THROW(Opm::NumericalProblem,"Well equation did not converge"); } - return iter_report; + return report; } @@ -439,7 +469,6 @@ namespace Opm { wellModel().applyScaleAdd(alpha, x, y); } - /// Solve the Jacobian system Jx = r where J is the Jacobian and /// r is the residual. void solveJacobianSystem(BVector& x, BVector& xw) const diff --git a/opm/autodiff/BlackoilMultiSegmentModel.hpp b/opm/autodiff/BlackoilMultiSegmentModel.hpp index d55eeea65..41b6187e7 100644 --- a/opm/autodiff/BlackoilMultiSegmentModel.hpp +++ b/opm/autodiff/BlackoilMultiSegmentModel.hpp @@ -106,7 +106,7 @@ namespace Opm { /// \param[in] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables /// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep - IterationReport + SimulatorReport assemble(const ReservoirState& reservoir_state, WellState& well_state, const bool initial_assembly); @@ -166,7 +166,7 @@ namespace Opm { const MultisegmentWells::MultisegmentWellOps& msWellOps() const { return well_model_.wellOps(); } - IterationReport + SimulatorReport solveWellEq(const std::vector& mob_perfcells, const std::vector& b_perfcells, const ReservoirState& reservoir_state, diff --git a/opm/autodiff/BlackoilMultiSegmentModel_impl.hpp b/opm/autodiff/BlackoilMultiSegmentModel_impl.hpp index f6558a7eb..f1ea6cb19 100644 --- a/opm/autodiff/BlackoilMultiSegmentModel_impl.hpp +++ b/opm/autodiff/BlackoilMultiSegmentModel_impl.hpp @@ -127,7 +127,7 @@ namespace Opm { template - IterationReport + SimulatorReport BlackoilMultiSegmentModel:: assemble(const ReservoirState& reservoir_state, WellState& well_state, @@ -188,9 +188,9 @@ namespace Opm { asImpl().assembleMassBalanceEq(state); // -------- Well equations ---------- - IterationReport iter_report = {false, false, 0, 0}; if ( ! wellsActive() ) { - return iter_report; + SimulatorReport report; + return report; } wellModel().computeSegmentFluidProperties(state); @@ -200,10 +200,11 @@ namespace Opm { std::vector mob_perfcells; std::vector b_perfcells; + SimulatorReport report; wellModel().extractWellPerfProperties(state, sd_.rq, mob_perfcells, b_perfcells); if (param_.solve_welleq_initially_ && initial_assembly) { // solve the well equations as a pre-processing step - iter_report = asImpl().solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state); + report = asImpl().solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state); } // the perforation flux here are different @@ -215,7 +216,7 @@ namespace Opm { wellModel().addWellFluxEq(cq_s, state, residual_); asImpl().addWellContributionToMassBalanceEq(cq_s, state, well_state); wellModel().addWellControlEq(state, well_state, aliveWells, residual_); - return iter_report; + return report; } @@ -223,16 +224,16 @@ namespace Opm { template - IterationReport + SimulatorReport BlackoilMultiSegmentModel::solveWellEq(const std::vector& mob_perfcells, const std::vector& b_perfcells, const ReservoirState& reservoir_state, SolutionState& state, WellState& well_state) { - IterationReport iter_report = Base::solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state); + SimulatorReport report = Base::solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state); - if (iter_report.converged) { + if (report.converged) { // We must now update the state.segp and state.segqs members, // that the base version does not know about. const int np = numPhases(); @@ -261,7 +262,7 @@ namespace Opm { asImpl().computeWellConnectionPressures(state, well_state); } - return iter_report; + return report; } diff --git a/opm/autodiff/BlackoilPressureModel.hpp b/opm/autodiff/BlackoilPressureModel.hpp index 3aab0a2a0..8ec744f7a 100644 --- a/opm/autodiff/BlackoilPressureModel.hpp +++ b/opm/autodiff/BlackoilPressureModel.hpp @@ -149,13 +149,14 @@ namespace Opm { - IterationReport + SimulatorReport assemble(const ReservoirState& reservoir_state, WellState& well_state, const bool initial_assembly) { - IterationReport iter_report = Base::assemble(reservoir_state, well_state, initial_assembly); + SimulatorReport report; + report += Base::assemble(reservoir_state, well_state, initial_assembly); if (initial_assembly) { } @@ -186,7 +187,8 @@ namespace Opm { assert(int(well_state.perfRates().size()) == wflux.size()); std::copy_n(wflux.data(), wflux.size(), well_state.perfRates().begin()); } - return iter_report; + + return report; } diff --git a/opm/autodiff/BlackoilSequentialModel.hpp b/opm/autodiff/BlackoilSequentialModel.hpp index 4c1ec6f8e..3aa23a3d9 100644 --- a/opm/autodiff/BlackoilSequentialModel.hpp +++ b/opm/autodiff/BlackoilSequentialModel.hpp @@ -123,7 +123,7 @@ namespace Opm { /// \param[in, out] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables template - IterationReport nonlinearIteration(const int iteration, + SimulatorReport nonlinearIteration(const int iteration, const SimulatorTimerInterface& timer, NonlinearSolverType& /* nonlinear_solver */, ReservoirState& reservoir_state, @@ -140,7 +140,8 @@ namespace Opm { OpmLog::info("Solving the pressure equation."); } ReservoirState initial_state = reservoir_state; - const int pressure_liniter = pressure_solver_.step(timer, reservoir_state, well_state); + const SimulatorReport pressure_report = pressure_solver_.step(timer, reservoir_state, well_state); + const int pressure_liniter = pressure_report.total_linear_iterations; if (pressure_liniter == -1) { OPM_THROW(std::runtime_error, "Pressure solver failed to converge."); } @@ -149,13 +150,17 @@ namespace Opm { if (terminalOutputEnabled()) { OpmLog::info("Solving the transport equations."); } - const int transport_liniter = transport_solver_.step(timer, initial_state, well_state, reservoir_state, well_state); + const SimulatorReport transport_report = transport_solver_.step(timer, initial_state, well_state, reservoir_state, well_state); + const int transport_liniter = transport_report.total_linear_iterations; if (transport_liniter == -1) { OPM_THROW(std::runtime_error, "Transport solver failed to converge."); } // Report and return. - return IterationReport { false, true, pressure_liniter + transport_liniter, 0 }; + SimulatorReport report; + report.converged = true; + report.total_linear_iterations = pressure_liniter + transport_liniter; + return report; } else { // Iterate to fully implicit solution. // This call is just for a single iteration (one pressure and one transport solve), @@ -168,7 +173,8 @@ namespace Opm { if (terminalOutputEnabled()) { OpmLog::info("Solving the pressure equation."); } - const int pressure_liniter = pressure_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state); + const SimulatorReport& pressure_report = pressure_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state); + const int pressure_liniter = pressure_report.total_linear_iterations; if (pressure_liniter == -1) { OPM_THROW(std::runtime_error, "Pressure solver failed to converge."); } @@ -177,14 +183,16 @@ namespace Opm { if (terminalOutputEnabled()) { OpmLog::info("Solving the transport equations."); } - const int transport_liniter = transport_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state); + const SimulatorReport& transport_report = transport_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state); + const int transport_liniter = transport_report.total_linear_iterations; if (transport_liniter == -1) { OPM_THROW(std::runtime_error, "Transport solver failed to converge."); } - // Report and return. - const bool converged = iteration >= 3; // TODO: replace this with a proper convergence check - return IterationReport { false, converged, pressure_liniter + transport_liniter, 0 }; + SimulatorReport report; + report.converged = iteration >= 3; // TODO: replace this with a proper convergence check + report.total_linear_iterations = pressure_liniter + transport_liniter; + return report; } } diff --git a/opm/autodiff/BlackoilTransportModel.hpp b/opm/autodiff/BlackoilTransportModel.hpp index fa9ad9a1a..11fca0914 100644 --- a/opm/autodiff/BlackoilTransportModel.hpp +++ b/opm/autodiff/BlackoilTransportModel.hpp @@ -83,14 +83,15 @@ namespace Opm { asImpl().makeConstantState(state0_); } - IterationReport + SimulatorReport assemble(const ReservoirState& reservoir_state, WellState& well_state, const bool initial_assembly) { - using namespace Opm::AutoDiffGrid; + SimulatorReport report; + // If we have VFP tables, we need the well connection // pressures for the "simple" hydrostatic correction // between well depth and vfp table depth. @@ -125,9 +126,8 @@ namespace Opm { asImpl().assembleMassBalanceEq(state); // -------- Well equations ---------- - IterationReport iter_report = {false, false, 0, 0}; if ( ! wellsActive() ) { - return iter_report; + return report; } std::vector mob_perfcells; @@ -135,7 +135,7 @@ namespace Opm { asImpl().wellModel().extractWellPerfProperties(state, sd_.rq, mob_perfcells, b_perfcells); if (param_.solve_welleq_initially_ && initial_assembly) { // solve the well equations as a pre-processing step - iter_report = asImpl().solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state); + report += asImpl().solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state); } V aliveWells; std::vector cq_s; @@ -154,7 +154,8 @@ namespace Opm { asImpl().makeConstantState(state0); asImpl().wellModel().computeWellPotentials(mob_perfcells, b_perfcells, state0, well_state); } - return iter_report; + + return report; } diff --git a/opm/autodiff/NonlinearSolver.hpp b/opm/autodiff/NonlinearSolver.hpp index fb414017b..2dfbab5f4 100644 --- a/opm/autodiff/NonlinearSolver.hpp +++ b/opm/autodiff/NonlinearSolver.hpp @@ -22,6 +22,7 @@ #define OPM_NONLINEARSOLVER_HEADER_INCLUDED #include +#include #include #include #include @@ -87,8 +88,7 @@ namespace Opm { /// \param[in] timer simulation timer /// \param[in, out] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables - /// \return number of linear iterations used - int + SimulatorReport step(const SimulatorTimerInterface& timer, ReservoirState& reservoir_state, WellState& well_state); @@ -103,7 +103,7 @@ namespace Opm { /// \param[in, out] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables /// \return number of linear iterations used - int + SimulatorReport step(const SimulatorTimerInterface& timer, const ReservoirState& initial_reservoir_state, const WellState& initial_well_state, @@ -178,7 +178,7 @@ namespace Opm { double relaxRelTol() const { return param_.relax_rel_tol_; } /// The maximum number of nonlinear iterations allowed. - double maxIter() const { return param_.max_iter_; } + int maxIter() const { return param_.max_iter_; } /// The minimum number of nonlinear iterations allowed. double minIter() const { return param_.min_iter_; } diff --git a/opm/autodiff/NonlinearSolver_impl.hpp b/opm/autodiff/NonlinearSolver_impl.hpp index 41d999591..3d92a8aac 100644 --- a/opm/autodiff/NonlinearSolver_impl.hpp +++ b/opm/autodiff/NonlinearSolver_impl.hpp @@ -102,7 +102,7 @@ namespace Opm } template - int + SimulatorReport NonlinearSolver:: step(const SimulatorTimerInterface& timer, ReservoirState& reservoir_state, @@ -114,7 +114,7 @@ namespace Opm template - int + SimulatorReport NonlinearSolver:: step(const SimulatorTimerInterface& timer, const ReservoirState& initial_reservoir_state, @@ -122,6 +122,9 @@ namespace Opm ReservoirState& reservoir_state, WellState& well_state) { + SimulatorReport iterReport; + SimulatorReport report; + // Do model-specific once-per-step calculations. model_->prepareStep(timer, initial_reservoir_state, initial_well_state); @@ -130,44 +133,31 @@ namespace Opm // Let the model do one nonlinear iteration. // Set up for main solver loop. - int linIters = 0; bool converged = false; - int wellIters = 0; // ---------- Main nonlinear solver loop ---------- do { // Do the nonlinear step. If we are in a converged state, the // model will usually do an early return without an expensive // solve, unless the minIter() count has not been reached yet. - IterationReport report = model_->nonlinearIteration(iteration, timer, *this, reservoir_state, well_state); - if (report.failed) { - OPM_THROW(Opm::NumericalProblem, "Failed to complete a nonlinear iteration."); - } + iterReport = model_->nonlinearIteration(iteration, timer, *this, reservoir_state, well_state); + + report += iterReport; + report.converged = iterReport.converged; + converged = report.converged; - linIters += report.linear_iterations; - wellIters += report.well_iterations; - ++iteration; - } while ( (!converged && (iteration <= maxIter())) || (iteration <= minIter())); + iteration += 1; + } while ( (!converged && (iteration <= maxIter())) || (iteration < minIter())); if (!converged) { - if (model_->terminalOutputEnabled()) { - std::cerr << "WARNING: Failed to compute converged solution in " << iteration - 1 << " iterations." << std::endl; - } - return -1; // -1 indicates that the solver has to be restarted + OPM_THROW(Opm::NumericalProblem, "Failed to complete a time step within "+std::to_string(maxIter())+" iterations."); } - linearIterations_ += linIters; - nonlinearIterations_ += iteration - 1; // Since the last one will always be trivial. - linearizations_ += iteration; - wellIterations_ += wellIters; - linearIterationsLast_ = linIters; - nonlinearIterationsLast_ = iteration; - wellIterationsLast_ = wellIters; - // Do model-specific post-step actions. model_->afterStep(timer, reservoir_state, well_state); + report.converged = true; - return linIters; + return report; } diff --git a/opm/autodiff/SimulatorBase_impl.hpp b/opm/autodiff/SimulatorBase_impl.hpp index fa1673e9e..de0eb8f06 100644 --- a/opm/autodiff/SimulatorBase_impl.hpp +++ b/opm/autodiff/SimulatorBase_impl.hpp @@ -98,7 +98,6 @@ namespace Opm // Create timers and file for writing timing info. Opm::time::StopWatch solver_timer; - double stime = 0.0; Opm::time::StopWatch step_timer; Opm::time::StopWatch total_timer; total_timer.start(); @@ -133,12 +132,11 @@ namespace Opm desiredRestoreStep ); } - unsigned int totalLinearizations = 0; - unsigned int totalNonlinearIterations = 0; - unsigned int totalLinearIterations = 0; bool is_well_potentials_computed = param_.getDefault("compute_well_potentials", false ); std::vector well_potentials; DynamicListEconLimited dynamic_list_econ_limited; + SimulatorReport report; + SimulatorReport stepReport; bool ooip_computed = false; std::vector fipnum_global = eclipse_state_->get3DProperties().getIntGridProperty("FIPNUM").getData(); @@ -186,9 +184,14 @@ namespace Opm // write the inital state at the report stage if (timer.initialStep()) { + Dune::Timer perfTimer; + perfTimer.start(); + // No per cell data is written for initial step, but will be // for subsequent steps, when we have started simulating output_writer_.writeTimeStepWithoutCellProperties( timer, state, well_state ); + + report.output_write_time += perfTimer.stop(); } // Max oil saturation (for VPPARS), hysteresis update. @@ -220,8 +223,7 @@ namespace Opm step_msg << "\nTime step " << std::setw(4) <nonlinearIterations(); - step_report.total_linear_iterations = solver->linearIterations(); - step_report.total_linearizations = solver->linearizations(); - if ( output_writer_.isIORank() ) { - step_report.reportParam(tstep_os); + stepReport.reportParam(tstep_os); } } @@ -318,8 +307,11 @@ namespace Opm ++timer; // write simulation state at the report stage + Dune::Timer perfTimer; + perfTimer.start(); const auto& physicalModel = solver->model(); output_writer_.writeTimeStep( timer, state, well_state, physicalModel ); + report.output_write_time += perfTimer.stop(); prev_well_state = well_state; // The well potentials are only computed if they are needed @@ -334,13 +326,8 @@ namespace Opm // Stop timer and create timing report total_timer.stop(); - SimulatorReport report; - report.pressure_time = stime; - report.transport_time = 0.0; report.total_time = total_timer.secsSinceStart(); - report.total_linearizations = totalLinearizations; - report.total_newton_iterations = totalNonlinearIterations; - report.total_linear_iterations = totalLinearIterations; + report.converged = true; return report; } diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp b/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp index 25bbd2ae9..fcb91801c 100644 --- a/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp +++ b/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp @@ -161,7 +161,6 @@ public: // Create timers and file for writing timing info. Opm::time::StopWatch solver_timer; - double stime = 0.0; Opm::time::StopWatch step_timer; Opm::time::StopWatch total_timer; total_timer.start(); @@ -188,12 +187,11 @@ public: // desiredRestoreStep ); } - unsigned int totalLinearizations = 0; - unsigned int totalNonlinearIterations = 0; - unsigned int totalLinearIterations = 0; bool is_well_potentials_computed = param_.getDefault("compute_well_potentials", false ); std::vector well_potentials; DynamicListEconLimited dynamic_list_econ_limited; + SimulatorReport report; + SimulatorReport stepReport; bool ooip_computed = false; std::vector fipnum_global = eclState().get3DProperties().getIntGridProperty("FIPNUM").getData(); @@ -253,22 +251,25 @@ public: // write the inital state at the report stage if (timer.initialStep()) { + Dune::Timer perfTimer; + perfTimer.start(); - // calculate Intensive Quantities + // make sure that the Intensive Quantities cache is up to date const auto& gridManager = ebosSimulator_.gridManager(); const auto& gridView = gridManager.gridView(); auto elemIt = gridView.template begin<0>(); auto elemEndIt = gridView.template end<0>(); ElementContext elemCtx(ebosSimulator_); for (; elemIt != elemEndIt; ++ elemIt) { - // this is convenient, but slightly inefficient: one only needs to update - // the primary intensive quantities - elemCtx.updateAll(*elemIt); + elemCtx.updatePrimaryStencil(*elemIt); + elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0); } // No per cell data is written for initial step, but will be // for subsequent steps, when we have started simulating output_writer_.writeTimeStepWithoutCellProperties( timer, state, well_state ); + + report.output_write_time += perfTimer.stop(); } // Compute orignal FIP; @@ -283,11 +284,12 @@ public: std::ostringstream step_msg; boost::posix_time::time_facet* facet = new boost::posix_time::time_facet("%d-%b-%Y"); step_msg.imbue(std::locale(std::locale::classic(), facet)); - step_msg << "\nTime step " << std::setw(4) <step( timer, *solver, state, well_state, output_writer_ ); + report += adaptiveTimeStepping->step( timer, *solver, state, well_state, output_writer_ ); } else { // solve for complete report step - solver->step(timer, state, well_state); + stepReport = solver->step(timer, state, well_state); + report += stepReport; if( terminal_output_ ) { + //stepReport.briefReport(); std::ostringstream iter_msg; iter_msg << "Stepsize " << (double)unit::convert::to(timer.currentStepLength(), unit::day); if (solver->wellIterations() != 0) { @@ -324,13 +328,8 @@ public: // take time that was used to solve system for this reportStep solver_timer.stop(); - // accumulate the number of nonlinear and linear Iterations - totalLinearizations += solver->linearizations(); - totalNonlinearIterations += solver->nonlinearIterations(); - totalLinearIterations += solver->linearIterations(); - - // Report timing. - const double st = solver_timer.secsSinceStart(); + // update timing. + report.solver_time += solver_timer.secsSinceStart(); // Compute current FIP. std::vector> COIP; @@ -345,30 +344,22 @@ public: for (size_t reg = 0; reg < OOIP.size(); ++reg) { outputFluidInPlace(OOIP[reg], COIP[reg], eclState().getUnits(), reg+1); } - } - // accumulate total time - stime += st; - - if ( terminal_output_ ) - { std::string msg; - msg = "Fully implicit solver took: " + std::to_string(st) + " seconds. Total solver time taken: " + std::to_string(stime) + " seconds."; + msg = + "Time step took " + std::to_string(stepReport.solver_time) + " seconds; " + "total solver time " + std::to_string(report.solver_time) + " seconds."; OpmLog::note(msg); } - if ( output_writer_.output() ) { - SimulatorReport step_report; - step_report.pressure_time = st; - step_report.total_time = step_timer.secsSinceStart(); - step_report.reportParam(tstep_os); - } - // Increment timer, remember well state. ++timer; // write simulation state at the report stage + Dune::Timer perfTimer; + perfTimer.start(); output_writer_.writeTimeStep( timer, state, well_state, solver->model() ); + report.output_write_time += perfTimer.stop(); prev_well_state = well_state; // The well potentials are only computed if they are needed @@ -383,13 +374,8 @@ public: // Stop timer and create timing report total_timer.stop(); - SimulatorReport report; - report.pressure_time = stime; - report.transport_time = 0.0; report.total_time = total_timer.secsSinceStart(); - report.total_linearizations = totalLinearizations; - report.total_newton_iterations = totalNonlinearIterations; - report.total_linear_iterations = totalLinearIterations; + report.converged = true; return report; } diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoilMultiSegment_impl.hpp b/opm/autodiff/SimulatorFullyImplicitBlackoilMultiSegment_impl.hpp index f29f105ce..354985c03 100644 --- a/opm/autodiff/SimulatorFullyImplicitBlackoilMultiSegment_impl.hpp +++ b/opm/autodiff/SimulatorFullyImplicitBlackoilMultiSegment_impl.hpp @@ -179,7 +179,7 @@ namespace Opm stime += st; if ( output_writer_.output() ) { SimulatorReport step_report; - step_report.pressure_time = st; + step_report.solver_time = st; step_report.total_time = step_timer.secsSinceStart(); step_report.reportParam(tstep_os); } @@ -198,9 +198,8 @@ namespace Opm // Stop timer and create timing report total_timer.stop(); SimulatorReport report; - report.pressure_time = stime; - report.transport_time = 0.0; report.total_time = total_timer.secsSinceStart(); + report.solver_time = stime; report.total_newton_iterations = totalNonlinearIterations; report.total_linear_iterations = totalLinearIterations; return report; diff --git a/opm/autodiff/StandardWellsDense.hpp b/opm/autodiff/StandardWellsDense.hpp index 0bebb8de6..0a7951085 100644 --- a/opm/autodiff/StandardWellsDense.hpp +++ b/opm/autodiff/StandardWellsDense.hpp @@ -177,14 +177,14 @@ enum WellVariablePositions { template - IterationReport assemble(Simulator& ebosSimulator, + SimulatorReport assemble(Simulator& ebosSimulator, const int iterationIdx, const double dt, WellState& well_state) { - IterationReport iter_report = {false, false, 0, 0}; + SimulatorReport report; if ( ! localWellsActive() ) { - return iter_report; + return report; } resetWellControlFromState(well_state); @@ -199,7 +199,7 @@ enum WellVariablePositions { if (param_.solve_welleq_initially_ && iterationIdx == 0) { // solve the well equations as a pre-processing step - iter_report = solveWellEq(ebosSimulator, dt, well_state); + report = solveWellEq(ebosSimulator, dt, well_state); } assembleWellEq(ebosSimulator, dt, well_state, false); @@ -207,7 +207,8 @@ enum WellVariablePositions { if (param_.compute_well_potentials_) { //wellModel().computeWellPotentials(mob_perfcells, b_perfcells, state0, well_state); } - return iter_report; + report.converged = true; + return report; } template @@ -660,9 +661,8 @@ enum WellVariablePositions { } } - template - IterationReport solveWellEq(Simulator& ebosSimulator, + SimulatorReport solveWellEq(Simulator& ebosSimulator, const double dt, WellState& well_state) { @@ -695,9 +695,10 @@ enum WellVariablePositions { well_state = well_state0; } - const bool failed = false; // Not needed in this method. - const int linear_iters = 0; // Not needed in this method - return IterationReport{failed, converged, linear_iters, it}; + SimulatorReport report; + report.converged = converged; + report.total_well_iterations = it; + return report; } void printIf(int c, double x, double y, double eps, std::string type) { @@ -1171,9 +1172,9 @@ enum WellVariablePositions { // We disregard terminal_ouput here as with it only messages // for wells on one process will be printed. std::ostringstream ss; - ss << "Switching control mode for well " << wells().name[w] - << " from " << modestring[well_controls_iget_type(wc, current)] - << " to " << modestring[well_controls_iget_type(wc, ctrl_index)] << std::endl; + ss << " Switching control mode for well " << wells().name[w] + << " from " << modestring[well_controls_iget_type(wc, current)] + << " to " << modestring[well_controls_iget_type(wc, ctrl_index)]; OpmLog::info(ss.str()); xw.currentControls()[w] = ctrl_index; current = xw.currentControls()[w]; diff --git a/opm/polymer/fullyimplicit/BlackoilPolymerModel.hpp b/opm/polymer/fullyimplicit/BlackoilPolymerModel.hpp index 4076c3d50..847ed810e 100644 --- a/opm/polymer/fullyimplicit/BlackoilPolymerModel.hpp +++ b/opm/polymer/fullyimplicit/BlackoilPolymerModel.hpp @@ -121,7 +121,7 @@ namespace Opm { /// \param[in] reservoir_state reservoir state variables /// \param[in, out] well_state well state variables /// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep - IterationReport + SimulatorReport assemble(const ReservoirState& reservoir_state, WellState& well_state, const bool initial_assembly); diff --git a/opm/polymer/fullyimplicit/BlackoilPolymerModel_impl.hpp b/opm/polymer/fullyimplicit/BlackoilPolymerModel_impl.hpp index 93fed4df4..c0eb6725f 100644 --- a/opm/polymer/fullyimplicit/BlackoilPolymerModel_impl.hpp +++ b/opm/polymer/fullyimplicit/BlackoilPolymerModel_impl.hpp @@ -494,13 +494,15 @@ namespace Opm { template - IterationReport + SimulatorReport BlackoilPolymerModel::assemble(const ReservoirState& reservoir_state, WellState& well_state, const bool initial_assembly) { using namespace Opm::AutoDiffGrid; + SimulatorReport report; + // Possibly switch well controls and updating well state to // get reasonable initial conditions for the wells // updateWellControls(well_state); @@ -531,10 +533,9 @@ namespace Opm { // -------- Mass balance equations -------- assembleMassBalanceEq(state); - IterationReport iter_report = {false, false, 0, 0}; // -------- Well equations ---------- if ( ! wellsActive() ) { - return iter_report; + return report; } std::vector mob_perfcells; @@ -572,7 +573,9 @@ namespace Opm { wellModel().addWellFluxEq(cq_s, state, residual_); addWellContributionToMassBalanceEq(cq_s, state, well_state); wellModel().addWellControlEq(state, well_state, aliveWells, residual_); - return iter_report; + + report.converged = true; + return report; } diff --git a/opm/simulators/WellSwitchingLogger.cpp b/opm/simulators/WellSwitchingLogger.cpp index d2e98ab5b..ea5c2974c 100644 --- a/opm/simulators/WellSwitchingLogger.cpp +++ b/opm/simulators/WellSwitchingLogger.cpp @@ -137,7 +137,7 @@ void WellSwitchingLogger::logSwitch(const char* name, std::array fromto, int rank) { std::ostringstream ss; - ss << "Switching control mode for well " << name + ss << " Switching control mode for well " << name << " from " << modestring[WellControlType(fromto[0])] << " to " << modestring[WellControlType(fromto[1])] << " on rank " << rank << std::endl; diff --git a/opm/simulators/WellSwitchingLogger.hpp b/opm/simulators/WellSwitchingLogger.hpp index 46fb40974..fe94e0602 100644 --- a/opm/simulators/WellSwitchingLogger.hpp +++ b/opm/simulators/WellSwitchingLogger.hpp @@ -74,7 +74,7 @@ public: else { std::ostringstream ss; - ss << "Switching control mode for well " << name + ss << " Switching control mode for well " << name << " from " << modestring[from] << " to " << modestring[to] << std::endl; OpmLog::info(ss.str());