OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2024-12-24 08:20:01 -06:00
Code Issues Packages Projects Releases Wiki Activity
df74d7e696
opm-simulators/opm/core/wells/WellsManager.cpp
Andreas Lauser 1ae94c8db3 do not explicitly pass the permeability to the well model anymore 
this information is already part of the EclipseState. The reason why
this should IMO be avoided is that this enforces an implementation
(ordering of the permeability matrices) the simulator on the well
model. If this needs to be done for performance reasons, IMO it would
be smarter to pass an array of matrices, instead of passing a raw
array of doubles.  I doubt that this is necessary, though: completing
the full Norne deck takes about 0.25 seconds longer on my machine,
that's substantially less than 0.1% of the total runtime.

in order to avoid code duplication, the permeability extraction
function of the RockFromDeck class is now made a public static
function and used as an implementation detail of the WellsManager.

finally, the permfield_valid_ attribute is removed from the
RockFromDeck class because this data was unused and not accessible via
the class' public API.
2017-01-27 12:51:12 +01:00

894 lines
42 KiB
C++
Raw Blame History

/*
Copyright 2012 SINTEF ICT, Applied Mathematics.
Copyright 2016 IRIS AS
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/>.
*/
#include "config.h"
#include <opm/core/wells/WellsManager.hpp>
#include <opm/core/grid.h>
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/common/ErrorMacros.hpp>
#include <opm/core/wells/WellCollection.hpp>
#include <opm/core/wells/WellsGroup.hpp>
#include <opm/core/props/phaseUsageFromDeck.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/ScheduleEnums.hpp>
#include <algorithm>
#include <cassert>
#include <cmath>
#include <cstddef>
#include <map>
#include <string>
#include <utility>
#include <iostream>
namespace
{
static double invalid_alq = -1e100;
static double invalid_vfp = -2147483647;
} //Namespace
// Helper structs and functions for the implementation.
namespace WellsManagerDetail
{
namespace ProductionControl
{
namespace Details {
std::map<std::string, Mode>
init_mode_map() {
std::map<std::string, Mode> m;
m.insert(std::make_pair("ORAT", ORAT));
m.insert(std::make_pair("WRAT", WRAT));
m.insert(std::make_pair("GRAT", GRAT));
m.insert(std::make_pair("LRAT", LRAT));
m.insert(std::make_pair("CRAT", CRAT));
m.insert(std::make_pair("RESV", RESV));
m.insert(std::make_pair("BHP" , BHP ));
m.insert(std::make_pair("THP" , THP ));
m.insert(std::make_pair("GRUP", GRUP));
return m;
}
} // namespace Details
Mode mode(const std::string& control)
{
static std::map<std::string, Mode>
mode_map = Details::init_mode_map();
std::map<std::string, Mode>::iterator
p = mode_map.find(control);
if (p != mode_map.end()) {
return p->second;
}
else {
OPM_THROW(std::runtime_error, "Unknown well control mode = "
<< control << " in input file");
}
}
Mode mode(Opm::WellProducer::ControlModeEnum controlMode)
{
switch( controlMode ) {
case Opm::WellProducer::ORAT:
return ORAT;
case Opm::WellProducer::WRAT:
return WRAT;
case Opm::WellProducer::GRAT:
return GRAT;
case Opm::WellProducer::LRAT:
return LRAT;
case Opm::WellProducer::CRAT:
return CRAT;
case Opm::WellProducer::RESV:
return RESV;
case Opm::WellProducer::BHP:
return BHP;
case Opm::WellProducer::THP:
return THP;
case Opm::WellProducer::GRUP:
return GRUP;
default:
throw std::invalid_argument("unhandled enum value");
}
}
} // namespace ProductionControl
namespace InjectionControl
{
namespace Details {
std::map<std::string, Mode>
init_mode_map() {
std::map<std::string, Mode> m;
m.insert(std::make_pair("RATE", RATE));
m.insert(std::make_pair("RESV", RESV));
m.insert(std::make_pair("BHP" , BHP ));
m.insert(std::make_pair("THP" , THP ));
m.insert(std::make_pair("GRUP", GRUP));
return m;
}
} // namespace Details
Mode mode(const std::string& control)
{
static std::map<std::string, Mode>
mode_map = Details::init_mode_map();
std::map<std::string, Mode>::iterator
p = mode_map.find(control);
if (p != mode_map.end()) {
return p->second;
}
else {
OPM_THROW(std::runtime_error, "Unknown well control mode = "
<< control << " in input file");
}
}
Mode mode(Opm::WellInjector::ControlModeEnum controlMode)
{
switch ( controlMode ) {
case Opm::WellInjector::GRUP:
return GRUP;
case Opm::WellInjector::RESV:
return RESV;
case Opm::WellInjector::RATE:
return RATE;
case Opm::WellInjector::THP:
return THP;
case Opm::WellInjector::BHP:
return BHP;
default:
throw std::invalid_argument("unhandled enum value");
}
}
} // namespace InjectionControl
// Compute direction permutation corresponding to completion's
// direction. First two elements of return value are directions
// perpendicular to completion while last element is direction
// along completion.
inline std::array< std::array<double,3>::size_type, 3 >
directionIndices(const Opm::WellCompletion::DirectionEnum direction)
{
typedef std::array<double,3>::size_type idx_t;
typedef std::array<idx_t,3> permutation;
switch (direction) {
case Opm::WellCompletion::DirectionEnum::X:
return permutation {{ idx_t(1), idx_t(2), idx_t(0) }};
case Opm::WellCompletion::DirectionEnum::Y:
return permutation {{ idx_t(2), idx_t(0), idx_t(1) }};
case Opm::WellCompletion::DirectionEnum::Z:
return permutation {{ idx_t(0), idx_t(1), idx_t(2) }};
}
// All enum values should be handled above. Therefore
// we should never reach this one. Anyway for the sake
// of reduced warnings we throw an exception.
throw std::invalid_argument("unhandled enum value");
}
// Permute (diagonal) permeability components according to
// completion's direction.
inline std::array<double,3>
permComponents(const Opm::WellCompletion::DirectionEnum direction,
const double* perm)
{
const auto p = directionIndices(direction);
const std::array<double,3>::size_type d = 3;
std::array<double,3>
K = {{ perm[ p[0]*(d + 1) ] ,
perm[ p[1]*(d + 1) ] ,
perm[ p[2]*(d + 1) ] }};
return K;
}
// Permute cell's geometric extent according to completion's
// direction. Honour net-to-gross ratio.
//
// Note: 'extent' is intentionally accepted by modifiable value
// rather than reference-to-const to support NTG manipulation.
inline std::array<double,3>
effectiveExtent(const Opm::WellCompletion::DirectionEnum direction,
const double ntg,
std::array<double,3> extent)
{
// Vertical extent affected by net-to-gross ratio.
extent[2] *= ntg;
const auto p = directionIndices(direction);
std::array<double,3>
D = {{ extent[ p[0] ] ,
extent[ p[1] ] ,
extent[ p[2] ] }};
return D;
}
// Compute Peaceman's effective radius of single completion.
inline double
effectiveRadius(const std::array<double,3>& K,
const std::array<double,3>& D)
{
const double K01 = K[0] / K[1];
const double K10 = K[1] / K[0];
const double D0_sq = D[0] * D[0];
const double D1_sq = D[1] * D[1];
const double num = std::sqrt((std::sqrt(K10) * D0_sq) +
(std::sqrt(K01) * D1_sq));
const double den = std::pow(K01, 0.25) + std::pow(K10, 0.25);
// Note: Analytic constant 0.28 derived for infintely sized
// formation with repeating well placement.
return 0.28 * (num / den);
}
// Use the Peaceman well model to compute well indices.
// radius is the radius of the well.
// cubical contains [dx, dy, dz] of the cell.
// (Note that the well model asumes that each cell is a cuboid).
// cell_permeability is the permeability tensor of the given cell.
// returns the well index of the cell.
double
computeWellIndex(const double radius,
const std::array<double, 3>& cubical,
const double* cell_permeability,
const double skin_factor,
const Opm::WellCompletion::DirectionEnum direction,
const double ntg)
{
const std::array<double,3>& K =
permComponents(direction, cell_permeability);
const std::array<double,3>& D =
effectiveExtent(direction, ntg, cubical);
const double r0 = effectiveRadius(K, D);
const double Kh = std::sqrt(K[0] * K[1]) * D[2];
// Angle of completion exposed to flow. We assume centre
// placement so there's complete exposure (= 2\pi).
const double angle =
6.2831853071795864769252867665590057683943387987502116419498;
double rw = radius;
if (r0 < rw) {
std::cerr << "Completion radius exceeds effective radius\n"
<< "Reset to effective\n";
rw = r0;
}
// NOTE: The formula is originally derived and valid for
// Cartesian grids only.
return (angle * Kh) / (std::log(r0 / rw) + skin_factor);
}
} // anonymous namespace
namespace Opm
{
/// Default constructor.
WellsManager::WellsManager()
: w_(0), is_parallel_run_(false)
{
}
/// Construct from existing wells object.
WellsManager::WellsManager(struct Wells* W)
: w_(clone_wells(W)), is_parallel_run_(false)
{
}
/// Construct wells from deck.
WellsManager::WellsManager(const Opm::EclipseState& eclipseState,
const size_t timeStep,
const UnstructuredGrid& grid)
: w_(0), is_parallel_run_(false)
{
std::vector<double> dummy_well_potentials;
// TODO: not sure about the usage of this WellsManager constructor
// TODO: not sure whether this is the correct thing to do here.
DynamicListEconLimited dummy_list_econ_limited;
init(eclipseState, timeStep, UgGridHelpers::numCells(grid),
UgGridHelpers::globalCell(grid), UgGridHelpers::cartDims(grid),
UgGridHelpers::dimensions(grid),
UgGridHelpers::cell2Faces(grid), UgGridHelpers::beginFaceCentroids(grid),
dummy_list_econ_limited, dummy_well_potentials,
std::unordered_set<std::string>());
}
/// Destructor.
WellsManager::~WellsManager()
{
destroy_wells(w_);
}
/// Does the "deck" define any wells?
bool WellsManager::empty() const
{
return (w_ == 0) || (w_->number_of_wells == 0);
}
/// Access the managed Wells.
/// The method is named similarly to c_str() in std::string,
/// to make it clear that we are returning a C-compatible struct.
const Wells* WellsManager::c_wells() const
{
return w_;
}
const WellCollection& WellsManager::wellCollection() const
{
return well_collection_;
}
WellCollection& WellsManager::wellCollection() {
return well_collection_;
}
bool WellsManager::conditionsMet(const std::vector<double>& well_bhp,
const std::vector<double>& well_reservoirrates_phase,
const std::vector<double>& well_surfacerates_phase)
{
return well_collection_.conditionsMet(well_bhp,
well_reservoirrates_phase,
well_surfacerates_phase);
}
/// Applies explicit reinjection controls. This must be called at each timestep to be correct.
/// \param[in] well_reservoirrates_phase
/// A vector containing reservoir rates by phase for each well.
/// Is assumed to be ordered the same way as the related Wells-struct,
/// with all phase rates of a single well adjacent in the array.
/// \param[in] well_surfacerates_phase
/// A vector containing surface rates by phase for each well.
/// Is assumed to be ordered the same way as the related Wells-struct,
/// with all phase rates of a single well adjacent in the array.
void WellsManager::applyExplicitReinjectionControls(const std::vector<double>& well_reservoirrates_phase,
const std::vector<double>& well_surfacerates_phase)
{
well_collection_.applyExplicitReinjectionControls(well_reservoirrates_phase, well_surfacerates_phase);
}
void WellsManager::setupCompressedToCartesian(const int* global_cell, int number_of_cells,
std::map<int,int>& cartesian_to_compressed ) {
// global_cell is a map from compressed cells to Cartesian grid cells.
// We must make the inverse lookup.
if (global_cell) {
for (int i = 0; i < number_of_cells; ++i) {
cartesian_to_compressed.insert(std::make_pair(global_cell[i], i));
}
}
else {
for (int i = 0; i < number_of_cells; ++i) {
cartesian_to_compressed.insert(std::make_pair(i, i));
}
}
}
void WellsManager::setupWellControls(std::vector< const Well* >& wells, size_t timeStep,
std::vector<std::string>& well_names, const PhaseUsage& phaseUsage,
const std::vector<int>& wells_on_proc,
const DynamicListEconLimited& list_econ_limited) {
int well_index = 0;
auto well_on_proc = wells_on_proc.begin();
for (auto wellIter= wells.begin(); wellIter != wells.end(); ++wellIter, ++well_on_proc) {
if( ! *well_on_proc )
{
// Wells not stored on the process are not in the list
continue;
}
const auto* well = (*wellIter);
if (well->getStatus(timeStep) == WellCommon::SHUT) {
//SHUT wells are not added to the well list
continue;
}
if (list_econ_limited.wellShutEconLimited(well->name())) {
continue;
}
if (well->getStatus(timeStep) == WellCommon::STOP || list_econ_limited.wellStoppedEconLimited(well->name())) {
// Stopped wells are kept in the well list but marked as stopped.
well_controls_stop_well(w_->ctrls[well_index]);
}
if (well->isInjector(timeStep)) {
const WellInjectionProperties& injectionProperties = well->getInjectionProperties(timeStep);
int ok = 1;
int control_pos[5] = { -1, -1, -1, -1, -1 };
clear_well_controls(well_index, w_);
if (injectionProperties.hasInjectionControl(WellInjector::RATE)) {
control_pos[WellsManagerDetail::InjectionControl::RATE] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 0.0, 0.0, 0.0 };
WellInjector::TypeEnum injectorType = injectionProperties.injectorType;
if (injectorType == WellInjector::TypeEnum::WATER) {
distr[phaseUsage.phase_pos[BlackoilPhases::Aqua]] = 1.0;
} else if (injectorType == WellInjector::TypeEnum::OIL) {
distr[phaseUsage.phase_pos[BlackoilPhases::Liquid]] = 1.0;
} else if (injectorType == WellInjector::TypeEnum::GAS) {
distr[phaseUsage.phase_pos[BlackoilPhases::Vapour]] = 1.0;
}
ok = append_well_controls(SURFACE_RATE,
injectionProperties.surfaceInjectionRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && injectionProperties.hasInjectionControl(WellInjector::RESV)) {
control_pos[WellsManagerDetail::InjectionControl::RESV] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 0.0, 0.0, 0.0 };
WellInjector::TypeEnum injectorType = injectionProperties.injectorType;
if (injectorType == WellInjector::TypeEnum::WATER) {
distr[phaseUsage.phase_pos[BlackoilPhases::Aqua]] = 1.0;
} else if (injectorType == WellInjector::TypeEnum::OIL) {
distr[phaseUsage.phase_pos[BlackoilPhases::Liquid]] = 1.0;
} else if (injectorType == WellInjector::TypeEnum::GAS) {
distr[phaseUsage.phase_pos[BlackoilPhases::Vapour]] = 1.0;
}
ok = append_well_controls(RESERVOIR_RATE,
injectionProperties.reservoirInjectionRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && injectionProperties.hasInjectionControl(WellInjector::BHP)) {
control_pos[WellsManagerDetail::InjectionControl::BHP] = well_controls_get_num(w_->ctrls[well_index]);
ok = append_well_controls(BHP,
injectionProperties.BHPLimit,
invalid_alq,
invalid_vfp,
NULL,
well_index,
w_);
}
if (ok && injectionProperties.hasInjectionControl(WellInjector::THP)) {
control_pos[WellsManagerDetail::InjectionControl::THP] = well_controls_get_num(w_->ctrls[well_index]);
const double thp_limit = injectionProperties.THPLimit;
const int vfp_number = injectionProperties.VFPTableNumber;
ok = append_well_controls(THP,
thp_limit,
invalid_alq,
vfp_number,
NULL,
well_index,
w_);
}
if (!ok) {
OPM_THROW(std::runtime_error, "Failure occured appending controls for well " << well_names[well_index]);
}
if (injectionProperties.controlMode != WellInjector::CMODE_UNDEFINED) {
WellsManagerDetail::InjectionControl::Mode mode = WellsManagerDetail::InjectionControl::mode( injectionProperties.controlMode );
int cpos = control_pos[mode];
if (cpos == -1 && mode != WellsManagerDetail::InjectionControl::GRUP) {
OPM_THROW(std::runtime_error, "Control not specified in well " << well_names[well_index]);
}
set_current_control(well_index, cpos, w_);
}
// Set well component fraction.
double cf[3] = { 0.0, 0.0, 0.0 };
{
WellInjector::TypeEnum injectorType = injectionProperties.injectorType;
if (injectorType == WellInjector::WATER) {
if (!phaseUsage.phase_used[BlackoilPhases::Aqua]) {
OPM_THROW(std::runtime_error, "Water phase not used, yet found water-injecting well.");
}
cf[phaseUsage.phase_pos[BlackoilPhases::Aqua]] = 1.0;
} else if (injectorType == WellInjector::OIL) {
if (!phaseUsage.phase_used[BlackoilPhases::Liquid]) {
OPM_THROW(std::runtime_error, "Oil phase not used, yet found oil-injecting well.");
}
cf[phaseUsage.phase_pos[BlackoilPhases::Liquid]] = 1.0;
} else if (injectorType == WellInjector::GAS) {
if (!phaseUsage.phase_used[BlackoilPhases::Vapour]) {
OPM_THROW(std::runtime_error, "Gas phase not used, yet found gas-injecting well.");
}
cf[phaseUsage.phase_pos[BlackoilPhases::Vapour]] = 1.0;
}
std::copy(cf, cf + phaseUsage.num_phases, w_->comp_frac + well_index*phaseUsage.num_phases);
}
}
if (well->isProducer(timeStep)) {
// Add all controls that are present in well.
// First we must clear existing controls, in case the
// current WCONPROD line is modifying earlier controls.
const WellProductionProperties& productionProperties = well->getProductionProperties(timeStep);
int control_pos[9] = { -1, -1, -1, -1, -1, -1, -1, -1, -1 };
int ok = 1;
clear_well_controls(well_index, w_);
if (ok && productionProperties.hasProductionControl(WellProducer::ORAT)) {
if (!phaseUsage.phase_used[BlackoilPhases::Liquid]) {
OPM_THROW(std::runtime_error, "Oil phase not active and ORAT control specified.");
}
control_pos[WellsManagerDetail::ProductionControl::ORAT] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 0.0, 0.0, 0.0 };
distr[phaseUsage.phase_pos[BlackoilPhases::Liquid]] = 1.0;
ok = append_well_controls(SURFACE_RATE,
-productionProperties.OilRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && productionProperties.hasProductionControl(WellProducer::WRAT)) {
if (!phaseUsage.phase_used[BlackoilPhases::Aqua]) {
OPM_THROW(std::runtime_error, "Water phase not active and WRAT control specified.");
}
control_pos[WellsManagerDetail::ProductionControl::WRAT] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 0.0, 0.0, 0.0 };
distr[phaseUsage.phase_pos[BlackoilPhases::Aqua]] = 1.0;
ok = append_well_controls(SURFACE_RATE,
-productionProperties.WaterRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && productionProperties.hasProductionControl(WellProducer::GRAT)) {
if (!phaseUsage.phase_used[BlackoilPhases::Vapour]) {
OPM_THROW(std::runtime_error, "Gas phase not active and GRAT control specified.");
}
control_pos[WellsManagerDetail::ProductionControl::GRAT] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 0.0, 0.0, 0.0 };
distr[phaseUsage.phase_pos[BlackoilPhases::Vapour]] = 1.0;
ok = append_well_controls(SURFACE_RATE,
-productionProperties.GasRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && productionProperties.hasProductionControl(WellProducer::LRAT)) {
if (!phaseUsage.phase_used[BlackoilPhases::Aqua]) {
OPM_THROW(std::runtime_error, "Water phase not active and LRAT control specified.");
}
if (!phaseUsage.phase_used[BlackoilPhases::Liquid]) {
OPM_THROW(std::runtime_error, "Oil phase not active and LRAT control specified.");
}
control_pos[WellsManagerDetail::ProductionControl::LRAT] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 0.0, 0.0, 0.0 };
distr[phaseUsage.phase_pos[BlackoilPhases::Aqua]] = 1.0;
distr[phaseUsage.phase_pos[BlackoilPhases::Liquid]] = 1.0;
ok = append_well_controls(SURFACE_RATE,
-productionProperties.LiquidRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && productionProperties.hasProductionControl(WellProducer::RESV)) {
control_pos[WellsManagerDetail::ProductionControl::RESV] = well_controls_get_num(w_->ctrls[well_index]);
double distr[3] = { 1.0, 1.0, 1.0 };
ok = append_well_controls(RESERVOIR_RATE,
-productionProperties.ResVRate,
invalid_alq,
invalid_vfp,
distr,
well_index,
w_);
}
if (ok && productionProperties.hasProductionControl(WellProducer::THP)) {
const double thp_limit = productionProperties.THPLimit;
const double alq_value = productionProperties.ALQValue;
const int vfp_number = productionProperties.VFPTableNumber;
control_pos[WellsManagerDetail::ProductionControl::THP] = well_controls_get_num(w_->ctrls[well_index]);
ok = append_well_controls(THP,
thp_limit,
alq_value,
vfp_number,
NULL,
well_index,
w_);
}
if (ok) {
// Always append a BHP control.
// If no explicit BHP control given, use a 1 atm control.
const bool has_explicit_limit = productionProperties.hasProductionControl(WellProducer::BHP);
const double bhp_limit = has_explicit_limit ? productionProperties.BHPLimit : unit::convert::from(1.0, unit::atm);
control_pos[WellsManagerDetail::ProductionControl::BHP] = well_controls_get_num(w_->ctrls[well_index]);
ok = append_well_controls(BHP,
bhp_limit,
invalid_alq,
invalid_vfp,
NULL,
well_index,
w_);
}
if (!ok) {
OPM_THROW(std::runtime_error, "Failure occured appending controls for well " << well_names[well_index]);
}
if (productionProperties.controlMode != WellProducer::CMODE_UNDEFINED) {
WellsManagerDetail::ProductionControl::Mode mode = WellsManagerDetail::ProductionControl::mode(productionProperties.controlMode);
int cpos = control_pos[mode];
if (cpos == -1 && mode != WellsManagerDetail::ProductionControl::GRUP) {
OPM_THROW(std::runtime_error, "Control mode type " << mode << " not present in well " << well_names[well_index]);
}
else {
set_current_control(well_index, cpos, w_);
}
}
// Set well component fraction to match preferred phase for the well.
double cf[3] = { 0.0, 0.0, 0.0 };
{
switch (well->getPreferredPhase()) {
case Phase::WATER:
if (!phaseUsage.phase_used[BlackoilPhases::Aqua]) {
OPM_THROW(std::runtime_error, "Water phase not used, yet found water-preferring well.");
}
cf[phaseUsage.phase_pos[BlackoilPhases::Aqua]] = 1.0;
break;
case Phase::OIL:
if (!phaseUsage.phase_used[BlackoilPhases::Liquid]) {
OPM_THROW(std::runtime_error, "Oil phase not used, yet found oil-preferring well.");
}
cf[phaseUsage.phase_pos[BlackoilPhases::Liquid]] = 1.0;
break;
case Phase::GAS:
if (!phaseUsage.phase_used[BlackoilPhases::Vapour]) {
OPM_THROW(std::runtime_error, "Gas phase not used, yet found gas-preferring well.");
}
cf[phaseUsage.phase_pos[BlackoilPhases::Vapour]] = 1.0;
break;
default:
OPM_THROW(std::logic_error, "Unknown preferred phase: " << well->getPreferredPhase());
}
std::copy(cf, cf + phaseUsage.num_phases, w_->comp_frac + well_index*phaseUsage.num_phases);
}
}
well_index++;
}
}
void WellsManager::setupGuideRates(std::vector< const Well* >& wells, const size_t timeStep, std::vector<WellData>& well_data, std::map<std::string, int>& well_names_to_index,
const PhaseUsage& phaseUsage, const std::vector<double>& well_potentials)
{
const int np = phaseUsage.num_phases;
for (auto wellIter = wells.begin(); wellIter != wells.end(); ++wellIter ) {
const auto* well = *wellIter;
if (well->getStatus(timeStep) == WellCommon::SHUT) {
//SHUT wells does not need guide rates
continue;
}
const int wix = well_names_to_index[well->name()];
WellNode& wellnode = *well_collection_.getLeafNodes()[wix];
if (well->getGuideRatePhase(timeStep) != GuideRate::UNDEFINED) {
if (well_data[wix].type == PRODUCER) {
wellnode.prodSpec().guide_rate_ = well->getGuideRate(timeStep);
if (well->getGuideRatePhase(timeStep) == GuideRate::OIL) {
wellnode.prodSpec().guide_rate_type_ = ProductionSpecification::OIL;
} else {
OPM_THROW(std::runtime_error, "Guide rate type " << GuideRate::GuideRatePhaseEnum2String(well->getGuideRatePhase(timeStep)) << " specified for producer "
<< well->name() << " in WGRUPCON, cannot handle.");
}
} else if (well_data[wix].type == INJECTOR) {
wellnode.injSpec().guide_rate_ = well->getGuideRate(timeStep);
if (well->getGuideRatePhase(timeStep) == GuideRate::RAT) {
wellnode.injSpec().guide_rate_type_ = InjectionSpecification::RAT;
} else {
OPM_THROW(std::runtime_error, "Guide rate type " << GuideRate::GuideRatePhaseEnum2String(well->getGuideRatePhase(timeStep)) << " specified for injector "
<< well->name() << " in WGRUPCON, cannot handle.");
}
} else {
OPM_THROW(std::runtime_error, "Unknown well type " << well_data[wix].type << " for well " << well->name());
}
} else if (well_potentials.size() > 0) { // default: calculate guide rates from well potentials
// Note: Modification of the guide rate using GUIDERAT is not supported
// only set the guide rates if there is a parent group with valied control
if (wellnode.getParent() != nullptr) {
const WellsGroupInterface& group = *wellnode.getParent();
if ( well->isProducer(timeStep) ) {
// The guide rates is calculated based on the group control
// Currently only supporting WRAT, ORAT and GRAT.
ProductionSpecification::ControlMode control_mode = group.prodSpec().control_mode_;
if (control_mode == ProductionSpecification::FLD) {
if (group.getParent() != nullptr) {
const WellsGroupInterface& higher_group = *(group.getParent());
control_mode = higher_group.prodSpec().control_mode_;
} else {
OPM_THROW(std::runtime_error, "Group " << group.name() << " is under FLD control while no higher level of group is specified.");
}
}
switch (control_mode) {
case ProductionSpecification::WRAT: {
if (!phaseUsage.phase_used[BlackoilPhases::Aqua]) {
OPM_THROW(std::runtime_error, "Water phase not used, yet found water rate controlled well.");
}
const int water_index = phaseUsage.phase_pos[BlackoilPhases::Aqua];
wellnode.prodSpec().guide_rate_ = well_potentials[np*wix + water_index];
wellnode.prodSpec().guide_rate_type_ = ProductionSpecification::WATER;
break;
}
case ProductionSpecification::ORAT: {
if (!phaseUsage.phase_used[BlackoilPhases::Liquid]) {
OPM_THROW(std::runtime_error, "Oil phase not used, yet found oil rate controlled well.");
}
const int oil_index = phaseUsage.phase_pos[BlackoilPhases::Liquid];
wellnode.prodSpec().guide_rate_ = well_potentials[np*wix + oil_index];
wellnode.prodSpec().guide_rate_type_ = ProductionSpecification::OIL;
break;
}
case ProductionSpecification::GRAT: {
if (!phaseUsage.phase_used[BlackoilPhases::Vapour]) {
OPM_THROW(std::runtime_error, "Gas phase not used, yet found gas rate controlled well.");
}
const int gas_index = phaseUsage.phase_pos[BlackoilPhases::Vapour];
wellnode.prodSpec().guide_rate_ = well_potentials[np*wix + gas_index];
wellnode.prodSpec().guide_rate_type_ = ProductionSpecification::GAS;
break;
}
case ProductionSpecification::FLD: {
OPM_THROW(std::logic_error, "Not support more than one continous level of FLD control");
}
case ProductionSpecification::LRAT: {
double guide_rate = 0;
if (phaseUsage.phase_used[BlackoilPhases::Liquid]) {
const int oil_index = phaseUsage.phase_pos[BlackoilPhases::Liquid];
const double potential_oil = well_potentials[np*wix + oil_index];
guide_rate += potential_oil;
}
if (phaseUsage.phase_used[BlackoilPhases::Aqua]) {
const int water_index = phaseUsage.phase_pos[BlackoilPhases::Aqua];
const double potential_water = well_potentials[np*wix + water_index];
guide_rate += potential_water;
}
// not sure if no water and no oil, what will happen here, zero guide_rate?
wellnode.prodSpec().guide_rate_ = guide_rate;
wellnode.prodSpec().guide_rate_type_ = ProductionSpecification::LIQ;
}
case ProductionSpecification::NONE: {
// Group control is not in use for this group.
break;
}
default:
OPM_THROW(std::logic_error, "Not supported control_mode for guide rate computed" <<
" from well potentials: " << group.prodSpec().control_mode_);
}
} else {
// The guide rates is calculated based on the group injector type
switch (group.injSpec().injector_type_) {
case InjectionSpecification::WATER: {
if (!phaseUsage.phase_used[BlackoilPhases::Aqua]) {
OPM_THROW(std::runtime_error, "Water phase not used, yet found water injecting well.");
}
const int water_index = phaseUsage.phase_pos[BlackoilPhases::Aqua];
wellnode.injSpec().guide_rate_ = well_potentials[np*wix + water_index];
// Guide rates applies to the phase that the well is injecting i.e water
wellnode.injSpec().guide_rate_type_ = InjectionSpecification::RAT;
break;
}
case InjectionSpecification::OIL: {
if (!phaseUsage.phase_used[BlackoilPhases::Liquid]) {
OPM_THROW(std::runtime_error, "Oil phase not used, yet found oil injecting well.");
}
const int oil_index = phaseUsage.phase_pos[BlackoilPhases::Liquid];
wellnode.injSpec().guide_rate_ = well_potentials[np*wix + oil_index];
// Guide rates applies to the phase that the well is injecting i.e. oil
wellnode.injSpec().guide_rate_type_ = InjectionSpecification::RAT;
break;
}
case InjectionSpecification::GAS: {
if (!phaseUsage.phase_used[BlackoilPhases::Vapour]) {
OPM_THROW(std::runtime_error, "Gas phase not used, yet found gas injecting well.");
}
const int gas_index = phaseUsage.phase_pos[BlackoilPhases::Vapour];
wellnode.injSpec().guide_rate_ = well_potentials[np*wix + gas_index];
// Guide rates applies to the phase that the well is injecting i.e gas
wellnode.injSpec().guide_rate_type_ = InjectionSpecification::RAT;
break;
}
default:
OPM_THROW(std::logic_error, "Not supported injector type for guide rate computed" <<
" from well potentials: " << group.injSpec().injector_type_);
}
}
}
} // if neither WGRUPCON nor well_potentials is given, distribute the flow equaly
}
}
} // namespace Opm
Reference in New Issue View Git Blame Copy Permalink