20 #ifndef OPM_WELLSTATEMULTISEGMENT_HEADER_INCLUDED 21 #define OPM_WELLSTATEMULTISEGMENT_HEADER_INCLUDED 24 #include <opm/core/wells.h> 25 #include <opm/core/well_controls.h> 26 #include <opm/common/ErrorMacros.hpp> 27 #include <opm/autodiff/AutoDiffBlock.hpp> 29 #include <opm/autodiff/MultisegmentWells.hpp> 30 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp> 55 int number_of_segments;
56 int start_perforation;
57 int number_of_perforations;
58 std::vector<int> start_perforation_segment;
59 std::vector<int> number_of_perforations_segment;
62 typedef std::map<std::string, SegmentedMapentryType> SegmentedWellMapType;
67 template <
class ReservoirState,
class PrevWellState>
68 void init(
const MultisegmentWells& ms_wells,
const ReservoirState& state,
const PrevWellState& prevState,
const Wells* legacy_wells_ptr)
71 this->wells_.reset( clone_wells( legacy_wells_ptr ) );
73 const std::vector<WellMultiSegmentConstPtr>& wells = ms_wells.msWells();
74 const int nw = wells.size();
80 segphaserates_.clear();
85 const int np = wells[0]->numberOfPhases();
87 for (
int iw = 0; iw < nw; ++iw) {
88 nperf_ += wells[iw]->numberOfPerforations();
89 nseg_ += wells[iw]->numberOfSegments();
94 top_segment_loc_.resize(nw);
95 temperature().resize(nw, 273.15 + 20);
97 wellRates().resize(nw * np, 0.0);
100 for(
int iw = 0; iw < nw; ++iw) {
101 currentControls()[iw] = well_controls_get_current(wells[iw]->wellControls());
104 for (
int iw = 0; iw < nw; ++iw) {
105 assert((wells[iw]->wellType() == INJECTOR) || (wells[iw]->wellType() == PRODUCER));
108 int start_segment = 0;
109 int start_perforation = 0;
112 perfPress().resize(nperf_, -1.0e100);
113 segphaserates_.resize(nseg_ * np, 0.0);
114 segpress_.resize(nseg_, -1.0e100);
117 for (
int w = 0; w < nw; ++w) {
118 assert((wells[w]->wellType() == INJECTOR) || (wells[w]->wellType() == PRODUCER));
119 const WellControls* ctrl = wells[w]->wellControls();
121 std::string well_name(wells[w]->name());
124 wellMapEntry.well_number = w;
125 wellMapEntry.start_segment = start_segment;
126 wellMapEntry.number_of_segments = wells[w]->numberOfSegments();
127 wellMapEntry.start_perforation = start_perforation;
128 wellMapEntry.number_of_perforations = wells[w]->numberOfPerforations();
130 top_segment_loc_[w] = start_segment;
132 int start_perforation_segment = 0;
133 wellMapEntry.start_perforation_segment.resize(wellMapEntry.number_of_segments);
134 wellMapEntry.number_of_perforations_segment.resize(wellMapEntry.number_of_segments);
136 for (
int i = 0; i < wellMapEntry.number_of_segments; ++i) {
137 wellMapEntry.start_perforation_segment[i] = start_perforation_segment;
138 wellMapEntry.number_of_perforations_segment[i] = wells[w]->segmentPerforations()[i].size();
139 start_perforation_segment += wellMapEntry.number_of_perforations_segment[i];
141 assert(start_perforation_segment == wellMapEntry.number_of_perforations);
144 if (well_controls_well_is_stopped(ctrl)) {
148 if (well_controls_get_current_type(ctrl) == BHP) {
149 bhp()[w] = well_controls_get_current_target(ctrl);
151 const int first_cell = wells[0]->wellCells()[0];
152 bhp()[w] = state.pressure()[first_cell];
156 if (well_controls_get_current_type(ctrl) == THP) {
157 thp()[w] = well_controls_get_current_target( ctrl );
169 if (well_controls_get_current_type(ctrl) == SURFACE_RATE) {
170 const double rate_target = well_controls_get_current_target(ctrl);
171 const double * distr = well_controls_get_current_distr( ctrl );
172 for (
int p = 0; p < np; ++p) {
173 wellRates()[np * w + p] = rate_target * distr[p];
176 const double small_rate = 1e-14;
177 const double sign = (wells[w]->wellType() == INJECTOR) ? 1.0 : -1.0;
178 for (
int p = 0; p < np; ++p) {
179 wellRates()[np * w + p] = small_rate *
sign;
184 if (well_controls_get_current_type(ctrl) == BHP) {
185 bhp()[w] = well_controls_get_current_target(ctrl);
187 const int first_cell = wells[w]->wellCells()[0];
188 const double safety_factor = (wells[w]->wellType() == INJECTOR) ? 1.01 : 0.99;
189 bhp()[w] = safety_factor* state.pressure()[first_cell];
192 if (well_controls_get_current_type(ctrl) == THP) {
193 thp()[w] = well_controls_get_current_target(ctrl);
199 int number_of_perforations = wellMapEntry.number_of_perforations;
200 for (
int i = 0; i < number_of_perforations; ++i) {
201 for (
int p = 0; p < np; ++p) {
202 perfPhaseRates()[np * (i + start_perforation) + p] = wellRates()[np * w + p] / double(number_of_perforations);
204 if (wells[w]->isMultiSegmented()) {
205 const double safety_factor = (wells[w]->wellType() == INJECTOR) ? 1.01 : 0.99;
206 perfPress()[i + start_perforation] = safety_factor * state.pressure()[wells[w]->wellCells()[i]];
208 perfPress()[i + start_perforation] = state.pressure()[wells[w]->wellCells()[i]];
213 int number_of_segments = wellMapEntry.number_of_segments;
218 segpress_[start_segment] = bhp()[w];
219 for (
int i = 1; i < number_of_segments; ++i) {
220 int first_perforation_segment = start_perforation + wellMapEntry.start_perforation_segment[i];
221 segpress_[i + start_segment] = perfPress()[first_perforation_segment];
225 for (
int p = 0; p < np; ++p) {
226 Eigen::VectorXd v_perf_rates(number_of_perforations);
227 for (
int i = 0; i < number_of_perforations; ++i) {
228 v_perf_rates[i] =
perfPhaseRates()[np * (i + start_perforation) + p];
231 Eigen::VectorXd v_segment_rates = wells[w]->wellOps().p2s_gather * v_perf_rates;
233 for (
int i = 0; i < number_of_segments; ++i) {
234 segphaserates_[np * (i + start_segment) + p] = v_segment_rates[i];
239 start_segment += wellMapEntry.number_of_segments;
240 start_perforation += wellMapEntry.number_of_perforations;
248 for (
int w = 0; w < nw; ++w) {
249 currentControls()[w] = well_controls_get_current(wells[w]->wellControls());
254 if ( !(prevState.segmentedWellMap().empty()) )
256 typedef typename SegmentedWellMapType::const_iterator const_iterator;
257 const_iterator end_old = prevState.segmentedWellMap().end();
259 for (
int w = 0; w < nw; ++w) {
260 std::string well_name(wells[w]->name());
261 const_iterator it_old = prevState.segmentedWellMap().find(well_name);
262 const_iterator it_this = segmentedWellMap().find(well_name);
264 assert(it_this != segmentedWellMap().end());
266 if (it_old != end_old) {
267 const int oldIndex = (*it_old).second.well_number;
268 const int newIndex = w;
271 bhp()[newIndex] = prevState.bhp()[oldIndex];
274 for(
int i=0, idx=newIndex*np, oldidx=oldIndex*np; i<np; ++i, ++idx, ++oldidx )
276 wellRates()[ idx ] = prevState.wellRates()[ oldidx ];
279 const int num_seg_old_well = (*it_old).second.number_of_segments;
280 const int num_perf_old_well = (*it_old).second.number_of_perforations;
282 const int num_seg_this_well = (*it_this).second.number_of_segments;
283 const int num_perf_this_well = (*it_this).second.number_of_perforations;
288 if ((num_perf_old_well == num_perf_this_well) && (num_seg_old_well == num_seg_this_well)) {
289 const int old_start_perforation = (*it_old).second.start_perforation;
290 const int old_start_segment = (*it_old).second.start_segment;
292 const int this_start_perforation = (*it_this).second.start_perforation;
293 const int this_start_segment = (*it_this).second.start_segment;
296 for (
int i = 0; i < num_seg_this_well * np; ++i) {
297 segphaserates_[this_start_segment * np + i] = prevState.segPhaseRates()[old_start_segment * np + i];
300 for (
int i = 0; i < num_perf_this_well * np; ++i) {
301 perfPhaseRates()[this_start_perforation * np + i] = prevState.perfPhaseRates()[old_start_perforation * np + i];
305 for (
int i = 0; i < num_perf_this_well; ++i) {
307 perfPress()[this_start_perforation + i] = prevState.perfPress()[old_start_perforation + i];
311 for (
int i = 0; i < num_seg_this_well; ++i) {
313 segpress_[this_start_segment + i] = prevState.segPress()[old_start_segment + i];
317 const int old_control_index = prevState.currentControls()[ oldIndex ];
318 if (old_control_index < well_controls_get_num(wells[w]->wellControls())) {
330 std::vector<double>& segPress() {
return segpress_; }
331 const std::vector<double>& segPress()
const {
return segpress_; }
333 std::vector<double>& segPhaseRates() {
return segphaserates_; }
334 const std::vector<double>& segPhaseRates()
const {
return segphaserates_; }
336 const std::vector<int>& topSegmentLoc()
const {
return top_segment_loc_; };
338 const SegmentedWellMapType& segmentedWellMap()
const {
return segmentedWellMap_; }
339 SegmentedWellMapType& segmentedWellMap() {
return segmentedWellMap_; }
341 int numSegments()
const {
return nseg_; }
342 int numPerforations()
const {
return nperf_; }
346 std::vector<double> segpress_;
348 std::vector<double> segphaserates_;
353 std::vector<int> top_segment_loc_;
355 SegmentedWellMapType segmentedWellMap_;
364 #endif // OPM_WELLSTATEMULTISEGMENT_HEADER_INCLUDE This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
void init(const MultisegmentWells &ms_wells, const ReservoirState &state, const PrevWellState &prevState, const Wells *legacy_wells_ptr)
Allocate and initialize if wells is non-null.
Definition: WellStateMultiSegment.hpp:68
std::vector< double > & perfPhaseRates()
One rate per phase and well connection.
Definition: WellStateFullyImplicitBlackoil.hpp:230
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
std::vector< int > & currentControls()
One current control per well.
Definition: WellStateFullyImplicitBlackoil.hpp:234
The state of a set of multi-sgemnet wells.
Definition: WellStateMultiSegment.hpp:45
Definition: WellStateMultiSegment.hpp:52
Class for handling the multi-segment well model.
Definition: MultisegmentWells.hpp:55
Eigen::ArrayXd sign(const Eigen::ArrayXd &x)
Return a vector of (-1.0, 0.0 or 1.0), depending on sign per element.
Definition: AutoDiffHelpers.hpp:719