20 #ifndef OPM_WELLSTATE_HEADER_INCLUDED
21 #define OPM_WELLSTATE_HEADER_INCLUDED
23 #include <opm/core/props/BlackoilPhases.hpp>
26 #include <opm/output/data/Wells.hpp>
42 typedef std::array< int, 3 > mapentry_t;
43 typedef std::map< std::string, mapentry_t > WellMapType;
45 template <
class State>
46 void init(
const Wells* wells,
const State& state)
48 init(wells, state.pressure());
56 void init(
const Wells* wells,
const std::vector<double>& cellPressures)
67 temperature_.resize(nw, 273.15 + 20);
68 wellrates_.resize(nw * np, 0.0);
69 for (
int w = 0; w < nw; ++w) {
76 assert( wells->
name[ w ] );
77 std::string name( wells->
name[ w ] );
78 assert( name.size() > 0 );
79 mapentry_t& wellMapEntry = wellMap_[name];
80 wellMapEntry[ 0 ] = w;
83 wellMapEntry[ 2 ] = num_perf_this_well;
86 if ( num_perf_this_well == 0 )
89 for (
int p = 0; p < np; ++p) {
90 wellrates_[np*w + p] = 0.0;
97 if (well_controls_well_is_stopped(ctrl)) {
100 for (
int p = 0; p < np; ++p) {
101 wellrates_[np*w + p] = 0.0;
106 if (well_controls_get_current_type(ctrl) ==
BHP) {
107 bhp_[w] = well_controls_get_current_target( ctrl );
110 bhp_[w] = cellPressures[first_cell];
120 if (well_controls_get_current_type(ctrl) ==
SURFACE_RATE) {
121 const double rate_target = well_controls_get_current_target(ctrl);
122 const double * distr = well_controls_get_current_distr( ctrl );
123 for (
int p = 0; p < np; ++p) {
124 wellrates_[np*w + p] = rate_target * distr[p];
127 const double small_rate = 1e-14;
128 const double sign = (wells->
type[w] ==
INJECTOR) ? 1.0 : -1.0;
129 for (
int p = 0; p < np; ++p) {
130 wellrates_[np*w + p] = small_rate * sign;
139 if (well_controls_get_current_type(ctrl) ==
BHP) {
140 bhp_[w] = well_controls_get_current_target( ctrl );
143 const double safety_factor = (wells->
type[w] ==
INJECTOR) ? 1.01 : 0.99;
144 bhp_[w] = safety_factor*cellPressures[first_cell];
151 const int nwc = well_controls_get_num(ctrl);
152 for (
int ctrl_index = 0; ctrl_index < nwc; ++ctrl_index) {
153 if (well_controls_iget_type(ctrl, ctrl_index) ==
THP) {
154 thp_[w] = well_controls_iget_target(ctrl, ctrl_index);
169 std::vector<double>&
bhp() {
return bhp_; }
170 const std::vector<double>&
bhp()
const {
return bhp_; }
173 std::vector<double>&
thp() {
return thp_; }
174 const std::vector<double>&
thp()
const {
return thp_; }
178 const std::vector<double>&
temperature()
const {
return temperature_; }
182 const std::vector<double>&
wellRates()
const {
return wellrates_; }
186 const std::vector<double>&
perfRates()
const {
return perfrates_; }
190 const std::vector<double>&
perfPress()
const {
return perfpress_; }
192 size_t getRestartBhpOffset()
const {
196 size_t getRestartPerfPressOffset()
const {
200 size_t getRestartPerfRatesOffset()
const {
201 return getRestartPerfPressOffset() + perfpress_.size();
204 size_t getRestartTemperatureOffset()
const {
205 return getRestartPerfRatesOffset() + perfrates_.size();
208 size_t getRestartWellRatesOffset()
const {
209 return getRestartTemperatureOffset() + temperature_.size();
212 const WellMapType& wellMap()
const {
return wellMap_; }
213 WellMapType& wellMap() {
return wellMap_; }
227 virtual data::Wells report(
const PhaseUsage& pu)
const
229 using rt = data::Rates::opt;
232 for(
const auto& itr : this->wellMap_ ) {
233 const auto well_index = itr.second[ 0 ];
235 auto& well = dw[ itr.first ];
236 well.bhp = this->
bhp().at( well_index );
237 well.thp = this->
thp().at( well_index );
238 well.temperature = this->
temperature().at( well_index );
240 const auto wellrate_index = well_index * pu.num_phases;
242 if( pu.phase_used[BlackoilPhases::Aqua] ) {
243 well.rates.set( rt::wat, wv[ wellrate_index + pu.phase_pos[BlackoilPhases::Aqua] ] );
246 if( pu.phase_used[BlackoilPhases::Liquid] ) {
247 well.rates.set( rt::oil, wv[ wellrate_index + pu.phase_pos[BlackoilPhases::Liquid] ] );
250 if( pu.phase_used[BlackoilPhases::Vapour] ) {
251 well.rates.set( rt::gas, wv[ wellrate_index + pu.phase_pos[BlackoilPhases::Vapour] ] );
254 const int num_perf_well = this->wells_->well_connpos[ well_index + 1 ]
255 - this->wells_->well_connpos[ well_index ];
256 well.completions.resize(num_perf_well);
258 for(
int i = 0; i < num_perf_well; ++i ) {
259 const auto wi = this->wells_->well_connpos[ well_index ] + i;
260 const auto active_index = this->wells_->well_cells[ wi ];
262 auto& completion = well.completions[ i ];
263 completion.index = active_index;
264 completion.pressure = this->
perfPress()[ itr.second[1] + i ];
265 completion.reservoir_rate = this->
perfRates()[ itr.second[1] + i ];
267 assert(num_perf_well ==
int(well.completions.size()));
274 virtual ~WellState() {}
276 WellState() =
default;
277 WellState(
const WellState& rhs ) :
280 temperature_( rhs.temperature_ ),
281 wellrates_( rhs.wellrates_ ),
282 perfrates_( rhs.perfrates_ ),
283 perfpress_( rhs.perfpress_ ),
284 wellMap_( rhs.wellMap_ ),
288 WellState& operator=(
const WellState& rhs ) {
289 this->bhp_ = rhs.bhp_;
290 this->thp_ = rhs.thp_;
291 this->temperature_ = rhs.temperature_;
292 this->wellrates_ = rhs.wellrates_;
293 this->perfrates_ = rhs.perfrates_;
294 this->perfpress_ = rhs.perfpress_;
295 this->wellMap_ = rhs.wellMap_;
296 this->wells_.reset(
clone_wells( rhs.wells_.get() ) );
302 std::vector<double> bhp_;
303 std::vector<double> thp_;
304 std::vector<double> temperature_;
305 std::vector<double> wellrates_;
306 std::vector<double> perfrates_;
307 std::vector<double> perfpress_;
309 WellMapType wellMap_;
315 std::unique_ptr< Wells, wdel > wells_;
320 #endif // OPM_WELLSTATE_HEADER_INCLUDED
Well constrained by BHP target.
Definition: well_controls.h:35
int numWells() const
The number of wells present.
Definition: WellState.hpp:216
API for managing sets of well controls.
char ** name
Well names.
Definition: wells.h:107
Controls for a single well.
Definition: well_controls.c:72
Well constrained by THP target.
Definition: well_controls.h:36
int * well_cells
Array of perforation cell indices.
Definition: wells.h:87
int number_of_wells
Number of wells.
Definition: wells.h:52
std::vector< double > & wellRates()
One rate per well and phase.
Definition: WellState.hpp:181
Data structure aggregating static information about all wells in a scenario.
Definition: wells.h:50
Well is an injector.
Definition: wells.h:42
std::vector< double > & temperature()
One temperature per well.
Definition: WellState.hpp:177
int * well_connpos
Array of indices into well_cells (and WI).
Definition: wells.h:81
Main OPM-Core well data structure along with functions to create, populate and destroy it...
Definition: WellState.hpp:312
enum WellType * type
Array of well types.
Definition: wells.h:58
int number_of_phases
Number of phases.
Definition: wells.h:53
std::vector< double > & bhp()
One bhp pressure per well.
Definition: WellState.hpp:169
Definition: BlackoilPhases.hpp:36
std::vector< double > & thp()
One thp pressure per well.
Definition: WellState.hpp:173
void init(const Wells *wells, const std::vector< double > &cellPressures)
Allocate and initialize if wells is non-null.
Definition: WellState.hpp:56
Well is a producer.
Definition: wells.h:43
struct Wells * clone_wells(const struct Wells *W)
Create a deep-copy (i.e., clone) of an existing Wells object, including its controls.
void destroy_wells(struct Wells *W)
Wells object destructor.
struct WellControls ** ctrls
Well controls, one set of controls for each well.
Definition: wells.h:102
std::vector< double > & perfRates()
One rate per well connection.
Definition: WellState.hpp:185
std::vector< double > & perfPress()
One pressure per well connection.
Definition: WellState.hpp:189
The state of a set of wells.
Definition: WellState.hpp:39
Well constrained by surface volume flow rate.
Definition: well_controls.h:38
int numPhases() const
The number of phases present.
Definition: WellState.hpp:222