Merge pull request #1501 from totto82/wtest

Implement WTEST support
This commit is contained in:
Kai Bao
2018-07-04 10:45:27 +02:00
committed by GitHub
8 changed files with 397 additions and 136 deletions

View File

@@ -259,6 +259,12 @@ add_test_compareECLFiles(CASENAME spe5
REL_TOL ${coarse_rel_tol}
TEST_ARGS max_iter=20)
add_test_compareECLFiles(CASENAME wecon_wtest
FILENAME 3D_WECON
SIMULATOR flow
ABS_TOL ${abs_tol}
REL_TOL ${coarse_rel_tol})
# Restart tests
opm_set_test_driver(${PROJECT_SOURCE_DIR}/tests/run-restart-regressionTest.sh "")

View File

@@ -206,7 +206,7 @@ namespace Opm {
wasSwitched_.resize(numDof);
std::fill(wasSwitched_.begin(), wasSwitched_.end(), false);
wellModel().beginTimeStep();
wellModel().beginTimeStep(timer.reportStepNum(), timer.simulationTimeElapsed());
if (param_.update_equations_scaling_) {
std::cout << "equation scaling not suported yet" << std::endl;
@@ -340,7 +340,7 @@ namespace Opm {
/// \param[in] timer simulation timer
void afterStep(const SimulatorTimerInterface& OPM_UNUSED timer)
{
wellModel().timeStepSucceeded();
wellModel().timeStepSucceeded(timer.simulationTimeElapsed());
aquiferModel().timeStepSucceeded(timer);
ebosSimulator_.problem().endTimeStep();

View File

@@ -2,7 +2,7 @@
Copyright 2016 SINTEF ICT, Applied Mathematics.
Copyright 2016 - 2017 Statoil ASA.
Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
Copyright 2016 - 2017 IRIS AS
Copyright 2016 - 2018 IRIS AS
This file is part of the Open Porous Media project (OPM).
@@ -33,6 +33,7 @@
#include <tuple>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/WellTestState.hpp>
#include <opm/core/wells.h>
#include <opm/core/wells/DynamicListEconLimited.hpp>
@@ -141,7 +142,7 @@ namespace Opm {
// compute the well fluxes and assemble them in to the reservoir equations as source terms
// and in the well equations.
void assemble(const int iterationIdx,
const double dt);
const double dt);
// substract Binv(D)rw from r;
void apply( BVector& r) const;
@@ -175,9 +176,9 @@ namespace Opm {
void setRestartWellState(const WellState& well_state);
// called at the beginning of a time step
void beginTimeStep();
void beginTimeStep(const int timeStepIdx,const double simulationTime);
// called at the end of a time step
void timeStepSucceeded();
void timeStepSucceeded(const double& simulationTime);
// called at the beginning of a report step
void beginReportStep(const int time_step);
@@ -236,7 +237,7 @@ namespace Opm {
using ConvergenceReport = typename WellInterface<TypeTag>::ConvergenceReport;
// create the well container
std::vector<WellInterfacePtr > createWellContainer(const int time_step) const;
std::vector<WellInterfacePtr > createWellContainer(const int time_step);
WellState well_state_;
WellState previous_well_state_;
@@ -254,12 +255,13 @@ namespace Opm {
std::vector<double> depth_;
bool initial_step_;
DynamicListEconLimited dynamic_list_econ_limited_;
std::unique_ptr<RateConverterType> rateConverter_;
std::unique_ptr<VFPProperties> vfp_properties_;
SimulatorReport last_report_;
WellTestState wellTestState_;
// used to better efficiency of calcuation
mutable BVector scaleAddRes_;
@@ -343,11 +345,13 @@ namespace Opm {
/// return true if wells are available on this process
bool localWellsActive() const;
/// upate the dynamic lists related to economic limits
void updateListEconLimited(DynamicListEconLimited& list_econ_limited) const;
/// upate the wellTestState related to economic limits
void updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const;
void updatePerforationIntensiveQuantities();
void wellTesting(const int timeStepIdx, const double simulationTime);
};

View File

@@ -45,6 +45,10 @@ namespace Opm {
wells_ecl_ = schedule().getWells(timeStepIdx);
// Create wells and well state.
// Pass empty dynamicListEconLimited class
// The closing of wells due to limites is
// handled by the wellTestState class
DynamicListEconLimited dynamic_list_econ_limited;
wells_manager_.reset( new WellsManager (eclState,
schedule(),
timeStepIdx,
@@ -54,7 +58,7 @@ namespace Opm {
Opm::UgGridHelpers::dimensions(grid),
Opm::UgGridHelpers::cell2Faces(grid),
Opm::UgGridHelpers::beginFaceCentroids(grid),
dynamic_list_econ_limited_,
dynamic_list_econ_limited,
grid.comm().size() > 1,
defunct_well_names) );
@@ -103,11 +107,35 @@ namespace Opm {
}
}
// update the previous well state. This is used to restart failed steps.
previous_well_state_ = well_state_;
// Compute reservoir volumes for RESV controls.
rateConverter_.reset(new RateConverterType (phase_usage_,
std::vector<int>(number_of_cells_, 0)));
std::vector<int>(number_of_cells_, 0)));
computeRESV(timeStepIdx);
// update VFP properties
vfp_properties_.reset (new VFPProperties (
schedule().getVFPInjTables(timeStepIdx),
schedule().getVFPProdTables(timeStepIdx)) );
}
// called at the beginning of a time step
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
beginTimeStep(const int timeStepIdx, const double simulationTime) {
well_state_ = previous_well_state_;
// test wells
wellTesting(timeStepIdx, simulationTime);
// create the well container
well_container_ = createWellContainer(timeStepIdx);
@@ -123,37 +151,124 @@ namespace Opm {
if (has_polymer_)
{
const Grid& grid = ebosSimulator_.vanguard().grid();
if (PolymerModule::hasPlyshlog()) {
computeRepRadiusPerfLength(grid);
}
}
// update VFP properties
vfp_properties_.reset (new VFPProperties (
schedule().getVFPInjTables(timeStepIdx),
schedule().getVFPProdTables(timeStepIdx)) );
for (auto& well : well_container_) {
well->setVFPProperties(vfp_properties_.get());
}
// update the previous well state. This is used to restart failed steps.
previous_well_state_ = well_state_;
// Close wells and completions due to economical reasons
for (auto& well : well_container_) {
well->closeCompletions(wellTestState_);
}
}
// called at the beginning of a time step
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
beginTimeStep() {
well_state_ = previous_well_state_;
BlackoilWellModel<TypeTag>::wellTesting(const int timeStepIdx, const double simulationTime) {
const auto& wtest_config = schedule().wtestConfig(timeStepIdx);
const auto& wellsForTesting = wellTestState_.updateWell(wtest_config, simulationTime);
if (wellCollection().havingVREPGroups() ) {
rateConverter_->template defineState<ElementContext>(ebosSimulator_);
// Do the well testing if enabled
if (wtest_config.size() > 0 && wellsForTesting.size() > 0) {
// solve the well equation isolated from the reservoir.
const int numComp = numComponents();
std::vector< Scalar > B_avg( numComp, Scalar() );
computeAverageFormationFactor(B_avg);
std::vector<WellInterfacePtr> well_container;
well_container.reserve(wellsForTesting.size());
for (auto& testWell : wellsForTesting) {
const std::string msg = std::string("well ") + testWell.first + std::string(" is tested");
OpmLog::info(msg);
// Finding the location of the well in wells_ecl
const int nw_wells_ecl = wells_ecl_.size();
int index_well = 0;
for (; index_well < nw_wells_ecl; ++index_well) {
if (testWell.first == wells_ecl_[index_well]->name()) {
break;
}
}
// It should be able to find in wells_ecl.
if (index_well == nw_wells_ecl) {
OPM_THROW(std::logic_error, "Could not find well " << testWell.first << " in wells_ecl ");
}
const Well* well_ecl = wells_ecl_[index_well];
// Finding the location of the well in wells struct.
const int nw = numWells();
int wellidx = -999;
for (int w = 0; w < nw; ++w) {
if (testWell.first == std::string(wells()->name[w])) {
wellidx = w;
break;
}
}
if (wellidx < 0) {
OPM_THROW(std::logic_error, "Could not find the well " << testWell.first << " in the well struct ");
}
// Use the pvtRegionIdx from the top cell
const int well_cell_top = wells()->well_cells[wells()->well_connpos[wellidx]];
const int pvtreg = pvt_region_idx_[well_cell_top];
if ( !well_ecl->isMultiSegment(timeStepIdx) || !param_.use_multisegment_well_) {
well_container.emplace_back(new StandardWell<TypeTag>(well_ecl, timeStepIdx, wells(),
param_, *rateConverter_, pvtreg, numComponents() ) );
} else {
well_container.emplace_back(new MultisegmentWell<TypeTag>(well_ecl, timeStepIdx, wells(),
param_, *rateConverter_, pvtreg, numComponents() ) );
}
}
for (auto& well : well_container) {
WellTestState wellTestStateForTheWellTest;
WellState wellStateCopy = well_state_;
well->init(&phase_usage_, depth_, gravity_, number_of_cells_);
const std::string& well_name = well->name();
const WellNode& well_node = wellCollection().findWellNode(well_name);
const double well_efficiency_factor = well_node.getAccumulativeEfficiencyFactor();
well->setWellEfficiencyFactor(well_efficiency_factor);
well->setVFPProperties(vfp_properties_.get());
well->updatePrimaryVariables(wellStateCopy);
well->initPrimaryVariablesEvaluation();
bool testWell = true;
// if a well is closed because all completions are closed, we need to check each completion
// individually. We first open all completions, then we close one by one by calling updateWellTestState
// untill the number of closed completions do not increase anymore.
while (testWell) {
const size_t numberOfClosedCompletions = wellTestStateForTheWellTest.sizeCompletions();
well->solveWellForTesting(ebosSimulator_, wellStateCopy, B_avg, terminal_output_);
well->updateWellTestState(wellStateCopy, simulationTime, wellTestStateForTheWellTest, /*writeMessageToOPMLog=*/ false);
well->closeCompletions(wellTestStateForTheWellTest);
// Stop testing if the well is closed or shut due to all completions shut
// Also check if number of completions has increased. If the number of closed completions do not increased
// we stop the testing.
if (wellTestStateForTheWellTest.sizeWells() > 0 || numberOfClosedCompletions == wellTestStateForTheWellTest.sizeCompletions())
testWell = false;
}
// update wellTestState if the well test succeeds
if (!wellTestStateForTheWellTest.hasWell(well->name(), WellTestConfig::Reason::ECONOMIC)) {
wellTestState_.openWell(well->name());
const std::string msg = std::string("well ") + well->name() + std::string(" is re-opened");
OpmLog::info(msg);
// also reopen completions
for (auto& completion : well->wellEcl()->getCompletions(timeStepIdx)) {
if (!wellTestStateForTheWellTest.hasCompletion(well->name(), completion.first))
wellTestState_.dropCompletion(well->name(), completion.first);
}
}
}
}
}
@@ -168,10 +283,6 @@ namespace Opm {
void
BlackoilWellModel<TypeTag>::
endReportStep() {
// update the list contanining information of closed wells
// and connections due to economical limits
// Used by the wellManager
updateListEconLimited(dynamic_list_econ_limited_);
}
// called at the end of a report step
@@ -184,20 +295,20 @@ namespace Opm {
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
timeStepSucceeded() {
timeStepSucceeded(const double& simulationTime) {
// TODO: when necessary
rateConverter_->template defineState<ElementContext>(ebosSimulator_);
for (const auto& well : well_container_) {
well->calculateReservoirRates(well_state_);
}
updateWellTestState(simulationTime, wellTestState_);
previous_well_state_ = well_state_;
}
template<typename TypeTag>
std::vector<typename BlackoilWellModel<TypeTag>::WellInterfacePtr >
BlackoilWellModel<TypeTag>::
createWellContainer(const int time_step) const
createWellContainer(const int time_step)
{
std::vector<WellInterfacePtr> well_container;
@@ -227,6 +338,24 @@ namespace Opm {
const Well* well_ecl = wells_ecl_[index_well];
// well is closed due to economical reasons
if (wellTestState_.hasWell(well_name, WellTestConfig::Reason::ECONOMIC)) {
if( well_ecl->getAutomaticShutIn() ) {
// shut wells are not added to the well container
well_state_.bhp()[w] = 0;
const int np = numPhases();
for (int p = 0; p < np; ++p) {
well_state_.wellRates()[np * w + p] = 0;
}
continue;
}
else {
// close wells are added to the container but marked as closed
struct WellControls* well_controls = wells()->ctrls[w];
well_controls_stop_well(well_controls);
}
}
// Use the pvtRegionIdx from the top cell
const int well_cell_top = wells()->well_cells[wells()->well_connpos[w]];
const int pvtreg = pvt_region_idx_[well_cell_top];
@@ -275,6 +404,7 @@ namespace Opm {
if (param_.solve_welleq_initially_ && iterationIdx == 0) {
// solve the well equations as a pre-processing step
last_report_ = solveWellEq(dt);
if (initial_step_) {
// update the explicit quantities to get the initial fluid distribution in the well correct.
calculateExplicitQuantities();
@@ -301,8 +431,8 @@ namespace Opm {
assembleWellEq(const double dt,
bool only_wells)
{
for (int w = 0; w < numWells(); ++w) {
well_container_[w]->assembleWellEq(ebosSimulator_, dt, well_state_, only_wells);
for (auto& well : well_container_) {
well->assembleWellEq(ebosSimulator_, dt, well_state_, only_wells);
}
}
@@ -392,13 +522,9 @@ namespace Opm {
BlackoilWellModel<TypeTag>::
resetWellControlFromState() const
{
const int nw = numWells();
assert(nw == int(well_container_.size()) );
for (int w = 0; w < nw; ++w) {
WellControls* wc = well_container_[w]->wellControls();
well_controls_set_current( wc, well_state_.currentControls()[w]);
for (auto& well : well_container_) {
WellControls* wc = well->wellControls();
well_controls_set_current( wc, well_state_.currentControls()[well->indexOfWell()]);
}
}
@@ -475,6 +601,7 @@ namespace Opm {
do {
assembleWellEq(dt, true);
//std::cout << "well convergence only wells " << std::endl;
converged = getWellConvergence(B_avg);
// checking whether the group targets are converged
@@ -632,10 +759,10 @@ namespace Opm {
template<typename TypeTag>
void
BlackoilWellModel<TypeTag>::
updateListEconLimited(DynamicListEconLimited& list_econ_limited) const
updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const
{
for (const auto& well : well_container_) {
well->updateListEconLimited(well_state_, list_econ_limited);
well->updateWellTestState(well_state_, simulationTime, wellTestState, /*writeMessageToOPMLog=*/ true);
}
}
@@ -651,13 +778,13 @@ namespace Opm {
const int np = numPhases();
well_potentials.resize(nw * np, 0.0);
for (int w = 0; w < nw; ++w) {
for (const auto& well : well_container_) {
std::vector<double> potentials;
well_container_[w]->computeWellPotentials(ebosSimulator_, well_state_, potentials);
well->computeWellPotentials(ebosSimulator_, well_state_, potentials);
// putting the sucessfully calculated potentials to the well_potentials
for (int p = 0; p < np; ++p) {
well_potentials[w * np + p] = std::abs(potentials[p]);
well_potentials[well->indexOfWell() * np + p] = std::abs(potentials[p]);
}
} // end of for (int w = 0; w < nw; ++w)
}
@@ -687,13 +814,14 @@ namespace Opm {
prepareGroupControl();
// since the controls are all updated, we should update well_state accordingly
for (int w = 0; w < numWells(); ++w) {
WellControls* wc = well_container_[w]->wellControls();
for (const auto& well : well_container_) {
const int w = well->indexOfWell();
WellControls* wc = well->wellControls();
const int control = well_controls_get_current(wc);
well_state_.currentControls()[w] = control;
// TODO: for VFP control, the perf_densities are still zero here, investigate better
// way to handle it later.
well_container_[w]->updateWellStateWithTarget(well_state_);
well->updateWellStateWithTarget(well_state_);
// The wells are not considered to be newly added
// for next time step
@@ -702,6 +830,7 @@ namespace Opm {
}
} // end of for (int w = 0; w < nw; ++w)
}
@@ -718,9 +847,9 @@ namespace Opm {
{
// group control related processing
if (wellCollection().groupControlActive()) {
for (int w = 0; w < numWells(); ++w) {
WellControls* wc = well_container_[w]->wellControls();
WellNode& well_node = wellCollection().findWellNode(well_container_[w]->name());
for (const auto& well : well_container_) {
WellControls* wc = well->wellControls();
WellNode& well_node = wellCollection().findWellNode(well->name());
// handling the situation that wells do not have a valid control
// it happens the well specified with GRUP and restarting due to non-convergencing
@@ -731,7 +860,7 @@ namespace Opm {
if (group_control_index >= 0 && ctrl_index < 0) {
// put well under group control
well_controls_set_current(wc, group_control_index);
well_state_.currentControls()[w] = group_control_index;
well_state_.currentControls()[well->indexOfWell()] = group_control_index;
}
// Final step, update whehter the well is under group control or individual control
@@ -804,15 +933,13 @@ namespace Opm {
return;
}
const int nw = numWells();
for (int w = 0; w < nw; ++w) {
const std::string well_name = well_container_[w]->name();
for (auto& well : well_container_) {
const std::string& well_name = well->name();
const WellNode& well_node = wellCollection().findWellNode(well_name);
const double well_efficiency_factor = well_node.getAccumulativeEfficiencyFactor();
well_container_[w]->setWellEfficiencyFactor(well_efficiency_factor);
well->setWellEfficiencyFactor(well_efficiency_factor);
}
}
@@ -846,9 +973,10 @@ namespace Opm {
std::vector<double> well_rates(np, 0.0);
std::vector<double> convert_coeff(np, 1.0);
for (int w = 0; w < nw; ++w) {
const bool is_producer = well_container_[w]->wellType() == PRODUCER;
const int well_cell_top = well_container_[w]->cells()[0];
for (auto& well : well_container_) {
const bool is_producer = well->wellType() == PRODUCER;
const int well_cell_top =well->cells()[0];
const int w = well->indexOfWell();
const int pvtRegionIdx = pvt_region_idx_[well_cell_top];
// not sure necessary to change all the value to be positive
@@ -917,13 +1045,13 @@ namespace Opm {
{
if (wellCollection().groupControlActive()) {
for (int w = 0; w < numWells(); ++w) {
for (auto& well : well_container_) {
// update whether well is under group control
// get well node in the well collection
WellNode& well_node = wellCollection().findWellNode(well_container_[w]->name());
WellNode& well_node = wellCollection().findWellNode(well->name());
// update whehter the well is under group control or individual control
const int current = well_state_.currentControls()[w];
const int current = well_state_.currentControls()[well->indexOfWell()];
if (well_node.groupControlIndex() >= 0 && current == well_node.groupControlIndex()) {
// under group control
well_node.setIndividualControl(false);
@@ -938,8 +1066,8 @@ namespace Opm {
// it will not change the control mode, only update the targets
wellCollection().updateWellTargets(well_state_.wellRates());
for (int w = 0; w < numWells(); ++w) {
well_container_[w]->updateWellStateWithTarget(well_state_);
for (auto& well : well_container_) {
well->updateWellStateWithTarget(well_state_);
}
}
}

View File

@@ -338,6 +338,7 @@ namespace Opm
// handle the non reasonable fractions due to numerical overshoot
void processFractions() const;
};
}

View File

@@ -28,10 +28,13 @@
#include <opm/common/Exceptions.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Well.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/WellTestState.hpp>
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/core/props/BlackoilPhases.hpp>
#include <opm/core/wells/WellsManager.hpp>
#include <opm/core/simulator/SimulatorReport.hpp>
#include <opm/autodiff/VFPProperties.hpp>
#include <opm/autodiff/VFPInjProperties.hpp>
@@ -113,6 +116,9 @@ namespace Opm
/// Well name.
const std::string& name() const;
/// Index of well in the wells struct and wellState
const int indexOfWell() const;
/// Well cells.
const std::vector<int>& cells() {return well_cells_; }
@@ -171,8 +177,10 @@ namespace Opm
WellState& well_state,
bool only_wells) = 0;
void updateListEconLimited(const WellState& well_state,
DynamicListEconLimited& list_econ_limited) const;
void updateWellTestState(const WellState& well_state,
const double& simulationTime,
WellTestState& wellTestState,
const bool& writeMessageToOPMLog) const;
void setWellEfficiencyFactor(const double efficiency_factor);
@@ -216,10 +224,18 @@ namespace Opm
// Add well contributions to matrix
virtual void addWellContributions(Mat&) const
{}
void solveWellForTesting(Simulator& ebosSimulator, WellState& well_state, const std::vector<double>& B_avg, bool terminal_output);
void closeCompletions(WellTestState& wellTestState);
const Well* wellEcl() const;
protected:
// to indicate a invalid connection
static const int INVALIDCONNECTION = -100000;
// to indicate a invalid completion
static const int INVALIDCOMPLETION = INT_MAX;
const Well* well_ecl_;
@@ -316,13 +332,12 @@ namespace Opm
double mostStrictBhpFromBhpLimits() const;
// a tuple type for ratio limit check.
// first value indicates whether ratio limit is violated, when the ratio limit is not violated, the following three
// first value indicates whether ratio limit is violated, when the ratio limit is not violated, the following two
// values should not be used.
// second value indicates whehter there is only one connection left.
// third value indicates the indx of the worst-offending connection.
// the last value indicates the extent of the violation for the worst-offending connection, which is defined by
// second value indicates the index of the worst-offending completion.
// the last value indicates the extent of the violation for the worst-offending completion, which is defined by
// the ratio of the actual value to the value of the violated limit.
using RatioCheckTuple = std::tuple<bool, bool, int, double>;
using RatioCheckTuple = std::tuple<bool, int, double>;
RatioCheckTuple checkMaxWaterCutLimit(const WellEconProductionLimits& econ_production_limits,
const WellState& well_state) const;
@@ -332,6 +347,7 @@ namespace Opm
double scalingFactor(const int comp_idx) const;
};
}

View File

@@ -1,6 +1,7 @@
/*
Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
Copyright 2017 Statoil ASA.
Copyright 2018 IRIS
This file is part of the Open Porous Media project (OPM).
@@ -101,7 +102,6 @@ namespace Opm
wells->sat_table_id + perf_index_end,
saturation_table_number_.begin() );
}
well_efficiency_factor_ = 1.0;
}
@@ -169,6 +169,15 @@ namespace Opm
return well_controls_;
}
template<typename TypeTag>
const int
WellInterface<TypeTag>::
indexOfWell() const
{
return index_of_well_;
}
@@ -194,6 +203,14 @@ namespace Opm
template<typename TypeTag>
const Well*
WellInterface<TypeTag>::
wellEcl() const
{
return well_ecl_;
}
template<typename TypeTag>
@@ -485,8 +502,7 @@ namespace Opm
const WellState& well_state) const
{
bool water_cut_limit_violated = false;
int worst_offending_connection = INVALIDCONNECTION;
bool last_connection = false;
int worst_offending_completion = INVALIDCOMPLETION;
double violation_extent = -1.0;
const int np = number_of_phases_;
@@ -528,29 +544,37 @@ namespace Opm
water_cut_perf[perf] = 0.;
}
}
const auto& completions = well_ecl_->getCompletions(current_step_);
const auto& connections = well_ecl_->getConnections(current_step_);
last_connection = (perf_number == 1);
if (last_connection) {
worst_offending_connection = 0;
violation_extent = water_cut_perf[0] / max_water_cut_limit;
return std::make_tuple(water_cut_limit_violated, last_connection, worst_offending_connection, violation_extent);
int complnumIdx = 0;
std::vector<double> water_cut_in_completions(completions.size(), 0.0);
for (const auto& completion : completions) {
int complnum = completion.first;
for (int perf = 0; perf < perf_number; ++perf) {
if (complnum == connections.get ( perf ).complnum) {
water_cut_in_completions[complnumIdx] += water_cut_perf[perf];
}
}
complnumIdx++;
}
double max_water_cut_perf = 0.;
for (int perf = 0; perf < perf_number; ++perf) {
if (water_cut_perf[perf] > max_water_cut_perf) {
worst_offending_connection = perf;
max_water_cut_perf = water_cut_perf[perf];
complnumIdx = 0;
for (const auto& completion : completions) {
if (water_cut_in_completions[complnumIdx] > max_water_cut_perf) {
worst_offending_completion = completion.first;
max_water_cut_perf = water_cut_in_completions[complnumIdx];
}
complnumIdx++;
}
assert(max_water_cut_perf != 0.);
assert((worst_offending_connection >= 0) && (worst_offending_connection < perf_number));
assert(max_water_cut_limit != 0.);
assert(worst_offending_completion != INVALIDCOMPLETION);
violation_extent = max_water_cut_perf / max_water_cut_limit;
}
return std::make_tuple(water_cut_limit_violated, last_connection, worst_offending_connection, violation_extent);
return std::make_tuple(water_cut_limit_violated, worst_offending_completion, violation_extent);
}
@@ -563,17 +587,16 @@ namespace Opm
checkRatioEconLimits(const WellEconProductionLimits& econ_production_limits,
const WellState& well_state) const
{
// TODO: not sure how to define the worst-offending connection when more than one
// TODO: not sure how to define the worst-offending completion when more than one
// ratio related limit is violated.
// The defintion used here is that we define the violation extent based on the
// ratio between the value and the corresponding limit.
// For each violated limit, we decide the worst-offending connection separately.
// Among the worst-offending connections, we use the one has the biggest violation
// For each violated limit, we decide the worst-offending completion separately.
// Among the worst-offending completions, we use the one has the biggest violation
// extent.
bool any_limit_violated = false;
bool last_connection = false;
int worst_offending_connection = INVALIDCONNECTION;
int worst_offending_completion = INVALIDCOMPLETION;
double violation_extent = -1.0;
if (econ_production_limits.onMaxWaterCut()) {
@@ -581,11 +604,10 @@ namespace Opm
bool water_cut_violated = std::get<0>(water_cut_return);
if (water_cut_violated) {
any_limit_violated = true;
const double violation_extent_water_cut = std::get<3>(water_cut_return);
const double violation_extent_water_cut = std::get<2>(water_cut_return);
if (violation_extent_water_cut > violation_extent) {
violation_extent = violation_extent_water_cut;
worst_offending_connection = std::get<2>(water_cut_return);
last_connection = std::get<1>(water_cut_return);
worst_offending_completion = std::get<1>(water_cut_return);
}
}
}
@@ -603,11 +625,11 @@ namespace Opm
}
if (any_limit_violated) {
assert(worst_offending_connection >=0);
assert(worst_offending_completion != INVALIDCOMPLETION);
assert(violation_extent > 1.);
}
return std::make_tuple(any_limit_violated, last_connection, worst_offending_connection, violation_extent);
return std::make_tuple(any_limit_violated, worst_offending_completion, violation_extent);
}
@@ -617,8 +639,10 @@ namespace Opm
template<typename TypeTag>
void
WellInterface<TypeTag>::
updateListEconLimited(const WellState& well_state,
DynamicListEconLimited& list_econ_limited) const
updateWellTestState(const WellState& well_state,
const double& simulationTime,
WellTestState& wellTestState,
const bool& writeMessageToOPMLog) const
{
// economic limits only apply for production wells.
if (wellType() != PRODUCER) {
@@ -660,14 +684,15 @@ namespace Opm
OpmLog::warning("NOT_SUPPORTING_FOLLOWONWELL", "opening following on well after well closed is not supported yet");
}
if (well_ecl_->getAutomaticShutIn()) {
list_econ_limited.addShutWell(well_name);
const std::string msg = std::string("well ") + well_name + std::string(" will be shut in due to rate economic limit");
wellTestState.addClosedWell(well_name, WellTestConfig::Reason::ECONOMIC, simulationTime);
if (writeMessageToOPMLog) {
if (well_ecl_->getAutomaticShutIn()) {
const std::string msg = std::string("well ") + well_name + std::string(" will be shut due to rate economic limit");
OpmLog::info(msg);
} else {
list_econ_limited.addStoppedWell(well_name);
const std::string msg = std::string("well ") + well_name + std::string(" will be stopped due to rate economic limit");
OpmLog::info(msg);
} else {
const std::string msg = std::string("well ") + well_name + std::string(" will be stopped due to rate economic limit");
OpmLog::info(msg);
}
}
// the well is closed, not need to check other limits
return;
@@ -687,36 +712,56 @@ namespace Opm
switch (workover) {
case WellEcon::CON:
{
const bool last_connection = std::get<1>(ratio_check_return);
const int worst_offending_connection = std::get<2>(ratio_check_return);
const int worst_offending_completion = std::get<1>(ratio_check_return);
assert((worst_offending_connection >= 0) && (worst_offending_connection < number_of_perforations_));
wellTestState.addClosedCompletion(well_name, worst_offending_completion, simulationTime);
if (writeMessageToOPMLog) {
if (worst_offending_completion < 0) {
const std::string msg = std::string("Connection ") + std::to_string(- worst_offending_completion) + std::string(" for well ")
+ well_name + std::string(" will be closed due to economic limit");
OpmLog::info(msg);
} else {
const std::string msg = std::string("Completion ") + std::to_string(worst_offending_completion) + std::string(" for well ")
+ well_name + std::string(" will be closed due to economic limit");
OpmLog::info(msg);
}
}
const int cell_worst_offending_connection = well_cells_[worst_offending_connection];
list_econ_limited.addClosedConnectionsForWell(well_name, cell_worst_offending_connection);
const std::string msg = std::string("Connection ") + std::to_string(worst_offending_connection) + std::string(" for well ")
+ well_name + std::string(" will be closed due to economic limit");
OpmLog::info(msg);
bool allCompletionsClosed = true;
const auto& connections = well_ecl_->getConnections(current_step_);
for (const auto& connection : connections) {
if (!wellTestState.hasCompletion(name(), connection.complnum)) {
allCompletionsClosed = false;
}
}
if (last_connection) {
// TODO: there is more things to check here
list_econ_limited.addShutWell(well_name);
const std::string msg2 = well_name + std::string(" will be shut due to the last connection closed");
OpmLog::info(msg2);
if (allCompletionsClosed) {
wellTestState.addClosedWell(well_name, WellTestConfig::Reason::ECONOMIC, simulationTime);
if (writeMessageToOPMLog) {
if (well_ecl_->getAutomaticShutIn()) {
const std::string msg = well_name + std::string(" will be shut due to last completion closed");
OpmLog::info(msg);
} else {
const std::string msg = well_name + std::string(" will be stopped due to last completion closed");
OpmLog::info(msg);
}
}
}
break;
}
case WellEcon::WELL:
{
wellTestState.addClosedWell(well_name, WellTestConfig::Reason::ECONOMIC, 0);
if (writeMessageToOPMLog) {
if (well_ecl_->getAutomaticShutIn()) {
list_econ_limited.addShutWell(well_name);
// tell the controll that the well is closed
const std::string msg = well_name + std::string(" will be shut due to ratio economic limit");
OpmLog::info(msg);
} else {
list_econ_limited.addStoppedWell(well_name);
const std::string msg = well_name + std::string(" will be stopped due to ratio economic limit");
OpmLog::info(msg);
}
}
break;
}
case WellEcon::NONE:
@@ -754,13 +799,13 @@ namespace Opm
bore_diameters_.reserve(nperf);
// COMPDAT handling
const auto& completionSet = well_ecl_->getConnections(current_step_);
for (size_t c=0; c<completionSet.size(); c++) {
const auto& completion = completionSet.get(c);
if (completion.state == WellCompletion::OPEN) {
const int i = completion.getI();
const int j = completion.getJ();
const int k = completion.getK();
const auto& connectionSet = well_ecl_->getConnections(current_step_);
for (size_t c=0; c<connectionSet.size(); c++) {
const auto& connection = connectionSet.get(c);
if (connection.state == WellCompletion::OPEN) {
const int i = connection.getI();
const int j = connection.getJ();
const int k = connection.getK();
const int* cpgdim = cart_dims;
const int cart_grid_indx = i + cpgdim[0]*(j + cpgdim[1]*k);
@@ -772,7 +817,7 @@ namespace Opm
const int cell = cgit->second;
{
double radius = 0.5*completion.getDiameter();
double radius = 0.5*connection.getDiameter();
if (radius <= 0.0) {
radius = 0.5*unit::feet;
OPM_MESSAGE("**** Warning: Well bore internal radius set to " << radius);
@@ -784,7 +829,7 @@ namespace Opm
double re; // area equivalent radius of the grid block
double perf_length; // the length of the well perforation
switch (completion.dir) {
switch (connection.dir) {
case Opm::WellCompletion::DirectionEnum::X:
re = std::sqrt(cubical[1] * cubical[2] / M_PI);
perf_length = cubical[0];
@@ -863,7 +908,59 @@ namespace Opm
for (int p = 0; p < np; ++p) {
well_state.wellReservoirRates()[well_rate_index + p] = voidage_rates[p];
}
}
template<typename TypeTag>
void
WellInterface<TypeTag>::closeCompletions(WellTestState& wellTestState)
{
const auto& connections = well_ecl_->getConnections(current_step_);
int perfIdx = 0;
for (const auto& connection : connections) {
if (wellTestState.hasCompletion(name(), connection.complnum)) {
well_index_[perfIdx] = 0.0;
}
perfIdx++;
}
}
template<typename TypeTag>
void
WellInterface<TypeTag>::solveWellForTesting(Simulator& ebosSimulator, WellState& well_state, const std::vector<double>& B_avg, bool terminal_output)
{
const int max_iter = param_.max_welleq_iter_;
int it = 0;
const double dt = 1.0; //not used for the well tests
bool converged;
WellState well_state0 = well_state;
do {
assembleWellEq(ebosSimulator, dt, well_state, true);
ConvergenceReport report;
report = getWellConvergence(B_avg);
converged = report.converged;
if (converged) {
break;
}
++it;
solveEqAndUpdateWellState(well_state);
wellhelpers::WellSwitchingLogger logger;
updateWellControl(well_state, logger);
initPrimaryVariablesEvaluation();
} while (it < max_iter);
if (converged) {
if ( terminal_output ) {
OpmLog::debug("WellTest: Well equation for well " + name() + " solution gets converged with " + std::to_string(it) + " iterations");
}
} else {
if ( terminal_output ) {
OpmLog::debug("WellTest: Well equation for well" +name() + " solution failed in getting converged with " + std::to_string(it) + " iterations");
}
}
}
}

View File

@@ -20,7 +20,7 @@ copyToReferenceDir () {
}
tests=${@:2}
test -z "$tests" && tests="spe11 spe12 spe12p spe1oilgas spe1nowells spe1thermal ctaquifer_2d_oilwater spe3 spe5 spe9 norne_init msw_2d_h msw_3d_hfa polymer2d spe9group polymer_oilwater"
test -z "$tests" && tests="spe11 spe12 spe12p spe1oilgas spe1nowells spe1thermal ctaquifer_2d_oilwater spe3 spe5 spe9 norne_init msw_2d_h msw_3d_hfa polymer2d spe9group polymer_oilwater wecon_wtest"
if grep -q -i "norne " <<< $ghprbCommentBody
then
if test -d $WORKSPACE/deps/opm-tests/norne/flow
@@ -218,6 +218,15 @@ for test_name in ${tests}; do
NORNE_ATW2013 \
UNSMRY
fi
if grep -q "wecon_wtest" <<< $test_name
then
copyToReferenceDir \
$configuration/build-opm-simulators/tests/results/flow+wecon_wtest/ \
$OPM_TESTS_ROOT/wecon_wtest/opm-simulation-reference/flow \
3D_WECON \
EGRID INIT SMSPEC UNRST UNSMRY
fi
done
if [ -z "${@:2}" ]