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 #2505 from atgeirr/group-controls-at-wells-cleaned

Browse Source 
Group controls at wells (cleaned)
This commit is contained in:
Tor Harald Sandve 2020-04-08 12:29:30 +02:00 committed by GitHub
commit 60bb9e4eaa
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
14 changed files with 2999 additions and 1751 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
CMakeLists_files.cmake
View File

@ -38,6 +38,7 @@ list (APPEND MAIN_SOURCE_FILES
opm/simulators/utils/ParallelRestart.cpp
opm/simulators/wells/VFPProdProperties.cpp
opm/simulators/wells/VFPInjProperties.cpp
opm/simulators/wells/WellGroupHelpers.cpp
)
if(CUDA_FOUND)
@ -188,6 +189,7 @@ list (APPEND PUBLIC_HEADER_FILES
opm/simulators/wells/PerforationData.hpp
opm/simulators/wells/RateConverter.hpp
opm/simulators/wells/SimFIBODetails.hpp
opm/simulators/wells/TargetCalculator.hpp
opm/simulators/wells/WellConnectionAuxiliaryModule.hpp
opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp
opm/simulators/wells/VFPProperties.hpp

12
opm/simulators/timestepping/ConvergenceReport.hpp
View File

@ -86,6 +86,7 @@ namespace Opm
: status_{AllGood}
, res_failures_{}
, well_failures_{}
, groupConverged_(true)
{
}
@ -94,6 +95,7 @@ namespace Opm
status_ = AllGood;
res_failures_.clear();
well_failures_.clear();
groupConverged_ = true;
}
void setReservoirFailed(const ReservoirFailure& rf)
@ -107,7 +109,12 @@ namespace Opm
status_ = static_cast<Status>(status_ | WellFailed);
well_failures_.push_back(wf);
}
void setGroupConverged(const bool groupConverged)
{
groupConverged_ = groupConverged;
}
ConvergenceReport& operator+=(const ConvergenceReport& other)
{
status_ = static_cast<Status>(status_ | other.status_);
@ -122,7 +129,7 @@ namespace Opm
bool converged() const
{
return status_ == AllGood;
return status_ == AllGood && groupConverged_;
}
bool reservoirFailed() const
@ -167,6 +174,7 @@ namespace Opm
Status status_;
std::vector<ReservoirFailure> res_failures_;
std::vector<WellFailure> well_failures_;
bool groupConverged_;
};
} // namespace Opm

19
opm/simulators/wells/BlackoilWellModel.hpp
View File

