20#ifndef OPM_WELLSTATE_HEADER_INCLUDED
21#define OPM_WELLSTATE_HEADER_INCLUDED
26#include <opm/output/data/Wells.hpp>
45 template <
class 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) {
71 const WellControls* ctrl = wells->
ctrls[w];
76 assert( wells->
name[ w ] );
77 std::string name( wells->
name[ w ] );
78 assert( name.size() > 0 );
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;
100 for (
int p = 0; p < np; ++p) {
101 wellrates_[np*w + p] = 0.0;
110 bhp_[w] = cellPressures[first_cell];
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;
143 const double safety_factor = (wells->
type[w] ==
INJECTOR) ? 1.01 : 0.99;
144 bhp_[w] = safety_factor*cellPressures[first_cell];
152 for (
int ctrl_index = 0; ctrl_index < nwc; ++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_; }
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;
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()));
280 temperature_( rhs.temperature_ ),
281 wellrates_( rhs.wellrates_ ),
282 perfrates_( rhs.perfrates_ ),
283 perfpress_( rhs.perfpress_ ),
284 wellMap_( rhs.wellMap_ ),
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_;
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_;
@ Liquid
Definition: BlackoilPhases.hpp:40
@ Aqua
Definition: BlackoilPhases.hpp:40
@ Vapour
Definition: BlackoilPhases.hpp:40
The state of a set of wells.
Definition: WellState.hpp:40
const std::vector< double > & perfRates() const
Definition: WellState.hpp:186
std::vector< double > & temperature()
One temperature per well.
Definition: WellState.hpp:177
virtual ~WellState()
Definition: WellState.hpp:274
WellState & operator=(const WellState &rhs)
Definition: WellState.hpp:288
size_t getRestartPerfRatesOffset() const
Definition: WellState.hpp:200
void init(const Wells *wells, const State &state)
Definition: WellState.hpp:46
int numWells() const
The number of wells present.
Definition: WellState.hpp:216
std::vector< double > & thp()
One thp pressure per well.
Definition: WellState.hpp:173
size_t getRestartBhpOffset() const
Definition: WellState.hpp:192
size_t getRestartTemperatureOffset() const
Definition: WellState.hpp:204
std::array< int, 3 > mapentry_t
Definition: WellState.hpp:42
const WellMapType & wellMap() const
Definition: WellState.hpp:212
const std::vector< double > & bhp() const
Definition: WellState.hpp:170
const std::vector< double > & perfPress() const
Definition: WellState.hpp:190
const std::vector< double > & thp() const
Definition: WellState.hpp:174
const std::vector< double > & temperature() const
Definition: WellState.hpp:178
std::vector< double > & perfPress()
One pressure per well connection.
Definition: WellState.hpp:189
std::vector< double > & perfRates()
One rate per well connection.
Definition: WellState.hpp:185
WellMapType & wellMap()
Definition: WellState.hpp:213
int numPhases() const
The number of phases present.
Definition: WellState.hpp:222
std::map< std::string, mapentry_t > WellMapType
Definition: WellState.hpp:43
std::vector< double > & bhp()
One bhp pressure per well.
Definition: WellState.hpp:169
size_t getRestartPerfPressOffset() const
Definition: WellState.hpp:196
const std::vector< double > & wellRates() const
Definition: WellState.hpp:182
void init(const Wells *wells, const std::vector< double > &cellPressures)
Definition: WellState.hpp:56
std::vector< double > & wellRates()
One rate per well and phase.
Definition: WellState.hpp:181
size_t getRestartWellRatesOffset() const
Definition: WellState.hpp:208
virtual data::Wells report(const PhaseUsage &pu) const
Definition: WellState.hpp:227
std::unique_ptr< Wells, wdel > wells_
Definition: WellState.hpp:315
WellState(const WellState &rhs)
Definition: WellState.hpp:277
@ BHP
Definition: legacy_well.h:42
@ INJECTOR
Definition: legacy_well.h:37
@ PRODUCER
Definition: legacy_well.h:37
void gas(const Grid &G, const Region ®, const std::array< double, 2 > &span, const double grav, double &po_goc, const CellRange &cells, std::vector< double > &press)
Definition: initStateEquil_impl.hpp:378
void oil(const Grid &G, const Region ®, const std::array< double, 2 > &span, const double grav, const CellRange &cells, std::vector< double > &press, double &po_woc, double &po_goc)
Definition: initStateEquil_impl.hpp:319
Definition: AnisotropicEikonal.hpp:44
Definition: BlackoilPhases.hpp:44
int phase_pos[MaxNumPhases+NumCryptoPhases]
Definition: BlackoilPhases.hpp:47
int phase_used[MaxNumPhases+NumCryptoPhases]
Definition: BlackoilPhases.hpp:46
int num_phases
Definition: BlackoilPhases.hpp:45
Definition: WellState.hpp:312
void operator()(Wells *w)
Definition: WellState.hpp:313
int * well_connpos
Definition: wells.h:81
enum WellType * type
Definition: wells.h:58
char ** name
Definition: wells.h:107
struct WellControls ** ctrls
Definition: wells.h:102
int number_of_wells
Definition: wells.h:52
int * well_cells
Definition: wells.h:87
int number_of_phases
Definition: wells.h:53
bool well_controls_well_is_stopped(const struct WellControls *ctrl)
double well_controls_iget_target(const struct WellControls *ctrl, int control_index)
enum WellControlType well_controls_get_current_type(const struct WellControls *ctrl)
double well_controls_get_current_target(const struct WellControls *ctrl)
int well_controls_get_num(const struct WellControls *ctrl)
const double * well_controls_get_current_distr(const struct WellControls *ctrl)
enum WellControlType well_controls_iget_type(const struct WellControls *ctrl, int control_index)
@ SURFACE_RATE
Definition: well_controls.h:38
@ THP
Definition: well_controls.h:36
struct Wells * clone_wells(const struct Wells *W)
void destroy_wells(struct Wells *W)