Split Restart File Loading Into Multiple Stages

Mostly for readability in preparation of initialising the guide rate
object with values from the restart file.  The new stages, each with
its own 'loadRestart*Data()' member function are

  - Connection
  - Segment
  - Well
  - Group

While here, also switch to using std::any_of() in place of a raw
loop.
This commit is contained in:
Bård Skaflestad 2021-10-07 13:21:35 +02:00
parent 37f589bff1
commit 8014442a1d
2 changed files with 237 additions and 146 deletions

View File

@ -21,6 +21,7 @@
*/
#include <config.h>
#include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
#include <opm/output/data/Groups.hpp>
@ -35,10 +36,13 @@
#include <opm/simulators/wells/GasLiftStage2.hpp>
#include <opm/simulators/wells/VFPProperties.hpp>
#include <opm/simulators/wells/WellGroupHelpers.hpp>
#include <opm/simulators/wells/WellInterfaceGeneric.hpp>
#include <opm/simulators/wells/WellState.hpp>
#include <algorithm>
#include <cassert>
#include <stdexcept>
#include <vector>
#include <fmt/format.h>
@ -95,9 +99,11 @@ bool
BlackoilWellModelGeneric::
hasWell(const std::string& wname)
{
auto iter = std::find_if(this->wells_ecl_.begin(), this->wells_ecl_.end(),
[&wname](const Well& well) { return well.name() == wname; });
return (iter != this->wells_ecl_.end());
return std::any_of(this->wells_ecl_.begin(), this->wells_ecl_.end(),
[&wname](const Well& well)
{
return well.name() == wname;
});
}
bool
@ -144,111 +150,160 @@ getWellEcl(const std::string& well_name) const
void
BlackoilWellModelGeneric::
loadRestartData(const data::Wells& rst_wells,
const data::GroupAndNetworkValues& grpNwrkValues,
const PhaseUsage& phases,
const bool handle_ms_well,
WellState& well_state)
loadRestartConnectionData(const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
SingleWellState& ws)
{
auto& perf_data = ws.perf_data;
auto perf_pressure = perf_data.pressure.begin();
auto perf_rates = perf_data.rates.begin();
auto perf_phase_rates = perf_data.phase_rates.begin();
for (const auto& pd : old_perf_data) {
const auto& rst_connection = rst_well.connections[pd.ecl_index];
*perf_pressure = rst_connection.pressure; ++perf_pressure;
*perf_rates = rst_connection.reservoir_rate; ++perf_rates;
for (const auto& phase : phs) {
*perf_phase_rates = rst_connection.rates.get(phase);
++perf_phase_rates;
}
}
}
void
BlackoilWellModelGeneric::
loadRestartSegmentData(const std::string& well_name,
const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
SingleWellState& ws)
{
const auto& segment_set = this->getWellEcl(well_name).getSegments();
const auto& rst_segments = rst_well.segments;
// \Note: Eventually we need to handle the situations that some segments are shut
assert(0u + segment_set.size() == rst_segments.size());
const auto np = phs.size();
const auto pres_idx = data::SegmentPressures::Value::Pressure;
auto& segments = ws.segments;
auto& segment_pressure = segments.pressure;
auto& segment_rates = segments.rates;
for (const auto& [segNum, rst_segment] : rst_segments) {
const int segment_index = segment_set.segmentNumberToIndex(segNum);
// Recovering segment rates and pressure from the restart values
segment_pressure[segment_index] = rst_segment.pressures[pres_idx];
const auto& rst_segment_rates = rst_segment.rates;
for (auto p = 0*np; p < np; ++p) {
segment_rates[segment_index*np + p] = rst_segment_rates.get(phs[p]);
}
}
}
void
BlackoilWellModelGeneric::
loadRestartWellData(const std::string& well_name,
const bool handle_ms_well,
const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
SingleWellState& ws)
{
const auto np = phs.size();
ws.bhp = rst_well.bhp;
ws.thp = rst_well.thp;
ws.temperature = rst_well.temperature;
if (rst_well.current_control.isProducer) {
ws.production_cmode = rst_well.current_control.prod;
}
else {
ws.injection_cmode = rst_well.current_control.inj;
}
for (auto i = 0*np; i < np; ++i) {
assert( rst_well.rates.has( phs[ i ] ) );
ws.surface_rates[i] = rst_well.rates.get(phs[i]);
}
this->loadRestartConnectionData(phs, rst_well, old_perf_data, ws);
if (handle_ms_well && !rst_well.segments.empty()) {
this->loadRestartSegmentData(well_name, phs, rst_well, ws);
}
}
void
BlackoilWellModelGeneric::
loadRestartGroupData(const std::string& group,
const data::GroupData& value)
{
using GPMode = Group::ProductionCMode;
using GIMode = Group::InjectionCMode;
const auto cpc = value.currentControl.currentProdConstraint;
const auto cgi = value.currentControl.currentGasInjectionConstraint;
const auto cwi = value.currentControl.currentWaterInjectionConstraint;
auto& grpState = this->groupState();
if (cpc != GPMode::NONE) {
grpState.production_control(group, cpc);
}
if (cgi != GIMode::NONE) {
grpState.injection_control(group, Phase::GAS, cgi);
}
if (cwi != GIMode::NONE) {
grpState.injection_control(group, Phase::WATER, cwi);
}
}
void
BlackoilWellModelGeneric::
loadRestartData(const data::Wells& rst_wells,
const data::GroupAndNetworkValues& grpNwrkValues,
const PhaseUsage& phases,
const bool handle_ms_well,
WellState& well_state)
{
using rt = data::Rates::opt;
const auto np = phases.num_phases;
std::vector< rt > phs( np );
if( phases.phase_used[BlackoilPhases::Aqua] ) {
phs.at( phases.phase_pos[BlackoilPhases::Aqua] ) = rt::wat;
std::vector<rt> phs(np);
if (phases.phase_used[BlackoilPhases::Aqua]) {
phs.at(phases.phase_pos[BlackoilPhases::Aqua]) = rt::wat;
}
if( phases.phase_used[BlackoilPhases::Liquid] ) {
if (phases.phase_used[BlackoilPhases::Liquid]) {
phs.at( phases.phase_pos[BlackoilPhases::Liquid] ) = rt::oil;
}
if( phases.phase_used[BlackoilPhases::Vapour] ) {
if (phases.phase_used[BlackoilPhases::Vapour]) {
phs.at( phases.phase_pos[BlackoilPhases::Vapour] ) = rt::gas;
}
for( std::size_t well_index = 0; well_index < well_state.size(); well_index++) {
for (auto well_index = 0*well_state.size();
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;
ws.temperature = rst_well.temperature;
if (rst_well.current_control.isProducer) {
ws.production_cmode = rst_well.current_control.prod;
}
else {
ws.injection_cmode = rst_well.current_control.inj;
}
for( size_t i = 0; i < phs.size(); ++i ) {
assert( rst_well.rates.has( phs[ i ] ) );
ws.surface_rates[i] = rst_well.rates.get(phs[i]);
}
auto& perf_data = ws.perf_data;
auto& perf_pressure = perf_data.pressure;
auto& perf_rates = perf_data.rates;
auto& perf_phase_rates = perf_data.phase_rates;
const auto& old_perf_data = this->well_perf_data_[well_index];
for (std::size_t perf_index = 0; perf_index < old_perf_data.size(); perf_index++) {
const auto& pd = old_perf_data[perf_index];
const auto& rst_connection = rst_well.connections[pd.ecl_index];
perf_pressure[perf_index] = rst_connection.pressure;
perf_rates[perf_index] = rst_connection.reservoir_rate;
for (int phase_index = 0; phase_index < np; ++phase_index)
perf_phase_rates[perf_index*np + phase_index] = rst_connection.rates.get(phs[phase_index]);
}
if (handle_ms_well && !rst_well.segments.empty()) {
// we need the well_ecl_ information
const Well& well_ecl = getWellEcl(well_name);
const WellSegments& segment_set = well_ecl.getSegments();
const auto& rst_segments = rst_well.segments;
// \Note: eventually we need to handle the situations that some segments are shut
assert(0u + segment_set.size() == rst_segments.size());
auto& segments = ws.segments;
auto& segment_pressure = segments.pressure;
auto& segment_rates = segments.rates;
for (const auto& rst_segment : rst_segments) {
const int segment_index = segment_set.segmentNumberToIndex(rst_segment.first);
// recovering segment rates and pressure from the restart values
const auto pres_idx = data::SegmentPressures::Value::Pressure;
segment_pressure[segment_index] = rst_segment.second.pressures[pres_idx];
const auto& rst_segment_rates = rst_segment.second.rates;
for (int p = 0; p < np; ++p) {
segment_rates[segment_index * np + p] = rst_segment_rates.get(phs[p]);
}
}
}
this->loadRestartWellData(well_name, handle_ms_well, phs,
rst_wells.at(well_name),
this->well_perf_data_[well_index],
well_state.well(well_index));
}
for (const auto& [group, value] : grpNwrkValues.groupData) {
const auto cpc = value.currentControl.currentProdConstraint;
const auto cgi = value.currentControl.currentGasInjectionConstraint;
const auto cwi = value.currentControl.currentWaterInjectionConstraint;
if (cpc != GPMode::NONE) {
this->groupState().production_control(group, cpc);
}
if (cgi != GIMode::NONE) {
this->groupState().injection_control(group, Phase::GAS, cgi);
}
if (cwi != GIMode::NONE) {
this->groupState().injection_control(group, Phase::WATER, cwi);
}
this->loadRestartGroupData(group, value);
}
}
@ -269,13 +324,17 @@ initFromRestartFile(const RestartValue& restartValues,
this->local_parallel_well_info_ = createLocalParallelWellInfo(wells_ecl_);
this->initializeWellProdIndCalculators();
initializeWellPerfData();
this->initializeWellPerfData();
const int nw = wells_ecl_.size();
if (nw > 0) {
if (! this->wells_ecl_.empty()) {
handle_ms_well &= anyMSWellOpenLocal();
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());
// Resize for restart step
this->wellState().resize(this->wells_ecl_, this->local_parallel_well_info_,
this->schedule(), handle_ms_well, numCells,
this->well_perf_data_, this->summaryState_);
loadRestartData(restartValues.wells, restartValues.grp_nwrk,
this->phase_usage_, handle_ms_well, this->wellState());
}
@ -1065,24 +1124,32 @@ updateEclWells(const int timeStepIdx,
const std::unordered_set<std::string>& wells)
{
for (const auto& wname : wells) {
auto well_iter = std::find_if( this->wells_ecl_.begin(), this->wells_ecl_.end(), [wname] (const auto& well) -> bool { return well.name() == wname;});
if (well_iter != this->wells_ecl_.end()) {
auto well_index = std::distance( this->wells_ecl_.begin(), well_iter );
this->wells_ecl_[well_index] = schedule_.getWell(wname, timeStepIdx);
auto well_iter = std::find_if(this->wells_ecl_.begin(), this->wells_ecl_.end(),
[&wname] (const auto& well) -> bool
{
return well.name() == wname;
});
const auto& well = this->wells_ecl_[well_index];
auto& pd = this->well_perf_data_[well_index];
auto pdIter = pd.begin();
for (const auto& conn : well.getConnections()) {
if (conn.state() != Connection::State::SHUT) {
pdIter->connection_transmissibility_factor = conn.CF();
++pdIter;
}
}
this->wellState().updateStatus(well_index, well.getStatus());
this->wellState().resetConnectionTransFactors(well_index, pd);
this->prod_index_calc_[well_index].reInit(well);
if (well_iter == this->wells_ecl_.end()) {
continue;
}
auto well_index = std::distance(this->wells_ecl_.begin(), well_iter);
this->wells_ecl_[well_index] = schedule_.getWell(wname, timeStepIdx);
const auto& well = this->wells_ecl_[well_index];
auto& pd = this->well_perf_data_[well_index];
auto pdIter = pd.begin();
for (const auto& conn : well.getConnections()) {
if (conn.state() != Connection::State::SHUT) {
pdIter->connection_transmissibility_factor = conn.CF();
++pdIter;
}
}
this->wellState().updateStatus(well_index, well.getStatus());
this->wellState().resetConnectionTransFactors(well_index, pd);
this->prod_index_calc_[well_index].reInit(well);
}
}
@ -1258,9 +1325,9 @@ getGuideRateValues(const Group& group) const
const auto& gname = group.name();
if (!this->groupState().has_production_rates(gname)) {
// No flow rates for production group 'gname' -- might be before group comes
// online (e.g., for the initial condition before simulation
// starts).
// No flow rates for production group 'gname' -- might be before
// group comes online (e.g., for the initial condition before
// simulation starts).
return grval;
}
@ -1269,7 +1336,8 @@ getGuideRateValues(const Group& group) const
return grval;
}
const auto qs = WellGroupHelpers::getProductionGroupRateVector(this->groupState(), this->phase_usage_, gname);
const auto qs = WellGroupHelpers::
getProductionGroupRateVector(this->groupState(), this->phase_usage_, gname);
const auto is_inj = false; // This procedure only applies to G*PGR.
this->getGuideRateValues(qs, is_inj, gname, grval);

View File

@ -20,10 +20,20 @@
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
#define OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
#include <opm/output/data/GuideRateValue.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/PerforationData.hpp>
#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
#include <opm/simulators/wells/WGState.hpp>
#include <functional>
#include <map>
#include <memory>
@ -32,46 +42,39 @@
#include <unordered_set>
#include <vector>
#include <dune/common/version.hh>
#include <opm/output/data/GuideRateValue.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
namespace Opm {
class DeferredLogger;
class EclipseState;
class GasLiftSingleWellGeneric;
class GasLiftWellState;
class Group;
class GuideRateConfig;
class ParallelWellInfo;
class RestartValue;
class Schedule;
class SummaryConfig;
class VFPProperties;
class WellInterfaceGeneric;
class WellState;
} // namespace Opm
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/ParallelWellInfo.hpp>
#include <opm/simulators/wells/PerforationData.hpp>
#include <opm/simulators/wells/WellInterfaceGeneric.hpp>
#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
#include <opm/simulators/wells/WGState.hpp>
namespace Opm { namespace data {
struct GroupData;
struct GroupGuideRates;
class GroupAndNetworkValues;
struct NodeData;
}} // namespace Opm::data
namespace Opm {
namespace data {
struct GroupData;
struct GroupGuideRates;
class GroupAndNetworkValues;
struct NodeData;
}
class DeferredLogger;
class EclipseState;
class GasLiftSingleWellGeneric;
class GasLiftWellState;
class Group;
class RestartValue;
class Schedule;
class SummaryConfig;
class VFPProperties;
class WellState;
/// Class for handling the blackoil well model.
class BlackoilWellModelGeneric
{
public:
// --------- Types ---------
using GLiftOptWells = std::map<std::string,std::unique_ptr<GasLiftSingleWellGeneric>>;
using GLiftProdWells = std::map<std::string,const WellInterfaceGeneric*>;
using GLiftWellStateMap = std::map<std::string,std::unique_ptr<GasLiftWellState>>;
using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric>>;
using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric*>;
using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState>>;
BlackoilWellModelGeneric(Schedule& schedule,
const SummaryState& summaryState,
@ -284,6 +287,26 @@ protected:
std::map<std::string, data::GroupData>& gvalues) const;
void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues) const;
void loadRestartConnectionData(const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
SingleWellState& ws);
void loadRestartSegmentData(const std::string& well_name,
const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
SingleWellState& ws);
void loadRestartWellData(const std::string& well_name,
const bool handle_ms_well,
const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
SingleWellState& ws);
void loadRestartGroupData(const std::string& group,
const data::GroupData& value);
std::unordered_map<std::string, data::GroupGuideRates>
calculateAllGroupGuiderates(const int reportStepIdx) const;
@ -426,7 +449,7 @@ protected:
bool glift_debug = false;
private:
private:
WellInterfaceGeneric* getGenWell(const std::string& well_name);
};