@ -243,7 +243,7 @@ namespace Opm {
void applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const;
// Check if well equations is converged.
ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg) const;
ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg, const bool checkGroupConvergence = false) const;
// return the internal well state, ignore the passed one.
// Used by the legacy code to make it compatible with the legacy well models.
@ -348,7 +348,9 @@ namespace Opm {
// xw to update Well State
void recoverWellSolutionAndUpdateWellState(const BVector& x);
void updateWellControls(Opm::DeferredLogger& deferred_logger, const bool checkGroupControl, const bool checkCurrentGroupControl);
void updateWellControls(Opm::DeferredLogger& deferred_logger, const bool checkGroupControls);
void updateAndCommunicateGroupData();
// setting the well_solutions_ based on well_state.
void updatePrimaryVariables(Opm::DeferredLogger& deferred_logger);
@ -421,11 +423,18 @@ namespace Opm {
const Well& getWellEcl(const std::string& well_name) const;
void checkGroupConstraints(const Group& group, const bool checkCurrentControl, Opm::DeferredLogger& deferred_logger);
void updateGroupIndividualControls(Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups);
void updateGroupIndividualControl(const Group& group, Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups);
bool checkGroupConstraints(const Group& group, Opm::DeferredLogger& deferred_logger) const;
Group::ProductionCMode checkGroupProductionConstraints(const Group& group, Opm::DeferredLogger& deferred_logger) const;
Group::InjectionCMode checkGroupInjectionConstraints(const Group& group, const Phase& phase) const;
void actionOnBrokenConstraints(const Group& group, const Group::ExceedAction& exceed_action, const Group::ProductionCMode& newControl, const int reportStepIdx, Opm::DeferredLogger& deferred_logger);
void updateGroupHigherControls(Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups);
void checkGroupHigherConstraints(const Group& group, Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups);
void actionOnBrokenConstraints(const Group& group, const Group::InjectionCMode& newControl, const Phase& topUpPhase, const int reportStepIdx, Opm::DeferredLogger& deferred_logger);
void actionOnBrokenConstraints(const Group& group, const Group::ExceedAction& exceed_action, const Group::ProductionCMode& newControl, Opm::DeferredLogger& deferred_logger);
void actionOnBrokenConstraints(const Group& group, const Group::InjectionCMode& newControl, const Phase& topUpPhase, Opm::DeferredLogger& deferred_logger);
WellInterfacePtr getWell(const std::string& well_name) const;

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

@ -292,7 +292,7 @@ namespace Opm {
}
}
const Group& fieldGroup = schedule().getGroup("FIELD", timeStepIdx);
wellGroupHelpers::setCmodeGroup(fieldGroup, schedule(), summaryState, timeStepIdx, well_state_);
WellGroupHelpers::setCmodeGroup(fieldGroup, schedule(), summaryState, timeStepIdx, well_state_);
// Compute reservoir volumes for RESV controls.
rateConverter_.reset(new RateConverterType (phase_usage_,
@ -381,7 +381,7 @@ namespace Opm {
//compute well guideRates
const auto& comm = ebosSimulator_.vanguard().grid().comm();
wellGroupHelpers::updateGuideRatesForWells(schedule(), phase_usage_, reportStepIdx, simulationTime, well_state_, comm, guideRate_.get());
WellGroupHelpers::updateGuideRatesForWells(schedule(), phase_usage_, reportStepIdx, simulationTime, well_state_, comm, guideRate_.get());
}
@ -408,7 +408,7 @@ namespace Opm {
well->init(&phase_usage_, depth_, gravity_, number_of_cells_);
const Well& wellEcl = schedule().getWell(well_name, timeStepIdx);
double well_efficiency_factor = wellEcl.getEfficiencyFactor();
wellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), timeStepIdx), schedule(), timeStepIdx, well_efficiency_factor);
WellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), timeStepIdx), schedule(), timeStepIdx, well_efficiency_factor);
well->setWellEfficiencyFactor(well_efficiency_factor);
well->setVFPProperties(vfp_properties_.get());
well->setGuideRate(guideRate_.get());
@ -820,11 +820,7 @@ namespace Opm {
calculateExplicitQuantities(local_deferredLogger);
prepareTimeStep(local_deferredLogger);
}
// check the current group control in the beginning of an episode for the first two iterations. The first itertion is needed for changes in group/well controls and closing of wells etc.
// a second check is needed for REIN and VREP controls since they depend on results from other wells.
// This check can probably be made more sofisticated, but this simple rule seems to work
bool checkCurrentGroupControls = (report_step_starts_ && iterationIdx < 2);
updateWellControls(local_deferredLogger, /*allow for switching to group controls*/true, checkCurrentGroupControls);
updateWellControls(local_deferredLogger, /* check group controls */ true);
// Set the well primary variables based on the value of well solutions
initPrimaryVariablesEvaluation();
@ -1065,7 +1061,7 @@ namespace Opm {
// are active wells anywhere in the global domain.
if( wellsActive() )
{
updateWellControls(deferred_logger, /*don't switch group controls*/false, /*don't check current group controls*/false);
updateWellControls(deferred_logger, /*don't switch group controls*/false);
initPrimaryVariablesEvaluation();
}
} catch (std::exception& e) {
@ -1107,7 +1103,7 @@ namespace Opm {
template<typename TypeTag>
ConvergenceReport
BlackoilWellModel<TypeTag>::
getWellConvergence(const std::vector<Scalar>& B_avg) const
getWellConvergence(const std::vector<Scalar>& B_avg, bool checkGroupConvergence) const
{
Opm::DeferredLogger local_deferredLogger;
@ -1118,7 +1114,6 @@ namespace Opm {
local_report += well->getWellConvergence(well_state_, B_avg, local_deferredLogger);
}
}
Opm::DeferredLogger global_deferredLogger = gatherDeferredLogger(local_deferredLogger);
if (terminal_output_) {
global_deferredLogger.logMessages();
@ -1137,6 +1132,12 @@ namespace Opm {
}
}
if (checkGroupConvergence) {
const int reportStepIdx = ebosSimulator_.episodeIndex();
const Group& fieldGroup = schedule().getGroup("FIELD", reportStepIdx);
bool violated = checkGroupConstraints(fieldGroup, global_deferredLogger);
report.setGroupConverged(!violated);
}
return report;
}
@ -1162,53 +1163,93 @@ namespace Opm {
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
updateWellControls(Opm::DeferredLogger& deferred_logger, const bool checkGroupControl, const bool checkCurrentGroupControl)
updateWellControls(Opm::DeferredLogger& deferred_logger, const bool checkGroupControls)
{
// Even if there are no wells active locally, we cannot
// return as the DeferredLogger uses global communication.
// For no well active globally we simply return.
if( !wellsActive() ) return ;
updateAndCommunicateGroupData();
std::set<std::string> switched_wells;
std::set<std::string> switched_groups;
if (checkGroupControls) {
// Check group individual constraints.
updateGroupIndividualControls(deferred_logger, switched_groups);
// Check group's constraints from higher levels.
updateGroupHigherControls(deferred_logger, switched_groups);
updateAndCommunicateGroupData();
// Check wells' group constraints and communicate.
for (const auto& well : well_container_) {
const auto mode = WellInterface<TypeTag>::IndividualOrGroup::Group;
const bool changed = well->updateWellControl(ebosSimulator_, mode, well_state_, deferred_logger);
if (changed) {
switched_wells.insert(well->name());
}
}
updateAndCommunicateGroupData();
}
// Check individual well constraints and communicate.
for (const auto& well : well_container_) {
if (switched_wells.count(well->name())) {
continue;
}
const auto mode = WellInterface<TypeTag>::IndividualOrGroup::Individual;
well->updateWellControl(ebosSimulator_, mode, well_state_, deferred_logger);
}
updateAndCommunicateGroupData();
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
updateAndCommunicateGroupData()
{
const int reportStepIdx = ebosSimulator_.episodeIndex();
const Group& fieldGroup = schedule().getGroup("FIELD", reportStepIdx);
// update group controls
if (checkGroupControl) {
checkGroupConstraints(fieldGroup, checkCurrentGroupControl, deferred_logger);
}
for (const auto& well : well_container_) {
well->updateWellControl(ebosSimulator_, well_state_, deferred_logger);
}
const int nupcol = schedule().getNupcol(reportStepIdx);
const int iterationIdx = ebosSimulator_.model().newtonMethod().numIterations();
// This builds some necessary lookup structures, so it must be called
// before we copy to well_state_nupcol_.
const auto& comm = ebosSimulator_.vanguard().grid().comm();
well_state_.updateGlobalIsGrup(schedule(), reportStepIdx, comm);
if (iterationIdx < nupcol) {
well_state_nupcol_ = well_state_;
}
// the group target reduction rates needs to be update since wells may have swicthed to/from GRUP control
// Currently the group targer reduction does not honor NUPCOL
// Currently the group target reduction does not honor NUPCOL. TODO: is that true?
std::vector<double> groupTargetReduction(numPhases(), 0.0);
wellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ false, phase_usage_, well_state_nupcol_, well_state_, groupTargetReduction);
WellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ false, phase_usage_, well_state_nupcol_, well_state_, groupTargetReduction);
std::vector<double> groupTargetReductionInj(numPhases(), 0.0);
wellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ true, phase_usage_, well_state_nupcol_, well_state_, groupTargetReductionInj);
const auto& comm = ebosSimulator_.vanguard().grid().comm();
well_state_.updateGlobalIsGrup(schedule(), reportStepIdx, comm);
WellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ true, phase_usage_, well_state_nupcol_, well_state_, groupTargetReductionInj);
const double simulationTime = ebosSimulator_.time();
std::vector<double> pot(numPhases(), 0.0);
wellGroupHelpers::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ false, well_state_, comm, guideRate_.get(), pot);
WellGroupHelpers::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ false, well_state_, comm, guideRate_.get(), pot);
std::vector<double> potInj(numPhases(), 0.0);
wellGroupHelpers::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ true, well_state_, comm, guideRate_.get(), potInj);
WellGroupHelpers::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ true, well_state_, comm, guideRate_.get(), potInj);
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
wellGroupHelpers::updateREINForGroups(fieldGroup, schedule(), reportStepIdx, phase_usage_, summaryState, well_state_nupcol_, well_state_);
wellGroupHelpers::updateVREPForGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_);
wellGroupHelpers::updateReservoirRatesInjectionGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_);
WellGroupHelpers::updateREINForGroups(fieldGroup, schedule(), reportStepIdx, phase_usage_, summaryState, well_state_nupcol_, well_state_);
WellGroupHelpers::updateVREPForGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_);
WellGroupHelpers::updateReservoirRatesInjectionGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_);
WellGroupHelpers::updateGroupProductionRates(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_);
WellGroupHelpers::updateWellRates(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_);
well_state_.communicateGroupRates(comm);
// compute wsolvent fraction for REIN wells
@ -1267,8 +1308,7 @@ namespace Opm {
summaryConfig.hasSummaryKey( "WOPP:" + well->name()) ||
summaryConfig.hasSummaryKey( "WGPP:" + well->name())) && well->isProducer());
const Well& eclWell = well->wellEcl();
bool needPotentialsForGuideRate = eclWell.getGuideRatePhase() == Well::GuideRateTarget::UNDEFINED;
bool needPotentialsForGuideRate = true;//eclWell.getGuideRatePhase() == Well::GuideRateTarget::UNDEFINED;
if (write_restart_file || needed_for_summary || needPotentialsForGuideRate)
{
std::vector<double> potentials;
@ -1357,7 +1397,7 @@ namespace Opm {
for (auto& well : well_container_) {
const Well& wellEcl = well->wellEcl();
double well_efficiency_factor = wellEcl.getEfficiencyFactor();
wellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), reportStepIdx), schedule(), reportStepIdx, well_efficiency_factor);
WellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), reportStepIdx), schedule(), reportStepIdx, well_efficiency_factor);
well->setWellEfficiencyFactor(well_efficiency_factor);
}
}
@ -1692,264 +1732,343 @@ namespace Opm {
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
checkGroupConstraints(const Group& group, const bool checkCurrentControl, Opm::DeferredLogger& deferred_logger) {
updateGroupIndividualControls(Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups)
{
const int reportStepIdx = ebosSimulator_.episodeIndex();
const Group& fieldGroup = schedule().getGroup("FIELD", reportStepIdx);
updateGroupIndividualControl(fieldGroup, deferred_logger, switched_groups);
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
updateGroupIndividualControl(const Group& group, Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups) {
const int reportStepIdx = ebosSimulator_.episodeIndex();
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
const auto& comm = ebosSimulator_.vanguard().grid().comm();
auto& well_state = well_state_;
if (group.isInjectionGroup())
const bool skip = switched_groups.count(group.name());
if (!skip && group.isInjectionGroup())
{
const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS};
for (Phase phase : all) {
if (!group.hasInjectionControl(phase)) {
continue;
}
int phasePos;
if (phase == Phase::GAS && phase_usage_.phase_used[BlackoilPhases::Vapour] )
phasePos = phase_usage_.phase_pos[BlackoilPhases::Vapour];
else if (phase == Phase::OIL && phase_usage_.phase_used[BlackoilPhases::Liquid])
phasePos = phase_usage_.phase_pos[BlackoilPhases::Liquid];
else if (phase == Phase::WATER && phase_usage_.phase_used[BlackoilPhases::Aqua] )
phasePos = phase_usage_.phase_pos[BlackoilPhases::Aqua];
else
continue;
const auto& controls = group.injectionControls(phase, summaryState);
const Group::InjectionCMode& currentControl = well_state.currentInjectionGroupControl(phase, group.name());
if (controls.has_control(Group::InjectionCMode::RATE))
Group::InjectionCMode newControl = checkGroupInjectionConstraints(group, phase);
if (newControl != Group::InjectionCMode::NONE)
{
if (checkCurrentControl || currentControl != Group::InjectionCMode::RATE)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.surface_max_rate < current_rate) {
actionOnBrokenConstraints(group, Group::InjectionCMode::RATE, phase, reportStepIdx, deferred_logger);
}
}
switched_groups.insert(group.name());
actionOnBrokenConstraints(group, newControl, phase, deferred_logger);
}
if (controls.has_control(Group::InjectionCMode::RESV))
{
if (checkCurrentControl || currentControl != Group::InjectionCMode::RESV)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.resv_max_rate < current_rate) {
actionOnBrokenConstraints(group, Group::InjectionCMode::RESV, phase, reportStepIdx, deferred_logger);
}
}
}
if (controls.has_control(Group::InjectionCMode::REIN))
{
if (checkCurrentControl || currentControl != Group::InjectionCMode::REIN)
{
double production_Rate = 0.0;
const Group& groupRein = schedule().getGroup(controls.reinj_group, reportStepIdx);
production_Rate += wellGroupHelpers::sumWellRates(groupRein, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/false);
// sum over all nodes
production_Rate = comm.sum(production_Rate);
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.target_reinj_fraction*production_Rate < current_rate) {
actionOnBrokenConstraints(group, Group::InjectionCMode::REIN, phase,reportStepIdx, deferred_logger);
}
}
}
if (controls.has_control(Group::InjectionCMode::VREP))
{
if (checkCurrentControl || currentControl != Group::InjectionCMode::VREP)
{
double voidage_rate = 0.0;
const Group& groupVoidage = schedule().getGroup(controls.voidage_group, reportStepIdx);
voidage_rate += wellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false);
voidage_rate += wellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false);
voidage_rate += wellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false);
// sum over all nodes
voidage_rate = comm.sum(voidage_rate);
double total_rate = 0.0;
total_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true);
total_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true);
total_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true);
// sum over all nodes
total_rate = comm.sum(total_rate);
if (controls.target_void_fraction*voidage_rate < total_rate) {
actionOnBrokenConstraints(group, Group::InjectionCMode::VREP, phase, reportStepIdx, deferred_logger);
}
}
}
// Handle GCONSALE
if (schedule().gConSale(reportStepIdx).has(group.name())) {
if (controls.phase != Phase::GAS)
OPM_THROW(std::runtime_error, "Group " + group.name() + " has GCONSALE control but is not a GAS group" );
const auto& gconsale = schedule().gConSale(reportStepIdx).get(group.name(), summaryState);
double sales_rate = 0.0;
int gasPos = phase_usage_.phase_pos[BlackoilPhases::Vapour];
sales_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/false);
sales_rate -= wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/true);
// sum over all nodes
sales_rate = comm.sum(sales_rate);
// add import rate and substract consumption rate for group for gas
if (schedule().gConSump(reportStepIdx).has(group.name())) {
const auto& gconsump = schedule().gConSump(reportStepIdx).get(group.name(), summaryState);
if (phase_usage_.phase_used[BlackoilPhases::Vapour]) {
sales_rate += gconsump.import_rate;
sales_rate -= gconsump.consumption_rate;
}
}
if (sales_rate > gconsale.max_sales_rate) {
OPM_THROW(std::runtime_error, "Group " + group.name() + " has sale rate more then the maximum permitted value. Not implemented in Flow" );
}
if (sales_rate < gconsale.min_sales_rate) {
OPM_THROW(std::runtime_error, "Group " + group.name() + " has sale rate less then minimum permitted value. Not implemented in Flow" );
}
if (gconsale.sales_target < 0.0) {
OPM_THROW(std::runtime_error, "Group " + group.name() + " has sale rate target less then zero. Not implemented in Flow" );
}
}
}
}
if (group.isProductionGroup()) {
if (!skip && group.isProductionGroup()) {
Group::ProductionCMode newControl = checkGroupProductionConstraints(group, deferred_logger);
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
const auto controls = group.productionControls(summaryState);
const Group::ProductionCMode& currentControl = well_state.currentProductionGroupControl(group.name());
if (group.has_control(Group::ProductionCMode::ORAT))
if (newControl != Group::ProductionCMode::NONE)
{
if (checkCurrentControl || currentControl != Group::ProductionCMode::ORAT)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.oil_target < current_rate ) {
actionOnBrokenConstraints(group, controls.exceed_action, Group::ProductionCMode::ORAT, reportStepIdx, deferred_logger);
}
}
switched_groups.insert(group.name());
actionOnBrokenConstraints(group, controls.exceed_action, newControl, deferred_logger);
}
if (group.has_control(Group::ProductionCMode::WRAT))
{
if (checkCurrentControl || currentControl != Group::ProductionCMode::WRAT)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.water_target < current_rate ) {
actionOnBrokenConstraints(group, controls.exceed_action, Group::ProductionCMode::WRAT, reportStepIdx, deferred_logger);
}
}
}
if (group.has_control(Group::ProductionCMode::GRAT))
{
if (checkCurrentControl || currentControl != Group::ProductionCMode::GRAT)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.gas_target < current_rate ) {
actionOnBrokenConstraints(group, controls.exceed_action, Group::ProductionCMode::GRAT, reportStepIdx, deferred_logger);
}
}
}
if (group.has_control(Group::ProductionCMode::LRAT))
{
if (checkCurrentControl || currentControl != Group::ProductionCMode::LRAT)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false);
current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.liquid_target < current_rate ) {
actionOnBrokenConstraints(group, controls.exceed_action, Group::ProductionCMode::LRAT, reportStepIdx, deferred_logger);
}
}
}
if (group.has_control(Group::ProductionCMode::CRAT))
{
OPM_DEFLOG_THROW(std::runtime_error, "Group " + group.name() + "CRAT control for production groups not implemented" , deferred_logger);
}
if (group.has_control(Group::ProductionCMode::RESV))
{
if (checkCurrentControl || currentControl != Group::ProductionCMode::RESV)
{
double current_rate = 0.0;
current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true);
current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true);
current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.resv_target < current_rate ) {
actionOnBrokenConstraints(group, controls.exceed_action, Group::ProductionCMode::RESV, reportStepIdx, deferred_logger);
}
}
}
if (group.has_control(Group::ProductionCMode::PRBL))
{
OPM_DEFLOG_THROW(std::runtime_error, "Group " + group.name() + "PRBL control for production groups not implemented", deferred_logger);
}
} else {
//neither production or injecting group FIELD?
}
// call recursively down the group hiearchy
for (const std::string& groupName : group.groups()) {
checkGroupConstraints( schedule().getGroup(groupName, reportStepIdx), checkCurrentControl, deferred_logger);
updateGroupIndividualControl( schedule().getGroup(groupName, reportStepIdx), deferred_logger, switched_groups);
}
}
template<typename TypeTag>
bool
BlackoilWellModel<TypeTag>::
checkGroupConstraints(const Group& group, Opm::DeferredLogger& deferred_logger) const {
const int reportStepIdx = ebosSimulator_.episodeIndex();
if (group.isInjectionGroup()) {
const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS};
for (Phase phase : all) {
if (!group.hasInjectionControl(phase)) {
continue;
}
Group::InjectionCMode newControl = checkGroupInjectionConstraints(group, phase);
if (newControl != Group::InjectionCMode::NONE) {
return true;
}
}
}
if (group.isProductionGroup()) {
Group::ProductionCMode newControl = checkGroupProductionConstraints(group, deferred_logger);
if (newControl != Group::ProductionCMode::NONE)
{
return true;
}
}
// call recursively down the group hiearchy
bool violated = false;
for (const std::string& groupName : group.groups()) {
violated = violated || checkGroupConstraints( schedule().getGroup(groupName, reportStepIdx), deferred_logger);
}
return violated;
}
template<typename TypeTag>
Group::ProductionCMode
BlackoilWellModel<TypeTag>::
checkGroupProductionConstraints(const Group& group, Opm::DeferredLogger& deferred_logger) const {
const int reportStepIdx = ebosSimulator_.episodeIndex();
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
const auto& comm = ebosSimulator_.vanguard().grid().comm();
const auto& well_state = well_state_;
const auto controls = group.productionControls(summaryState);
const Group::ProductionCMode& currentControl = well_state.currentProductionGroupControl(group.name());
if (group.has_control(Group::ProductionCMode::ORAT))
{
if (currentControl != Group::ProductionCMode::ORAT)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.oil_target < current_rate ) {
return Group::ProductionCMode::ORAT;
}
}
}
if (group.has_control(Group::ProductionCMode::WRAT))
{
if (currentControl != Group::ProductionCMode::WRAT)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.water_target < current_rate ) {
return Group::ProductionCMode::WRAT;
}
}
}
if (group.has_control(Group::ProductionCMode::GRAT))
{
if (currentControl != Group::ProductionCMode::GRAT)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.gas_target < current_rate ) {
return Group::ProductionCMode::GRAT;
}
}
}
if (group.has_control(Group::ProductionCMode::LRAT))
{
if (currentControl != Group::ProductionCMode::LRAT)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false);
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.liquid_target < current_rate ) {
return Group::ProductionCMode::LRAT;
}
}
}
if (group.has_control(Group::ProductionCMode::CRAT))
{
OPM_DEFLOG_THROW(std::runtime_error, "Group " + group.name() + "CRAT control for production groups not implemented" , deferred_logger);
}
if (group.has_control(Group::ProductionCMode::RESV))
{
if (currentControl != Group::ProductionCMode::RESV)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true);
current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true);
current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.resv_target < current_rate ) {
return Group::ProductionCMode::RESV;
}
}
}
if (group.has_control(Group::ProductionCMode::PRBL))
{
OPM_DEFLOG_THROW(std::runtime_error, "Group " + group.name() + "PRBL control for production groups not implemented", deferred_logger);
}
return Group::ProductionCMode::NONE;
}
template<typename TypeTag>
Group::InjectionCMode
BlackoilWellModel<TypeTag>::
checkGroupInjectionConstraints(const Group& group, const Phase& phase) const {
const int reportStepIdx = ebosSimulator_.episodeIndex();
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
const auto& comm = ebosSimulator_.vanguard().grid().comm();
const auto& well_state = well_state_;
int phasePos;
if (phase == Phase::GAS && phase_usage_.phase_used[BlackoilPhases::Vapour] )
phasePos = phase_usage_.phase_pos[BlackoilPhases::Vapour];
else if (phase == Phase::OIL && phase_usage_.phase_used[BlackoilPhases::Liquid])
phasePos = phase_usage_.phase_pos[BlackoilPhases::Liquid];
else if (phase == Phase::WATER && phase_usage_.phase_used[BlackoilPhases::Aqua] )
phasePos = phase_usage_.phase_pos[BlackoilPhases::Aqua];
else
OPM_THROW(std::runtime_error, "Unknown phase" );
const auto& controls = group.injectionControls(phase, summaryState);
const Group::InjectionCMode& currentControl = well_state.currentInjectionGroupControl(phase, group.name());
if (controls.has_control(Group::InjectionCMode::RATE))
{
if (currentControl != Group::InjectionCMode::RATE)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.surface_max_rate < current_rate) {
return Group::InjectionCMode::RATE;
}
}
}
if (controls.has_control(Group::InjectionCMode::RESV))
{
if (currentControl != Group::InjectionCMode::RESV)
{
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.resv_max_rate < current_rate) {
return Group::InjectionCMode::RESV;
}
}
}
if (controls.has_control(Group::InjectionCMode::REIN))
{
if (currentControl != Group::InjectionCMode::REIN)
{
double production_Rate = 0.0;
const Group& groupRein = schedule().getGroup(controls.reinj_group, reportStepIdx);
production_Rate += WellGroupHelpers::sumWellRates(groupRein, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/false);
// sum over all nodes
production_Rate = comm.sum(production_Rate);
double current_rate = 0.0;
current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true);
// sum over all nodes
current_rate = comm.sum(current_rate);
if (controls.target_reinj_fraction*production_Rate < current_rate) {
return Group::InjectionCMode::REIN;
}
}
}
if (controls.has_control(Group::InjectionCMode::VREP))
{
if (currentControl != Group::InjectionCMode::VREP)
{
double voidage_rate = 0.0;
const Group& groupVoidage = schedule().getGroup(controls.voidage_group, reportStepIdx);
voidage_rate += WellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false);
voidage_rate += WellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false);
voidage_rate += WellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false);
// sum over all nodes
voidage_rate = comm.sum(voidage_rate);
double total_rate = 0.0;
total_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true);
total_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true);
total_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true);
// sum over all nodes
total_rate = comm.sum(total_rate);
if (controls.target_void_fraction*voidage_rate < total_rate) {
return Group::InjectionCMode::VREP;
}
}
}
// Handle GCONSALE
if (schedule().gConSale(reportStepIdx).has(group.name())) {
if (controls.phase != Phase::GAS)
OPM_THROW(std::runtime_error, "Group " + group.name() + " has GCONSALE control but is not a GAS group" );
const auto& gconsale = schedule().gConSale(reportStepIdx).get(group.name(), summaryState);
double sales_rate = 0.0;
int gasPos = phase_usage_.phase_pos[BlackoilPhases::Vapour];
sales_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/false);
sales_rate -= WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/true);
// sum over all nodes
sales_rate = comm.sum(sales_rate);
// add import rate and substract consumption rate for group for gas
if (schedule().gConSump(reportStepIdx).has(group.name())) {
const auto& gconsump = schedule().gConSump(reportStepIdx).get(group.name(), summaryState);
if (phase_usage_.phase_used[BlackoilPhases::Vapour]) {
sales_rate += gconsump.import_rate;
sales_rate -= gconsump.consumption_rate;
}
}
if (sales_rate > gconsale.max_sales_rate) {
OPM_THROW(std::runtime_error, "Group " + group.name() + " has sale rate more then the maximum permitted value. Not implemented in Flow" );
}
if (sales_rate < gconsale.min_sales_rate) {
OPM_THROW(std::runtime_error, "Group " + group.name() + " has sale rate less then minimum permitted value. Not implemented in Flow" );
}
if (gconsale.sales_target < 0.0) {
OPM_THROW(std::runtime_error, "Group " + group.name() + " has sale rate target less then zero. Not implemented in Flow" );
}
}
return Group::InjectionCMode::NONE;
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
actionOnBrokenConstraints(const Group& group, const Group::ExceedAction& exceed_action, const Group::ProductionCMode& newControl, const int reportStepIdx, Opm::DeferredLogger& deferred_logger) {
actionOnBrokenConstraints(const Group& group, const Group::ExceedAction& exceed_action, const Group::ProductionCMode& newControl, Opm::DeferredLogger& deferred_logger) {
auto& well_state = well_state_;
const Group::ProductionCMode& oldControl = well_state.currentProductionGroupControl(group.name());
@ -1984,9 +2103,6 @@ namespace Opm {
well_state.setCurrentProductionGroupControl(group.name(), newControl);
ss << "Switching control mode for group "<< group.name() << " to " << Group::ProductionCMode2String(newControl);
}
// Pass a dummy phase for producer groups. The topUpPhase is only relevant for injector groups
const Phase topUpPhase = Phase::WATER;
wellGroupHelpers::setGroupControl(group, schedule(), topUpPhase, reportStepIdx, false, well_state, ss);
break;
}
default:
@ -2007,7 +2123,7 @@ namespace Opm {
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
actionOnBrokenConstraints(const Group& group, const Group::InjectionCMode& newControl, const Phase& controlPhase, const int reportStepIdx, Opm::DeferredLogger& deferred_logger) {
actionOnBrokenConstraints(const Group& group, const Group::InjectionCMode& newControl, const Phase& controlPhase, Opm::DeferredLogger& deferred_logger) {
auto& well_state = well_state_;
const Group::InjectionCMode& oldControl = well_state.currentInjectionGroupControl(controlPhase, group.name());
@ -2023,13 +2139,130 @@ namespace Opm {
}
well_state.setCurrentInjectionGroupControl(controlPhase, group.name(), newControl);
}
wellGroupHelpers::setGroupControl(group, schedule(), controlPhase, reportStepIdx, /*isInjector*/true, well_state, ss);
if (!ss.str().empty())
deferred_logger.info(ss.str());
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
updateGroupHigherControls(Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups)
{
const int reportStepIdx = ebosSimulator_.episodeIndex();
const Group& fieldGroup = schedule().getGroup("FIELD", reportStepIdx);
checkGroupHigherConstraints(fieldGroup, deferred_logger, switched_groups);
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
checkGroupHigherConstraints(const Group& group, Opm::DeferredLogger& deferred_logger, std::set<std::string>& switched_groups)
{
// Set up coefficients for RESV <-> surface rate conversion.
// Use the pvtRegionIdx from the top cell of the first well.
// TODO fix this!
// This is only used for converting RESV rates.
// What is the proper approach?
const int fipnum = 0;
const int pvtreg = well_perf_data_.empty()
? pvt_region_idx_[0]
: pvt_region_idx_[well_perf_data_[0][0].cell_index];
std::vector<double> resv_coeff(phase_usage_.num_phases, 0.0);
rateConverter_->calcCoeff(fipnum, pvtreg, resv_coeff);
const int reportStepIdx = ebosSimulator_.episodeIndex();
const auto& summaryState = ebosSimulator_.vanguard().summaryState();
std::vector<double> rates(phase_usage_.num_phases, 0.0);
const auto& comm = ebosSimulator_.vanguard().grid().comm();
const bool skip = switched_groups.count(group.name()) || group.name() == "FIELD";
if (!skip && group.isInjectionGroup()) {
// Obtain rates for group.
for (int phasePos = 0; phasePos < phase_usage_.num_phases; ++phasePos) {
const double local_current_rate = WellGroupHelpers::sumWellRates(
group, schedule(), well_state_, reportStepIdx, phasePos, /* isInjector */ true);
// Sum over all processes
rates[phasePos] = comm.sum(local_current_rate);
}
const Phase all[] = { Phase::WATER, Phase::OIL, Phase::GAS };
for (Phase phase : all) {
// Check higher up only if under individual (not FLD) control.
const Group::InjectionCMode& currentControl = well_state_.currentInjectionGroupControl(phase, group.name());
if (currentControl != Group::InjectionCMode::FLD) {
const Group& parentGroup = schedule().getGroup(group.parent(), reportStepIdx);
const std::pair<bool, double> changed = WellGroupHelpers::checkGroupConstraintsInj(
group.name(),
group.parent(),
parentGroup,
well_state_,
reportStepIdx,
guideRate_.get(),
rates.data(),
phase,
phase_usage_,
group.getGroupEfficiencyFactor(),
schedule(),
summaryState,
resv_coeff,
deferred_logger);
if (changed.first) {
switched_groups.insert(group.name());
actionOnBrokenConstraints(group, Group::InjectionCMode::FLD, phase, deferred_logger);
}
}
}
}
if (!skip && group.isProductionGroup()) {
// Obtain rates for group.
for (int phasePos = 0; phasePos < phase_usage_.num_phases; ++phasePos) {
const double local_current_rate = WellGroupHelpers::sumWellRates(
group, schedule(), well_state_, reportStepIdx, phasePos, /* isInjector */ false);
// Sum over all processes
rates[phasePos] = -comm.sum(local_current_rate);
}
// Check higher up only if under individual (not FLD) control.
const Group::ProductionCMode& currentControl = well_state_.currentProductionGroupControl(group.name());
if (currentControl != Group::ProductionCMode::FLD) {
const Group& parentGroup = schedule().getGroup(group.parent(), reportStepIdx);
const std::pair<bool, double> changed = WellGroupHelpers::checkGroupConstraintsProd(
group.name(),
group.parent(),
parentGroup,
well_state_,
reportStepIdx,
guideRate_.get(),
rates.data(),
phase_usage_,
group.getGroupEfficiencyFactor(),
schedule(),
summaryState,
resv_coeff,
deferred_logger);
if (changed.first) {
switched_groups.insert(group.name());
const auto exceed_action = group.productionControls(summaryState).exceed_action;
actionOnBrokenConstraints(group, exceed_action, Group::ProductionCMode::FLD, deferred_logger);
}
}
}
// call recursively down the group hiearchy
for (const std::string& groupName : group.groups()) {
checkGroupHigherConstraints( schedule().getGroup(groupName, reportStepIdx), deferred_logger, switched_groups);
}
}
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
@ -2045,8 +2278,8 @@ namespace Opm {
const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(Phase::GAS, group.name());
if( currentGroupControl == Group::InjectionCMode::REIN ) {
int gasPos = phase_usage_.phase_pos[BlackoilPhases::Vapour];
double gasProductionRate = wellGroupHelpers::sumWellRates(group, schedule, wellState, reportStepIdx, gasPos, /*isInjector*/false);
double solventProductionRate = wellGroupHelpers::sumSolventRates(group, schedule, wellState, reportStepIdx, /*isInjector*/false);
double gasProductionRate = WellGroupHelpers::sumWellRates(group, schedule, wellState, reportStepIdx, gasPos, /*isInjector*/false);
double solventProductionRate = WellGroupHelpers::sumSolventRates(group, schedule, wellState, reportStepIdx, /*isInjector*/false);
const auto& comm = ebosSimulator_.vanguard().grid().comm();
solventProductionRate = comm.sum(solventProductionRate);

12
opm/simulators/wells/MultisegmentWell.hpp
View File

@ -283,11 +283,10 @@ namespace Opm
void initMatrixAndVectors(const int num_cells) const;
// protected functions
// EvalWell getBhp(); this one should be something similar to getSegmentPressure();
// EvalWell getQs(); this one should be something similar to getSegmentRates()
// EValWell wellVolumeFractionScaled, wellVolumeFraction, wellSurfaceVolumeFraction ... these should have different names, and probably will be needed.
// bool crossFlowAllowed(const Simulator& ebosSimulator) const; probably will be needed
EvalWell getBhp() const;
EvalWell getQs(const int comp_idx) const;
EvalWell getWQTotal() const;
// xw = inv(D)*(rw - C*x)
void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
@ -381,9 +380,6 @@ namespace Opm
const Well::ProductionControls& prod_controls,
Opm::DeferredLogger& deferred_logger);
void assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger);
void assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger);
void assemblePressureEq(const int seg) const;
// hytrostatic pressure loss

559
opm/simulators/wells/MultisegmentWell_impl.hpp
View File

@ -675,7 +675,7 @@ namespace Opm
}
// If the well is pressure controlled the potential equals the rate.
{
/* {
bool pressure_controlled_well = false;
if (this->isInjector()) {
const Opm::Well::InjectorCMode& current = well_state.currentInjectionControls()[index_of_well_];
@ -690,12 +690,12 @@ namespace Opm
}
if (pressure_controlled_well) {
for (int compIdx = 0; compIdx < num_components_; ++compIdx) {
const EvalWell rate = this->getSegmentRate(0, compIdx);
const EvalWell rate = this->getQs(compIdx);
well_potentials[ebosCompIdxToFlowCompIdx(compIdx)] = rate.value();
}
return;
}
}
} */
// creating a copy of the well itself, to avoid messing up the explicit informations
// during this copy, the only information not copied properly is the well controls
@ -790,7 +790,7 @@ namespace Opm
const int np = number_of_phases_;
well_flux.resize(np, 0.0);
for (int compIdx = 0; compIdx < num_components_; ++compIdx) {
const EvalWell rate = well_copy.getSegmentRate(0, compIdx);
const EvalWell rate = well_copy.getQs(compIdx);
well_flux[ebosCompIdxToFlowCompIdx(compIdx)] = rate.value();
}
debug_cost_counter_ += well_copy.debug_cost_counter_;
@ -1621,6 +1621,18 @@ namespace Opm
template <typename TypeTag>
typename MultisegmentWell<TypeTag>::EvalWell
MultisegmentWell<TypeTag>::
getBhp() const
{
return getSegmentPressure(0);
}
template <typename TypeTag>
typename MultisegmentWell<TypeTag>::EvalWell
MultisegmentWell<TypeTag>::
@ -1634,6 +1646,18 @@ namespace Opm
template <typename TypeTag>
typename MultisegmentWell<TypeTag>::EvalWell
MultisegmentWell<TypeTag>::
getQs(const int comp_idx) const
{
return getSegmentRate(0, comp_idx);
}
template <typename TypeTag>
typename MultisegmentWell<TypeTag>::EvalWell
MultisegmentWell<TypeTag>::
@ -1690,6 +1714,18 @@ namespace Opm
template <typename TypeTag>
typename MultisegmentWell<TypeTag>::EvalWell
MultisegmentWell<TypeTag>::
getWQTotal() const
{
return getSegmentGTotal(0);
}
template <typename TypeTag>
void
MultisegmentWell<TypeTag>::
@ -1759,20 +1795,28 @@ namespace Opm
EvalWell control_eq(0.0);
const auto& well = well_ecl_;
const int well_index = index_of_well_;
double efficiencyFactor = well.getEfficiencyFactor();
auto getRates = [&]() {
std::vector<EvalWell> rates(3, 0.0);
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
rates[Water] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
rates[Oil] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
rates[Gas] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
}
return rates;
};
if (wellIsStopped_) {
control_eq = getSegmentGTotal(0);
control_eq = getWQTotal();
} else if (this->isInjector() ) {
const Opm::Well::InjectorCMode& current = well_state.currentInjectionControls()[well_index];
const auto& controls = inj_controls;
InjectorType injectorType = controls.injector_type;
// Find scaling factor to get injection rate,
const InjectorType injectorType = inj_controls.injector_type;
double scaling = 1.0;
const auto& pu = phaseUsage();
switch (injectorType) {
case InjectorType::WATER:
{
@ -1792,195 +1836,23 @@ namespace Opm
default:
throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
}
switch(current) {
case Well::InjectorCMode::RATE:
{
control_eq = getSegmentGTotal(0) / scaling - controls.surface_rate;
break;
}
case Well::InjectorCMode::RESV:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
double coeff = 1.0;
switch (injectorType) {
case InjectorType::WATER:
{
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Aqua]];
break;
}
case InjectorType::OIL:
{
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Liquid]];
break;
}
case InjectorType::GAS:
{
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Vapour]];
break;
}
default:
throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
}
control_eq = coeff*getSegmentGTotal(0) / scaling - controls.reservoir_rate;
break;
}
case Well::InjectorCMode::THP:
{
std::vector<EvalWell> rates(3, 0.);
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
rates[ Water ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
rates[ Oil ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
rates[ Gas ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
}
control_eq = getSegmentPressure(0) - calculateBhpFromThp(rates, well, summaryState, deferred_logger);
break;
}
case Well::InjectorCMode::BHP:
{
const auto& bhp = controls.bhp_limit;
control_eq = getSegmentPressure(0) - bhp;
break;
}
case Well::InjectorCMode::GRUP:
{
assert(well.isAvailableForGroupControl());
const auto& group = schedule.getGroup( well.groupName(), current_step_ );
assembleGroupInjectionControl(group, well_state, schedule, summaryState, controls.injector_type, control_eq, efficiencyFactor, deferred_logger);
break;
}
case Well::InjectorCMode::CMODE_UNDEFINED:
{
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger);
}
}
}
//Producer
else
{
const Well::ProducerCMode& current = well_state.currentProductionControls()[well_index];
const auto& controls = prod_controls;
switch (current) {
case Well::ProducerCMode::ORAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - controls.oil_rate;
break;
}
case Well::ProducerCMode::WRAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
control_eq = rate - controls.water_rate;
break;
}
case Well::ProducerCMode::GRAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx));
const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
control_eq = rate - controls.gas_rate;
break;
}
case Well::ProducerCMode::LRAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx))
-getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - controls.liquid_rate;
break;
}
case Well::ProducerCMode::CRAT:
{
OPM_DEFLOG_THROW(std::runtime_error, "CRAT control not supported " << name(), deferred_logger);
}
case Well::ProducerCMode::RESV:
{
EvalWell total_rate(0.); // reservoir rate
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
for (int phase = 0; phase < number_of_phases_; ++phase) {
total_rate += getSegmentRate(0, flowPhaseToEbosCompIdx(phase) ) * convert_coeff[phase];
}
if (controls.prediction_mode) {
control_eq = total_rate - controls.resv_rate;
} else {
std::vector<double> hrates(number_of_phases_,0.);
const PhaseUsage& pu = phaseUsage();
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
hrates[pu.phase_pos[Water]] = controls.water_rate;
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
hrates[pu.phase_pos[Oil]] = controls.oil_rate;
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
hrates[pu.phase_pos[Gas]] = controls.gas_rate;
}
std::vector<double> hrates_resv(number_of_phases_,0.);
Base::rateConverter_.calcReservoirVoidageRates(/*fipreg*/ 0, Base::pvtRegionIdx_, hrates, hrates_resv);
double target = -std::accumulate(hrates_resv.begin(), hrates_resv.end(), 0.0);
control_eq = total_rate - target;
}
break;
}
case Well::ProducerCMode::BHP:
{
control_eq = getSegmentPressure(0) - controls.bhp_limit;
break;
}
case Well::ProducerCMode::THP:
{
std::vector<EvalWell> rates(3, 0.);
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
rates[ Water ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
rates[ Oil ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
rates[ Gas ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
}
control_eq = getSegmentPressure(0) - calculateBhpFromThp(rates, well, summaryState, deferred_logger);
break;
}
case Well::ProducerCMode::GRUP:
{
assert(well.isAvailableForGroupControl());
const auto& group = schedule.getGroup( well.groupName(), current_step_ );
assembleGroupProductionControl(group, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
break;
}
case Well::ProducerCMode::CMODE_UNDEFINED:
{
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger );
}
case Well::ProducerCMode::NONE:
{
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger );
}
}
const EvalWell injection_rate = getWQTotal() / scaling;
// Setup function for evaluation of BHP from THP (used only if needed).
auto bhp_from_thp = [&]() {
const auto rates = getRates();
return calculateBhpFromThp(rates, well, summaryState, deferred_logger);
};
// Call generic implementation.
Base::assembleControlEqInj(well_state, schedule, summaryState, inj_controls, getBhp(), injection_rate, bhp_from_thp, control_eq, deferred_logger);
} else {
// Find rates.
const auto rates = getRates();
// Setup function for evaluation of BHP from THP (used only if needed).
auto bhp_from_thp = [&]() {
return calculateBhpFromThp(rates, well, summaryState, deferred_logger);
};
// Call generic implementation.
Base::assembleControlEqProd(well_state, schedule, summaryState, prod_controls, getBhp(), rates, bhp_from_thp, control_eq, deferred_logger);
}
// using control_eq to update the matrix and residuals
@ -2105,288 +1977,8 @@ namespace Opm
else {
OPM_DEFLOG_THROW(std::logic_error, "Expected INJECTOR or PRODUCER well", deferred_logger);
}
}
template <typename TypeTag>
void
MultisegmentWell<TypeTag>::
assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger)
{
const auto& well = well_ecl_;
const auto& pu = phaseUsage();
int phasePos;
Well::GuideRateTarget wellTarget;
Phase injectionPhase;
double scaling = 1.0;
switch (injectorType) {
case InjectorType::WATER:
{
phasePos = pu.phase_pos[BlackoilPhases::Aqua];
wellTarget = Well::GuideRateTarget::WAT;
injectionPhase = Phase::WATER;
scaling = scalingFactor(pu.phase_pos[BlackoilPhases::Aqua]);
break;
}
case InjectorType::OIL:
{
phasePos = pu.phase_pos[BlackoilPhases::Liquid];
wellTarget = Well::GuideRateTarget::OIL;
injectionPhase = Phase::OIL;
scaling = scalingFactor(pu.phase_pos[BlackoilPhases::Liquid]);
break;
}
case InjectorType::GAS:
{
phasePos = pu.phase_pos[BlackoilPhases::Vapour];
wellTarget = Well::GuideRateTarget::GAS;
injectionPhase = Phase::GAS;
scaling = scalingFactor(pu.phase_pos[BlackoilPhases::Vapour]);
break;
}
default:
throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
}
const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name());
if (currentGroupControl == Group::InjectionCMode::FLD) {
// Inject share of parents control
const auto& parent = schedule.getGroup( group.parent(), current_step_ );
if (group.getTransferGroupEfficiencyFactor())
efficiencyFactor *= group.getGroupEfficiencyFactor();
assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger);
return;
}
if (!group.isInjectionGroup() || currentGroupControl == Group::InjectionCMode::NONE) {
// use bhp as control eq and let the updateControl code find a valid control
const auto& controls = well.injectionControls(summaryState);
control_eq = getSegmentPressure(0) - controls.bhp_limit;
return;
}
assert(group.hasInjectionControl(injectionPhase));
const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState);
const std::vector<double>& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name());
double groupTargetReduction = groupInjectionReductions[phasePos];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, wellTarget, /*isInjector*/true);
wellGroupHelpers::accumulateGroupInjectionPotentialFractions(well.groupName(), group.name(), schedule, well_state, pu, current_step_, injectionPhase, fraction);
switch(currentGroupControl) {
case Group::InjectionCMode::NONE:
{
// The NONE case is handled earlier
assert(false);
break;
}
case Group::InjectionCMode::RATE:
{
control_eq = getSegmentGTotal(0) / scaling - fraction * (groupcontrols.surface_max_rate / efficiencyFactor - groupTargetReduction);
break;
}
case Group::InjectionCMode::RESV:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
double coeff = convert_coeff[phasePos];
double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff/efficiencyFactor - groupTargetReduction));
control_eq = getSegmentGTotal(0) / scaling - fraction * target;
break;
}
case Group::InjectionCMode::REIN:
{
double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos];
productionRate /= efficiencyFactor;
double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction));
control_eq = getSegmentGTotal(0) / scaling - fraction * target;
break;
}
case Group::InjectionCMode::VREP:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
double coeff = convert_coeff[phasePos];
double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction;
double injReduction = 0.0;
std::vector<double> groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name());
if (groupcontrols.phase != Phase::WATER)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]];
if (groupcontrols.phase != Phase::OIL)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]];
if (groupcontrols.phase != Phase::GAS)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]];
voidageRate -= injReduction;
voidageRate /= efficiencyFactor;
double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction));
control_eq = getSegmentGTotal(0) / scaling - fraction * target;
break;
}
case Group::InjectionCMode::FLD:
{
// The FLD case is handled earlier
assert(false);
break;
}
case Group::InjectionCMode::SALE:
{
// only for gas injectors
assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]);
// Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate;
double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos];
if (schedule.gConSump(current_step_).has(group.name())) {
const auto& gconsump = schedule.gConSump(current_step_).get(group.name(), summaryState);
if (pu.phase_used[BlackoilPhases::Vapour]) {
inj_rate += gconsump.import_rate;
inj_rate -= gconsump.consumption_rate;
}
}
const auto& gconsale = schedule.gConSale(current_step_).get(group.name(), summaryState);
inj_rate -= gconsale.sales_target;
inj_rate /= efficiencyFactor;
double target = std::max(0.0, (inj_rate - groupTargetReduction));
control_eq = getSegmentGTotal(0) /scaling - fraction * target;
break;
}
default:
OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger );
}
}
template <typename TypeTag>
void
MultisegmentWell<TypeTag>::
assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger)
{
const auto& well = well_ecl_;
const auto pu = phaseUsage();
const Group::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name());
if (currentGroupControl == Group::ProductionCMode::FLD ) {
// Produce share of parents control
const auto& parent = schedule.getGroup( group.parent(), current_step_ );
if (group.getTransferGroupEfficiencyFactor())
efficiencyFactor *= group.getGroupEfficiencyFactor();
assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
return;
}
if (!group.isProductionGroup() || currentGroupControl == Group::ProductionCMode::NONE) {
// use bhp as control eq and let the updateControl code find a vallied control
const auto& controls = well.productionControls(summaryState);
control_eq = getSegmentPressure(0) - controls.bhp_limit;
return;
}
const auto& groupcontrols = group.productionControls(summaryState);
const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name());
switch(currentGroupControl) {
case Group::ProductionCMode::NONE:
{
// The NONE case is handled earlier
assert(false);
break;
}
case Group::ProductionCMode::ORAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::OIL, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::OIL, fraction);
const double rate_target = std::max(0.0, groupcontrols.oil_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::WRAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Water]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::WAT, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::WAT, fraction);
const double rate_target = std::max(0.0, groupcontrols.water_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::GRAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Gas]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::GAS, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::GAS, fraction);
const double rate_target = std::max(0.0, groupcontrols.gas_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::gasCompIdx));
const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::LRAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] + groupTargetReductions[pu.phase_pos[Water]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::LIQ, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::LIQ, fraction);
const double rate_target = std::max(0.0, groupcontrols.liquid_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx))
-getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::CRAT:
{
OPM_DEFLOG_THROW(std::runtime_error, "CRAT group control not implemented for producers", deferred_logger );
break;
}
case Group::ProductionCMode::RESV:
{
OPM_DEFLOG_THROW(std::runtime_error, "RESV group control not implemented for producers", deferred_logger );
break;
}
case Group::ProductionCMode::PRBL:
{
OPM_DEFLOG_THROW(std::runtime_error, "PRBL group control not implemented for producers", deferred_logger );
break;
}
case Group::ProductionCMode::FLD:
{
// The FLD case is handled earlier
assert(false);
break;
}
default:
OPM_DEFLOG_THROW(std::runtime_error, "Unvallied group control specified for group " + well.groupName(), deferred_logger );
}
}
template <typename TypeTag>
@ -4096,5 +3688,4 @@ namespace Opm
const double sign = mass_rate <= 0. ? 1.0 : -1.0;
return sign * (friction_pressure_loss + constriction_pressure_loss);
}
}
}

