OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #3551 from joakim-hove/pinfo-ws

Browse Source 
Use SingleWellState to manage parallel well info
This commit is contained in:
Markus Blatt 2021-09-29 09:59:00 +02:00 committed by GitHub
commit f6ecacf1b7
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
11 changed files with 96 additions and 101 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

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

@ -265,7 +265,7 @@ initFromRestartFile(const RestartValue& restartValues,
const int report_step = std::max(eclState_.getInitConfig().getRestartStep() - 1, 0);
// wells_ecl_ should only contain wells on this processor.
wells_ecl_ = getLocalWells(report_step);
local_parallel_well_info_ = createLocalParallelWellInfo(wells_ecl_);
this->local_parallel_well_info_ = createLocalParallelWellInfo(wells_ecl_);
this->initializeWellProdIndCalculators();
initializeWellPerfData();
@ -273,7 +273,7 @@ initFromRestartFile(const RestartValue& restartValues,
const int nw = wells_ecl_.size();
if (nw > 0) {
handle_ms_well &= anyMSWellOpenLocal();
this->wellState().resize(wells_ecl_, local_parallel_well_info_, schedule(), handle_ms_well, numCells, well_perf_data_, summaryState_); // Resize for restart step
this->wellState().resize(wells_ecl_, this->local_parallel_well_info_, schedule(), handle_ms_well, numCells, well_perf_data_, summaryState_); // Resize for restart step
loadRestartData(restartValues.wells, restartValues.grp_nwrk, phase_usage_, handle_ms_well, this->wellState());
}
@ -297,11 +297,11 @@ getLocalWells(const int timeStepIdx) const
return w;
}
std::vector<ParallelWellInfo*>
std::vector<std::reference_wrapper<ParallelWellInfo>>
BlackoilWellModelGeneric::
createLocalParallelWellInfo(const std::vector<Well>& wells)
{
std::vector<ParallelWellInfo*> local_parallel_well_info;
std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info;
local_parallel_well_info.reserve(wells.size());
for (const auto& well : wells)
{
@ -311,7 +311,7 @@ createLocalParallelWellInfo(const std::vector<Well>& wells)
wellPair);
assert(pwell != parallel_well_info_.end() &&
*pwell == wellPair);
local_parallel_well_info.push_back(&(*pwell));
local_parallel_well_info.push_back(std::ref(*pwell));
}
return local_parallel_well_info;
}
@ -339,10 +339,10 @@ initializeWellPerfData()
int completion_index_above = ParallelWellInfo::INVALID_ECL_INDEX;
well_perf_data_[well_index].clear();
well_perf_data_[well_index].reserve(well.getConnections().size());
CheckDistributedWellConnections checker(well, *local_parallel_well_info_[well_index]);
CheckDistributedWellConnections checker(well, local_parallel_well_info_[well_index].get());
bool hasFirstPerforation = false;
bool firstOpenCompletion = true;
auto& parallelWellInfo = *local_parallel_well_info_[well_index];
auto& parallelWellInfo = this->local_parallel_well_info_[well_index].get();
parallelWellInfo.beginReset();
for (const auto& completion : well.getConnections()) {

4
opm/simulators/wells/BlackoilWellModelGeneric.hpp
View File

@ -241,7 +241,7 @@ protected:
/// \brief Create the parallel well information
/// \param localWells The local wells from ECL schedule
std::vector<ParallelWellInfo*> createLocalParallelWellInfo(const std::vector<Well>& wells);
std::vector<std::reference_wrapper<ParallelWellInfo>> createLocalParallelWellInfo(const std::vector<Well>& wells);
void initializeWellProdIndCalculators();
void initializeWellPerfData();
@ -396,7 +396,7 @@ protected:
std::vector<int> local_shut_wells_{};
std::vector<ParallelWellInfo> parallel_well_info_;
std::vector<ParallelWellInfo*> local_parallel_well_info_;
std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
std::vector<WellProdIndexCalculator> prod_index_calc_;

10
opm/simulators/wells/BlackoilWellModel_impl.hpp
View File

@ -183,7 +183,7 @@ namespace Opm {
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
// Make wells_ecl_ contain only this partition's wells.
wells_ecl_ = getLocalWells(timeStepIdx);
local_parallel_well_info_ = createLocalParallelWellInfo(wells_ecl_);
this->local_parallel_well_info_ = createLocalParallelWellInfo(wells_ecl_);
// at least initializeWellState might be throw
// exception in opm-material (UniformTabulated2DFunction.hpp)
@ -416,9 +416,9 @@ namespace Opm {
endReportStep()
{
// Clear the communication data structures for above values.
for (auto&& pinfo : local_parallel_well_info_)
for (auto&& pinfo : this->local_parallel_well_info_)
{
pinfo->clear();
pinfo.get().clear();
}
}
@ -736,7 +736,7 @@ namespace Opm {
const auto pvtreg = perf_data.empty()
? 0 : pvt_region_idx_[perf_data.front().cell_index];
const auto& parallel_well_info = *local_parallel_well_info_[wellID];
const auto& parallel_well_info = this->local_parallel_well_info_[wellID].get();
const auto global_pvtreg = parallel_well_info.broadcastFirstPerforationValue(pvtreg);
return std::make_unique<WellType>(this->wells_ecl_[wellID],
@ -1685,7 +1685,7 @@ namespace Opm {
double weighted_temperature = 0.0;
double total_weight = 0.0;
auto& well_info = *local_parallel_well_info_[wellID];
auto& well_info = local_parallel_well_info_[wellID].get();
const int num_perf_this_well = well_info.communication().sum(well_perf_data_[wellID].size());
auto& ws = this->wellState().well(wellID);
auto& perf_data = ws.perf_data;

24
opm/simulators/wells/SingleWellState.cpp
View File

@ -21,8 +21,9 @@
namespace Opm {
SingleWellState::SingleWellState(bool is_producer, std::size_t num_perf, std::size_t num_phases, double temp)
: producer(is_producer)
SingleWellState::SingleWellState(const ParallelWellInfo& pinfo, bool is_producer, std::size_t num_perf, std::size_t num_phases, double temp)
: parallel_info(pinfo)
, producer(is_producer)
, temperature(temp)
, well_potentials(num_phases)
, productivity_index(num_phases)
@ -66,6 +67,25 @@ void SingleWellState::open() {
this->status = Well::Status::OPEN;
}
double SingleWellState::sum_connection_rates(const std::vector<double> connection_rates) const {
return this->parallel_info.get().sumPerfValues(connection_rates.begin(), connection_rates.end());
}
double SingleWellState::sum_brine_rates() const {
return this->sum_connection_rates(this->perf_data.brine_rates);
}
double SingleWellState::sum_polymer_rates() const {
return this->sum_connection_rates(this->perf_data.polymer_rates);
}
double SingleWellState::sum_solvent_rates() const {
return this->sum_connection_rates(this->perf_data.solvent_rates);
}
}

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

@ -20,18 +20,22 @@
#ifndef OPM_SINGLE_WELL_STATE_HEADER_INCLUDED
#define OPM_SINGLE_WELL_STATE_HEADER_INCLUDED
#include <functional>
#include <vector>
#include <opm/simulators/wells/SegmentState.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Well/Well.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Events.hpp>
#include <opm/simulators/wells/PerfData.hpp>
#include <opm/simulators/wells/ParallelWellInfo.hpp>
namespace Opm {
class SingleWellState {
public:
SingleWellState(bool is_producer, std::size_t num_perf, std::size_t num_phases, double temp);
SingleWellState(const ParallelWellInfo& pinfo, bool is_producer, std::size_t num_perf, std::size_t num_phases, double temp);
std::reference_wrapper<const ParallelWellInfo> parallel_info;
Well::Status status{Well::Status::OPEN};
bool producer;
@ -55,6 +59,15 @@ public:
void shut();
void stop();
void open();
// The sum_xxx_rates() functions sum over all connection rates of pertinent
// types. In the case of distributed wells this involves an MPI
// communication.
double sum_solvent_rates() const;
double sum_polymer_rates() const;
double sum_brine_rates() const;
private:
double sum_connection_rates(const std::vector<double> connection_rates) const;
};

4
opm/simulators/wells/StandardWellEval.cpp
View File

@ -291,10 +291,10 @@ updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_log
}
if constexpr (has_gfrac_variable) {
primary_variables_[GFrac] = baseif_.scalingFactor(pu.phase_pos[Gas]) * (ws.surface_rates[pu.phase_pos[Gas]]
- (Indices::enableSolvent ? well_state.solventWellRate(well_index) : 0.0) ) / total_well_rate ;
- (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0) ) / total_well_rate ;
}
if constexpr (Indices::enableSolvent) {
primary_variables_[SFrac] = baseif_.scalingFactor(pu.phase_pos[Gas]) * well_state.solventWellRate(well_index) / total_well_rate ;
primary_variables_[SFrac] = baseif_.scalingFactor(pu.phase_pos[Gas]) * ws.sum_solvent_rates() / total_well_rate ;
}
} else { // total_well_rate == 0
if (baseif_.isInjector()) {

7
opm/simulators/wells/StandardWell_impl.hpp
View File

@ -1210,7 +1210,6 @@ namespace Opm
const PhaseUsage& pu = phaseUsage();
b_perf.resize(nperf * this->num_components_);
surf_dens_perf.resize(nperf * this->num_components_);
const int w = this->index_of_well_;
const auto& ws = well_state.well(this->index_of_well_);
const bool waterPresent = FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx);
@ -1234,8 +1233,6 @@ namespace Opm
const auto& intQuants = *(ebosSimulator.model().cachedIntensiveQuantities(cell_idx, /*timeIdx=*/0));
const auto& fs = intQuants.fluidState();
// TODO: this is another place to show why WellState need to be a vector of WellState.
// TODO: to check why should be perf - 1
const double p_avg = (perf_press[perf] + p_above[perf])/2;
const double temperature = fs.temperature(FluidSystem::oilPhaseIdx).value();
const double saltConcentration = fs.saltConcentration().value();
@ -1254,7 +1251,7 @@ namespace Opm
const double oilrate = std::abs(ws.surface_rates[pu.phase_pos[Oil]]); //in order to handle negative rates in producers
rvmax_perf[perf] = FluidSystem::gasPvt().saturatedOilVaporizationFactor(fs.pvtRegionIndex(), temperature, p_avg);
if (oilrate > 0) {
const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? well_state.solventWellRate(w) : 0.0);
const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? ws.sum_solvent_rates() : 0.0);
double rv = 0.0;
if (gasrate > 0) {
rv = oilrate / gasrate;
@ -1277,7 +1274,7 @@ namespace Opm
const int oilpos = oilCompIdx + perf * this->num_components_;
if (gasPresent) {
rsmax_perf[perf] = FluidSystem::oilPvt().saturatedGasDissolutionFactor(fs.pvtRegionIndex(), temperature, p_avg);
const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? well_state.solventWellRate(w) : 0.0);
const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? ws.sum_solvent_rates() : 0.0);
if (gasrate > 0) {
const double oilrate = std::abs(ws.surface_rates[pu.phase_pos[Oil]]);
double rs = 0.0;

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

@ -240,11 +240,12 @@ namespace WellGroupHelpers
if (wellEcl.getStatus() == Well::Status::SHUT)
continue;
const auto& ws = wellState.well(well_index.value());
double factor = wellEcl.getEfficiencyFactor();
if (injector)
rate += factor * wellState.solventWellRate(well_index.value());
rate += factor * ws.sum_solvent_rates();
else
rate -= factor * wellState.solventWellRate(well_index.value());
rate -= factor * ws.sum_solvent_rates();
}
return rate;
}

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

@ -33,13 +33,12 @@ namespace Opm
{
void WellState::base_init(const std::vector<double>& cellPressures,
const std::vector<Well>& wells_ecl,
const std::vector<ParallelWellInfo*>& parallel_well_info,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state)
const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state)
{
// clear old name mapping
this->parallel_well_info_.clear();
this->wells_.clear();
{
// const int nw = wells->number_of_wells;
@ -61,7 +60,7 @@ void WellState::base_init(const std::vector<double>& cellPressures,
void WellState::initSingleWell(const std::vector<double>& cellPressures,
const Well& well,
const std::vector<PerforationData>& well_perf_data,
const ParallelWellInfo* well_info,
const ParallelWellInfo& well_info,
const SummaryState& summary_state)
{
assert(well.isInjector() || well.isProducer());
@ -72,8 +71,7 @@ void WellState::initSingleWell(const std::vector<double>& cellPressures,
const int np = pu.num_phases;
double temp = well.isInjector() ? well.injectionControls(summary_state).temperature : 273.15 + 15.56;
this->parallel_well_info_.add(well.name(), well_info);
auto& ws = this->wells_.add(well.name(), SingleWellState{well.isProducer(), well_perf_data.size(), static_cast<std::size_t>(np), temp});
auto& ws = this->wells_.add(well.name(), SingleWellState{well_info, well.isProducer(), well_perf_data.size(), static_cast<std::size_t>(np), temp});
if ( ws.perf_data.empty())
return;
@ -91,7 +89,7 @@ void WellState::initSingleWell(const std::vector<double>& cellPressures,
const double local_pressure = well_perf_data.empty() ?
0 : cellPressures[well_perf_data[0].cell_index];
const double global_pressure = well_info->broadcastFirstPerforationValue(local_pressure);
const double global_pressure = well_info.broadcastFirstPerforationValue(local_pressure);
if (well.getStatus() == Well::Status::OPEN) {
ws.status = Well::Status::OPEN;
@ -194,13 +192,13 @@ void WellState::initSingleWell(const std::vector<double>& cellPressures,
void WellState::init(const std::vector<double>& cellPressures,
const Schedule& schedule,
const std::vector<Well>& wells_ecl,
const std::vector<ParallelWellInfo*>& parallel_well_info,
const int report_step,
const WellState* prevState,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state)
const Schedule& schedule,
const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const int report_step,
const WellState* prevState,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state)
{
// call init on base class
this->base_init(cellPressures, wells_ecl, parallel_well_info, well_perf_data, summary_state);
@ -211,7 +209,7 @@ void WellState::init(const std::vector<double>& cellPressures,
}
for (const auto& winfo: parallel_well_info)
{
well_rates[winfo->name()].first = winfo->isOwner();
well_rates[winfo.get().name()].first = winfo.get().isOwner();
}
const int nw = wells_ecl.size();
@ -371,12 +369,12 @@ void WellState::init(const std::vector<double>& cellPressures,
}
void WellState::resize(const std::vector<Well>& wells_ecl,
const std::vector<ParallelWellInfo*>& parallel_well_info,
const Schedule& schedule,
const bool handle_ms_well,
const size_t numCells,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state)
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const Schedule& schedule,
const bool handle_ms_well,
const size_t numCells,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state)
{
const std::vector<double> tmp(numCells, 0.0); // <- UGLY HACK to pass the size
init(tmp, schedule, wells_ecl, parallel_well_info, 0, nullptr, well_perf_data, summary_state);
@ -471,15 +469,15 @@ WellState::report(const int* globalCellIdxMap,
}
if (pu.has_solvent || pu.has_zFraction) {
well.rates.set(rt::solvent, solventWellRate(well_index));
well.rates.set(rt::solvent, ws.sum_solvent_rates());
}
if (pu.has_polymer) {
well.rates.set(rt::polymer, polymerWellRate(well_index));
well.rates.set(rt::polymer, ws.sum_polymer_rates());
}
if (pu.has_brine) {
well.rates.set(rt::brine, brineWellRate(well_index));
well.rates.set(rt::brine, ws.sum_brine_rates());
}
if (ws.producer) {
@ -500,7 +498,7 @@ WellState::report(const int* globalCellIdxMap,
curr.inj = ws.injection_cmode;
}
const auto& pwinfo = *this->parallel_well_info_[well_index];
const auto& pwinfo = ws.parallel_info.get();
if (pwinfo.communication().size()==1)
{
reportConnections(well.connections, pu, well_index, globalCellIdxMap);
@ -737,31 +735,6 @@ WellState::calculateSegmentRates(const std::vector<std::vector<int>>& segment_in
}
}
double WellState::solventWellRate(const int w) const
{
auto& ws = this->well(w);
const auto& perf_data = ws.perf_data;
const auto& perf_rates_solvent = perf_data.solvent_rates;
return parallel_well_info_[w]->sumPerfValues(perf_rates_solvent.begin(), perf_rates_solvent.end());
}
double WellState::polymerWellRate(const int w) const
{
auto& ws = this->well(w);
const auto& perf_data = ws.perf_data;
const auto& perf_rates_polymer = perf_data.polymer_rates;
return parallel_well_info_[w]->sumPerfValues(perf_rates_polymer.begin(), perf_rates_polymer.end());
}
double WellState::brineWellRate(const int w) const
{
auto& ws = this->well(w);
const auto& perf_data = ws.perf_data;
const auto& perf_rates_brine = perf_data.brine_rates;
return parallel_well_info_[w]->sumPerfValues(perf_rates_brine.begin(), perf_rates_brine.end());
}
void WellState::stopWell(int well_index)
{
auto& ws = this->well(well_index);
@ -908,7 +881,7 @@ WellState::reportSegmentResults(const PhaseUsage& pu,
bool WellState::wellIsOwned(std::size_t well_index,
[[maybe_unused]] const std::string& wellName) const
{
const auto& well_info = parallelWellInfo(well_index);
const auto& well_info = this->parallelWellInfo(well_index);
assert(well_info.name() == wellName);
return well_info.isOwner();
@ -986,7 +959,8 @@ void WellState::resetConnectionTransFactors(const int well_index,
const ParallelWellInfo&
WellState::parallelWellInfo(std::size_t well_index) const
{
return *parallel_well_info_[well_index];
const auto& ws = this->well(well_index);
return ws.parallel_info;
}
template void WellState::updateGlobalIsGrup<ParallelWellInfo::Communication>(const ParallelWellInfo::Communication& comm);

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

@ -86,14 +86,14 @@ public:
void init(const std::vector<double>& cellPressures,
const Schedule& schedule,
const std::vector<Well>& wells_ecl,
const std::vector<ParallelWellInfo*>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const int report_step,
const WellState* prevState,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state);
void resize(const std::vector<Well>& wells_ecl,
const std::vector<ParallelWellInfo*>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const Schedule& schedule,
const bool handle_ms_well,
const size_t numCells,
@ -132,15 +132,6 @@ public:
static void calculateSegmentRates(const std::vector<std::vector<int>>& segment_inlets, const std::vector<std::vector<int>>&segment_perforations,
const std::vector<double>& perforation_rates, const int np, const int segment, std::vector<double>& segment_rates);
/// One rate pr well
double solventWellRate(const int w) const;
/// One rate pr well
double polymerWellRate(const int w) const;
/// One rate pr well
double brineWellRate(const int w) const;
template<class Comm>
void communicateGroupRates(const Comm& comm);
@ -268,7 +259,6 @@ private:
PhaseUsage phase_usage_;
WellContainer<SingleWellState> wells_;
WellContainer<const ParallelWellInfo*> parallel_well_info_;
// The well_rates variable is defined for all wells on all processors. The
// bool in the value pair is whether the current process owns the well or
// not.
@ -303,14 +293,14 @@ private:
/// with -1e100.
void base_init(const std::vector<double>& cellPressures,
const std::vector<Well>& wells_ecl,
const std::vector<ParallelWellInfo*>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const SummaryState& summary_state);
void initSingleWell(const std::vector<double>& cellPressures,
const Well& well,
const std::vector<PerforationData>& well_perf_data,
const ParallelWellInfo* well_info,
const ParallelWellInfo& well_info,
const SummaryState& summary_state);

14
tests/test_wellstate.cpp
View File

@ -18,6 +18,7 @@
*/
#include <config.h>
#include <functional>
#define BOOST_TEST_MODULE WellStateFIBOTest
@ -138,17 +139,15 @@ namespace {
auto wells = setup.sched.getWells(timeStep);
pinfos.resize(wells.size());
std::vector<Opm::ParallelWellInfo*> ppinfos(wells.size());
std::vector<std::reference_wrapper<Opm::ParallelWellInfo>> ppinfos;
auto pw = pinfos.begin();
auto ppw = ppinfos.begin();
for (const auto& well : wells)
{
*pw = {well.name()};
*ppw = &(*pw);
ppinfos.push_back(std::ref(*pw));
pw->communicateFirstPerforation(true);
++pw;
++ppw;
}
state.init(cpress, setup.sched,
@ -569,9 +568,10 @@ BOOST_AUTO_TEST_CASE(TESTPerfData) {
BOOST_AUTO_TEST_CASE(TestSingleWellState) {
Opm::SingleWellState ws1(true, 10, 3, 1);
Opm::SingleWellState ws2(true, 10, 3, 2);
Opm::SingleWellState ws3(false, 10, 3, 3);
Opm::ParallelWellInfo pinfo;
Opm::SingleWellState ws1(pinfo, true, 10, 3, 1);
Opm::SingleWellState ws2(pinfo, true, 10, 3, 2);
Opm::SingleWellState ws3(pinfo, false, 10, 3, 3);
ws1.bhp = 100;
ws1.thp = 200;