From 6606bb75b2a51c8a460ee376d78ba73254d391ae Mon Sep 17 00:00:00 2001 From: Tor Harald Sandve Date: Fri, 29 Jun 2018 09:12:14 +0200 Subject: [PATCH] Minor fix based on code review --- opm/autodiff/BlackoilWellModel.hpp | 3 +-- opm/autodiff/BlackoilWellModel_impl.hpp | 27 ++++++++--------------- opm/autodiff/StandardWell.hpp | 1 - opm/autodiff/WellInterface.hpp | 4 +--- opm/autodiff/WellInterface_impl.hpp | 19 +++++++++++----- opm/core/wells/DynamicListEconLimited.hpp | 16 -------------- 6 files changed, 24 insertions(+), 46 deletions(-) diff --git a/opm/autodiff/BlackoilWellModel.hpp b/opm/autodiff/BlackoilWellModel.hpp index 7b754d272..034c5c9c1 100644 --- a/opm/autodiff/BlackoilWellModel.hpp +++ b/opm/autodiff/BlackoilWellModel.hpp @@ -142,8 +142,7 @@ namespace Opm { // compute the well fluxes and assemble them in to the reservoir equations as source terms // and in the well equations. void assemble(const int iterationIdx, - const double dt, - bool wtest = false); + const double dt); // substract Binv(D)rw from r; void apply( BVector& r) const; diff --git a/opm/autodiff/BlackoilWellModel_impl.hpp b/opm/autodiff/BlackoilWellModel_impl.hpp index 581a10d41..780e68019 100644 --- a/opm/autodiff/BlackoilWellModel_impl.hpp +++ b/opm/autodiff/BlackoilWellModel_impl.hpp @@ -110,14 +110,9 @@ namespace Opm { // update the previous well state. This is used to restart failed steps. previous_well_state_ = well_state_; - - if (wellCollection().havingVREPGroups() ) { - rateConverter_->template defineState(ebosSimulator_); - } - // Compute reservoir volumes for RESV controls. rateConverter_.reset(new RateConverterType (phase_usage_, - std::vector(number_of_cells_, 0))); + std::vector(number_of_cells_, 0))); computeRESV(timeStepIdx); // update VFP properties @@ -154,10 +149,9 @@ namespace Opm { // calculate the efficiency factors for each well calculateEfficiencyFactors(); - const Grid& grid = ebosSimulator_.vanguard().grid(); - if (has_polymer_) { + const Grid& grid = ebosSimulator_.vanguard().grid(); if (PolymerModule::hasPlyshlog()) { computeRepRadiusPerfLength(grid); } @@ -194,7 +188,7 @@ namespace Opm { const std::string msg = std::string("well ") + testWell.first + std::string(" is tested"); OpmLog::info(msg); - // finding the location of the well in wells_ecl + // Finding the location of the well in wells_ecl const int nw_wells_ecl = wells_ecl_.size(); int index_well = 0; for (; index_well < nw_wells_ecl; ++index_well) { @@ -202,15 +196,13 @@ namespace Opm { break; } } - // It should be able to find in wells_ecl. if (index_well == nw_wells_ecl) { OPM_THROW(std::logic_error, "Could not find well " << testWell.first << " in wells_ecl "); } const Well* well_ecl = wells_ecl_[index_well]; - - // Find the index in the wells() struct + // Finding the location of the well in wells struct. const int nw = numWells(); int wellidx = -999; for (int w = 0; w < nw; ++w) { @@ -218,9 +210,11 @@ namespace Opm { wellidx = w; break; } + } + if (wellidx < 0) { + OPM_THROW(std::logic_error, "Could not find the well " << testWell.first << " in the well struct "); } - // Use the pvtRegionIdx from the top cell const int well_cell_top = wells()->well_cells[wells()->well_connpos[wellidx]]; const int pvtreg = pvt_region_idx_[well_cell_top]; @@ -373,8 +367,7 @@ namespace Opm { void BlackoilWellModel:: assemble(const int iterationIdx, - const double dt, - bool onlyDoTheWellTest) + const double dt) { @@ -398,8 +391,6 @@ namespace Opm { if (param_.solve_welleq_initially_ && iterationIdx == 0) { // solve the well equations as a pre-processing step last_report_ = solveWellEq(dt); - if (onlyDoTheWellTest) - return; if (initial_step_) { // update the explicit quantities to get the initial fluid distribution in the well correct. @@ -930,7 +921,7 @@ namespace Opm { } for (auto& well : well_container_) { - const std::string well_name = well->name(); + const std::string& well_name = well->name(); const WellNode& well_node = wellCollection().findWellNode(well_name); const double well_efficiency_factor = well_node.getAccumulativeEfficiencyFactor(); diff --git a/opm/autodiff/StandardWell.hpp b/opm/autodiff/StandardWell.hpp index 78f893890..ff6dac53c 100644 --- a/opm/autodiff/StandardWell.hpp +++ b/opm/autodiff/StandardWell.hpp @@ -188,7 +188,6 @@ namespace Opm using Base::wellHasTHPConstraints; using Base::mostStrictBhpFromBhpLimits; using Base::scalingFactor; - using Base::updateWellControl; // protected member variables from the Base class using Base::current_step_; diff --git a/opm/autodiff/WellInterface.hpp b/opm/autodiff/WellInterface.hpp index 1c15d0fd4..844d4240e 100644 --- a/opm/autodiff/WellInterface.hpp +++ b/opm/autodiff/WellInterface.hpp @@ -229,7 +229,7 @@ namespace Opm void closeWellsAndCompletions(WellTestState& wellTestState); - const Well* wellEcl() const { return well_ecl_;} + const Well* wellEcl() const; protected: @@ -347,8 +347,6 @@ namespace Opm double scalingFactor(const int comp_idx) const; - int numberOfCompletions() const { return well_ecl_->getCompletions(current_step_).size();} - }; diff --git a/opm/autodiff/WellInterface_impl.hpp b/opm/autodiff/WellInterface_impl.hpp index 7926c329f..87cff3dbe 100644 --- a/opm/autodiff/WellInterface_impl.hpp +++ b/opm/autodiff/WellInterface_impl.hpp @@ -203,6 +203,14 @@ namespace Opm + template + const Well* + WellInterface:: + wellEcl() const + { + return well_ecl_; + } + template @@ -540,7 +548,7 @@ namespace Opm const auto& connections = well_ecl_->getConnections(current_step_); int complnumIdx = 0; - std::vector water_cut_in_completions(numberOfCompletions(), 0.0); + std::vector water_cut_in_completions(completions.size(), 0.0); for (const auto& completion : completions) { int complnum = completion.first; for (int perf = 0; perf < perf_number; ++perf) { @@ -731,10 +739,10 @@ namespace Opm wellTestState.addClosedWell(well_name, WellTestConfig::Reason::ECONOMIC, simulationTime); if (writeMessageToOPMLog) { if (well_ecl_->getAutomaticShutIn()) { - const std::string msg = well_name + std::string(" will be shut due to last compleation closed"); + const std::string msg = well_name + std::string(" will be shut due to last completion closed"); OpmLog::info(msg); } else { - const std::string msg = well_name + std::string(" will be stopped due to last compleation closed"); + const std::string msg = well_name + std::string(" will be stopped due to last completion closed"); OpmLog::info(msg); } } @@ -908,7 +916,6 @@ namespace Opm WellInterface::closeWellsAndCompletions(WellTestState& wellTestState) { if (wellTestState.hasWell(name(), WellTestConfig::Reason::ECONOMIC)) { - assert(!well_ecl_->getAutomaticShutIn()); well_controls_stop_well(wellControls()); } @@ -951,11 +958,11 @@ namespace Opm if (converged) { if ( terminal_output ) { - OpmLog::debug("Well equation solution gets converged with " + std::to_string(it) + " iterations"); + OpmLog::debug("Well equation for well " + name() + " solution gets converged with " + std::to_string(it) + " iterations"); } } else { if ( terminal_output ) { - OpmLog::debug("Well equation solution failed in getting converged with " + std::to_string(it) + " iterations"); + OpmLog::debug("Well equation for well" +name() + " solution failed in getting converged with " + std::to_string(it) + " iterations"); well_state = well_state0; updatePrimaryVariables(well_state); // also recover the old well controls diff --git a/opm/core/wells/DynamicListEconLimited.hpp b/opm/core/wells/DynamicListEconLimited.hpp index a428ebaf0..22bf2355d 100644 --- a/opm/core/wells/DynamicListEconLimited.hpp +++ b/opm/core/wells/DynamicListEconLimited.hpp @@ -80,22 +80,6 @@ namespace Opm } } - void removeShutWell(const std::string& well_name) { - auto itr = std::find(m_shut_wells.begin(), m_shut_wells.end(), well_name); - if (itr != m_shut_wells.end()) - m_shut_wells.erase(itr); - } - - void removeStoppedWell(const std::string& well_name) { - auto itr = std::find(m_stopped_wells.begin(), m_stopped_wells.end(), well_name); - if (itr != m_stopped_wells.end()) - m_stopped_wells.erase(itr); - } - - void removeClosedConnectionsForWell(const std::string& well_name) { - m_cells_closed_connections.erase(well_name); - } - private: std::vector m_shut_wells; std::vector m_stopped_wells;