8
opm/simulators/wells/StandardWell.hpp
View File

@ -410,10 +410,10 @@ namespace Opm
void updateThp(WellState& well_state, Opm::DeferredLogger& deferred_logger) const;
void assembleControlEq(const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, Opm::DeferredLogger& deferred_logger);
void assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger);
void assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger);
void assembleControlEq(const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
Opm::DeferredLogger& deferred_logger);
// handle the non reasonable fractions due to numerical overshoot
void processFractions() const;

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

@ -764,209 +764,53 @@ namespace Opm
template <typename TypeTag>
void
StandardWell<TypeTag>::
assembleControlEq(const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, Opm::DeferredLogger& deferred_logger)
StandardWell<TypeTag>::assembleControlEq(const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
Opm::DeferredLogger& deferred_logger)
{
EvalWell control_eq(numWellEq_ + numEq, 0.);
EvalWell control_eq(numWellEq_ + numEq, 0.0);
const auto& well = well_ecl_;
const int well_index = index_of_well_;
double efficiencyFactor = well.getEfficiencyFactor();
auto getRates = [&]() {
std::vector<EvalWell> rates(3, EvalWell(numWellEq_ + numEq, 0.0));
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
rates[Water] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
rates[Oil] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
rates[Gas] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
}
return rates;
};
if (wellIsStopped_) {
control_eq = getWQTotal();
} else if (this->isInjector()) {
const Opm::Well::InjectorCMode& current = well_state.currentInjectionControls()[well_index];
const auto& controls = well.injectionControls(summaryState);
switch(current) {
case Well::InjectorCMode::RATE:
{
control_eq = getWQTotal() - controls.surface_rate;
break;
}
case Well::InjectorCMode::RESV:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
const auto& pu = phaseUsage();
InjectorType injectorType = controls.injector_type;
double coeff = 1.0;
switch (injectorType) {
case InjectorType::WATER:
{
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Aqua]];
break;
}
case InjectorType::OIL:
{
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Liquid]];
break;
}
case InjectorType::GAS:
{
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Vapour]];
break;
}
default:
throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
}
control_eq = coeff*getWQTotal() - controls.reservoir_rate;
break;
}
case Well::InjectorCMode::THP:
{
std::vector<EvalWell> rates(3, {numWellEq_ + numEq, 0.});
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
rates[ Water ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
rates[ Oil ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
rates[ Gas ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
}
control_eq = getBhp() - calculateBhpFromThp(rates, well, summaryState, deferred_logger);
break;
}
case Well::InjectorCMode::BHP:
{
const auto& bhp = controls.bhp_limit;
control_eq = getBhp() - bhp;
break;
}
case Well::InjectorCMode::GRUP:
{
assert(well.isAvailableForGroupControl());
const auto& group = schedule.getGroup( well.groupName(), current_step_ );
assembleGroupInjectionControl(group, well_state, schedule, summaryState, controls.injector_type, control_eq, efficiencyFactor, deferred_logger);
break;
}
case Well::InjectorCMode::CMODE_UNDEFINED:
{
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name() , deferred_logger);
}
}
// Find injection rate.
const EvalWell injection_rate = getWQTotal();
// Setup function for evaluation of BHP from THP (used only if needed).
auto bhp_from_thp = [&]() {
const auto rates = getRates();
return calculateBhpFromThp(rates, well, summaryState, deferred_logger);
};
// Call generic implementation.
const auto& inj_controls = well.injectionControls(summaryState);
Base::assembleControlEqInj(well_state, schedule, summaryState, inj_controls, getBhp(), injection_rate, bhp_from_thp, control_eq, deferred_logger);
} else {
// Find rates.
const auto rates = getRates();
// Setup function for evaluation of BHP from THP (used only if needed).
auto bhp_from_thp = [&]() {
return calculateBhpFromThp(rates, well, summaryState, deferred_logger);
};
// Call generic implementation.
const auto& prod_controls = well.productionControls(summaryState);
Base::assembleControlEqProd(well_state, schedule, summaryState, prod_controls, getBhp(), rates, bhp_from_thp, control_eq, deferred_logger);
}
//Producer
else
{
const Well::ProducerCMode& current = well_state.currentProductionControls()[well_index];
const auto& controls = well.productionControls(summaryState);
switch (current) {
case Well::ProducerCMode::ORAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - controls.oil_rate;
break;
}
case Well::ProducerCMode::WRAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
control_eq = rate - controls.water_rate;
break;
}
case Well::ProducerCMode::GRAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
control_eq = rate - controls.gas_rate;
break;
}
case Well::ProducerCMode::LRAT:
{
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx))
- getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - controls.liquid_rate;
break;
}
case Well::ProducerCMode::CRAT:
{
OPM_DEFLOG_THROW(std::runtime_error, "CRAT control not supported " << name(), deferred_logger);
}
case Well::ProducerCMode::RESV:
{
EvalWell total_rate(numWellEq_ + numEq, 0.); // reservoir rate
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
for (int phase = 0; phase < number_of_phases_; ++phase) {
total_rate -= getQs( flowPhaseToEbosCompIdx(phase) ) * convert_coeff[phase];
}
if (controls.prediction_mode) {
control_eq = total_rate - controls.resv_rate;
} else {
std::vector<double> hrates(number_of_phases_,0.);
const PhaseUsage& pu = phaseUsage();
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
hrates[pu.phase_pos[Water]] = controls.water_rate;
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
hrates[pu.phase_pos[Oil]] = controls.oil_rate;
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
hrates[pu.phase_pos[Gas]] = controls.gas_rate;
}
std::vector<double> hrates_resv(number_of_phases_,0.);
Base::rateConverter_.calcReservoirVoidageRates(/*fipreg*/ 0, Base::pvtRegionIdx_, hrates, hrates_resv);
double target = std::accumulate(hrates_resv.begin(), hrates_resv.end(), 0.0);
control_eq = total_rate - target;
}
break;
}
case Well::ProducerCMode::BHP:
{
control_eq = getBhp() - controls.bhp_limit;
break;
}
case Well::ProducerCMode::THP:
{
std::vector<EvalWell> rates(3, {numWellEq_ + numEq, 0.});
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
rates[ Water ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
rates[ Oil ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
rates[ Gas ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
}
control_eq = getBhp() - calculateBhpFromThp(rates, well, summaryState, deferred_logger);
break;
}
case Well::ProducerCMode::GRUP:
{
assert(well.isAvailableForGroupControl());
const auto& group = schedule.getGroup( well.groupName(), current_step_ );
assembleGroupProductionControl(group, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
break;
}
case Well::ProducerCMode::CMODE_UNDEFINED:
{
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(),deferred_logger);
}
case Well::ProducerCMode::NONE:
{
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger);
}
}
}
// using control_eq to update the matrix and residuals
// TODO: we should use a different index system for the well equations
@ -976,293 +820,6 @@ namespace Opm
}
}
template <typename TypeTag>
void
StandardWell<TypeTag>::
assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger)
{
const auto& well = well_ecl_;
const auto pu = phaseUsage();
int phasePos;
Well::GuideRateTarget wellTarget;
Phase injectionPhase;
switch (injectorType) {
case InjectorType::WATER:
{
phasePos = pu.phase_pos[BlackoilPhases::Aqua];
wellTarget = Well::GuideRateTarget::WAT;
injectionPhase = Phase::WATER;
break;
}
case InjectorType::OIL:
{
phasePos = pu.phase_pos[BlackoilPhases::Liquid];
wellTarget = Well::GuideRateTarget::OIL;
injectionPhase = Phase::OIL;
break;
}
case InjectorType::GAS:
{
phasePos = pu.phase_pos[BlackoilPhases::Vapour];
wellTarget = Well::GuideRateTarget::GAS;
injectionPhase = Phase::GAS;
break;
}
default:
throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
}
const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name());
if (currentGroupControl == Group::InjectionCMode::FLD) {
// Inject share of parents control
const auto& parent = schedule.getGroup( group.parent(), current_step_ );
if (group.getTransferGroupEfficiencyFactor())
efficiencyFactor *= group.getGroupEfficiencyFactor();
assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger);
return;
}
if (!group.isInjectionGroup() || currentGroupControl == Group::InjectionCMode::NONE) {
// use bhp as control eq and let the updateControl code find a valid control
const auto& controls = well.injectionControls(summaryState);
control_eq = getBhp() - controls.bhp_limit;
return;
}
assert(group.hasInjectionControl(injectionPhase));
const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState);
const std::vector<double>& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name());
double groupTargetReduction = groupInjectionReductions[phasePos];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, wellTarget, /*isInjector*/true);
wellGroupHelpers::accumulateGroupInjectionPotentialFractions(well.groupName(), group.name(), schedule, well_state, pu, current_step_, injectionPhase, fraction);
switch(currentGroupControl) {
case Group::InjectionCMode::NONE:
{
// The NONE case is handled earlier
assert(false);
break;
}
case Group::InjectionCMode::RATE:
{
double target = std::max(0.0, (groupcontrols.surface_max_rate / efficiencyFactor - groupTargetReduction) );
control_eq = getWQTotal() - fraction * target;
break;
}
case Group::InjectionCMode::RESV:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
double coeff = convert_coeff[phasePos];
double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff/efficiencyFactor - groupTargetReduction));
control_eq = getWQTotal() - fraction * target;
break;
}
case Group::InjectionCMode::REIN:
{
double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos];
productionRate /= efficiencyFactor;
double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction));
control_eq = getWQTotal() - fraction * target;
break;
}
case Group::InjectionCMode::VREP:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
double coeff = convert_coeff[phasePos];
double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction;
double injReduction = 0.0;
std::vector<double> groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name());
if (groupcontrols.phase != Phase::WATER)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]];
if (groupcontrols.phase != Phase::OIL)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]];
if (groupcontrols.phase != Phase::GAS)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]];
voidageRate -= injReduction;
voidageRate /= efficiencyFactor;
double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction));
control_eq = getWQTotal() - fraction * target;
break;
}
case Group::InjectionCMode::FLD:
{
// The FLD case is handled earlier
assert(false);
break;
}
case Group::InjectionCMode::SALE:
{
// only for gas injectors
assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]);
// Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate;
double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos];
if (schedule.gConSump(current_step_).has(group.name())) {
const auto& gconsump = schedule.gConSump(current_step_).get(group.name(), summaryState);
if (pu.phase_used[BlackoilPhases::Vapour]) {
inj_rate += gconsump.import_rate;
inj_rate -= gconsump.consumption_rate;
}
}
const auto& gconsale = schedule.gConSale(current_step_).get(group.name(), summaryState);
inj_rate -= gconsale.sales_target;
inj_rate /= efficiencyFactor;
double target = std::max(0.0, (inj_rate - groupTargetReduction));
control_eq = getWQTotal() - fraction * target;
break;
}
default:
OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger );
}
}
template <typename TypeTag>
void
StandardWell<TypeTag>::
assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger)
{
const auto& well = well_ecl_;
const auto pu = phaseUsage();
const Group::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name());
if (currentGroupControl == Group::ProductionCMode::FLD) {
// Produce share of parents control
const auto& parent = schedule.getGroup( group.parent(), current_step_ );
if (group.getTransferGroupEfficiencyFactor())
efficiencyFactor *= group.getGroupEfficiencyFactor();
assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
return;
}
if (!group.isProductionGroup() || currentGroupControl == Group::ProductionCMode::NONE) {
// use bhp as control eq and let the updateControl code find a vallied control
const auto& controls = well.productionControls(summaryState);
control_eq = getBhp() - controls.bhp_limit;
return;
}
const auto& groupcontrols = group.productionControls(summaryState);
const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name());
switch(currentGroupControl) {
case Group::ProductionCMode::NONE:
{
// The NONE case is handled earlier
assert(false);
}
case Group::ProductionCMode::ORAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::OIL, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::OIL, fraction);
const double rate_target = std::max(0.0, groupcontrols.oil_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::WRAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Water]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::WAT, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::WAT, fraction);
const double rate_target = std::max(0.0, groupcontrols.water_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::GRAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Gas]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::GAS, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::GAS, fraction);
const double rate_target = std::max(0.0, groupcontrols.gas_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::gasCompIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::LRAT:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] + groupTargetReductions[pu.phase_pos[Water]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::LIQ, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::LIQ, fraction);
const double rate_target = std::max(0.0, groupcontrols.liquid_target / efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx))
- getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx));
control_eq = rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::CRAT:
{
OPM_DEFLOG_THROW(std::runtime_error, "CRAT group control not implemented for producers", deferred_logger );
break;
}
case Group::ProductionCMode::RESV:
{
double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]]
+ groupTargetReductions[pu.phase_pos[Gas]]
+ groupTargetReductions[pu.phase_pos[Water]];
double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well::GuideRateTarget::RES, /*isInjector*/false);
wellGroupHelpers::accumulateGroupFractionsFromGuideRates(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group::GuideRateTarget::RES, fraction);
EvalWell total_rate(numWellEq_ + numEq, 0.); // reservoir rate
std::vector<double> convert_coeff(number_of_phases_, 1.0);
Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
for (int phase = 0; phase < number_of_phases_; ++phase) {
total_rate -= getQs( flowPhaseToEbosCompIdx(phase) ) * convert_coeff[phase];
}
const double rate_target = std::max(0.0, groupcontrols.resv_target/ efficiencyFactor - groupTargetReduction);
assert(FluidSystem::phaseIsActive(FluidSystem::gasCompIdx));
control_eq = total_rate - fraction * rate_target;
break;
}
case Group::ProductionCMode::PRBL:
{
OPM_DEFLOG_THROW(std::runtime_error, "PRBL group control not implemented for producers", deferred_logger );
break;
}
case Group::ProductionCMode::FLD:
{
// The FLD case is handled earlier
assert(false);
break;
}
default:
OPM_DEFLOG_THROW(std::runtime_error, "Unvallied group control specified for group " + well.groupName(), deferred_logger );
}
}

