Use well index for perforation rates

This commit is contained in:
Joakim Hove
2021-06-03 11:08:55 +02:00
parent 739f6ad28a
commit 1a6737aff3
6 changed files with 22 additions and 23 deletions

View File

@@ -1997,7 +1997,7 @@ namespace Opm {
auto& perf_pressure = well_state.perfPress(well_index);
auto& perf_rates = well_state.perfRates(well_index);
auto * perf_phase_rates = &well_state.mutable_perfPhaseRates()[wm.second[1]*np];
auto * perf_phase_rates = well_state.perfPhaseRates(well_index);
const auto& perf_data = this->well_perf_data_[well_index];
for (std::size_t perf_index = 0; perf_index < perf_data.size(); perf_index++) {
@@ -3374,7 +3374,7 @@ namespace Opm {
auto& well_info = *local_parallel_well_info_[wellID];
const int num_perf_this_well = well_info.communication().sum(well_perf_data_[wellID].size());
auto * perf_phase_rate = &this->wellState().perfPhaseRates()[connpos];
auto * perf_phase_rate = this->wellState().perfPhaseRates(wellID);
for (int perf = 0; perf < num_perf_this_well; ++perf) {
const int cell_idx = well_perf_data_[wellID][perf].cell_index;

View File

@@ -2637,8 +2637,7 @@ namespace Opm
// calculating the perforation rate for each perforation that belongs to this segment
const EvalWell seg_pressure = getSegmentPressure(seg);
const int rate_start_offset = first_perf_ * number_of_phases_;
auto * perf_rates = &well_state.mutable_perfPhaseRates()[rate_start_offset];
auto * perf_rates = well_state.perfPhaseRates(this->index_of_well_);
auto& perf_press_state = well_state.perfPress(this->index_of_well_);
for (const int perf : segment_perforations_[seg]) {
const int cell_idx = well_cells_[perf];

View File

@@ -591,8 +591,7 @@ namespace Opm
const int np = number_of_phases_;
std::vector<RateVector> connectionRates = connectionRates_; // Copy to get right size.
const int rate_start_offset = first_perf_ * number_of_phases_;
auto * perf_rates = &well_state.mutable_perfPhaseRates()[rate_start_offset];
auto * perf_rates = well_state.perfPhaseRates(this->index_of_well_);
for (int perf = 0; perf < number_of_perforations_; ++perf) {
// Calculate perforation quantities.
std::vector<EvalWell> cq_s(num_components_, {numWellEq_ + numEq, 0.0});
@@ -2105,7 +2104,7 @@ namespace Opm
const int nperf = number_of_perforations_;
const int np = number_of_phases_;
std::vector<double> perfRates(b_perf.size(),0.0);
const auto * perf_rates_state = &well_state.perfPhaseRates()[first_perf_ * np];
const auto * perf_rates_state = well_state.perfPhaseRates(this->index_of_well_);
for (int perf = 0; perf < nperf; ++perf) {
for (int comp = 0; comp < np; ++comp) {

View File

@@ -846,7 +846,7 @@ checkMaxRatioLimitCompletions(const WellState& well_state,
double max_ratio_completion = 0;
const int np = number_of_phases_;
const auto * perf_phase_rates = &well_state.perfPhaseRates()[first_perf_ * np];
const auto * perf_phase_rates = well_state.perfPhaseRates(this->index_of_well_);
// look for the worst_offending_completion
for (const auto& completion : completions_) {
std::vector<double> completion_rates(np, 0.0);

View File

@@ -298,7 +298,9 @@ void WellState::init(const std::vector<double>& cellPressures,
const int num_perf_this_well = well_info[2];
const int global_num_perf_this_well = parallel_well_info[w]->communication().sum(num_perf_this_well);
auto& perf_press = this->perfPress(w);
auto * phase_rates = &this->mutable_perfPhaseRates()[connpos * this->numPhases()];
first_perf_index_[w] = connpos;
auto phase_rates = this->perfPhaseRates(w);
for (int perf = 0; perf < num_perf_this_well; ++perf) {
if (wells_ecl[w].getStatus() == Well::Status::OPEN) {
@@ -308,7 +310,6 @@ void WellState::init(const std::vector<double>& cellPressures,
}
perf_press[perf] = cellPressures[well_perf_data[w][perf].cell_index];
}
first_perf_index_[w] = connpos;
this->well_reservoir_rates_.add(wname, std::vector<double>(np, 0));
this->well_dissolved_gas_rates_.add(wname, 0);
@@ -408,7 +409,6 @@ void WellState::init(const std::vector<double>& cellPressures,
}
// perfPhaseRates
const int oldPerf_idx_beg = (*it).second[ 1 ];
const int num_perf_old_well = (*it).second[ 2 ];
const auto new_iter = this->wellMap().find(well.name());
if (new_iter == this->wellMap().end()) {
@@ -418,7 +418,6 @@ void WellState::init(const std::vector<double>& cellPressures,
};
}
const int connpos = new_iter->second[1];
const int num_perf_this_well = new_iter->second[2];
const int num_perf_changed = parallel_well_info[w]->communication()
@@ -432,8 +431,8 @@ void WellState::init(const std::vector<double>& cellPressures,
// number of perforations.
if (global_num_perf_same)
{
const auto * src_rates = &prevState->perfPhaseRates()[oldPerf_idx_beg* np];
auto * target_rates = &this->mutable_perfPhaseRates()[connpos*np];
const auto * src_rates = prevState->perfPhaseRates(oldIndex);
auto * target_rates = this->perfPhaseRates(newIndex);
for (int perf_index = 0; perf_index < num_perf_this_well; perf_index++) {
for (int p = 0; p < np; p++) {
target_rates[perf_index*np + p] = src_rates[perf_index*np + p];
@@ -441,7 +440,7 @@ void WellState::init(const std::vector<double>& cellPressures,
}
} else {
const int global_num_perf_this_well = parallel_well_info[w]->communication().sum(num_perf_this_well);
auto * target_rates = &this->mutable_perfPhaseRates()[connpos*np];
auto * target_rates = this->perfPhaseRates(newIndex);
for (int perf_index = 0; perf_index < num_perf_this_well; perf_index++) {
for (int p = 0; p < np; ++p) {
target_rates[perf_index*np + p] = wellRates(w)[p] / double(global_num_perf_this_well);
@@ -771,7 +770,7 @@ void WellState::reportConnections(data::Well& well,
for( auto& comp : well.connections) {
const auto connPhaseOffset = np * (wt.second[1] + local_comp_index);
const auto * rates = &this->perfPhaseRates()[connPhaseOffset];
const auto * rates = &this->perfPhaseRates(well_index)[np*local_comp_index];
const auto connPI = this->connectionProductivityIndex().begin() + connPhaseOffset;
for( int i = 0; i < np; ++i ) {
@@ -813,7 +812,6 @@ void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
const auto& well_ecl = wells_ecl[w];
const auto& wname = wells_ecl[w].name();
const auto& well_info = this->wellMap().at(wname);
const int connpos = well_info[1];
const int num_perf_this_well = well_info[2];
if ( well_ecl.isMultiSegment() ) {
@@ -852,14 +850,12 @@ void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
// for the seg_rates_, now it becomes a recursive solution procedure.
{
const int start_perf = connpos;
// make sure the information from wells_ecl consistent with wells
assert((n_activeperf == num_perf_this_well) &&
"Inconsistent number of reservoir connections in well");
if (pu.phase_used[Gas]) {
auto * perf_rates = &this->mutable_perfPhaseRates()[np * start_perf];
auto * perf_rates = this->perfPhaseRates(w);
const int gaspos = pu.phase_pos[Gas];
// scale the phase rates for Gas to avoid too bad initial guess for gas fraction
// it will probably benefit the standard well too, while it needs to be justified
@@ -870,7 +866,7 @@ void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
perf_rates[perf*np + gaspos] *= 100;
}
const auto * perf_rates = &perfPhaseRates()[np*start_perf];
const auto * perf_rates = this->perfPhaseRates(w);
std::vector<double> perforation_rates(perf_rates, perf_rates + num_perf_this_well*np);
auto& segments = this->segments(w);

View File

@@ -100,8 +100,13 @@ public:
const SummaryState& summary_state);
/// One rate per phase and well connection.
std::vector<double>& mutable_perfPhaseRates() { return perfphaserates_; }
const std::vector<double>& perfPhaseRates() const { return perfphaserates_; }
double * perfPhaseRates(std::size_t well_index) {
return &this->perfphaserates_[this->first_perf_index_[well_index] * this->numPhases()];
}
const double * perfPhaseRates(std::size_t well_index) const {
return &this->perfphaserates_[this->first_perf_index_[well_index] * this->numPhases()];
}
/// One current control per injecting well.
Well::InjectorCMode currentInjectionControl(std::size_t well_index) const { return current_injection_controls_[well_index]; }