OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-03 16:40:29 -06:00
Code Issues Packages Projects Releases Wiki Activity

Remove WellMap from WellState

Browse Source
This commit is contained in:
Joakim Hove 2021-08-06 13:18:18 +02:00
parent 8f9e7d8e96
commit 5a721b8cd2
6 changed files with 71 additions and 139 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
ebos/ecltracermodel.hh
View File

@ -433,7 +433,7 @@ protected:
// TODO: Some inconsistencies here that perhaps should be clarified. The "offical" rate as reported below is
// occasionally significant different from the sum over connections (as calculated above). Only observed
// for small values, neglible for the rate itself, but matters when used to calculate tracer concentrations.
std::size_t well_index = simulator_.problem().wellModel().wellState().wellIndex(well.name());
std::size_t well_index = simulator_.problem().wellModel().wellState().index(well.name()).value();
Scalar official_well_rate_total = simulator_.problem().wellModel().wellState().well(well_index).surface_rates[tr.phaseIdx_];
rateWellTotal = official_well_rate_total;

9
opm/simulators/wells/BlackoilWellModelGeneric.cpp
View File

@ -169,9 +169,9 @@ loadRestartData(const data::Wells& rst_wells,
phs.at( phases.phase_pos[BlackoilPhases::Vapour] ) = rt::gas;
}
for( const auto& wm : well_state.wellMap() ) {
const auto well_index = wm.second[ 0 ];
const auto& rst_well = rst_wells.at( wm.first );
for( std::size_t well_index = 0; well_index < well_state.size(); well_index++) {
const auto& well_name = well_state.name(well_index);
const auto& rst_well = rst_wells.at(well_name);
auto& ws = well_state.well(well_index);
ws.bhp = rst_well.bhp;
ws.thp = rst_well.thp;
@ -206,7 +206,6 @@ loadRestartData(const data::Wells& rst_wells,
if (handle_ms_well && !rst_well.segments.empty()) {
// we need the well_ecl_ information
const std::string& well_name = wm.first;
const Well& well_ecl = getWellEcl(well_name);
const WellSegments& segment_set = well_ecl.getSegments();
@ -1528,7 +1527,7 @@ updateAndCommunicateGroupData(const int reportStepIdx,
// Set ALQ for off-process wells to zero
for (const auto& wname : schedule().wellNames(reportStepIdx)) {
const bool is_producer = schedule().getWell(wname, reportStepIdx).isProducer();
const bool not_on_this_process = well_state.wellMap().count(wname) == 0;
const bool not_on_this_process = !well_state.has(wname);
if (is_producer && not_on_this_process) {
well_state.setALQ(wname, 0.0);
}

2
opm/simulators/wells/GasLiftSingleWellGeneric.cpp
View File

@ -1235,7 +1235,7 @@ bool
GasLiftSingleWellGeneric::OptimizeState::
checkThpControl()
{
const int well_index = this->parent.well_state_.wellIndex(this->parent.well_name_);
const int well_index = this->parent.well_state_.index(this->parent.well_name_).value();
const Well::ProducerCMode& control_mode = this->parent.well_state_.well(well_index).production_cmode;
return control_mode == Well::ProducerCMode::THP;
}

76
opm/simulators/wells/WellGroupHelpers.cpp
View File

@ -74,16 +74,13 @@ namespace {
const auto& gefac = groupTmp.getGroupEfficiencyFactor();
rate += gefac * sumWellPhaseRates(res_rates, groupTmp, schedule, wellState, reportStepIdx, phasePos, injector);
}
const auto& end = wellState.wellMap().end();
for (const std::string& wellName : group.wells()) {
const auto& it = wellState.wellMap().find(wellName);
if (it == end) // the well is not found
const auto& well_index = wellState.index(wellName);
if (!well_index.has_value())
continue;
int well_index = it->second[0];
if (! wellState.wellIsOwned(well_index, wellName) ) // Only sum once
if (! wellState.wellIsOwned(well_index.value(), wellName) ) // Only sum once
{
continue;
}
@ -97,7 +94,7 @@ namespace {
continue;
double factor = wellEcl.getEfficiencyFactor();
const auto& ws = wellState.well(well_index);
const auto& ws = wellState.well(well_index.value());
if (res_rates) {
const auto& well_rates = ws.reservoir_rates;
if (injector)
@ -216,15 +213,13 @@ namespace WellGroupHelpers
const auto& gefac = groupTmp.getGroupEfficiencyFactor();
rate += gefac * sumSolventRates(groupTmp, schedule, wellState, reportStepIdx, injector);
}
const auto& end = wellState.wellMap().end();
for (const std::string& wellName : group.wells()) {
const auto& it = wellState.wellMap().find(wellName);
if (it == end) // the well is not found
const auto& well_index = wellState.index(wellName);
if (!well_index.has_value())
continue;
int well_index = it->second[0];
if (! wellState.wellIsOwned(well_index, wellName) ) // Only sum once
if (! wellState.wellIsOwned(well_index.value(), wellName) ) // Only sum once
{
continue;
}
@ -239,9 +234,9 @@ namespace WellGroupHelpers
double factor = wellEcl.getEfficiencyFactor();
if (injector)
rate += factor * wellState.solventWellRate(well_index);
rate += factor * wellState.solventWellRate(well_index.value());
else
rate -= factor * wellState.solventWellRate(well_index);
rate -= factor * wellState.solventWellRate(well_index.value());
}
return rate;
}
@ -401,22 +396,19 @@ namespace WellGroupHelpers
if (wellTmp.getStatus() == Well::Status::SHUT)
continue;
const auto& end = wellState.wellMap().end();
const auto& it = wellState.wellMap().find(wellName);
if (it == end) // the well is not found
const auto& well_index = wellState.index(wellName);
if (!well_index.has_value())
continue;
int well_index = it->second[0];
if (! wellState.wellIsOwned(well_index, wellName) ) // Only sum once
if (! wellState.wellIsOwned(well_index.value(), wellName) ) // Only sum once
{
continue;
}
const double efficiency = wellTmp.getEfficiencyFactor();
// add contributino from wells not under group control
const auto& ws_nupcol = wellStateNupcol.well(well_index);
const auto& ws = wellState.well(well_index);
const auto& ws_nupcol = wellStateNupcol.well(well_index.value());
const auto& ws = wellState.well(well_index.value());
if (isInjector) {
if (ws.injection_cmode != Well::InjectorCMode::GRUP)
for (int phase = 0; phase < np; phase++) {
@ -476,20 +468,17 @@ namespace WellGroupHelpers
if (wellTmp.getStatus() == Well::Status::SHUT)
continue;
const auto& end = wellState.wellMap().end();
const auto& it = wellState.wellMap().find(wellName);
if (it == end) // the well is not found
const auto& well_index = wellState.index(wellName);
if (!well_index.has_value())
continue;
int well_index = it->second[0];
if (! wellState.wellIsOwned(well_index, wellName) ) // Only sum once
if (! wellState.wellIsOwned(well_index.value(), wellName) ) // Only sum once
{
continue;
}
// scale rates
auto& ws = wellState.well(well_index);
auto& ws = wellState.well(well_index.value());
if (isInjector) {
if (ws.injection_cmode == Well::InjectorCMode::GRUP)
for (int phase = 0; phase < np; phase++) {
@ -568,19 +557,17 @@ namespace WellGroupHelpers
updateWellRates(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState);
}
const int np = wellState.numPhases();
const auto& end = wellState.wellMap().end();
for (const std::string& wellName : group.wells()) {
std::vector<double> rates(np, 0.0);
const auto& it = wellState.wellMap().find(wellName);
if (it != end) { // the well is found on this node
int well_index = it->second[0];
const auto& well_index = wellState.index(wellName);
if (well_index.has_value()) { // the well is found on this node
const auto& wellTmp = schedule.getWell(wellName, reportStepIdx);
int sign = 1;
// production wellRates are negative. The users of currentWellRates uses the convention in
// opm-common that production and injection rates are positive.
if (!wellTmp.isInjector())
sign = -1;
const auto& ws = wellStateNupcol.well(well_index);
const auto& ws = wellStateNupcol.well(well_index.value());
for (int phase = 0; phase < np; ++phase) {
rates[phase] = sign * ws.surface_rates[phase];
}
@ -1380,19 +1367,16 @@ namespace WellGroupHelpers
if (wellTmp.getStatus() == Well::Status::SHUT)
continue;
const auto& end = wellState.wellMap().end();
const auto& it = wellState.wellMap().find(wellName);
if (it == end) // the well is not found
const auto& well_index = wellState.index(wellName);
if (!well_index.has_value()) // the well is not found
continue;
int well_index = it->second[0];
if (! wellState.wellIsOwned(well_index, wellName) ) // Only sum once
if (! wellState.wellIsOwned(well_index.value(), wellName) ) // Only sum once
{
continue;
}
const auto& ws = wellState.well(well_index);
const auto& ws = wellState.well(well_index.value());
// add contribution from wells unconditionally
for (int phase = 0; phase < np; phase++) {
pot[phase] += wefac * ws.well_potentials[phase];
@ -1430,18 +1414,16 @@ namespace WellGroupHelpers
const Comm& comm,
GuideRate* guideRate)
{
const auto& end = wellState.wellMap().end();
for (const auto& well : schedule.getWells(reportStepIdx)) {
double oilpot = 0.0;
double gaspot = 0.0;
double waterpot = 0.0;
const auto& it = wellState.wellMap().find(well.name());
if (it != end && wellState.wellIsOwned(it->second[0], well.name()))
const auto& well_index = wellState.index(well.name());
if (well_index.has_value() && wellState.wellIsOwned(well_index.value(), well.name()))
{
// the well is found and owned
int well_index = it->second[0];
const auto& ws = wellState.well(well_index);
const auto& ws = wellState.well(well_index.value());
const auto& wpot = ws.well_potentials;
if (pu.phase_used[BlackoilPhases::Liquid] > 0)
oilpot = wpot[pu.phase_pos[BlackoilPhases::Liquid]];

106
opm/simulators/wells/WellState.cpp
View File

@ -39,29 +39,17 @@ void WellState::base_init(const std::vector<double>& cellPressures,
const SummaryState& summary_state)
{
// clear old name mapping
this->wellMap_.clear();
this->parallel_well_info_.clear();
this->wells_.clear();
{
// const int nw = wells->number_of_wells;
const int nw = wells_ecl.size();
// const int np = wells->number_of_phases;
int connpos = 0;
for (int w = 0; w < nw; ++w) {
const Well& well = wells_ecl[w];
// Initialize bhp(), thp(), wellRates(), temperature().
initSingleWell(cellPressures, well, well_perf_data[w], parallel_well_info[w], summary_state);
// Setup wellname -> well index mapping.
const int num_perf_this_well = well_perf_data[w].size();
std::string name = well.name();
assert( !name.empty() );
mapentry_t& wellMapEntry = wellMap_[name];
wellMapEntry[ 0 ] = w;
wellMapEntry[ 1 ] = connpos;
wellMapEntry[ 2 ] = num_perf_this_well;
connpos += num_perf_this_well;
}
}
}
@ -85,10 +73,9 @@ void WellState::initSingleWell(const std::vector<double>& cellPressures,
double temp = well.isInjector() ? well.injectionControls(summary_state).temperature : 273.15 + 15.56;
this->parallel_well_info_.add(well.name(), well_info);
const int num_perf_this_well = well_info->communication().sum(well_perf_data.size());
auto& ws = this->wells_.add(well.name(), SingleWellState{well.isProducer(), static_cast<std::size_t>(num_perf_this_well), static_cast<std::size_t>(np), temp});
auto& ws = this->wells_.add(well.name(), SingleWellState{well.isProducer(), well_perf_data.size(), static_cast<std::size_t>(np), temp});
if ( num_perf_this_well == 0 )
if ( ws.perf_data.empty())
return;
const auto inj_controls = well.isInjector() ? well.injectionControls(summary_state) : Well::InjectionControls(0);
@ -256,11 +243,10 @@ void WellState::init(const std::vector<double>& cellPressures,
// rates divided by the number of perforations.
const auto& ecl_well = wells_ecl[w];
const auto& wname = ecl_well.name();
const auto& well_info = this->wellMap().at(wname);
const int num_perf_this_well = well_info[2];
const int global_num_perf_this_well = ecl_well.getConnections().num_open();
auto& ws = this->well(w);
auto& perf_data = ws.perf_data;
const int num_perf_this_well = perf_data.size();
const int global_num_perf_this_well = ecl_well.getConnections().num_open();
const auto& perf_input = well_perf_data[w];
for (int perf = 0; perf < num_perf_this_well; ++perf) {
@ -309,20 +295,16 @@ void WellState::init(const std::vector<double>& cellPressures,
// intialize wells that have been there before
// order may change so the mapping is based on the well name
if (prevState && !prevState->wellMap().empty()) {
auto end = prevState->wellMap().end();
if (prevState && prevState->size() > 0) {
for (int w = 0; w < nw; ++w) {
const Well& well = wells_ecl[w];
if (well.getStatus() == Well::Status::SHUT) {
continue;
}
auto& new_well = this->well(w);
auto it = prevState->wellMap().find(well.name());
if (it != end)
{
const int oldIndex = it->second[ 0 ];
const auto& prev_well = prevState->well(oldIndex);
const auto& old_index = prevState->index(well.name());
if (old_index.has_value()) {
const auto& prev_well = prevState->well(old_index.value());
new_well.init_timestep(prev_well);
@ -346,16 +328,8 @@ void WellState::init(const std::vector<double>& cellPressures,
new_well.well_potentials = prev_well.well_potentials;
// perfPhaseRates
const int num_perf_old_well = (*it).second[ 2 ];
const auto new_iter = this->wellMap().find(well.name());
if (new_iter == this->wellMap().end()) {
throw std::logic_error {
well.name() + " is not in internal well map - "
"Bug in WellState"
};
}
const int num_perf_this_well = new_iter->second[2];
const int num_perf_old_well = prev_well.perf_data.size();
const int num_perf_this_well = new_well.perf_data.size();
const bool global_num_perf_same = (num_perf_this_well == num_perf_old_well);
// copy perforation rates when the number of
@ -459,8 +433,7 @@ WellState::report(const int* globalCellIdxMap,
const auto& pu = this->phaseUsage();
data::Wells res;
for( const auto& [wname, winfo]: this->wellMap() ) {
const auto well_index = winfo[ 0 ];
for( std::size_t well_index = 0; well_index < this->size(); well_index++) {
const auto& ws = this->well(well_index);
if ((ws.status == Well::Status::SHUT) && !wasDynamicallyClosed(well_index))
{
@ -471,6 +444,8 @@ WellState::report(const int* globalCellIdxMap,
const auto& well_potentials = ws.well_potentials;
const auto& wpi = ws.productivity_index;
const auto& wv = ws.surface_rates;
const auto& wname = this->name(well_index);
data::Well well;
well.bhp = ws.bhp;
@ -572,8 +547,6 @@ void WellState::reportConnections(std::vector<data::Connection>& connections,
connection.reservoir_rate = perf_rates[i];
connection.trans_factor = perf_data.connection_transmissibility_factor[i];
}
assert(num_perf_well == int(connections.size()));
const int np = pu.num_phases;
size_t local_comp_index = 0;
@ -616,7 +589,6 @@ void WellState::reportConnections(std::vector<data::Connection>& connections,
++local_comp_index;
}
assert(local_comp_index == this->perfdata[well_index].size());
}
void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
@ -672,28 +644,22 @@ void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
auto& perf_data = ws.perf_data;
// for the seg_rates_, now it becomes a recursive solution procedure.
{
// make sure the information from wells_ecl consistent with wells
assert((n_activeperf == this->wellMap().at(well_ecl.name())[2]) &&
"Inconsistent number of reservoir connections in well");
if (pu.phase_used[Gas]) {
auto& perf_rates = perf_data.phase_rates;
const int gaspos = pu.phase_pos[Gas];
// scale the phase rates for Gas to avoid too bad initial guess for gas fraction
// it will probably benefit the standard well too, while it needs to be justified
// TODO: to see if this strategy can benefit StandardWell too
// TODO: it might cause big problem for gas rate control or if there is a gas rate limit
// maybe the best way is to initialize the fractions first then get the rates
for (int perf = 0; perf < n_activeperf; perf++)
perf_rates[perf*np + gaspos] *= 100;
}
const auto& perf_rates = perf_data.phase_rates;
std::vector<double> perforation_rates(perf_rates.begin(), perf_rates.end());
calculateSegmentRates(segment_inlets, segment_perforations, perforation_rates, np, 0 /* top segment */, ws.segments.rates);
if (pu.phase_used[Gas]) {
auto& perf_rates = perf_data.phase_rates;
const int gaspos = pu.phase_pos[Gas];
// scale the phase rates for Gas to avoid too bad initial guess for gas fraction
// it will probably benefit the standard well too, while it needs to be justified
// TODO: to see if this strategy can benefit StandardWell too
// TODO: it might cause big problem for gas rate control or if there is a gas rate limit
// maybe the best way is to initialize the fractions first then get the rates
for (int perf = 0; perf < n_activeperf; perf++)
perf_rates[perf*np + gaspos] *= 100;
}
const auto& perf_rates = perf_data.phase_rates;
std::vector<double> perforation_rates(perf_rates.begin(), perf_rates.end());
calculateSegmentRates(segment_inlets, segment_perforations, perforation_rates, np, 0 /* top segment */, ws.segments.rates);
// for the segment pressure, the segment pressure is the same with the first perforation belongs to the segment
// if there is no perforation associated with this segment, it uses the pressure from the outlet segment
// which requres the ordering is successful
@ -947,21 +913,11 @@ bool WellState::wellIsOwned(std::size_t well_index,
bool WellState::wellIsOwned(const std::string& wellName) const
{
const auto& it = this->wellMap_.find( wellName );
if (it == this->wellMap_.end()) {
const auto& well_index = this->index(wellName);
if (!well_index.has_value())
OPM_THROW(std::logic_error, "Could not find well " << wellName << " in well map");
}
const int well_index = it->second[0];
return wellIsOwned(well_index, wellName);
}
int WellState::wellIndex(const std::string& wellName) const
{
const auto& it = this->wellMap_.find( wellName );
if (it == this->wellMap_.end()) {
OPM_THROW(std::logic_error, "Could not find well " << wellName << " in well map");
}
return it->second[0];
return wellIsOwned(well_index.value(), wellName);
}
int WellState::numSegments(const int well_id) const

15
opm/simulators/wells/WellState.hpp
View File

@ -52,9 +52,6 @@ class Schedule;
class WellState
{
public:
using mapentry_t = std::array<int, 3>;
using WellMapType = std::map<std::string, mapentry_t>;
static const uint64_t event_mask = ScheduleEvents::WELL_STATUS_CHANGE + ScheduleEvents::PRODUCTION_UPDATE + ScheduleEvents::INJECTION_UPDATE;
virtual ~WellState() = default;
@ -69,11 +66,8 @@ public:
this->phase_usage_ = pu;
}
const WellMapType& wellMap() const { return wellMap_; }
WellMapType& wellMap() { return wellMap_; }
std::size_t size() const {
return this->wellMap_.size();
return this->wells_.size();
}
@ -82,8 +76,6 @@ public:
return this->size();
}
int wellIndex(const std::string& wellName) const;
const ParallelWellInfo& parallelWellInfo(std::size_t well_index) const;
@ -250,6 +242,10 @@ public:
return this->wells_.well_name(well_index);
}
std::optional<std::size_t> index(const std::string& well_name) const {
return this->wells_.well_index(well_name);
}
const SingleWellState& well(std::size_t well_index) const {
return this->wells_[well_index];
}
@ -271,7 +267,6 @@ public:
}
private:
WellMapType wellMap_;
// Use of std::optional<> here is a technical crutch, the
// WellStateFullyImplicitBlackoil class should be default constructible,
// whereas the GlobalWellInfo is not.