160
opm/simulators/wells/TargetCalculator.hpp Normal file
View File

@ -0,0 +1,160 @@
/*
Copyright 2020 Equinor ASA.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_TARGETCALCULATOR_HEADER_INCLUDED
#define OPM_TARGETCALCULATOR_HEADER_INCLUDED
#include <opm/parser/eclipse/EclipseState/Schedule/Group/Group.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/simulators/utils/DeferredLogger.hpp>
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp>
#include <algorithm>
#include <cassert>
#include <type_traits>
#include <vector>
namespace Opm
{
namespace WellGroupHelpers
{
/// Based on a group control mode, extract or calculate rates, and
/// provide other conveniences.
class TargetCalculator
{
public:
TargetCalculator(const Group::ProductionCMode cmode,
const PhaseUsage& pu,
const std::vector<double>& resv_coeff)
: cmode_(cmode)
, pu_(pu)
, resv_coeff_(resv_coeff)
{
}
template <typename RateVec>
auto calcModeRateFromRates(const RateVec& rates) const
{
// ElemType is just the plain element type of the rates container,
// without any reference, const or volatile modifiers.
using ElemType = std::remove_cv_t<std::remove_reference_t<decltype(rates[0])>>;
switch (cmode_) {
case Group::ProductionCMode::ORAT: {
assert(pu_.phase_used[BlackoilPhases::Liquid]);
const int pos = pu_.phase_pos[BlackoilPhases::Liquid];
return rates[pos];
}
case Group::ProductionCMode::WRAT: {
assert(pu_.phase_used[BlackoilPhases::Aqua]);
const int pos = pu_.phase_pos[BlackoilPhases::Aqua];
return rates[pos];
}
case Group::ProductionCMode::GRAT: {
assert(pu_.phase_used[BlackoilPhases::Vapour]);
const int pos = pu_.phase_pos[BlackoilPhases::Vapour];
return rates[pos];
}
case Group::ProductionCMode::LRAT: {
assert(pu_.phase_used[BlackoilPhases::Liquid]);
assert(pu_.phase_used[BlackoilPhases::Aqua]);
const int opos = pu_.phase_pos[BlackoilPhases::Liquid];
const int wpos = pu_.phase_pos[BlackoilPhases::Aqua];
return rates[opos] + rates[wpos];
}
case Group::ProductionCMode::RESV: {
assert(pu_.phase_used[BlackoilPhases::Liquid]);
assert(pu_.phase_used[BlackoilPhases::Aqua]);
assert(pu_.phase_used[BlackoilPhases::Vapour]);
ElemType mode_rate = zero<ElemType>();
for (int phase = 0; phase < pu_.num_phases; ++phase) {
mode_rate += rates[phase] * resv_coeff_[phase];
}
return mode_rate;
}
default:
// Should never be here.
assert(false);
return zero<ElemType>();
}
}
double groupTarget(const Group::ProductionControls ctrl) const
{
switch (cmode_) {
case Group::ProductionCMode::ORAT:
return ctrl.oil_target;
case Group::ProductionCMode::WRAT:
return ctrl.water_target;
case Group::ProductionCMode::GRAT:
return ctrl.gas_target;
case Group::ProductionCMode::LRAT:
return ctrl.liquid_target;
case Group::ProductionCMode::RESV:
return ctrl.resv_target;
default:
// Should never be here.
assert(false);
return 0.0;
}
}
GuideRateModel::Target guideTargetMode() const
{
switch (cmode_) {
case Group::ProductionCMode::ORAT:
return GuideRateModel::Target::OIL;
case Group::ProductionCMode::WRAT:
return GuideRateModel::Target::WAT;
case Group::ProductionCMode::GRAT:
return GuideRateModel::Target::GAS;
case Group::ProductionCMode::LRAT:
return GuideRateModel::Target::LIQ;
case Group::ProductionCMode::RESV:
return GuideRateModel::Target::RES;
default:
// Should never be here.
assert(false);
return GuideRateModel::Target::NONE;
}
}
private:
template <typename ElemType>
static ElemType zero()
{
// This is for Evaluation types.
ElemType x;
x = 0.0;
return x;
}
Group::ProductionCMode cmode_;
const PhaseUsage& pu_;
const std::vector<double>& resv_coeff_;
};
} // namespace WellGroupHelpers
} // namespace Opm
#endif

1164
opm/simulators/wells/WellGroupHelpers.cpp Normal file
View File

File diff suppressed because it is too large Load Diff

700
opm/simulators/wells/WellGroupHelpers.hpp
View File

@ -21,339 +21,104 @@
#ifndef OPM_WELLGROUPHELPERS_HEADER_INCLUDED
#define OPM_WELLGROUPHELPERS_HEADER_INCLUDED
#include <opm/output/data/Groups.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Group/Group.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/simulators/utils/DeferredLogger.hpp>
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp>
#include <algorithm>
#include <cassert>
#include <type_traits>
#include <vector>
#include <opm/parser/eclipse/EclipseState/Schedule/ScheduleTypes.hpp>
namespace Opm {
namespace Opm
{
namespace wellGroupHelpers
{
inline void setGroupControl(const Group& group, const Schedule& schedule, const Phase& groupInjectionPhase, const int reportStepIdx, const bool injector, WellStateFullyImplicitBlackoil& wellState, std::ostringstream& ss) {
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
setGroupControl(groupTmp, schedule, groupInjectionPhase, reportStepIdx, injector, wellState, ss);
if (injector) {
wellState.setCurrentInjectionGroupControl(groupInjectionPhase, groupName, Group::InjectionCMode::FLD);
} else {
wellState.setCurrentProductionGroupControl(groupName, Group::ProductionCMode::FLD);
}
}
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
continue;
int well_index = it->second[0];
const auto& wellEcl = schedule.getWell(wellName, reportStepIdx);
if (wellEcl.getStatus() == Well::Status::SHUT)
continue;
if (!wellEcl.isAvailableForGroupControl())
continue;
if (wellEcl.isProducer() && !injector) {
if (wellState.currentProductionControls()[well_index] != Well::ProducerCMode::GRUP) {
wellState.currentProductionControls()[well_index] = Well::ProducerCMode::GRUP;
ss <<"\n Producer " << wellName << " switches to GRUP control limit";
}
}
if (wellEcl.isInjector() && injector) {
// only switch if the well phase is the same as the group phase
// Get the current controls.
const InjectorType& injectorType = wellEcl.getInjectionProperties().injectorType;
if (injectorType == InjectorType::WATER && groupInjectionPhase != Phase::WATER)
continue;
if (injectorType == InjectorType::OIL && groupInjectionPhase != Phase::OIL)
continue;
if (injectorType == InjectorType::GAS && groupInjectionPhase != Phase::GAS)
continue;
if (injectorType == InjectorType::MULTI)
throw("Expected WATER, OIL or GAS as type for injectors " + wellEcl.name());
if (wellState.currentInjectionControls()[well_index] != Well::InjectorCMode::GRUP) {
wellState.currentInjectionControls()[well_index] = Well::InjectorCMode::GRUP;
ss <<"\n Injector " << wellName << " switches to GRUP control limit";
}
}
}
}
inline void setCmodeGroup(const Group& group, const Schedule& schedule, const SummaryState& summaryState, const int reportStepIdx, WellStateFullyImplicitBlackoil& wellState) {
for (const std::string& groupName : group.groups()) {
setCmodeGroup( schedule.getGroup(groupName, reportStepIdx), schedule, summaryState, reportStepIdx, wellState);
}
// use NONE as default control
const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS};
for (Phase phase : all) {
if (!wellState.hasInjectionGroupControl(phase, group.name())) {
wellState.setCurrentInjectionGroupControl(phase, group.name(), Group::InjectionCMode::NONE);
}
}
if (!wellState.hasProductionGroupControl(group.name())) {
wellState.setCurrentProductionGroupControl(group.name(), Group::ProductionCMode::NONE);
}
if (group.isInjectionGroup() && schedule.hasWellGroupEvent(group.name(), ScheduleEvents::GROUP_INJECTION_UPDATE, reportStepIdx)) {
for (Phase phase : all) {
if (!group.hasInjectionControl(phase))
continue;
const auto& controls = group.injectionControls(phase, summaryState);
wellState.setCurrentInjectionGroupControl(phase, group.name(), controls.cmode);
}
}
if (group.isProductionGroup() && schedule.hasWellGroupEvent(group.name(), ScheduleEvents::GROUP_PRODUCTION_UPDATE, reportStepIdx)) {
const auto controls = group.productionControls(summaryState);
wellState.setCurrentProductionGroupControl(group.name(), controls.cmode);
}
if (schedule.gConSale(reportStepIdx).has(group.name())) {
wellState.setCurrentInjectionGroupControl(Phase::GAS, group.name(), Group::InjectionCMode::SALE);
std::ostringstream ss;
setGroupControl(group, schedule, Phase::GAS, reportStepIdx, /*injector*/true, wellState, ss);
}
}
inline void accumulateGroupEfficiencyFactor(const Group& group, const Schedule& schedule, const int reportStepIdx, double& factor) {
factor *= group.getGroupEfficiencyFactor();
if (group.parent() != "FIELD")
accumulateGroupEfficiencyFactor(schedule.getGroup(group.parent(), reportStepIdx), schedule, reportStepIdx, factor);
}
namespace WellGroupHelpers
{
inline void computeGroupTargetReduction(const Group& group, const WellStateFullyImplicitBlackoil& wellState, const Schedule& schedule, const int reportStepIdx, const int phasePos, const bool isInjector, double& groupTargetReduction )
{
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
computeGroupTargetReduction(groupTmp, wellState, schedule, reportStepIdx, phasePos, isInjector, groupTargetReduction);
}
for (const std::string& wellName : group.wells()) {
const auto& wellTmp = schedule.getWell(wellName, reportStepIdx);
void setCmodeGroup(const Group& group,
const Schedule& schedule,
const SummaryState& summaryState,
const int reportStepIdx,
WellStateFullyImplicitBlackoil& wellState);
if (wellTmp.isProducer() && isInjector)
continue;
void accumulateGroupEfficiencyFactor(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
double& factor);
if (wellTmp.isInjector() && !isInjector)
continue;
double sumWellPhaseRates(const std::vector<double>& rates,
const Group& group,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const int phasePos,
const bool injector);
if (wellTmp.getStatus() == Well::Status::SHUT)
continue;
double sumWellRates(const Group& group,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const int phasePos,
const bool injector);
const auto& end = wellState.wellMap().end();
const auto& it = wellState.wellMap().find( wellName );
if (it == end) // the well is not found
continue;
double sumWellResRates(const Group& group,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const int phasePos,
const bool injector);
int well_index = it->second[0];
const auto wellrate_index = well_index * wellState.numPhases();
// add contributino from wells not under group control
if (isInjector) {
if (wellState.currentInjectionControls()[well_index] != Well::InjectorCMode::GRUP)
groupTargetReduction += wellState.wellRates()[wellrate_index + phasePos];
} else {
if (wellState.currentProductionControls()[well_index] != Well::ProducerCMode::GRUP)
groupTargetReduction -= wellState.wellRates()[wellrate_index + phasePos];
}
}
}
double sumSolventRates(const Group& group,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const bool injector);
inline double sumWellPhaseRates(const std::vector<double>& rates, const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const int phasePos,
const bool injector) {
double rate = 0.0;
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
rate += groupTmp.getGroupEfficiencyFactor()*sumWellPhaseRates(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
continue;
int well_index = it->second[0];
const auto& wellEcl = schedule.getWell(wellName, reportStepIdx);
//only count producers or injectors
if ( (wellEcl.isProducer() && injector) || (wellEcl.isInjector() && !injector))
continue;
if (wellEcl.getStatus() == Well::Status::SHUT)
continue;
double factor = wellEcl.getEfficiencyFactor();
const auto wellrate_index = well_index * wellState.numPhases();
if (injector)
rate += factor * rates[ wellrate_index + phasePos];
else
rate -= factor * rates[ wellrate_index + phasePos];
}
return rate;
}
inline double sumWellRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const int phasePos, const bool injector) {
return sumWellPhaseRates(wellState.wellRates(), group, schedule, wellState, reportStepIdx, phasePos, injector);
}
inline double sumWellResRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const int phasePos, const bool injector) {
return sumWellPhaseRates(wellState.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phasePos, injector);
}
inline double sumSolventRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const bool injector) {
double rate = 0.0;
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
rate += groupTmp.getGroupEfficiencyFactor()*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
continue;
int well_index = it->second[0];
const auto& wellEcl = schedule.getWell(wellName, reportStepIdx);
//only count producers or injectors
if ( (wellEcl.isProducer() && injector) || (wellEcl.isInjector() && !injector))
continue;
if (wellEcl.getStatus() == Well::Status::SHUT)
continue;
double factor = wellEcl.getEfficiencyFactor();
if (injector)
rate += factor * wellState.solventWellRate(well_index);
else
rate -= factor * wellState.solventWellRate(well_index);
}
return rate;
}
inline void updateGroupTargetReduction(const Group& group, const Schedule& schedule, const int reportStepIdx, const bool isInjector, const PhaseUsage& pu, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState, std::vector<double>& groupTargetReduction)
{
const int np = wellState.numPhases();
for (const std::string& groupName : group.groups()) {
std::vector<double> thisGroupTargetReduction(np, 0.0);
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
updateGroupTargetReduction(groupTmp, schedule, reportStepIdx, isInjector, pu, wellStateNupcol, wellState, thisGroupTargetReduction);
// accumulate group contribution from sub group
if (isInjector) {
const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS};
for (Phase phase : all) {
const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(phase, groupName);
int phasePos;
if (phase == Phase::GAS && pu.phase_used[BlackoilPhases::Vapour] )
phasePos = pu.phase_pos[BlackoilPhases::Vapour];
else if (phase == Phase::OIL && pu.phase_used[BlackoilPhases::Liquid])
phasePos = pu.phase_pos[BlackoilPhases::Liquid];
else if (phase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua] )
phasePos = pu.phase_pos[BlackoilPhases::Aqua];
else
continue;
if (currentGroupControl != Group::InjectionCMode::FLD) {
groupTargetReduction[phasePos] += sumWellRates(groupTmp, schedule, wellStateNupcol, reportStepIdx, phasePos, isInjector);
} else {
groupTargetReduction[phasePos] += thisGroupTargetReduction[phasePos];
}
}
} else {
const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName);
if (currentGroupControl != Group::ProductionCMode::FLD) {
for (int phase = 0; phase < np; phase++) {
groupTargetReduction[phase] += sumWellRates(groupTmp, schedule, wellStateNupcol, reportStepIdx, phase, isInjector);
}
} else {
// or accumulate directly from the wells if controled from its parents
for (int phase = 0; phase < np; phase++) {
groupTargetReduction[phase] += thisGroupTargetReduction[phase];
}
}
}
}
for (const std::string& wellName : group.wells()) {
const auto& wellTmp = schedule.getWell(wellName, reportStepIdx);
if (wellTmp.isProducer() && isInjector)
continue;
if (wellTmp.isInjector() && !isInjector)
continue;
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
continue;
int well_index = it->second[0];
const auto wellrate_index = well_index * wellState.numPhases();
// add contributino from wells not under group control
if (isInjector) {
if (wellState.currentInjectionControls()[well_index] != Well::InjectorCMode::GRUP)
for (int phase = 0; phase < np; phase++) {
groupTargetReduction[phase] += wellStateNupcol.wellRates()[wellrate_index + phase];
}
} else {
if (wellState.currentProductionControls()[well_index] != Well::ProducerCMode::GRUP)
for (int phase = 0; phase < np; phase++) {
groupTargetReduction[phase] -= wellStateNupcol.wellRates()[wellrate_index + phase];
}
}
}
if (isInjector)
wellState.setCurrentInjectionGroupReductionRates(group.name(), groupTargetReduction);
else
wellState.setCurrentProductionGroupReductionRates(group.name(), groupTargetReduction);
}
void updateGroupTargetReduction(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
const bool isInjector,
const PhaseUsage& pu,
const WellStateFullyImplicitBlackoil& wellStateNupcol,
WellStateFullyImplicitBlackoil& wellState,
std::vector<double>& groupTargetReduction);
template <class Comm>
inline void updateGuideRateForGroups(const Group& group, const Schedule& schedule, const PhaseUsage& pu, const int reportStepIdx, const double& simTime, const bool isInjector, WellStateFullyImplicitBlackoil& wellState, const Comm& comm, GuideRate* guideRate, std::vector<double>& pot)
void updateGuideRateForGroups(const Group& group,
const Schedule& schedule,
const PhaseUsage& pu,
const int reportStepIdx,
const double& simTime,
const bool isInjector,
WellStateFullyImplicitBlackoil& wellState,
const Comm& comm,
GuideRate* guideRate,
std::vector<double>& pot)
{
const int np = pu.num_phases;
for (const std::string& groupName : group.groups()) {
std::vector<double> thisPot(np, 0.0);
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
updateGuideRateForGroups(groupTmp, schedule, pu, reportStepIdx, simTime, isInjector, wellState, comm, guideRate, thisPot);
updateGuideRateForGroups(
groupTmp, schedule, pu, reportStepIdx, simTime, isInjector, wellState, comm, guideRate, thisPot);
// accumulate group contribution from sub group if FLD
// accumulate group contribution from sub group unconditionally
if (isInjector) {
const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS};
for (Phase phase : all) {
const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(phase, groupName);
if (currentGroupControl != Group::InjectionCMode::FLD) {
continue;
}
int phasePos;
if (phase == Phase::GAS && pu.phase_used[BlackoilPhases::Vapour] )
if (phase == Phase::GAS && pu.phase_used[BlackoilPhases::Vapour])
phasePos = pu.phase_pos[BlackoilPhases::Vapour];
else if (phase == Phase::OIL && pu.phase_used[BlackoilPhases::Liquid])
phasePos = pu.phase_pos[BlackoilPhases::Liquid];
else if (phase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua] )
else if (phase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua])
phasePos = pu.phase_pos[BlackoilPhases::Aqua];
else
continue;
@ -362,14 +127,14 @@ namespace Opm {
}
} else {
const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName);
if (currentGroupControl != Group::ProductionCMode::FLD) {
if (currentGroupControl != Group::ProductionCMode::FLD
&& currentGroupControl != Group::ProductionCMode::NONE) {
continue;
}
for (int phase = 0; phase < np; phase++) {
pot[phase] += thisPot[phase];
}
}
}
for (const std::string& wellName : group.wells()) {
const auto& wellTmp = schedule.getWell(wellName, reportStepIdx);
@ -383,41 +148,35 @@ namespace Opm {
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& it = wellState.wellMap().find(wellName);
if (it == end) // the well is not found
continue;
int well_index = it->second[0];
const auto wellrate_index = well_index * wellState.numPhases();
// add contribution from wells not under group control
if (isInjector) {
if (wellState.currentInjectionControls()[well_index] == Well::InjectorCMode::GRUP)
for (int phase = 0; phase < np; phase++) {
pot[phase] += wellState.wellPotentials()[wellrate_index + phase];
}
} else {
if (wellState.currentProductionControls()[well_index] == Well::ProducerCMode::GRUP)
for (int phase = 0; phase < np; phase++) {
pot[phase] -= wellState.wellPotentials()[wellrate_index + phase];
}
// add contribution from wells unconditionally
for (int phase = 0; phase < np; phase++) {
pot[phase] += wellState.wellPotentials()[wellrate_index + phase];
}
}
double oilPot = 0.0;
if (pu.phase_used[BlackoilPhases::Liquid])
oilPot = pot [ pu.phase_pos[BlackoilPhases::Liquid]];
oilPot = pot[pu.phase_pos[BlackoilPhases::Liquid]];
double gasPot = 0.0;
if (pu.phase_used[BlackoilPhases::Vapour])
gasPot = pot [ pu.phase_pos[BlackoilPhases::Vapour]];
gasPot = pot[pu.phase_pos[BlackoilPhases::Vapour]];
double waterPot = 0.0;
if (pu.phase_used[BlackoilPhases::Aqua])
waterPot = pot [pu.phase_pos[BlackoilPhases::Aqua]];
waterPot = pot[pu.phase_pos[BlackoilPhases::Aqua]];
oilPot = comm.sum(oilPot);
gasPot = comm.sum(gasPot);
waterPot = comm.sum(waterPot);
const double gefac = group.getGroupEfficiencyFactor();
oilPot = comm.sum(oilPot) * gefac;
gasPot = comm.sum(gasPot) * gefac;
waterPot = comm.sum(waterPot) * gefac;
if (isInjector) {
wellState.setCurrentGroupInjectionPotentials(group.name(), pot);
@ -427,7 +186,14 @@ namespace Opm {
}
template <class Comm>
inline void updateGuideRatesForWells(const Schedule& schedule, const PhaseUsage& pu, const int reportStepIdx, const double& simTime, const WellStateFullyImplicitBlackoil& wellState, const Comm& comm, GuideRate* guideRate) {
void updateGuideRatesForWells(const Schedule& schedule,
const PhaseUsage& pu,
const int reportStepIdx,
const double& simTime,
const WellStateFullyImplicitBlackoil& wellState,
const Comm& comm,
GuideRate* guideRate)
{
const auto& end = wellState.wellMap().end();
for (const auto& well : schedule.getWells(reportStepIdx)) {
@ -435,12 +201,12 @@ namespace Opm {
double gaspot = 0.0;
double waterpot = 0.0;
const auto& it = wellState.wellMap().find( well.name());
if (it != end) { // the well is found
const auto& it = wellState.wellMap().find(well.name());
if (it != end) { // the well is found
int well_index = it->second[0];
const auto wpot = wellState.wellPotentials().data() + well_index*wellState.numPhases();
const auto wpot = wellState.wellPotentials().data() + well_index * wellState.numPhases();
if (pu.phase_used[BlackoilPhases::Liquid] > 0)
oilpot = wpot[pu.phase_pos[BlackoilPhases::Liquid]];
@ -450,170 +216,172 @@ namespace Opm {
if (pu.phase_used[BlackoilPhases::Aqua] > 0)
waterpot = wpot[pu.phase_pos[BlackoilPhases::Aqua]];
}
oilpot = comm.sum(oilpot);
gaspot = comm.sum(gaspot);
waterpot = comm.sum(waterpot);
const double wefac = well.getEfficiencyFactor();
oilpot = comm.sum(oilpot) * wefac;
gaspot = comm.sum(gaspot) * wefac;
waterpot = comm.sum(waterpot) * wefac;
guideRate->compute(well.name(), reportStepIdx, simTime, oilpot, gaspot, waterpot);
}
}
inline void updateVREPForGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) {
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
updateVREPForGroups(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState);
}
const int np = wellState.numPhases();
double resv = 0.0;
for (int phase = 0; phase < np; ++phase) {
resv += sumWellPhaseRates(wellStateNupcol.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false);
}
wellState.setCurrentInjectionVREPRates(group.name(), resv);
}
void updateVREPForGroups(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
const WellStateFullyImplicitBlackoil& wellStateNupcol,
WellStateFullyImplicitBlackoil& wellState);
inline void updateReservoirRatesInjectionGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) {
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
updateReservoirRatesInjectionGroups(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState);
}
const int np = wellState.numPhases();
std::vector<double> resv(np, 0.0);
for (int phase = 0; phase < np; ++phase) {
resv[phase] = sumWellPhaseRates(wellStateNupcol.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ true);
}
wellState.setCurrentInjectionGroupReservoirRates(group.name(), resv);
}
void updateReservoirRatesInjectionGroups(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
const WellStateFullyImplicitBlackoil& wellStateNupcol,
WellStateFullyImplicitBlackoil& wellState);
inline void updateREINForGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, const PhaseUsage& pu, const SummaryState& st, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) {
const int np = wellState.numPhases();
for (const std::string& groupName : group.groups()) {
const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
updateREINForGroups(groupTmp, schedule, reportStepIdx, pu, st, wellStateNupcol, wellState);
}
void updateWellRates(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
const WellStateFullyImplicitBlackoil& wellStateNupcol,
WellStateFullyImplicitBlackoil& wellState);
std::vector<double> rein(np, 0.0);
for (int phase = 0; phase < np; ++phase) {
rein[phase] = sumWellPhaseRates(wellStateNupcol.wellRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false);
}
void updateGroupProductionRates(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
const WellStateFullyImplicitBlackoil& wellStateNupcol,
WellStateFullyImplicitBlackoil& wellState);
// add import rate and substract consumption rate for group for gas
if (schedule.gConSump(reportStepIdx).has(group.name())) {
const auto& gconsump = schedule.gConSump(reportStepIdx).get(group.name(), st);
if (pu.phase_used[BlackoilPhases::Vapour]) {
rein[pu.phase_pos[BlackoilPhases::Vapour]] += gconsump.import_rate;
rein[pu.phase_pos[BlackoilPhases::Vapour]] -= gconsump.consumption_rate;
}
}
void updateREINForGroups(const Group& group,
const Schedule& schedule,
const int reportStepIdx,
const PhaseUsage& pu,
const SummaryState& st,
const WellStateFullyImplicitBlackoil& wellStateNupcol,
WellStateFullyImplicitBlackoil& wellState);
wellState.setCurrentInjectionREINRates(group.name(), rein);
}
GuideRate::RateVector
getRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name);
inline double wellFractionFromGuideRates(const Well& well, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const GuideRate* guideRate, const Well::GuideRateTarget& wellTarget, const bool isInjector) {
double groupTotalGuideRate = 0.0;
const Group& groupTmp = schedule.getGroup(well.groupName(), reportStepIdx);
int global_well_index = -1;
for (const std::string& wellName : groupTmp.wells()) {
const auto& wellTmp = schedule.getWell(wellName, reportStepIdx);
global_well_index++;
double getGuideRate(const std::string& name,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const GuideRate* guideRate,
const GuideRateModel::Target target,
const PhaseUsage& pu);
if (wellTmp.isProducer() && isInjector)
continue;
if (wellTmp.isInjector() && !isInjector)
continue;
double getGuideRateInj(const std::string& name,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const GuideRate* guideRate,
const GuideRateModel::Target target,
const Phase& injectionPhase,
const PhaseUsage& pu);
if (wellTmp.getStatus() == Well::Status::SHUT)
continue;
// only count wells under group control
if (isInjector) {
if (!wellState.isInjectionGrup(wellName))
continue;
} else {
if (!wellState.isProductionGrup(wellName))
continue;
}
int groupControlledWells(const Schedule& schedule,
const WellStateFullyImplicitBlackoil& well_state,
const int report_step,
const std::string& group_name,
const std::string& always_included_child);
groupTotalGuideRate += guideRate->get(wellName, wellTarget);
}
if (groupTotalGuideRate == 0.0)
return 0.0;
class FractionCalculator
{
public:
FractionCalculator(const Schedule& schedule,
const WellStateFullyImplicitBlackoil& well_state,
const int report_step,
const GuideRate* guide_rate,
const GuideRateModel::Target target,
const PhaseUsage& pu);
double fraction(const std::string& name, const std::string& control_group_name, const bool always_include_this);
double localFraction(const std::string& name, const std::string& always_included_child);
double wellGuideRate = guideRate->get(well.name(), wellTarget);
return wellGuideRate / groupTotalGuideRate;
}
private:
std::string parent(const std::string& name);
double guideRateSum(const Group& group, const std::string& always_included_child, const bool include_all);
double guideRate(const std::string& name, const std::string& always_included_child);
int groupControlledWells(const std::string& group_name, const std::string& always_included_child);
GuideRate::RateVector getGroupRateVector(const std::string& group_name);
const Schedule& schedule_;
const WellStateFullyImplicitBlackoil& well_state_;
int report_step_;
const GuideRate* guide_rate_;
GuideRateModel::Target target_;
PhaseUsage pu_;
};
inline double groupFractionFromGuideRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const GuideRate* guideRate, const Group::GuideRateTarget& groupTarget) {
double groupTotalGuideRate = 0.0;
const Group& groupParent = schedule.getGroup(group.parent(), reportStepIdx);
for (const std::string& groupName : groupParent.groups()) {
// only count group under group control from its parent
const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName);
if (currentGroupControl != Group::ProductionCMode::FLD)
continue;
GuideRate::RateVector getGroupRateVector(const std::string& group_name);
groupTotalGuideRate += guideRate->get(groupName, groupTarget);
}
if (groupTotalGuideRate == 0.0)
return 1.0;
double groupGuideRate = guideRate->get(group.name(), groupTarget);
return groupGuideRate / groupTotalGuideRate;
}
double fractionFromGuideRates(const std::string& name,
const std::string& controlGroupName,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const GuideRate* guideRate,
const GuideRateModel::Target target,
const PhaseUsage& pu,
const bool alwaysIncludeThis = false);
inline void accumulateGroupFractionsFromGuideRates(const std::string& groupName, const std::string& controlGroupName, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState,const int reportStepIdx, const GuideRate* guideRate, const Group::GuideRateTarget& groupTarget, double& fraction) {
const Group& group = schedule.getGroup(groupName, reportStepIdx);
if (groupName != controlGroupName) {
fraction *= groupFractionFromGuideRates(group, schedule, wellState, reportStepIdx, guideRate, groupTarget);
accumulateGroupFractionsFromGuideRates(group.parent(), controlGroupName, schedule, wellState, reportStepIdx, guideRate, groupTarget, fraction);
}
return;
}
double fractionFromInjectionPotentials(const std::string& name,
const std::string& controlGroupName,
const Schedule& schedule,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const GuideRate* guideRate,
const GuideRateModel::Target target,
const PhaseUsage& pu,
const Phase& injectionPhase,
const bool alwaysIncludeThis = false);
inline double groupFractionFromInjectionPotentials(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const PhaseUsage& pu, const int reportStepIdx, const Phase& injectionPhase) {
double groupTotalGuideRate = 0.0;
const Group& groupParent = schedule.getGroup(group.parent(), reportStepIdx);
int phasePos;
if (injectionPhase == Phase::GAS && pu.phase_used[BlackoilPhases::Vapour] )
phasePos = pu.phase_pos[ pu.phase_pos[BlackoilPhases::Vapour] ];
else if (injectionPhase == Phase::OIL && pu.phase_used[BlackoilPhases::Liquid])
phasePos = pu.phase_pos[ pu.phase_pos[BlackoilPhases::Liquid] ];
else if (injectionPhase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua] )
phasePos = pu.phase_pos[ pu.phase_pos[BlackoilPhases::Aqua] ];
else
throw("this should not happen");
for (const std::string& groupName : groupParent.groups()) {
// only count group under group control from its parent
const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(injectionPhase, groupName);
if (currentGroupControl != Group::InjectionCMode::FLD)
continue;
groupTotalGuideRate += wellState.currentGroupInjectionPotentials(groupName)[phasePos];
}
if (groupTotalGuideRate == 0.0)
return 1.0;
double groupGuideRate = wellState.currentGroupInjectionPotentials(group.name())[phasePos];
return groupGuideRate / groupTotalGuideRate;
}
inline void accumulateGroupInjectionPotentialFractions(const std::string& groupName, const std::string& controlGroupName, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const PhaseUsage& pu, const int reportStepIdx, const Phase& injectionPhase, double& fraction) {
const Group& group = schedule.getGroup(groupName, reportStepIdx);
if (groupName != controlGroupName) {
fraction *= groupFractionFromInjectionPotentials(group, schedule, wellState, pu, reportStepIdx, injectionPhase);
accumulateGroupInjectionPotentialFractions(group.parent(), controlGroupName, schedule, wellState, pu, reportStepIdx, injectionPhase, fraction);
}
return;
}
std::pair<bool, double> checkGroupConstraintsInj(const std::string& name,
const std::string& parent,
const Group& group,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const GuideRate* guideRate,
const double* rates,
Phase injectionPhase,
const PhaseUsage& pu,
const double efficiencyFactor,
const Schedule& schedule,
const SummaryState& summaryState,
const std::vector<double>& resv_coeff,
DeferredLogger& deferred_logger);
} // namespace wellGroupHelpers
}
std::vector<std::string> groupChainTopBot(const std::string& bottom,
const std::string& top,
const Schedule& schedule,
const int report_step);
std::pair<bool, double> checkGroupConstraintsProd(const std::string& name,
const std::string& parent,
const Group& group,
const WellStateFullyImplicitBlackoil& wellState,
const int reportStepIdx,
const GuideRate* guideRate,
const double* rates,
const PhaseUsage& pu,
const double efficiencyFactor,
const Schedule& schedule,
const SummaryState& summaryState,
const std::vector<double>& resv_coeff,
DeferredLogger& deferred_logger);
} // namespace WellGroupHelpers
} // namespace Opm
#endif

71
opm/simulators/wells/WellInterface.hpp
View File

@ -201,7 +201,9 @@ namespace Opm
WellState& well_state,
Opm::DeferredLogger& deferred_logger) const = 0;
void updateWellControl(const Simulator& ebos_simulator,
enum class IndividualOrGroup { Individual, Group, Both };
bool updateWellControl(const Simulator& ebos_simulator,
const IndividualOrGroup iog,
WellState& well_state,
Opm::DeferredLogger& deferred_logger) /* const */;
@ -478,10 +480,75 @@ namespace Opm
// index calculations
int well_productivity_index_logger_counter_;
bool checkConstraints(WellState& well_state, const SummaryState& summaryState);
bool checkConstraints(WellState& well_state,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const;
bool checkIndividualConstraints(WellState& well_state,
const SummaryState& summaryState) const;
bool checkGroupConstraints(WellState& well_state,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const;
std::pair<bool, double> checkGroupConstraintsProd(const Group& group,
const WellState& well_state,
const double efficiencyFactor,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const;
std::pair<bool, double> checkGroupConstraintsInj(const Group& group,
const WellState& well_state,
const double efficiencyFactor,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const;
template <class EvalWell>
void getGroupInjectionControl(const Group& group,
const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const InjectorType& injectorType,
const EvalWell& bhp,
const EvalWell& injection_rate,
EvalWell& control_eq,
double efficiencyFactor);
template <class EvalWell>
void getGroupProductionControl(const Group& group,
const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const EvalWell& bhp,
const std::vector<EvalWell>& rates,
EvalWell& control_eq,
double efficiencyFactor);
template <class EvalWell, class BhpFromThpFunc>
void assembleControlEqInj(const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const Well::InjectionControls& controls,
const EvalWell& bhp,
const EvalWell& injection_rate,
BhpFromThpFunc bhp_from_thp,
EvalWell& control_eq,
Opm::DeferredLogger& deferred_logger);
template <class EvalWell, class BhpFromThpFunc>
void assembleControlEqProd(const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const Well::ProductionControls& controls,
const EvalWell& bhp,
const std::vector<EvalWell>& rates, // Always 3 canonical rates.
BhpFromThpFunc bhp_from_thp,
EvalWell& control_eq,
Opm::DeferredLogger& deferred_logger);
};

673
opm/simulators/wells/WellInterface_impl.hpp
View File

@ -21,6 +21,7 @@
#include <opm/parser/eclipse/EclipseState/Schedule/ScheduleTypes.hpp>
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/TargetCalculator.hpp>
namespace Opm
{
@ -461,17 +462,19 @@ namespace Opm
template<typename TypeTag>
void
bool
WellInterface<TypeTag>::
updateWellControl(const Simulator& ebos_simulator,
const IndividualOrGroup iog,
WellState& well_state,
Opm::DeferredLogger& deferred_logger) /* const */
{
if (this->wellIsStopped()) {
return;
return false;
}
const auto& summaryState = ebos_simulator.vanguard().summaryState();
const auto& schedule = ebos_simulator.vanguard().schedule();
const auto& well = well_ecl_;
std::string from;
if (well.isInjector()) {
@ -480,7 +483,15 @@ namespace Opm
from = Well::ProducerCMode2String(well_state.currentProductionControls()[index_of_well_]);
}
bool changed = checkConstraints(well_state, summaryState);
bool changed = false;
if (iog == IndividualOrGroup::Individual) {
changed = checkIndividualConstraints(well_state, summaryState);
} else if (iog == IndividualOrGroup::Group) {
changed = checkGroupConstraints(well_state, schedule, summaryState, deferred_logger);
} else {
assert(iog == IndividualOrGroup::Both);
changed = checkConstraints(well_state, schedule, summaryState, deferred_logger);
}
auto cc = Dune::MPIHelper::getCollectiveCommunication();
@ -503,6 +514,8 @@ namespace Opm
updateWellStateWithTarget(ebos_simulator, well_state, deferred_logger);
updatePrimaryVariables(well_state, deferred_logger);
}
return changed;
}
@ -1436,20 +1449,38 @@ namespace Opm
return operability_status_.isOperable();
}
template<typename TypeTag>
bool
WellInterface<TypeTag>::
checkConstraints(WellState& well_state, const SummaryState& summaryState) {
template <typename TypeTag>
bool
WellInterface<TypeTag>::checkConstraints(WellState& well_state,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const
{
const bool ind_broken = checkIndividualConstraints(well_state, summaryState);
if (ind_broken) {
return true;
} else {
return checkGroupConstraints(well_state, schedule, summaryState, deferred_logger);
}
}
template <typename TypeTag>
bool
WellInterface<TypeTag>::checkIndividualConstraints(WellState& well_state,
const SummaryState& summaryState) const
{
const auto& well = well_ecl_;
const PhaseUsage& pu = phaseUsage();
const int well_index = index_of_well_;
const auto wellrate_index = well_index * pu.num_phases;
// bool changed = false;
// // Stopped wells can not change control
// if (currentControl == "STOP")
// return newControl;
if (well.isInjector()) {
const auto controls = well.injectionControls(summaryState);
@ -1633,4 +1664,620 @@ namespace Opm
}
}
template <typename TypeTag>
bool
WellInterface<TypeTag>::checkGroupConstraints(WellState& well_state,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const
{
const auto& well = well_ecl_;
const int well_index = index_of_well_;
if (well.isInjector()) {
Opm::Well::InjectorCMode& currentControl = well_state.currentInjectionControls()[well_index];
if (currentControl != Well::InjectorCMode::GRUP) {
// This checks only the first encountered group limit,
// in theory there could be several, and then we should
// test all but the one currently applied. At that point,
// this if-statement should be removed and we should always
// check, skipping over only the single group parent whose
// control is the active one for the well (if any).
const auto& group = schedule.getGroup( well.groupName(), current_step_ );
const double efficiencyFactor = well.getEfficiencyFactor();
const std::pair<bool, double> group_constraint = checkGroupConstraintsInj(
group, well_state, efficiencyFactor, schedule, summaryState, deferred_logger);
// If a group constraint was broken, we set the current well control to
// be GRUP.
if (group_constraint.first) {
well_state.currentInjectionControls()[index_of_well_] = Well::InjectorCMode::GRUP;
const int np = well_state.numPhases();
for (int p = 0; p<np; ++p) {
well_state.wellRates()[index_of_well_*np + p] *= group_constraint.second;
}
}
return group_constraint.first;
}
}
if (well.isProducer( )) {
Well::ProducerCMode& currentControl = well_state.currentProductionControls()[well_index];
if (currentControl != Well::ProducerCMode::GRUP) {
// This checks only the first encountered group limit,
// in theory there could be several, and then we should
// test all but the one currently applied. At that point,
// this if-statement should be removed and we should always
// check, skipping over only the single group parent whose
// control is the active one for the well (if any).
const auto& group = schedule.getGroup( well.groupName(), current_step_ );
const double efficiencyFactor = well.getEfficiencyFactor();
const std::pair<bool, double> group_constraint = checkGroupConstraintsProd(
group, well_state, efficiencyFactor, schedule, summaryState, deferred_logger);
// If a group constraint was broken, we set the current well control to
// be GRUP.
if (group_constraint.first) {
well_state.currentProductionControls()[index_of_well_] = Well::ProducerCMode::GRUP;
const int np = well_state.numPhases();
for (int p = 0; p<np; ++p) {
well_state.wellRates()[index_of_well_*np + p] *= group_constraint.second;
}
}
return group_constraint.first;
}
}
return false;
}
template <typename TypeTag>
std::pair<bool, double>
WellInterface<TypeTag>::checkGroupConstraintsInj(const Group& group,
const WellState& well_state,
const double efficiencyFactor,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const
{
// Translate injector type from control to Phase.
const auto& well_controls = well_ecl_.injectionControls(summaryState);
auto injectorType = well_controls.injector_type;
Phase injectionPhase;
switch (injectorType) {
case InjectorType::WATER:
{
injectionPhase = Phase::WATER;
break;
}
case InjectorType::OIL:
{
injectionPhase = Phase::OIL;
break;
}
case InjectorType::GAS:
{
injectionPhase = Phase::GAS;
break;
}
default:
throw("Expected WATER, OIL or GAS as type for injector " + name());
}
// Make conversion factors for RESV <-> surface rates.
std::vector<double> resv_coeff(phaseUsage().num_phases, 1.0);
rateConverter_.calcCoeff(0, pvtRegionIdx_, resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS.
// Call check for the well's injection phase.
return WellGroupHelpers::checkGroupConstraintsInj(name(),
well_ecl_.groupName(),
group,
well_state,
current_step_,
guide_rate_,
well_state.wellRates().data() + index_of_well_ * phaseUsage().num_phases,
injectionPhase,
phaseUsage(),
efficiencyFactor,
schedule,
summaryState,
resv_coeff,
deferred_logger);
}
template <typename TypeTag>
std::pair<bool, double>
WellInterface<TypeTag>::checkGroupConstraintsProd(const Group& group,
const WellState& well_state,
const double efficiencyFactor,
const Schedule& schedule,
const SummaryState& summaryState,
DeferredLogger& deferred_logger) const
{
// Make conversion factors for RESV <-> surface rates.
std::vector<double> resv_coeff(phaseUsage().num_phases, 1.0);
rateConverter_.calcCoeff(0, pvtRegionIdx_, resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS.
return WellGroupHelpers::checkGroupConstraintsProd(name(),
well_ecl_.groupName(),
group,
well_state,
current_step_,
guide_rate_,
well_state.wellRates().data() + index_of_well_ * phaseUsage().num_phases,
phaseUsage(),
efficiencyFactor,
schedule,
summaryState,
resv_coeff,
deferred_logger);
}
template <typename TypeTag>
template <class EvalWell, class BhpFromThpFunc>
void
WellInterface<TypeTag>::assembleControlEqInj(const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const Well::InjectionControls& controls,
const EvalWell& bhp,
const EvalWell& injection_rate,
BhpFromThpFunc bhp_from_thp,
EvalWell& control_eq,
Opm::DeferredLogger& deferred_logger)
{
const Opm::Well::InjectorCMode& current = well_state.currentInjectionControls()[index_of_well_];
const InjectorType injectorType = controls.injector_type;
const auto& pu = phaseUsage();
const double efficiencyFactor = well_ecl_.getEfficiencyFactor();
switch (current) {
case Well::InjectorCMode::RATE: {
control_eq = injection_rate - controls.surface_rate;
break;
}
case Well::InjectorCMode::RESV: {
std::vector<double> convert_coeff(number_of_phases_, 1.0);
rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff);
double coeff = 1.0;
switch (injectorType) {
case InjectorType::WATER: {
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Aqua]];
break;
}
case InjectorType::OIL: {
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Liquid]];
break;
}
case InjectorType::GAS: {
coeff = convert_coeff[pu.phase_pos[BlackoilPhases::Vapour]];
break;
}
default:
throw("Expected WATER, OIL or GAS as type for injectors " + well_ecl_.name());
}
control_eq = coeff * injection_rate - controls.reservoir_rate;
break;
}
case Well::InjectorCMode::THP: {
control_eq = bhp - bhp_from_thp();
break;
}
case Well::InjectorCMode::BHP: {
control_eq = bhp - controls.bhp_limit;
break;
}
case Well::InjectorCMode::GRUP: {
assert(well_ecl_.isAvailableForGroupControl());
const auto& group = schedule.getGroup(well_ecl_.groupName(), current_step_);
getGroupInjectionControl(group,
well_state,
schedule,
summaryState,
injectorType,
bhp,
injection_rate,
control_eq,
efficiencyFactor);
break;
}
case Well::InjectorCMode::CMODE_UNDEFINED: {
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger);
}
}
}
template <typename TypeTag>
template <class EvalWell, class BhpFromThpFunc>
void
WellInterface<TypeTag>::assembleControlEqProd(const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const Well::ProductionControls& controls,
const EvalWell& bhp,
const std::vector<EvalWell>& rates, // Always 3 canonical rates.
BhpFromThpFunc bhp_from_thp,
EvalWell& control_eq,
Opm::DeferredLogger& deferred_logger)
{
const Well::ProducerCMode& current = well_state.currentProductionControls()[index_of_well_];
const auto& pu = phaseUsage();
const double efficiencyFactor = well_ecl_.getEfficiencyFactor();
switch (current) {
case Well::ProducerCMode::ORAT: {
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
const EvalWell rate = -rates[BlackoilPhases::Liquid];
control_eq = rate - controls.oil_rate;
break;
}
case Well::ProducerCMode::WRAT: {
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
const EvalWell rate = -rates[BlackoilPhases::Aqua];
control_eq = rate - controls.water_rate;
break;
}
case Well::ProducerCMode::GRAT: {
assert(FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx));
const EvalWell rate = -rates[BlackoilPhases::Vapour];
control_eq = rate - controls.gas_rate;
break;
}
case Well::ProducerCMode::LRAT: {
assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx));
EvalWell rate = -rates[BlackoilPhases::Aqua] - rates[BlackoilPhases::Liquid];
control_eq = rate - controls.liquid_rate;
break;
}
case Well::ProducerCMode::CRAT: {
OPM_DEFLOG_THROW(std::runtime_error, "CRAT control not supported " << name(), deferred_logger);
}
case Well::ProducerCMode::RESV: {
auto total_rate = rates[0]; // To get the correct type only.
total_rate = 0.0;
std::vector<double> convert_coeff(number_of_phases_, 1.0);
rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff);
for (int phase = 0; phase < 3; ++phase) {
if (pu.phase_used[phase]) {
const int pos = pu.phase_pos[phase];
total_rate -= rates[phase] * convert_coeff[pos]; // Note different indices.
}
}
if (controls.prediction_mode) {
control_eq = total_rate - controls.resv_rate;
} else {
std::vector<double> hrates(number_of_phases_, 0.);
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
hrates[pu.phase_pos[Water]] = controls.water_rate;
}
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
hrates[pu.phase_pos[Oil]] = controls.oil_rate;
}
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
hrates[pu.phase_pos[Gas]] = controls.gas_rate;
}
std::vector<double> hrates_resv(number_of_phases_, 0.);
rateConverter_.calcReservoirVoidageRates(/*fipreg*/ 0, pvtRegionIdx_, hrates, hrates_resv);
double target = std::accumulate(hrates_resv.begin(), hrates_resv.end(), 0.0);
control_eq = total_rate - target;
}
break;
}
case Well::ProducerCMode::BHP: {
control_eq = bhp - controls.bhp_limit;
break;
}
case Well::ProducerCMode::THP: {
control_eq = bhp - bhp_from_thp();
break;
}
case Well::ProducerCMode::GRUP: {
assert(well_ecl_.isAvailableForGroupControl());
const auto& group = schedule.getGroup(well_ecl_.groupName(), current_step_);
// Annoying thing: the rates passed to this function are
// always of size 3 and in canonical (for PhaseUsage)
// order. This is what is needed for VFP calculations if
// they are required (THP controlled well). But for the
// group production control things we must pass only the
// active phases' rates.
std::vector<EvalWell> active_rates(pu.num_phases);
for (int canonical_phase = 0; canonical_phase < 3; ++canonical_phase) {
if (pu.phase_used[canonical_phase]) {
active_rates[pu.phase_pos[canonical_phase]] = rates[canonical_phase];
}
}
getGroupProductionControl(group, well_state, schedule, summaryState, bhp, active_rates, control_eq, efficiencyFactor);
break;
}
case Well::ProducerCMode::CMODE_UNDEFINED: {
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger);
}
case Well::ProducerCMode::NONE: {
OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + name(), deferred_logger);
}
}
}
template <typename TypeTag>
template <class EvalWell>
void
WellInterface<TypeTag>::getGroupInjectionControl(const Group& group,
const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const InjectorType& injectorType,
const EvalWell& bhp,
const EvalWell& injection_rate,
EvalWell& control_eq,
double efficiencyFactor)
{
const auto& well = well_ecl_;
const auto pu = phaseUsage();
// Setting some defaults to silence warnings below.
// Will be overwritten in the switch statement.
int phasePos = -1;
Well::GuideRateTarget wellTarget = Well::GuideRateTarget::UNDEFINED;
Phase injectionPhase = Phase::WATER;
switch (injectorType) {
case InjectorType::WATER:
{
phasePos = pu.phase_pos[BlackoilPhases::Aqua];
wellTarget = Well::GuideRateTarget::WAT;
injectionPhase = Phase::WATER;
break;
}
case InjectorType::OIL:
{
phasePos = pu.phase_pos[BlackoilPhases::Liquid];
wellTarget = Well::GuideRateTarget::OIL;
injectionPhase = Phase::OIL;
break;
}
case InjectorType::GAS:
{
phasePos = pu.phase_pos[BlackoilPhases::Vapour];
wellTarget = Well::GuideRateTarget::GAS;
injectionPhase = Phase::GAS;
break;
}
default:
// Should not be here.
assert(false);
}
const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name());
if (currentGroupControl == Group::InjectionCMode::FLD ||
currentGroupControl == Group::InjectionCMode::NONE) {
if (!group.isAvailableForGroupControl()) {
// We cannot go any further up the hierarchy. This could
// be the FIELD group, or any group for which this has
// been set in GCONINJE or GCONPROD. If we are here
// anyway, it is likely that the deck set inconsistent
// requirements, such as GRUP control mode on a well with
// no appropriate controls defined on any of its
// containing groups. We will therefore use the wells' bhp
// limit equation as a fallback.
const auto& controls = well_ecl_.injectionControls(summaryState);
control_eq = bhp - controls.bhp_limit;
return;
} else {
// Inject share of parents control
const auto& parent = schedule.getGroup( group.parent(), current_step_ );
efficiencyFactor *= group.getGroupEfficiencyFactor();
getGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, bhp, injection_rate, control_eq, efficiencyFactor);
return;
}
}
assert(group.hasInjectionControl(injectionPhase));
const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState);
const std::vector<double>& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name());
double groupTargetReduction = groupInjectionReductions[phasePos];
double fraction = WellGroupHelpers::fractionFromInjectionPotentials(well.name(),
group.name(),
schedule,
well_state,
current_step_,
guide_rate_,
GuideRateModel::convert_target(wellTarget),
pu,
injectionPhase,
false);
switch (currentGroupControl) {
case Group::InjectionCMode::NONE:
{
// The NONE case is handled earlier
assert(false);
break;
}
case Group::InjectionCMode::RATE:
{
double target = std::max(0.0, (groupcontrols.surface_max_rate - groupTargetReduction)) / efficiencyFactor;
control_eq = injection_rate - fraction * target;
break;
}
case Group::InjectionCMode::RESV:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff);
double coeff = convert_coeff[phasePos];
double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff - groupTargetReduction)) / efficiencyFactor;
control_eq = injection_rate - fraction * target;
break;
}
case Group::InjectionCMode::REIN:
{
double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos];
double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction)) / efficiencyFactor;
control_eq = injection_rate - fraction * target;
break;
}
case Group::InjectionCMode::VREP:
{
std::vector<double> convert_coeff(number_of_phases_, 1.0);
rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff);
double coeff = convert_coeff[phasePos];
double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction;
double injReduction = 0.0;
std::vector<double> groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name());
if (groupcontrols.phase != Phase::WATER)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]];
if (groupcontrols.phase != Phase::OIL)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]];
if (groupcontrols.phase != Phase::GAS)
injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]];
voidageRate -= injReduction;
double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction)) / efficiencyFactor;
control_eq = injection_rate - fraction * target;
break;
}
case Group::InjectionCMode::FLD:
{
// The FLD case is handled earlier
assert(false);
break;
}
case Group::InjectionCMode::SALE:
{
// only for gas injectors
assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]);
// Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate;
double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos];
if (schedule.gConSump(current_step_).has(group.name())) {
const auto& gconsump = schedule.gConSump(current_step_).get(group.name(), summaryState);
if (pu.phase_used[BlackoilPhases::Vapour]) {
inj_rate += gconsump.import_rate;
inj_rate -= gconsump.consumption_rate;
}
}
const auto& gconsale = schedule.gConSale(current_step_).get(group.name(), summaryState);
inj_rate -= gconsale.sales_target;
double target = std::max(0.0, (inj_rate - groupTargetReduction)) / efficiencyFactor;
control_eq = injection_rate - fraction * target;
break;
}
// default:
// OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger );
}
}
template <typename TypeTag>
template <class EvalWell>
void
WellInterface<TypeTag>::getGroupProductionControl(const Group& group,
const WellState& well_state,
const Opm::Schedule& schedule,
const SummaryState& summaryState,
const EvalWell& bhp,
const std::vector<EvalWell>& rates,
EvalWell& control_eq,
double efficiencyFactor)
{
const Group::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name());
if (currentGroupControl == Group::ProductionCMode::FLD ||
currentGroupControl == Group::ProductionCMode::NONE) {
if (!group.isAvailableForGroupControl()) {
// We cannot go any further up the hierarchy. This could
// be the FIELD group, or any group for which this has
// been set in GCONINJE or GCONPROD. If we are here
// anyway, it is likely that the deck set inconsistent
// requirements, such as GRUP control mode on a well with
// no appropriate controls defined on any of its
// containing groups. We will therefore use the wells' bhp
// limit equation as a fallback.
const auto& controls = well_ecl_.productionControls(summaryState);
control_eq = bhp - controls.bhp_limit;
return;
} else {
// Produce share of parents control
const auto& parent = schedule.getGroup( group.parent(), current_step_ );
efficiencyFactor *= group.getGroupEfficiencyFactor();
getGroupProductionControl(parent, well_state, schedule, summaryState, bhp, rates, control_eq, efficiencyFactor);
return;
}
}
const auto& well = well_ecl_;
const auto pu = phaseUsage();
if (!group.isProductionGroup()) {
// use bhp as control eq and let the updateControl code find a valid control
const auto& controls = well.productionControls(summaryState);
control_eq = bhp - controls.bhp_limit;
return;
}
// If we are here, we are at the topmost group to be visited in the recursion.
// This is the group containing the control we will check against.
// Make conversion factors for RESV <-> surface rates.
std::vector<double> resv_coeff(phaseUsage().num_phases, 1.0);
rateConverter_.calcCoeff(0, pvtRegionIdx_, resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS.
WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff);
WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu);
auto localFraction = [&](const std::string& child) {
return fcalc.localFraction(child, "");
};
auto localReduction = [&](const std::string& group_name) {
const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group_name);
return tcalc.calcModeRateFromRates(groupTargetReductions);
};
const double orig_target = tcalc.groupTarget(group.productionControls(summaryState));
const auto chain = WellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, current_step_);
// Because 'name' is the last of the elements, and not an ancestor, we subtract one below.
const size_t num_ancestors = chain.size() - 1;
double target = orig_target;
for (size_t ii = 0; ii < num_ancestors; ++ii) {
target -= localReduction(chain[ii]);
target *= localFraction(chain[ii+1]);
}
// Avoid negative target rates comming from too large local reductions.
const double target_rate = std::max(0.0, target / efficiencyFactor);
const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers.
control_eq = current_rate - target_rate;
}
} // namespace Opm

72
opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp
View File

@ -350,8 +350,33 @@ namespace Opm
return it->second;
}
void setCurrentWellRates(const std::string& wellName, const std::vector<double>& rates ) {
well_rates[wellName] = rates;
}
const std::vector<double>& currentWellRates(const std::string& wellName) const {
auto it = well_rates.find(wellName);
if (it == well_rates.end())
OPM_THROW(std::logic_error, "Could not find any rates for well " << wellName);
return it->second;
}
void setCurrentProductionGroupRates(const std::string& groupName, const std::vector<double>& rates ) {
production_group_rates[groupName] = rates;
}
const std::vector<double>& currentProductionGroupRates(const std::string& groupName) const {
auto it = production_group_rates.find(groupName);
if (it == production_group_rates.end())
OPM_THROW(std::logic_error, "Could not find any rates for productino group " << groupName);
return it->second;
}
void setCurrentProductionGroupReductionRates(const std::string& groupName, const std::vector<double>& target ) {
production_group_reduction_rates[groupName] = target;
}
@ -880,27 +905,46 @@ namespace Opm
for (auto& x : injection_group_reservoir_rates) {
comm.sum(x.second.data(), x.second.size());
}
for (auto& x : production_group_rates) {
comm.sum(x.second.data(), x.second.size());
}
for (auto& x : well_rates) {
comm.sum(x.second.data(), x.second.size());
}
}
template<class Comm>
void updateGlobalIsGrup(const Schedule& schedule, const int reportStepIdx, const Comm& comm) {
int global_well_index = -1;
void updateGlobalIsGrup(const Schedule& schedule, const int reportStepIdx, const Comm& comm)
{
std::fill(globalIsInjectionGrup_.begin(), globalIsInjectionGrup_.end(), 0);
std::fill(globalIsProductionGrup_.begin(), globalIsProductionGrup_.end(), 0);
int global_well_index = 0;
const auto& end = wellMap().end();
for (const auto& well : schedule.getWells(reportStepIdx)) {
global_well_index ++;
// Build global name->index map.
wellNameToGlobalIdx_[well.name()] = global_well_index;
// For wells on this process...
const auto& it = wellMap().find( well.name());
if (it == end) // the well is not found
continue;
int well_index = it->second[0];
if (well.isInjector())
globalIsInjectionGrup_[global_well_index] = (currentInjectionControls()[well_index] == Well::InjectorCMode::GRUP);
else
globalIsProductionGrup_[global_well_index] = (currentProductionControls()[well_index] == Well::ProducerCMode::GRUP);
if (it != end) {
// ... set the GRUP/not GRUP states.
const int well_index = it->second[0];
if (!this->open_for_output_[well_index]) {
// Well is shut.
if (well.isInjector()) {
globalIsInjectionGrup_[global_well_index] = 0;
} else {
globalIsProductionGrup_[global_well_index] = 0;
}
} else {
if (well.isInjector()) {
globalIsInjectionGrup_[global_well_index] = (currentInjectionControls()[well_index] == Well::InjectorCMode::GRUP);
} else {
globalIsProductionGrup_[global_well_index] = (currentProductionControls()[well_index] == Well::ProducerCMode::GRUP);
}
}
}
++global_well_index;
}
comm.sum(globalIsInjectionGrup_.data(), globalIsInjectionGrup_.size());
comm.sum(globalIsProductionGrup_.data(), globalIsProductionGrup_.size());
@ -939,6 +983,8 @@ namespace Opm
std::map<std::string, Group::ProductionCMode> current_production_group_controls_;
std::map<std::pair<Opm::Phase, std::string>, Group::InjectionCMode> current_injection_group_controls_;
std::map<std::string, std::vector<double>> well_rates;
std::map<std::string, std::vector<double>> production_group_rates;
std::map<std::string, std::vector<double>> production_group_reduction_rates;
std::map<std::string, std::vector<double>> injection_group_reduction_rates;
std::map<std::string, std::vector<double>> injection_group_reservoir_rates;