WellInterface_impl.hpp
Go to the documentation of this file.
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2018 IRIS
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_WELLINTERFACE_IMPL_HEADER_INCLUDED
23#define OPM_WELLINTERFACE_IMPL_HEADER_INCLUDED
24
25// Improve IDE experience
26#ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
27#include <config.h>
29#endif
30
31#include <opm/common/Exceptions.hpp>
32
33#include <opm/input/eclipse/Schedule/ScheduleTypes.hpp>
34#include <opm/input/eclipse/Schedule/Well/WDFAC.hpp>
35
37
42
43#include <dune/common/version.hh>
44
45#include <algorithm>
46#include <cassert>
47#include <cstddef>
48#include <utility>
49
50#include <fmt/format.h>
51
52namespace Opm
53{
54
55
56 template<typename TypeTag>
58 WellInterface(const Well& well,
59 const ParallelWellInfo<Scalar>& pw_info,
60 const int time_step,
61 const ModelParameters& param,
62 const RateConverterType& rate_converter,
63 const int pvtRegionIdx,
64 const int num_conservation_quantities,
65 const int num_phases,
66 const int index_of_well,
67 const std::vector<PerforationData<Scalar>>& perf_data)
69 pw_info,
70 time_step,
71 param,
72 rate_converter,
73 pvtRegionIdx,
74 num_conservation_quantities,
75 num_phases,
76 index_of_well,
77 perf_data)
78 {
80
81 if constexpr (has_solvent || has_zFraction) {
82 if (well.isInjector()) {
83 auto injectorType = this->well_ecl_.injectorType();
84 if (injectorType == InjectorType::GAS) {
85 this->wsolvent_ = this->well_ecl_.getSolventFraction();
86 }
87 }
88 }
89 }
90
91
92 template<typename TypeTag>
93 void
95 init(const std::vector<Scalar>& /* depth_arg */,
96 const Scalar gravity_arg,
97 const std::vector<Scalar>& B_avg,
98 const bool changed_to_open_this_step)
99 {
100 this->gravity_ = gravity_arg;
101 B_avg_ = B_avg;
102 this->changed_to_open_this_step_ = changed_to_open_this_step;
103 }
104
105
106
107
108 template<typename TypeTag>
111 wpolymer() const
112 {
113 if constexpr (has_polymer) {
114 return this->wpolymer_();
115 }
116
117 return 0.0;
118 }
119
120
121
122
123
124 template<typename TypeTag>
127 wfoam() const
128 {
129 if constexpr (has_foam) {
130 return this->wfoam_();
131 }
132
133 return 0.0;
134 }
135
136
137
138 template<typename TypeTag>
141 wsalt() const
142 {
143 if constexpr (has_brine) {
144 return this->wsalt_();
145 }
146
147 return 0.0;
148 }
149
150 template<typename TypeTag>
153 wmicrobes() const
154 {
155 if constexpr (has_micp) {
156 return this->wmicrobes_();
157 }
158
159 return 0.0;
160 }
161
162 template<typename TypeTag>
165 woxygen() const
166 {
167 if constexpr (has_micp) {
168 return this->woxygen_();
169 }
170
171 return 0.0;
172 }
173
174 template<typename TypeTag>
177 wurea() const
178 {
179 if constexpr (has_micp) {
180 return this->wurea_();
181 }
182
183 return 0.0;
184 }
185
186 template<typename TypeTag>
187 bool
189 updateWellControl(const Simulator& simulator,
190 const IndividualOrGroup iog,
191 WellStateType& well_state,
192 const GroupState<Scalar>& group_state,
193 DeferredLogger& deferred_logger) /* const */
194 {
195 OPM_TIMEFUNCTION();
196 if (stoppedOrZeroRateTarget(simulator, well_state, deferred_logger)) {
197 return false;
198 }
199
200 const auto& summaryState = simulator.vanguard().summaryState();
201 const auto& schedule = simulator.vanguard().schedule();
202 const auto& well = this->well_ecl_;
203 auto& ws = well_state.well(this->index_of_well_);
204 std::string from;
205 bool is_grup = false;
206 if (well.isInjector()) {
207 from = WellInjectorCMode2String(ws.injection_cmode);
208 is_grup = ws.injection_cmode == Well::InjectorCMode::GRUP;
209 } else {
210 from = WellProducerCMode2String(ws.production_cmode);
211 is_grup = ws.production_cmode == Well::ProducerCMode::GRUP;
212 }
213
214 const int episodeIdx = simulator.episodeIndex();
215 const int iterationIdx = simulator.model().newtonMethod().numIterations();
216 const int nupcol = schedule[episodeIdx].nupcol();
217 const bool oscillating = std::count(this->well_control_log_.begin(), this->well_control_log_.end(), from) >= this->param_.max_number_of_well_switches_;
218 if (oscillating && !is_grup) { // we would like to avoid ending up as GRUP
219 // only output first time
220 const bool output = std::count(this->well_control_log_.begin(), this->well_control_log_.end(), from) == this->param_.max_number_of_well_switches_;
221 if (output) {
222 const auto msg = fmt::format(" The control mode for well {} is oscillating. \n"
223 "We don't allow for more than {} switches after NUPCOL iterations. (NUPCOL = {}) \n"
224 "The control is kept at {}.",
225 this->name(), this->param_.max_number_of_well_switches_, nupcol, from);
226 deferred_logger.info(msg);
227 // add one more to avoid outputting the same info again
228 this->well_control_log_.push_back(from);
229 }
230 return false;
231 }
232 bool changed = false;
233 if (iog == IndividualOrGroup::Individual) {
234 changed = this->checkIndividualConstraints(ws, summaryState, deferred_logger);
235 } else if (iog == IndividualOrGroup::Group) {
236 changed = this->checkGroupConstraints(well_state, group_state, schedule, summaryState, true, deferred_logger);
237 } else {
238 assert(iog == IndividualOrGroup::Both);
239 changed = this->checkConstraints(well_state, group_state, schedule, summaryState, deferred_logger);
240 }
241 Parallel::Communication cc = simulator.vanguard().grid().comm();
242 // checking whether control changed
243 if (changed) {
244 std::string to;
245 if (well.isInjector()) {
246 to = WellInjectorCMode2String(ws.injection_cmode);
247 } else {
248 to = WellProducerCMode2String(ws.production_cmode);
249 }
250 std::ostringstream ss;
251 ss << " Switching control mode for well " << this->name()
252 << " from " << from
253 << " to " << to;
254 if (cc.size() > 1) {
255 ss << " on rank " << cc.rank();
256 }
257 deferred_logger.debug(ss.str());
258
259 // We always store the current control as it is used for output
260 // and only after iteration >= nupcol
261 // we log all switches to check if the well controls oscillates
262 if (iterationIdx >= nupcol || this->well_control_log_.empty()) {
263 this->well_control_log_.push_back(from);
264 }
265 updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger);
266 updatePrimaryVariables(simulator, well_state, deferred_logger);
267 }
268
269 return changed;
270 }
271
272 template<typename TypeTag>
273 bool
276 WellStateType& well_state,
277 const GroupState<Scalar>& group_state,
278 const Well::InjectionControls& inj_controls,
279 const Well::ProductionControls& prod_controls,
280 const Scalar wqTotal,
281 DeferredLogger& deferred_logger,
282 const bool fixed_control,
283 const bool fixed_status)
284 {
285 OPM_TIMEFUNCTION();
286 const auto& summary_state = simulator.vanguard().summaryState();
287 const auto& schedule = simulator.vanguard().schedule();
288 auto& ws = well_state.well(this->index_of_well_);
289 std::string from;
290 if (this->isInjector()) {
291 from = WellInjectorCMode2String(ws.injection_cmode);
292 } else {
293 from = WellProducerCMode2String(ws.production_cmode);
294 }
295 const bool oscillating = std::count(this->well_control_log_.begin(), this->well_control_log_.end(), from) >= this->param_.max_number_of_well_switches_;
296
297 if (oscillating || this->wellUnderZeroRateTarget(simulator, well_state, deferred_logger) || !(well_state.well(this->index_of_well_).status == WellStatus::OPEN)) {
298 return false;
299 }
300
301 const Scalar sgn = this->isInjector() ? 1.0 : -1.0;
302 if (!this->wellIsStopped()){
303 if (wqTotal*sgn <= 0.0 && !fixed_status){
304 this->stopWell();
305 return true;
306 } else {
307 bool changed = false;
308 if (!fixed_control) {
309 // Changing to group controls here may lead to inconsistencies in the group handling which in turn
310 // may result in excessive back and forth switching. However, we currently allow this by default.
311 // The switch check_group_constraints_inner_well_iterations_ is a temporary solution.
312
313 const bool hasGroupControl = this->isInjector() ? inj_controls.hasControl(Well::InjectorCMode::GRUP) :
314 prod_controls.hasControl(Well::ProducerCMode::GRUP);
315 bool isGroupControl = ws.production_cmode == Well::ProducerCMode::GRUP || ws.injection_cmode == Well::InjectorCMode::GRUP;
316 if (! (isGroupControl && !this->param_.check_group_constraints_inner_well_iterations_)) {
317 changed = this->checkIndividualConstraints(ws, summary_state, deferred_logger, inj_controls, prod_controls);
318 }
319 if (hasGroupControl && this->param_.check_group_constraints_inner_well_iterations_) {
320 changed = changed || this->checkGroupConstraints(well_state, group_state, schedule, summary_state, false, deferred_logger);
321 }
322
323 if (changed) {
324 const bool thp_controlled = this->isInjector() ? ws.injection_cmode == Well::InjectorCMode::THP :
325 ws.production_cmode == Well::ProducerCMode::THP;
326 if (thp_controlled){
327 ws.thp = this->getTHPConstraint(summary_state);
328 } else {
329 // don't call for thp since this might trigger additional local solve
330 updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger);
331 }
332 updatePrimaryVariables(simulator, well_state, deferred_logger);
333 }
334 }
335 return changed;
336 }
337 } else if (!fixed_status){
338 // well is stopped, check if current bhp allows reopening
339 const Scalar bhp = well_state.well(this->index_of_well_).bhp;
340 Scalar prod_limit = prod_controls.bhp_limit;
341 Scalar inj_limit = inj_controls.bhp_limit;
342 const bool has_thp = this->wellHasTHPConstraints(summary_state);
343 if (has_thp){
344 std::vector<Scalar> rates(this->num_conservation_quantities_);
345 if (this->isInjector()){
346 const Scalar bhp_thp = WellBhpThpCalculator(*this).
347 calculateBhpFromThp(well_state, rates,
348 this->well_ecl_,
349 summary_state,
350 this->getRefDensity(),
351 deferred_logger);
352 inj_limit = std::min(bhp_thp, static_cast<Scalar>(inj_controls.bhp_limit));
353 } else {
354 // if the well can operate, it must at least be able to produce
355 // at the lowest bhp of the bhp-curve (explicit fractions)
356 const Scalar bhp_min = WellBhpThpCalculator(*this).
357 calculateMinimumBhpFromThp(well_state,
358 this->well_ecl_,
359 summary_state,
360 this->getRefDensity());
361 prod_limit = std::max(bhp_min, static_cast<Scalar>(prod_controls.bhp_limit));
362 }
363 }
364 const Scalar bhp_diff = (this->isInjector())? inj_limit - bhp: bhp - prod_limit;
365 if (bhp_diff > 0){
366 this->openWell();
367 well_state.well(this->index_of_well_).bhp = (this->isInjector())? inj_limit : prod_limit;
368 if (has_thp) {
369 well_state.well(this->index_of_well_).thp = this->getTHPConstraint(summary_state);
370 }
371 return true;
372 } else {
373 return false;
374 }
375 } else {
376 return false;
377 }
378 }
379
380 template<typename TypeTag>
381 void
383 wellTesting(const Simulator& simulator,
384 const double simulation_time,
385 /* const */ WellStateType& well_state,
386 const GroupState<Scalar>& group_state,
387 WellTestState& well_test_state,
388 GLiftEclWells& ecl_well_map,
389 std::map<std::string, double>& open_times,
390 DeferredLogger& deferred_logger)
391 {
392 OPM_TIMEFUNCTION();
393 deferred_logger.info(" well " + this->name() + " is being tested");
394
395 WellStateType well_state_copy = well_state;
396 auto& ws = well_state_copy.well(this->indexOfWell());
397
398 const auto& summary_state = simulator.vanguard().summaryState();
399 const bool has_thp_limit = this->wellHasTHPConstraints(summary_state);
400 if (this->isProducer()) {
401 ws.production_cmode = has_thp_limit ? Well::ProducerCMode::THP : Well::ProducerCMode::BHP;
402 } else {
403 ws.injection_cmode = has_thp_limit ? Well::InjectorCMode::THP : Well::InjectorCMode::BHP;
404 }
405 // We test the well as an open well during the well testing
406 ws.open();
407
408 scaleSegmentRatesAndPressure(well_state_copy);
409 calculateExplicitQuantities(simulator, well_state_copy, deferred_logger);
410 updatePrimaryVariables(simulator, well_state_copy, deferred_logger);
411
412 if (this->isProducer()) {
413 const auto& schedule = simulator.vanguard().schedule();
414 const auto report_step = simulator.episodeIndex();
415 const auto& glo = schedule.glo(report_step);
416 if (glo.active()) {
417 gliftBeginTimeStepWellTestUpdateALQ(simulator,
418 well_state_copy,
419 group_state,
420 ecl_well_map,
421 deferred_logger);
422 }
423 }
424
425 WellTestState welltest_state_temp;
426
427 bool testWell = true;
428 // if a well is closed because all completions are closed, we need to check each completion
429 // individually. We first open all completions, then we close one by one by calling updateWellTestState
430 // untill the number of closed completions do not increase anymore.
431 while (testWell) {
432 const std::size_t original_number_closed_completions = welltest_state_temp.num_closed_completions();
433 bool converged = solveWellForTesting(simulator, well_state_copy, group_state, deferred_logger);
434 if (!converged) {
435 const auto msg = fmt::format("WTEST: Well {} is not solvable (physical)", this->name());
436 deferred_logger.debug(msg);
437 return;
438 }
439
440
441 updateWellOperability(simulator, well_state_copy, deferred_logger);
442 if ( !this->isOperableAndSolvable() ) {
443 const auto msg = fmt::format("WTEST: Well {} is not operable (physical)", this->name());
444 deferred_logger.debug(msg);
445 return;
446 }
447 std::vector<Scalar> potentials;
448 try {
449 computeWellPotentials(simulator, well_state_copy, potentials, deferred_logger);
450 } catch (const std::exception& e) {
451 const std::string msg = fmt::format("well {}: computeWellPotentials() "
452 "failed during testing for re-opening: ",
453 this->name(), e.what());
454 deferred_logger.info(msg);
455 return;
456 }
457 const int np = well_state_copy.numPhases();
458 for (int p = 0; p < np; ++p) {
459 ws.well_potentials[p] = std::max(Scalar{0.0}, potentials[p]);
460 }
461 const bool under_zero_target = this->wellUnderZeroGroupRateTarget(simulator, well_state_copy, deferred_logger);
462 this->updateWellTestState(well_state_copy.well(this->indexOfWell()),
463 simulation_time,
464 /*writeMessageToOPMLog=*/ false,
465 under_zero_target,
466 welltest_state_temp,
467 deferred_logger);
468 this->closeCompletions(welltest_state_temp);
469
470 // Stop testing if the well is closed or shut due to all completions shut
471 // Also check if number of completions has increased. If the number of closed completions do not increased
472 // we stop the testing.
473 // TODO: it can be tricky here, if the well is shut/closed due to other reasons
474 if ( welltest_state_temp.num_closed_wells() > 0 ||
475 (original_number_closed_completions == welltest_state_temp.num_closed_completions()) ) {
476 testWell = false; // this terminates the while loop
477 }
478 }
479
480 // update wellTestState if the well test succeeds
481 if (!welltest_state_temp.well_is_closed(this->name())) {
482 well_test_state.open_well(this->name());
483
484 std::string msg = std::string("well ") + this->name() + std::string(" is re-opened");
485 deferred_logger.info(msg);
486
487 // also reopen completions
488 for (const auto& completion : this->well_ecl_.getCompletions()) {
489 if (!welltest_state_temp.completion_is_closed(this->name(), completion.first))
490 well_test_state.open_completion(this->name(), completion.first);
491 }
492 well_state = well_state_copy;
493 open_times.try_emplace(this->name(), well_test_state.lastTestTime(this->name()));
494 }
495 }
496
497
498
499
500 template<typename TypeTag>
501 bool
503 iterateWellEquations(const Simulator& simulator,
504 const double dt,
505 WellStateType& well_state,
506 const GroupState<Scalar>& group_state,
507 DeferredLogger& deferred_logger)
508 {
509 OPM_TIMEFUNCTION();
510 const auto& summary_state = simulator.vanguard().summaryState();
511 const auto inj_controls = this->well_ecl_.isInjector() ? this->well_ecl_.injectionControls(summary_state) : Well::InjectionControls(0);
512 const auto prod_controls = this->well_ecl_.isProducer() ? this->well_ecl_.productionControls(summary_state) : Well::ProductionControls(0);
513 const auto& ws = well_state.well(this->indexOfWell());
514 const auto pmode_orig = ws.production_cmode;
515 const auto imode_orig = ws.injection_cmode;
516 bool converged = false;
517 try {
518 // TODO: the following two functions will be refactored to be one to reduce the code duplication
519 if (!this->param_.local_well_solver_control_switching_){
520 converged = this->iterateWellEqWithControl(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
521 } else {
522 if (this->param_.use_implicit_ipr_ && this->well_ecl_.isProducer() && (well_state.well(this->index_of_well_).status == WellStatus::OPEN)) {
523 converged = solveWellWithOperabilityCheck(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
524 } else {
525 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
526 }
527 }
528
529 } catch (NumericalProblem& e ) {
530 const std::string msg = "Inner well iterations failed for well " + this->name() + " Treat the well as unconverged. ";
531 deferred_logger.warning("INNER_ITERATION_FAILED", msg);
532 converged = false;
533 }
534 if (converged) {
535 // Add debug info for switched controls
536 if (ws.production_cmode != pmode_orig || ws.injection_cmode != imode_orig) {
537 std::string from,to;
538 if (this->isInjector()) {
539 from = WellInjectorCMode2String(imode_orig);
540 to = WellInjectorCMode2String(ws.injection_cmode);
541 } else {
542 from = WellProducerCMode2String(pmode_orig);
543 to = WellProducerCMode2String(ws.production_cmode);
544 }
545 const auto msg = fmt::format(" Well {} switched from {} to {} during local solve", this->name(), from, to);
546 deferred_logger.debug(msg);
547 const int episodeIdx = simulator.episodeIndex();
548 const int iterationIdx = simulator.model().newtonMethod().numIterations();
549 const auto& schedule = simulator.vanguard().schedule();
550 const int nupcol = schedule[episodeIdx].nupcol();
551 // We always store the current control as it is used for output
552 // and only after iteration >= nupcol
553 // we log all switches to check if the well controls oscillates
554 if (iterationIdx >= nupcol || this->well_control_log_.empty()) {
555 this->well_control_log_.push_back(from);
556 }
557 }
558 }
559
560 return converged;
561 }
562
563 template<typename TypeTag>
564 bool
567 const double dt,
568 const Well::InjectionControls& inj_controls,
569 const Well::ProductionControls& prod_controls,
570 WellStateType& well_state,
571 const GroupState<Scalar>& group_state,
572 DeferredLogger& deferred_logger)
573 {
574 OPM_TIMEFUNCTION();
575 const auto& summary_state = simulator.vanguard().summaryState();
576 bool converged = true;
577 auto& ws = well_state.well(this->index_of_well_);
578 // if well is stopped, check if we can reopen
579 if (this->wellIsStopped()) {
580 this->openWell();
581 auto bhp_target = estimateOperableBhp(simulator, dt, well_state, summary_state, deferred_logger);
582 if (!bhp_target.has_value()) {
583 // no intersection with ipr
584 const auto msg = fmt::format("estimateOperableBhp: Did not find operable BHP for well {}", this->name());
585 deferred_logger.debug(msg);
586 // well can't operate using explicit fractions stop the well
587 // solve with zero rates
588 converged = solveWellWithZeroRate(simulator, dt, well_state, deferred_logger);
589 this->stopWell();
590 this->operability_status_.can_obtain_bhp_with_thp_limit = false;
591 this->operability_status_.obey_thp_limit_under_bhp_limit = false;
592 return converged;
593 } else {
594 // solve well with the estimated target bhp (or limit)
595 ws.thp = this->getTHPConstraint(summary_state);
596 const Scalar bhp = std::max(bhp_target.value(),
597 static_cast<Scalar>(prod_controls.bhp_limit));
598 solveWellWithBhp(simulator, dt, bhp, well_state, deferred_logger);
599 }
600 }
601 // solve well-equation
602 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
603
604
605 const bool isThp = ws.production_cmode == Well::ProducerCMode::THP;
606 // check stability of solution under thp-control
607 if (converged && !stoppedOrZeroRateTarget(simulator, well_state, deferred_logger) && isThp) {
608 auto rates = well_state.well(this->index_of_well_).surface_rates;
609 this->adaptRatesForVFP(rates);
610 this->updateIPRImplicit(simulator, well_state, deferred_logger);
611 bool is_stable = WellBhpThpCalculator(*this).isStableSolution(well_state, this->well_ecl_, rates, summary_state);
612 if (!is_stable) {
613 // solution converged to an unstable point!
614 this->operability_status_.use_vfpexplicit = true;
615 auto bhp_stable = WellBhpThpCalculator(*this).estimateStableBhp(well_state, this->well_ecl_, rates, this->getRefDensity(), summary_state);
616 // if we find an intersection with a sufficiently lower bhp, re-solve equations
617 const Scalar reltol = 1e-3;
618 const Scalar cur_bhp = ws.bhp;
619 if (bhp_stable.has_value() && cur_bhp - bhp_stable.value() > cur_bhp*reltol){
620 const auto msg = fmt::format("Well {} converged to an unstable solution, re-solving", this->name());
621 deferred_logger.debug(msg);
622 solveWellWithBhp(simulator, dt, bhp_stable.value(), well_state, deferred_logger);
623 // re-solve with hopefully good initial guess
624 ws.thp = this->getTHPConstraint(summary_state);
625 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
626 }
627 }
628 }
629
630 if (!converged) {
631 // Well did not converge, switch to explicit fractions
632 this->operability_status_.use_vfpexplicit = true;
633 this->openWell();
634 auto bhp_target = estimateOperableBhp(simulator, dt, well_state, summary_state, deferred_logger);
635 if (!bhp_target.has_value()) {
636 // solve with zero rate
637 // well can't operate using explicit fractions stop the well
638 converged = solveWellWithZeroRate(simulator, dt, well_state, deferred_logger);
639 this->stopWell();
640 this->operability_status_.can_obtain_bhp_with_thp_limit = false;
641 this->operability_status_.obey_thp_limit_under_bhp_limit = false;
642 return converged;
643 } else {
644 // solve well with the estimated target bhp (or limit)
645 const Scalar bhp = std::max(bhp_target.value(),
646 static_cast<Scalar>(prod_controls.bhp_limit));
647 solveWellWithBhp(simulator, dt, bhp, well_state, deferred_logger);
648 ws.thp = this->getTHPConstraint(summary_state);
649 const auto msg = fmt::format("Well {} did not converge, re-solving with explicit fractions for VFP caculations.", this->name());
650 deferred_logger.debug(msg);
651 converged = this->iterateWellEqWithSwitching(simulator, dt,
652 inj_controls,
653 prod_controls,
654 well_state,
655 group_state,
656 deferred_logger);
657 }
658 }
659 // update operability
660 this->operability_status_.can_obtain_bhp_with_thp_limit = !this->wellIsStopped();
661 this->operability_status_.obey_thp_limit_under_bhp_limit = !this->wellIsStopped();
662 return converged;
663 }
664
665 template<typename TypeTag>
666 std::optional<typename WellInterface<TypeTag>::Scalar>
668 estimateOperableBhp(const Simulator& simulator,
669 const double dt,
670 WellStateType& well_state,
671 const SummaryState& summary_state,
672 DeferredLogger& deferred_logger)
673 {
674 if (!this->wellHasTHPConstraints(summary_state)) {
675 const Scalar bhp_limit = WellBhpThpCalculator(*this).mostStrictBhpFromBhpLimits(summary_state);
676 const bool converged = solveWellWithBhp(simulator, dt, bhp_limit, well_state, deferred_logger);
677 if (!converged || this->wellIsStopped()) {
678 return std::nullopt;
679 }
680
681 return bhp_limit;
682 }
683 OPM_TIMEFUNCTION();
684 // Given an unconverged well or closed well, estimate an operable bhp (if any)
685 // Get minimal bhp from vfp-curve
686 Scalar bhp_min = WellBhpThpCalculator(*this).calculateMinimumBhpFromThp(well_state, this->well_ecl_, summary_state, this->getRefDensity());
687 // Solve
688 const bool converged = solveWellWithBhp(simulator, dt, bhp_min, well_state, deferred_logger);
689 if (!converged || this->wellIsStopped()) {
690 return std::nullopt;
691 }
692 this->updateIPRImplicit(simulator, well_state, deferred_logger);
693 auto rates = well_state.well(this->index_of_well_).surface_rates;
694 this->adaptRatesForVFP(rates);
695 return WellBhpThpCalculator(*this).estimateStableBhp(well_state, this->well_ecl_, rates, this->getRefDensity(), summary_state);
696 }
697
698 template<typename TypeTag>
699 bool
701 solveWellWithBhp(const Simulator& simulator,
702 const double dt,
703 const Scalar bhp,
704 WellStateType& well_state,
705 DeferredLogger& deferred_logger)
706 {
707 OPM_TIMEFUNCTION();
708 // Solve a well using single bhp-constraint (but close if not operable under this)
709 auto group_state = GroupState<Scalar>(); // empty group
710 auto inj_controls = Well::InjectionControls(0);
711 auto prod_controls = Well::ProductionControls(0);
712 auto& ws = well_state.well(this->index_of_well_);
713 auto cmode_inj = ws.injection_cmode;
714 auto cmode_prod = ws.production_cmode;
715 if (this->isInjector()) {
716 inj_controls.addControl(Well::InjectorCMode::BHP);
717 inj_controls.bhp_limit = bhp;
718 inj_controls.cmode = Well::InjectorCMode::BHP;
719 ws.injection_cmode = Well::InjectorCMode::BHP;
720 } else {
721 prod_controls.addControl(Well::ProducerCMode::BHP);
722 prod_controls.bhp_limit = bhp;
723 prod_controls.cmode = Well::ProducerCMode::BHP;
724 ws.production_cmode = Well::ProducerCMode::BHP;
725 }
726 // update well-state
727 ws.bhp = bhp;
728 // solve
729 const bool converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger, /*fixed_control*/true);
730 ws.injection_cmode = cmode_inj;
731 ws.production_cmode = cmode_prod;
732 return converged;
733 }
734
735 template<typename TypeTag>
736 bool
738 solveWellWithZeroRate(const Simulator& simulator,
739 const double dt,
740 WellStateType& well_state,
741 DeferredLogger& deferred_logger)
742 {
743 OPM_TIMEFUNCTION();
744 // Solve a well as stopped
745 const auto well_status_orig = this->wellStatus_;
746 this->stopWell();
747
748 auto group_state = GroupState<Scalar>(); // empty group
749 auto inj_controls = Well::InjectionControls(0);
750 auto prod_controls = Well::ProductionControls(0);
751 const bool converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger, /*fixed_control*/true, /*fixed_status*/ true);
752 this->wellStatus_ = well_status_orig;
753 return converged;
754 }
755
756 template<typename TypeTag>
757 bool
759 solveWellForTesting(const Simulator& simulator,
760 WellStateType& well_state,
761 const GroupState<Scalar>& group_state,
762 DeferredLogger& deferred_logger)
763 {
764 OPM_TIMEFUNCTION();
765 const double dt = simulator.timeStepSize();
766
767 const auto& summary_state = simulator.vanguard().summaryState();
768 auto inj_controls = this->well_ecl_.isInjector() ? this->well_ecl_.injectionControls(summary_state) : Well::InjectionControls(0);
769 auto prod_controls = this->well_ecl_.isProducer() ? this->well_ecl_.productionControls(summary_state) : Well::ProductionControls(0);
770 this->onlyKeepBHPandTHPcontrols(summary_state, well_state, inj_controls, prod_controls);
771
772 bool converged = false;
773 try {
774 // TODO: the following two functions will be refactored to be one to reduce the code duplication
775 if (!this->param_.local_well_solver_control_switching_){
776 converged = this->iterateWellEqWithControl(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
777 } else {
778 if (this->param_.use_implicit_ipr_ && this->well_ecl_.isProducer() && (well_state.well(this->index_of_well_).status == WellStatus::OPEN)) {
779 converged = solveWellWithOperabilityCheck(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
780 } else {
781 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
782 }
783 }
784
785 } catch (NumericalProblem& e ) {
786 const std::string msg = "Inner well iterations failed for well " + this->name() + " Treat the well as unconverged. ";
787 deferred_logger.warning("INNER_ITERATION_FAILED", msg);
788 converged = false;
789 }
790
791 if (converged) {
792 deferred_logger.debug("WellTest: Well equation for well " + this->name() + " converged");
793 return true;
794 }
795 const int max_iter = this->param_.max_welleq_iter_;
796 deferred_logger.debug("WellTest: Well equation for well " + this->name() + " failed converging in "
797 + std::to_string(max_iter) + " iterations");
798 return false;
799 }
800
801
802 template<typename TypeTag>
803 void
805 solveWellEquation(const Simulator& simulator,
806 WellStateType& well_state,
807 const GroupState<Scalar>& group_state,
808 DeferredLogger& deferred_logger)
809 {
810 OPM_TIMEFUNCTION();
811 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
812 return;
813
814 // keep a copy of the original well state
815 const WellStateType well_state0 = well_state;
816 const double dt = simulator.timeStepSize();
817 bool converged = iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
818
819 // Newly opened wells with THP control sometimes struggles to
820 // converge due to bad initial guess. Or due to the simple fact
821 // that the well needs to change to another control.
822 // We therefore try to solve the well with BHP control to get
823 // an better initial guess.
824 // If the well is supposed to operate under THP control
825 // "updateWellControl" will switch it back to THP later.
826 if (!converged) {
827 auto& ws = well_state.well(this->indexOfWell());
828 bool thp_control = false;
829 if (this->well_ecl_.isInjector()) {
830 thp_control = ws.injection_cmode == Well::InjectorCMode::THP;
831 if (thp_control) {
832 ws.injection_cmode = Well::InjectorCMode::BHP;
833 if (this->well_control_log_.empty()) { // only log the first control
834 this->well_control_log_.push_back(WellInjectorCMode2String(Well::InjectorCMode::THP));
835 }
836 }
837 } else {
838 thp_control = ws.production_cmode == Well::ProducerCMode::THP;
839 if (thp_control) {
840 ws.production_cmode = Well::ProducerCMode::BHP;
841 if (this->well_control_log_.empty()) { // only log the first control
842 this->well_control_log_.push_back(WellProducerCMode2String(Well::ProducerCMode::THP));
843 }
844 }
845 }
846 if (thp_control) {
847 const std::string msg = std::string("The newly opened well ") + this->name()
848 + std::string(" with THP control did not converge during inner iterations, we try again with bhp control");
849 deferred_logger.debug(msg);
850 converged = this->iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
851 }
852 }
853
854 if (!converged) {
855 const int max_iter = this->param_.max_welleq_iter_;
856 deferred_logger.debug("Compute initial well solution for well " + this->name() + ". Failed to converge in "
857 + std::to_string(max_iter) + " iterations");
858 well_state = well_state0;
859 }
860 }
861
862
863
864 template <typename TypeTag>
865 void
867 assembleWellEq(const Simulator& simulator,
868 const double dt,
869 WellStateType& well_state,
870 const GroupState<Scalar>& group_state,
871 DeferredLogger& deferred_logger)
872 {
873 OPM_TIMEFUNCTION();
874 prepareWellBeforeAssembling(simulator, dt, well_state, group_state, deferred_logger);
875 assembleWellEqWithoutIteration(simulator, dt, well_state, group_state, deferred_logger);
876 }
877
878
879
880 template <typename TypeTag>
881 void
884 const double dt,
885 WellStateType& well_state,
886 const GroupState<Scalar>& group_state,
887 DeferredLogger& deferred_logger)
888 {
889 OPM_TIMEFUNCTION();
890 const auto& summary_state = simulator.vanguard().summaryState();
891 const auto inj_controls = this->well_ecl_.isInjector() ? this->well_ecl_.injectionControls(summary_state) : Well::InjectionControls(0);
892 const auto prod_controls = this->well_ecl_.isProducer() ? this->well_ecl_.productionControls(summary_state) : Well::ProductionControls(0);
893 // TODO: the reason to have inj_controls and prod_controls in the arguments, is that we want to change the control used for the well functions
894 // TODO: maybe we can use std::optional or pointers to simplify here
895 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
896 }
897
898
899
900 template<typename TypeTag>
901 void
904 const double dt,
905 WellStateType& well_state,
906 const GroupState<Scalar>& group_state,
907 DeferredLogger& deferred_logger)
908 {
909 OPM_TIMEFUNCTION();
910 const bool old_well_operable = this->operability_status_.isOperableAndSolvable();
911
912 if (this->param_.check_well_operability_iter_)
913 checkWellOperability(simulator, well_state, deferred_logger);
914
915 // only use inner well iterations for the first newton iterations.
916 const int iteration_idx = simulator.model().newtonMethod().numIterations();
917 if (iteration_idx < this->param_.max_niter_inner_well_iter_ || this->well_ecl_.isMultiSegment()) {
918 const auto& ws = well_state.well(this->indexOfWell());
919 const bool nonzero_rate_original =
920 std::any_of(ws.surface_rates.begin(),
921 ws.surface_rates.begin() + well_state.numPhases(),
922 [](Scalar rate) { return rate != Scalar(0.0); });
923
924 this->operability_status_.solvable = true;
925 if (number_of_well_reopenings_ >= this->param_.max_well_status_switch_) {
926 // only output the first time
927 if (number_of_well_reopenings_ == this->param_.max_well_status_switch_) {
928 const std::string msg = fmt::format("well {} is oscillating between open and stop. \n"
929 "We don't allow for more than {} re-openings "
930 "and the well is therefore kept stopped.",
931 this->name(), number_of_well_reopenings_);
932 deferred_logger.debug(msg);
933 }
934 this->stopWell();
935 changed_to_stopped_this_step_ = true;
936 bool converged_zero_rate = this->solveWellWithZeroRate(simulator, dt, well_state, deferred_logger);
937 if (this->param_.shut_unsolvable_wells_ && !converged_zero_rate ) {
938 this->operability_status_.solvable = false;
939 }
940 // we increse the number of reopenings to avoid output in the next iteration
941 number_of_well_reopenings_++;
942 return;
943 }
944 bool converged = this->iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
945
946 if (converged) {
947 const bool zero_target = this->wellUnderZeroRateTarget(simulator, well_state, deferred_logger);
948 if (this->wellIsStopped() && !zero_target && nonzero_rate_original) {
949 // Well had non-zero rate, but was stopped during local well-solve. We re-open the well
950 // for the next global iteration, but if the zero rate persists, it will be stopped.
951 // This logic is introduced to prevent/ameliorate stopped/revived oscillations
952 this->operability_status_.resetOperability();
953 this->openWell();
954 deferred_logger.debug(" " + this->name() + " is re-opened after being stopped during local solve");
955 number_of_well_reopenings_++;
956 }
957 } else {
958 // unsolvable wells are treated as not operable and will not be solved for in this iteration.
959 if (this->param_.shut_unsolvable_wells_) {
960 this->operability_status_.solvable = false;
961 }
962 }
963 }
964 if (this->operability_status_.has_negative_potentials) {
965 auto well_state_copy = well_state;
966 std::vector<Scalar> potentials;
967 try {
968 computeWellPotentials(simulator, well_state_copy, potentials, deferred_logger);
969 } catch (const std::exception& e) {
970 const std::string msg = fmt::format("well {}: computeWellPotentials() failed "
971 "during attempt to recompute potentials for well: ",
972 this->name(), e.what());
973 deferred_logger.info(msg);
974 this->operability_status_.has_negative_potentials = true;
975 }
976 auto& ws = well_state.well(this->indexOfWell());
977 const int np = well_state.numPhases();
978 for (int p = 0; p < np; ++p) {
979 ws.well_potentials[p] = std::max(Scalar{0.0}, potentials[p]);
980 }
981 }
982 this->changed_to_open_this_step_ = false;
983 changed_to_stopped_this_step_ = false;
984
985 const bool well_operable = this->operability_status_.isOperableAndSolvable();
986 if (!well_operable) {
987 this->stopWell();
988 if (old_well_operable) {
989 deferred_logger.debug(" well " + this->name() + " gets STOPPED during iteration ");
990 changed_to_stopped_this_step_ = true;
991 }
992 } else if (well_state.isOpen(this->name())) {
993 this->openWell();
994 if (!old_well_operable) {
995 deferred_logger.debug(" well " + this->name() + " gets REVIVED during iteration ");
996 this->changed_to_open_this_step_ = true;
997 }
998 }
999 }
1000
1001 template<typename TypeTag>
1002 void
1004 {
1005 if(!this->isOperableAndSolvable() && !this->wellIsStopped())
1006 return;
1007
1008 for (int perfIdx = 0; perfIdx < this->number_of_local_perforations_; ++perfIdx) {
1009 if (this->cells()[perfIdx] == cellIdx) {
1010 for (int i = 0; i < RateVector::dimension; ++i) {
1011 rates[i] += connectionRates_[perfIdx][i];
1012 }
1013 }
1014 }
1015 }
1016
1017 template<typename TypeTag>
1020 {
1021 for (int perfIdx = 0; perfIdx < this->number_of_local_perforations_; ++perfIdx) {
1022 if (this->cells()[perfIdx] == cellIdx) {
1023 const unsigned activeCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::solventComponentIndex(phaseIdx));
1024 return connectionRates_[perfIdx][activeCompIdx].value();
1025 }
1026 }
1027 // this is not thread safe
1028 OPM_THROW(std::invalid_argument, "The well with name " + this->name()
1029 + " does not perforate cell " + std::to_string(cellIdx));
1030 return 0.0;
1031 }
1032
1033
1034
1035
1036 template<typename TypeTag>
1037 void
1039 checkWellOperability(const Simulator& simulator,
1040 const WellStateType& well_state,
1041 DeferredLogger& deferred_logger)
1042 {
1043 OPM_TIMEFUNCTION();
1044 if (!this->param_.check_well_operability_) {
1045 return;
1046 }
1047
1048 if (this->wellIsStopped() && !changed_to_stopped_this_step_) {
1049 return;
1050 }
1051
1052 updateWellOperability(simulator, well_state, deferred_logger);
1053 if (!this->operability_status_.isOperableAndSolvable()) {
1054 this->operability_status_.use_vfpexplicit = true;
1055 deferred_logger.debug("EXPLICIT_LOOKUP_VFP",
1056 "well not operable, trying with explicit vfp lookup: " + this->name());
1057 updateWellOperability(simulator, well_state, deferred_logger);
1058 }
1059 }
1060
1061
1062
1063 template<typename TypeTag>
1064 void
1067 WellStateType& well_state,
1068 const GroupState<Scalar>& group_state,
1069 GLiftEclWells& ecl_well_map,
1070 DeferredLogger& deferred_logger)
1071 {
1072 OPM_TIMEFUNCTION();
1073 const auto& summary_state = simulator.vanguard().summaryState();
1074 const auto& well_name = this->name();
1075 if (!this->wellHasTHPConstraints(summary_state)) {
1076 const std::string msg = fmt::format("GLIFT WTEST: Well {} does not have THP constraints", well_name);
1077 deferred_logger.info(msg);
1078 return;
1079 }
1080 const auto& schedule = simulator.vanguard().schedule();
1081 const auto report_step_idx = simulator.episodeIndex();
1082 const auto& glo = schedule.glo(report_step_idx);
1083 if (!glo.has_well(well_name)) {
1084 const std::string msg = fmt::format(
1085 "GLIFT WTEST: Well {} : Gas lift not activated: "
1086 "WLIFTOPT is probably missing. Skipping.", well_name);
1087 deferred_logger.info(msg);
1088 return;
1089 }
1090 const auto& gl_well = glo.well(well_name);
1091
1092 // Use gas lift optimization to get ALQ for well test
1093 std::unique_ptr<GasLiftSingleWell> glift =
1094 initializeGliftWellTest_<GasLiftSingleWell>(simulator,
1095 well_state,
1096 group_state,
1097 ecl_well_map,
1098 deferred_logger);
1099 auto [wtest_alq, success] = glift->wellTestALQ();
1100 std::string msg;
1101 const auto& unit_system = schedule.getUnits();
1102 if (success) {
1103 well_state.well(well_name).alq_state.set(wtest_alq);
1104 msg = fmt::format(
1105 "GLIFT WTEST: Well {} : Setting ALQ to optimized value = {}",
1106 well_name, unit_system.from_si(UnitSystem::measure::gas_surface_rate, wtest_alq));
1107 }
1108 else {
1109 if (!gl_well.use_glo()) {
1110 msg = fmt::format(
1111 "GLIFT WTEST: Well {} : Gas lift optimization deactivated. Setting ALQ to WLIFTOPT item 3 = {}",
1112 well_name,
1113 unit_system.from_si(UnitSystem::measure::gas_surface_rate, well_state.well(well_name).alq_state.get()));
1114
1115 }
1116 else {
1117 msg = fmt::format(
1118 "GLIFT WTEST: Well {} : Gas lift optimization failed, no ALQ set.",
1119 well_name);
1120 }
1121 }
1122 deferred_logger.info(msg);
1123 }
1124
1125 template<typename TypeTag>
1126 void
1128 updateWellOperability(const Simulator& simulator,
1129 const WellStateType& well_state,
1130 DeferredLogger& deferred_logger)
1131 {
1132 OPM_TIMEFUNCTION();
1133 if (this->param_.local_well_solver_control_switching_) {
1134 const bool success = updateWellOperabilityFromWellEq(simulator, well_state, deferred_logger);
1135 if (!success) {
1136 this->operability_status_.solvable = false;
1137 deferred_logger.debug("Operability check using well equations did not converge for well "
1138 + this->name() + ". Mark the well as unsolvable." );
1139 }
1140 return;
1141 }
1142 this->operability_status_.resetOperability();
1143
1144 bool thp_controlled = this->isInjector() ? well_state.well(this->index_of_well_).injection_cmode == Well::InjectorCMode::THP:
1145 well_state.well(this->index_of_well_).production_cmode == Well::ProducerCMode::THP;
1146 bool bhp_controlled = this->isInjector() ? well_state.well(this->index_of_well_).injection_cmode == Well::InjectorCMode::BHP:
1147 well_state.well(this->index_of_well_).production_cmode == Well::ProducerCMode::BHP;
1148
1149 // Operability checking is not free
1150 // Only check wells under BHP and THP control
1151 bool check_thp = thp_controlled || this->operability_status_.thp_limit_violated_but_not_switched;
1152 if (check_thp || bhp_controlled) {
1153 updateIPR(simulator, deferred_logger);
1154 checkOperabilityUnderBHPLimit(well_state, simulator, deferred_logger);
1155 }
1156 // we do some extra checking for wells under THP control.
1157 if (check_thp) {
1158 checkOperabilityUnderTHPLimit(simulator, well_state, deferred_logger);
1159 }
1160 }
1161
1162 template<typename TypeTag>
1163 bool
1166 const WellStateType& well_state,
1167 DeferredLogger& deferred_logger)
1168 {
1169 OPM_TIMEFUNCTION();
1170 // only makes sense if we're using this parameter is true
1171 assert(this->param_.local_well_solver_control_switching_);
1172 this->operability_status_.resetOperability();
1173 WellStateType well_state_copy = well_state;
1174 const auto& group_state = simulator.problem().wellModel().groupState();
1175 const double dt = simulator.timeStepSize();
1176 // equations should be converged at this stage, so only one it is needed
1177 bool converged = iterateWellEquations(simulator, dt, well_state_copy, group_state, deferred_logger);
1178 return converged;
1179 }
1180
1181 template<typename TypeTag>
1182 void
1184 scaleSegmentRatesAndPressure([[maybe_unused]] WellStateType& well_state) const
1185 {
1186 // only relevant for MSW
1187 }
1188
1189 template<typename TypeTag>
1190 void
1192 updateWellStateWithTarget(const Simulator& simulator,
1193 const GroupState<Scalar>& group_state,
1194 WellStateType& well_state,
1195 DeferredLogger& deferred_logger) const
1196 {
1197 OPM_TIMEFUNCTION();
1198 // only bhp and wellRates are used to initilize the primaryvariables for standard wells
1199 const auto& well = this->well_ecl_;
1200 const int well_index = this->index_of_well_;
1201 auto& ws = well_state.well(well_index);
1202 const int np = well_state.numPhases();
1203 const auto& summaryState = simulator.vanguard().summaryState();
1204 const auto& schedule = simulator.vanguard().schedule();
1205
1206 // Discard old primary variables, the new well state
1207 // may not be anywhere near the old one.
1208 ws.primaryvar.resize(0);
1209
1210 if (this->wellIsStopped()) {
1211 for (int p = 0; p<np; ++p) {
1212 ws.surface_rates[p] = 0;
1213 }
1214 ws.thp = 0;
1215 return;
1216 }
1217
1218 if (this->isInjector() )
1219 {
1220 const auto& controls = well.injectionControls(summaryState);
1221
1222 InjectorType injectorType = controls.injector_type;
1223 int phasePos;
1224 switch (injectorType) {
1225 case InjectorType::WATER:
1226 {
1227 phasePos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::waterPhaseIdx);
1228 break;
1229 }
1230 case InjectorType::OIL:
1231 {
1232 phasePos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
1233 break;
1234 }
1235 case InjectorType::GAS:
1236 {
1237 phasePos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
1238 break;
1239 }
1240 default:
1241 OPM_DEFLOG_THROW(std::runtime_error, "Expected WATER, OIL or GAS as type for injectors " + this->name(), deferred_logger );
1242 }
1243
1244 const auto current = ws.injection_cmode;
1245
1246 switch (current) {
1247 case Well::InjectorCMode::RATE:
1248 {
1249 ws.surface_rates[phasePos] = (1.0 - this->rsRvInj()) * controls.surface_rate;
1250 if(this->rsRvInj() > 0) {
1251 if (injectorType == InjectorType::OIL && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1252 const int gas_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
1253 ws.surface_rates[gas_pos] = controls.surface_rate * this->rsRvInj();
1254 } else if (injectorType == InjectorType::GAS && FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
1255 const int oil_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
1256 ws.surface_rates[oil_pos] = controls.surface_rate * this->rsRvInj();
1257 } else {
1258 OPM_DEFLOG_THROW(std::runtime_error, "Expected OIL or GAS as type for injectors when RS/RV (item 10) is non-zero " + this->name(), deferred_logger );
1259 }
1260 }
1261 break;
1262 }
1263
1264 case Well::InjectorCMode::RESV:
1265 {
1266 std::vector<Scalar> convert_coeff(this->number_of_phases_, 1.0);
1267 this->rateConverter_.calcCoeff(/*fipreg*/ 0, this->pvtRegionIdx_, convert_coeff);
1268 const Scalar coeff = convert_coeff[phasePos];
1269 ws.surface_rates[phasePos] = controls.reservoir_rate/coeff;
1270 break;
1271 }
1272
1273 case Well::InjectorCMode::THP:
1274 {
1275 auto rates = ws.surface_rates;
1276 Scalar bhp = WellBhpThpCalculator(*this).calculateBhpFromThp(well_state,
1277 rates,
1278 well,
1279 summaryState,
1280 this->getRefDensity(),
1281 deferred_logger);
1282 ws.bhp = bhp;
1283 ws.thp = this->getTHPConstraint(summaryState);
1284
1285 // if the total rates are negative or zero
1286 // we try to provide a better intial well rate
1287 // using the well potentials
1288 Scalar total_rate = std::accumulate(rates.begin(), rates.end(), 0.0);
1289 if (total_rate <= 0.0)
1290 ws.surface_rates = ws.well_potentials;
1291
1292 break;
1293 }
1294 case Well::InjectorCMode::BHP:
1295 {
1296 ws.bhp = controls.bhp_limit;
1297 Scalar total_rate = 0.0;
1298 for (int p = 0; p<np; ++p) {
1299 total_rate += ws.surface_rates[p];
1300 }
1301 // if the total rates are negative or zero
1302 // we try to provide a better intial well rate
1303 // using the well potentials
1304 if (total_rate <= 0.0)
1305 ws.surface_rates = ws.well_potentials;
1306
1307 break;
1308 }
1309 case Well::InjectorCMode::GRUP:
1310 {
1311 assert(well.isAvailableForGroupControl());
1312 const auto& group = schedule.getGroup(well.groupName(), this->currentStep());
1313 const Scalar efficiencyFactor = well.getEfficiencyFactor() *
1314 well_state[well.name()].efficiency_scaling_factor;
1315 std::optional<Scalar> target =
1316 this->getGroupInjectionTargetRate(group,
1317 well_state,
1318 group_state,
1319 schedule,
1320 summaryState,
1321 injectorType,
1322 efficiencyFactor,
1323 deferred_logger);
1324 if (target)
1325 ws.surface_rates[phasePos] = *target;
1326 break;
1327 }
1328 case Well::InjectorCMode::CMODE_UNDEFINED:
1329 {
1330 OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + this->name(), deferred_logger );
1331 }
1332
1333 }
1334 // for wells with zero injection rate, if we assign exactly zero rate,
1335 // we will have to assume some trivial composition in the wellbore.
1336 // here, we use some small value (about 0.01 m^3/day ~= 1.e-7) to initialize
1337 // the zero rate target, then we can use to retain the composition information
1338 // within the wellbore from the previous result, and hopefully it is a good
1339 // initial guess for the zero rate target.
1340 ws.surface_rates[phasePos] = std::max(Scalar{1.e-7}, ws.surface_rates[phasePos]);
1341
1342 if (ws.bhp == 0.) {
1343 ws.bhp = controls.bhp_limit;
1344 }
1345 }
1346 //Producer
1347 else
1348 {
1349 const auto current = ws.production_cmode;
1350 const auto& controls = well.productionControls(summaryState);
1351 switch (current) {
1352 case Well::ProducerCMode::ORAT:
1353 {
1354 const int oil_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
1355 Scalar current_rate = -ws.surface_rates[oil_pos];
1356 // for trivial rates or opposite direction we don't just scale the rates
1357 // but use either the potentials or the mobility ratio to initial the well rates
1358 if (current_rate > 0.0) {
1359 for (int p = 0; p<np; ++p) {
1360 ws.surface_rates[p] *= controls.oil_rate/current_rate;
1361 }
1362 } else {
1363 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1364 double control_fraction = fractions[oil_pos];
1365 if (control_fraction != 0.0) {
1366 for (int p = 0; p<np; ++p) {
1367 ws.surface_rates[p] = - fractions[p] * controls.oil_rate/control_fraction;
1368 }
1369 }
1370 }
1371 break;
1372 }
1373 case Well::ProducerCMode::WRAT:
1374 {
1375 const int water_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::waterPhaseIdx);
1376 Scalar current_rate = -ws.surface_rates[water_pos];
1377 // for trivial rates or opposite direction we don't just scale the rates
1378 // but use either the potentials or the mobility ratio to initial the well rates
1379 if (current_rate > 0.0) {
1380 for (int p = 0; p<np; ++p) {
1381 ws.surface_rates[p] *= controls.water_rate/current_rate;
1382 }
1383 } else {
1384 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1385 const Scalar control_fraction = fractions[water_pos];
1386 if (control_fraction != 0.0) {
1387 for (int p = 0; p<np; ++p) {
1388 ws.surface_rates[p] = - fractions[p] * controls.water_rate / control_fraction;
1389 }
1390 }
1391 }
1392 break;
1393 }
1394 case Well::ProducerCMode::GRAT:
1395 {
1396 const int gas_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
1397 Scalar current_rate = -ws.surface_rates[gas_pos];
1398 // or trivial rates or opposite direction we don't just scale the rates
1399 // but use either the potentials or the mobility ratio to initial the well rates
1400 if (current_rate > 0.0) {
1401 for (int p = 0; p<np; ++p) {
1402 ws.surface_rates[p] *= controls.gas_rate/current_rate;
1403 }
1404 } else {
1405 const std::vector<Scalar > fractions = initialWellRateFractions(simulator, well_state);
1406 const Scalar control_fraction = fractions[gas_pos];
1407 if (control_fraction != 0.0) {
1408 for (int p = 0; p<np; ++p) {
1409 ws.surface_rates[p] = - fractions[p] * controls.gas_rate / control_fraction;
1410 }
1411 }
1412 }
1413
1414 break;
1415
1416 }
1417 case Well::ProducerCMode::LRAT:
1418 {
1419 const int water_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::waterPhaseIdx);
1420 const int oil_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
1421 Scalar current_rate = - ws.surface_rates[water_pos]
1422 - ws.surface_rates[oil_pos];
1423 // or trivial rates or opposite direction we don't just scale the rates
1424 // but use either the potentials or the mobility ratio to initial the well rates
1425 if (current_rate > 0.0) {
1426 for (int p = 0; p<np; ++p) {
1427 ws.surface_rates[p] *= controls.liquid_rate/current_rate;
1428 }
1429 } else {
1430 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1431 const Scalar control_fraction = fractions[water_pos] + fractions[oil_pos];
1432 if (control_fraction != 0.0) {
1433 for (int p = 0; p<np; ++p) {
1434 ws.surface_rates[p] = - fractions[p] * controls.liquid_rate / control_fraction;
1435 }
1436 }
1437 }
1438 break;
1439 }
1440 case Well::ProducerCMode::CRAT:
1441 {
1442 OPM_DEFLOG_THROW(std::runtime_error,
1443 fmt::format("CRAT control not supported, well {}", this->name()),
1444 deferred_logger);
1445 }
1446 case Well::ProducerCMode::RESV:
1447 {
1448 std::vector<Scalar> convert_coeff(this->number_of_phases_, 1.0);
1449 this->rateConverter_.calcCoeff(/*fipreg*/ 0, this->pvtRegionIdx_, ws.surface_rates, convert_coeff);
1450 Scalar total_res_rate = 0.0;
1451 for (int p = 0; p<np; ++p) {
1452 total_res_rate -= ws.surface_rates[p] * convert_coeff[p];
1453 }
1454 if (controls.prediction_mode) {
1455 // or trivial rates or opposite direction we don't just scale the rates
1456 // but use either the potentials or the mobility ratio to initial the well rates
1457 if (total_res_rate > 0.0) {
1458 for (int p = 0; p<np; ++p) {
1459 ws.surface_rates[p] *= controls.resv_rate/total_res_rate;
1460 }
1461 } else {
1462 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1463 for (int p = 0; p<np; ++p) {
1464 ws.surface_rates[p] = - fractions[p] * controls.resv_rate / convert_coeff[p];
1465 }
1466 }
1467 } else {
1468 std::vector<Scalar> hrates(this->number_of_phases_,0.);
1469 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
1470 const int phase_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::waterPhaseIdx);
1471 hrates[phase_pos] = controls.water_rate;
1472 }
1473 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
1474 const int phase_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
1475 hrates[phase_pos] = controls.oil_rate;
1476 }
1477 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1478 const int phase_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
1479 hrates[phase_pos] = controls.gas_rate;
1480 }
1481 std::vector<Scalar> hrates_resv(this->number_of_phases_,0.);
1482 this->rateConverter_.calcReservoirVoidageRates(/*fipreg*/ 0, this->pvtRegionIdx_, hrates, hrates_resv);
1483 Scalar target = std::accumulate(hrates_resv.begin(), hrates_resv.end(), 0.0);
1484 // or trivial rates or opposite direction we don't just scale the rates
1485 // but use either the potentials or the mobility ratio to initial the well rates
1486 if (total_res_rate > 0.0) {
1487 for (int p = 0; p<np; ++p) {
1488 ws.surface_rates[p] *= target/total_res_rate;
1489 }
1490 } else {
1491 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1492 for (int p = 0; p<np; ++p) {
1493 ws.surface_rates[p] = - fractions[p] * target / convert_coeff[p];
1494 }
1495 }
1496 }
1497 break;
1498 }
1499 case Well::ProducerCMode::BHP:
1500 {
1501 ws.bhp = controls.bhp_limit;
1502 Scalar total_rate = 0.0;
1503 for (int p = 0; p<np; ++p) {
1504 total_rate -= ws.surface_rates[p];
1505 }
1506 // if the total rates are negative or zero
1507 // we try to provide a better intial well rate
1508 // using the well potentials
1509 if (total_rate <= 0.0){
1510 for (int p = 0; p<np; ++p) {
1511 ws.surface_rates[p] = -ws.well_potentials[p];
1512 }
1513 }
1514 break;
1515 }
1516 case Well::ProducerCMode::THP:
1517 {
1518 const bool update_success = updateWellStateWithTHPTargetProd(simulator, well_state, deferred_logger);
1519
1520 if (!update_success) {
1521 // the following is the original way of initializing well state with THP constraint
1522 // keeping it for robust reason in case that it fails to get a bhp value with THP constraint
1523 // more sophisticated design might be needed in the future
1524 auto rates = ws.surface_rates;
1525 this->adaptRatesForVFP(rates);
1527 well_state, rates, well, summaryState, this->getRefDensity(), deferred_logger);
1528 ws.bhp = bhp;
1529 ws.thp = this->getTHPConstraint(summaryState);
1530 // if the total rates are negative or zero
1531 // we try to provide a better initial well rate
1532 // using the well potentials
1533 const Scalar total_rate = -std::accumulate(rates.begin(), rates.end(), 0.0);
1534 if (total_rate <= 0.0) {
1535 for (int p = 0; p < this->number_of_phases_; ++p) {
1536 ws.surface_rates[p] = -ws.well_potentials[p];
1537 }
1538 }
1539 }
1540 break;
1541 }
1542 case Well::ProducerCMode::GRUP:
1543 {
1544 assert(well.isAvailableForGroupControl());
1545 const auto& group = schedule.getGroup(well.groupName(), this->currentStep());
1546 const Scalar efficiencyFactor = well.getEfficiencyFactor() *
1547 well_state[well.name()].efficiency_scaling_factor;
1548 Scalar scale = this->getGroupProductionTargetRate(group,
1549 well_state,
1550 group_state,
1551 schedule,
1552 summaryState,
1553 efficiencyFactor,
1554 deferred_logger);
1555
1556 // we don't want to scale with zero and get zero rates.
1557 if (scale > 0) {
1558 for (int p = 0; p<np; ++p) {
1559 ws.surface_rates[p] *= scale;
1560 }
1561 ws.trivial_group_target = false;
1562 } else {
1563 // If group target is trivial we dont want to flip to other controls. To avoid oscillation we store
1564 // this information in the well state and explicitly check for this condition when evaluating well controls.
1565 ws.trivial_group_target = true;
1566 }
1567 break;
1568 }
1569 case Well::ProducerCMode::CMODE_UNDEFINED:
1571 {
1572 OPM_DEFLOG_THROW(std::runtime_error, "Well control must be specified for well " + this->name() , deferred_logger);
1573 break;
1574 }
1575 } // end of switch
1576
1577 if (ws.bhp == 0.) {
1578 ws.bhp = controls.bhp_limit;
1579 }
1580 }
1581 }
1582
1583 template<typename TypeTag>
1584 bool
1586 wellUnderZeroRateTarget(const Simulator& simulator,
1587 const WellStateType& well_state,
1588 DeferredLogger& deferred_logger) const
1589 {
1590 OPM_TIMEFUNCTION();
1591 // Check if well is under zero rate control, either directly or from group
1592 const bool isGroupControlled = this->wellUnderGroupControl(well_state.well(this->index_of_well_));
1593 if (!isGroupControlled) {
1594 // well is not under group control, check "individual" version
1595 const auto& summaryState = simulator.vanguard().summaryState();
1596 return this->wellUnderZeroRateTargetIndividual(summaryState, well_state);
1597 } else {
1598 return this->wellUnderZeroGroupRateTarget(simulator, well_state, deferred_logger, isGroupControlled);
1599 }
1600 }
1601
1602 template <typename TypeTag>
1603 bool
1605 const WellStateType& well_state,
1606 DeferredLogger& deferred_logger,
1607 const std::optional<bool> group_control) const
1608 {
1609 // Check if well is under zero rate target from group
1610 const bool isGroupControlled = group_control.value_or(this->wellUnderGroupControl(well_state.well(this->index_of_well_)));
1611 if (isGroupControlled) {
1612 const auto& summaryState = simulator.vanguard().summaryState();
1613 const auto& group_state = simulator.problem().wellModel().groupState();
1614 const auto& schedule = simulator.vanguard().schedule();
1615 return this->zeroGroupRateTarget(summaryState, schedule, well_state, group_state, deferred_logger);
1616 }
1617 return false;
1618 }
1619
1620 template<typename TypeTag>
1621 bool
1623 stoppedOrZeroRateTarget(const Simulator& simulator,
1624 const WellStateType& well_state,
1625 DeferredLogger& deferred_logger) const
1626 {
1627 // Check if well is stopped or under zero rate control, either
1628 // directly or from group.
1629 return this->wellIsStopped()
1630 || this->wellUnderZeroRateTarget(simulator, well_state, deferred_logger);
1631 }
1632
1633 template<typename TypeTag>
1634 std::vector<typename WellInterface<TypeTag>::Scalar>
1636 initialWellRateFractions(const Simulator& simulator,
1637 const WellStateType& well_state) const
1638 {
1639 OPM_TIMEFUNCTION();
1640 const int np = this->number_of_phases_;
1641 std::vector<Scalar> scaling_factor(np);
1642 const auto& ws = well_state.well(this->index_of_well_);
1643
1644 Scalar total_potentials = 0.0;
1645 for (int p = 0; p<np; ++p) {
1646 total_potentials += ws.well_potentials[p];
1647 }
1648 if (total_potentials > 0) {
1649 for (int p = 0; p<np; ++p) {
1650 scaling_factor[p] = ws.well_potentials[p] / total_potentials;
1651 }
1652 return scaling_factor;
1653 }
1654 // if we don't have any potentials we weight it using the mobilites
1655 // We only need approximation so we don't bother with the vapporized oil and dissolved gas
1656 Scalar total_tw = 0;
1657 const int nperf = this->number_of_local_perforations_;
1658 for (int perf = 0; perf < nperf; ++perf) {
1659 total_tw += this->well_index_[perf];
1660 }
1661 total_tw = this->parallelWellInfo().communication().sum(total_tw);
1662
1663 for (int perf = 0; perf < nperf; ++perf) {
1664 const int cell_idx = this->well_cells_[perf];
1665 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/0);
1666 const auto& fs = intQuants.fluidState();
1667 const Scalar well_tw_fraction = this->well_index_[perf] / total_tw;
1668 Scalar total_mobility = 0.0;
1669 for (int p = 0; p < np; ++p) {
1670 const int canonical_phase_idx = FluidSystem::activeToCanonicalPhaseIdx(p);
1671 total_mobility += fs.invB(canonical_phase_idx).value() * intQuants.mobility(canonical_phase_idx).value();
1672 }
1673 for (int p = 0; p < np; ++p) {
1674 const int canonical_phase_idx = FluidSystem::activeToCanonicalPhaseIdx(p);
1675 scaling_factor[p] += well_tw_fraction * fs.invB(canonical_phase_idx).value() * intQuants.mobility(canonical_phase_idx).value() / total_mobility;
1676 }
1677 }
1678 return scaling_factor;
1679 }
1680
1681
1682
1683 template <typename TypeTag>
1684 void
1687 WellStateType& well_state,
1688 DeferredLogger& deferred_logger) const
1689 {
1690 assert(this->isProducer());
1691 OPM_TIMEFUNCTION();
1692 // Check if the rates of this well only are single-phase, do nothing
1693 // if more than one nonzero rate.
1694 auto& ws = well_state.well(this->index_of_well_);
1695 int nonzero_rate_index = -1;
1696 const Scalar floating_point_error_epsilon = 1e-14;
1697 for (int p = 0; p < this->number_of_phases_; ++p) {
1698 if (std::abs(ws.surface_rates[p]) > floating_point_error_epsilon) {
1699 if (nonzero_rate_index == -1) {
1700 nonzero_rate_index = p;
1701 } else {
1702 // More than one nonzero rate.
1703 return;
1704 }
1705 }
1706 }
1707
1708 // Calculate rates at bhp limit, or 1 bar if no limit.
1709 std::vector<Scalar> well_q_s(this->number_of_phases_, 0.0);
1710 bool rates_evaluated_at_1bar = false;
1711 {
1712 const auto& summary_state = simulator.vanguard().summaryState();
1713 const auto& prod_controls = this->well_ecl_.productionControls(summary_state);
1714 const double bhp_limit = std::max(prod_controls.bhp_limit, 1.0 * unit::barsa);
1715 this->computeWellRatesWithBhp(simulator, bhp_limit, well_q_s, deferred_logger);
1716 // Remember of we evaluated the rates at (approx.) 1 bar or not.
1717 rates_evaluated_at_1bar = (bhp_limit < 1.1 * unit::barsa);
1718 // Check that no rates are positive.
1719 if (std::any_of(well_q_s.begin(), well_q_s.end(), [](Scalar q) { return q > 0.0; })) {
1720 // Did we evaluate at 1 bar? If not, then we can try again at 1 bar.
1721 if (!rates_evaluated_at_1bar) {
1722 this->computeWellRatesWithBhp(simulator, 1.0 * unit::barsa, well_q_s, deferred_logger);
1723 rates_evaluated_at_1bar = true;
1724 }
1725 // At this point we can only set the wrong-direction (if any) values to zero.
1726 for (auto& q : well_q_s) {
1727 q = std::min(q, Scalar{0.0});
1728 }
1729 }
1730 }
1731
1732 if (nonzero_rate_index == -1) {
1733 // No nonzero rates on input.
1734 // Use the computed rate directly, or scaled by a factor
1735 // 0.5 (to avoid too high values) if it was evaluated at 1 bar.
1736 const Scalar factor = rates_evaluated_at_1bar ? 0.5 : 1.0;
1737 for (int p = 0; p < this->number_of_phases_; ++p) {
1738 ws.surface_rates[p] = factor * well_q_s[p];
1739 }
1740 return;
1741 }
1742
1743 // If we are here, we had a single nonzero rate for the well,
1744 // typically from a rate constraint. We must make sure it is
1745 // respected, so if it was lower than the calculated rate for
1746 // the same phase we scale all rates to match.
1747 const Scalar initial_nonzero_rate = ws.surface_rates[nonzero_rate_index];
1748 const Scalar computed_rate = well_q_s[nonzero_rate_index];
1749 if (std::abs(initial_nonzero_rate) < std::abs(computed_rate)) {
1750 // Note that both rates below are negative. The factor should be < 1.0.
1751 const Scalar factor = initial_nonzero_rate / computed_rate;
1752 assert(factor < 1.0);
1753 for (int p = 0; p < this->number_of_phases_; ++p) {
1754 // We skip the nonzero_rate_index, as that should remain as it was.
1755 if (p != nonzero_rate_index) {
1756 ws.surface_rates[p] = factor * well_q_s[p];
1757 }
1758 }
1759 return;
1760 }
1761
1762 // If we are here, we had a single nonzero rate, but it was
1763 // higher than the one calculated from the bhp limit, so we
1764 // use the calculated rates.
1765 for (int p = 0; p < this->number_of_phases_; ++p) {
1766 ws.surface_rates[p] = well_q_s[p];
1767 }
1768 }
1769
1770 template <typename TypeTag>
1771 std::vector<typename WellInterface<TypeTag>::Scalar>
1773 wellIndex(const int perf,
1774 const IntensiveQuantities& intQuants,
1775 const Scalar trans_mult,
1776 const SingleWellStateType& ws) const
1777 {
1778 OPM_TIMEFUNCTION_LOCAL(Subsystem::Wells);
1779 // Add a Forchheimer term to the gas phase CTF if the run uses
1780 // either of the WDFAC or the WDFACCOR keywords.
1781 if (static_cast<std::size_t>(perf) >= this->well_cells_.size()) {
1782 OPM_THROW(std::invalid_argument,"The perforation index exceeds the size of the local containers - possibly wellIndex was called with a global instead of a local perforation index!");
1783 }
1784 auto wi = std::vector<Scalar>
1785 (this->num_conservation_quantities_, this->well_index_[perf] * trans_mult);
1786
1787 if constexpr (! Indices::gasEnabled) {
1788 return wi;
1789 }
1790
1791 const auto& wdfac = this->well_ecl_.getWDFAC();
1792
1793 if (! wdfac.useDFactor() || (this->well_index_[perf] == 0.0)) {
1794 return wi;
1795 }
1796
1797 const Scalar d = this->computeConnectionDFactor(perf, intQuants, ws);
1798 if (d < 1.0e-15) {
1799 return wi;
1800 }
1801
1802 // Solve quadratic equations for connection rates satisfying the ipr and the flow-dependent skin.
1803 // If more than one solution, pick the one corresponding to lowest absolute rate (smallest skin).
1804 const auto& connection = this->well_ecl_.getConnections()[ws.perf_data.ecl_index[perf]];
1805 const Scalar Kh = connection.Kh();
1806 const Scalar scaling = 3.141592653589 * Kh * connection.wpimult();
1807 const unsigned gas_comp_idx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx);
1808
1809 const Scalar connection_pressure = ws.perf_data.pressure[perf];
1810 const Scalar cell_pressure = getValue(intQuants.fluidState().pressure(FluidSystem::gasPhaseIdx));
1811 const Scalar drawdown = cell_pressure - connection_pressure;
1812 const Scalar invB = getValue(intQuants.fluidState().invB(FluidSystem::gasPhaseIdx));
1813 const Scalar mob_g = getValue(intQuants.mobility(FluidSystem::gasPhaseIdx)) * invB;
1814 const Scalar a = d;
1815 const Scalar b = 2*scaling/wi[gas_comp_idx];
1816 const Scalar c = -2*scaling*mob_g*drawdown;
1817
1818 Scalar consistent_Q = -1.0e20;
1819 // Find and check negative solutions (a --> -a)
1820 const Scalar r2n = b*b + 4*a*c;
1821 if (r2n >= 0) {
1822 const Scalar rn = std::sqrt(r2n);
1823 const Scalar xn1 = (b-rn)*0.5/a;
1824 if (xn1 <= 0) {
1825 consistent_Q = xn1;
1826 }
1827 const Scalar xn2 = (b+rn)*0.5/a;
1828 if (xn2 <= 0 && xn2 > consistent_Q) {
1829 consistent_Q = xn2;
1830 }
1831 }
1832 // Find and check positive solutions
1833 consistent_Q *= -1;
1834 const Scalar r2p = b*b - 4*a*c;
1835 if (r2p >= 0) {
1836 const Scalar rp = std::sqrt(r2p);
1837 const Scalar xp1 = (rp-b)*0.5/a;
1838 if (xp1 > 0 && xp1 < consistent_Q) {
1839 consistent_Q = xp1;
1840 }
1841 const Scalar xp2 = -(rp+b)*0.5/a;
1842 if (xp2 > 0 && xp2 < consistent_Q) {
1843 consistent_Q = xp2;
1844 }
1845 }
1846 wi[gas_comp_idx] = 1.0/(1.0/(trans_mult * this->well_index_[perf]) + (consistent_Q/2 * d / scaling));
1847
1848 return wi;
1849 }
1850
1851 template <typename TypeTag>
1852 void
1854 updateConnectionDFactor(const Simulator& simulator,
1855 SingleWellStateType& ws) const
1856 {
1857 if (! this->well_ecl_.getWDFAC().useDFactor()) {
1858 return;
1859 }
1860
1861 auto& d_factor = ws.perf_data.connection_d_factor;
1862
1863 for (int perf = 0; perf < this->number_of_local_perforations_; ++perf) {
1864 const int cell_idx = this->well_cells_[perf];
1865 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/ 0);
1866
1867 d_factor[perf] = this->computeConnectionDFactor(perf, intQuants, ws);
1868 }
1869 }
1870
1871 template <typename TypeTag>
1874 computeConnectionDFactor(const int perf,
1875 const IntensiveQuantities& intQuants,
1876 const SingleWellStateType& ws) const
1877 {
1878 auto rhoGS = [regIdx = this->pvtRegionIdx()]() {
1879 return FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, regIdx);
1880 };
1881
1882 // Viscosity is evaluated at connection pressure.
1883 auto gas_visc = [connection_pressure = ws.perf_data.pressure[perf],
1884 temperature = ws.temperature,
1885 regIdx = this->pvtRegionIdx(), &intQuants]()
1886 {
1887 const auto rv = getValue(intQuants.fluidState().Rv());
1888
1889 const auto& gasPvt = FluidSystem::gasPvt();
1890
1891 // Note that rv here is from grid block with typically
1892 // p_block > connection_pressure
1893 // so we may very well have rv > rv_sat
1894 const Scalar rv_sat = gasPvt.saturatedOilVaporizationFactor
1895 (regIdx, temperature, connection_pressure);
1896
1897 if (! (rv < rv_sat)) {
1898 return gasPvt.saturatedViscosity(regIdx, temperature,
1899 connection_pressure);
1900 }
1901
1902 return gasPvt.viscosity(regIdx, temperature, connection_pressure,
1903 rv, getValue(intQuants.fluidState().Rvw()));
1904 };
1905
1906 const auto& connection = this->well_ecl_.getConnections()
1907 [ws.perf_data.ecl_index[perf]];
1908
1909 return this->well_ecl_.getWDFAC().getDFactor(rhoGS, gas_visc, connection);
1910 }
1911
1912
1913 template <typename TypeTag>
1914 void
1917 SingleWellStateType& ws) const
1918 {
1919 auto connCF = [&connIx = std::as_const(ws.perf_data.ecl_index),
1920 &conns = this->well_ecl_.getConnections()]
1921 (const int perf)
1922 {
1923 return conns[connIx[perf]].CF();
1924 };
1925
1926 auto& tmult = ws.perf_data.connection_compaction_tmult;
1927 auto& ctf = ws.perf_data.connection_transmissibility_factor;
1928
1929 for (int perf = 0; perf < this->number_of_local_perforations_; ++perf) {
1930 const int cell_idx = this->well_cells_[perf];
1931
1932 const auto& intQuants = simulator.model()
1933 .intensiveQuantities(cell_idx, /*timeIdx=*/ 0);
1934
1935 tmult[perf] = simulator.problem()
1936 .template wellTransMultiplier<double>(intQuants, cell_idx);
1937
1938 ctf[perf] = connCF(perf) * tmult[perf];
1939 }
1940 }
1941
1942
1943 template<typename TypeTag>
1946 {
1947 if constexpr (Indices::oilEnabled) {
1948 return fs.pressure(FluidSystem::oilPhaseIdx);
1949 } else if constexpr (Indices::gasEnabled) {
1950 return fs.pressure(FluidSystem::gasPhaseIdx);
1951 } else {
1952 return fs.pressure(FluidSystem::waterPhaseIdx);
1953 }
1954 }
1955
1956 template <typename TypeTag>
1957 template<class Value, class Callback>
1958 void
1960 getMobility(const Simulator& simulator,
1961 const int local_perf_index,
1962 std::vector<Value>& mob,
1963 Callback& extendEval,
1964 [[maybe_unused]] DeferredLogger& deferred_logger) const
1965 {
1966 auto relpermArray = []()
1967 {
1968 if constexpr (std::is_same_v<Value, Scalar>) {
1969 return std::array<Scalar,3>{};
1970 } else {
1971 return std::array<Eval,3>{};
1972 }
1973 };
1974 if (static_cast<std::size_t>(local_perf_index) >= this->well_cells_.size()) {
1975 OPM_THROW(std::invalid_argument,"The perforation index exceeds the size of the local containers - possibly getMobility was called with a global instead of a local perforation index!");
1976 }
1977 const int cell_idx = this->well_cells_[local_perf_index];
1978 assert (int(mob.size()) == this->num_conservation_quantities_);
1979 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/0);
1980 const auto& materialLawManager = simulator.problem().materialLawManager();
1981
1982 // either use mobility of the perforation cell or calculate its own
1983 // based on passing the saturation table index
1984 const int satid = this->saturation_table_number_[local_perf_index] - 1;
1985 const int satid_elem = materialLawManager->satnumRegionIdx(cell_idx);
1986 if (satid == satid_elem) { // the same saturation number is used. i.e. just use the mobilty from the cell
1987 for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
1988 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1989 continue;
1990 }
1991
1992 const unsigned activeCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::solventComponentIndex(phaseIdx));
1993 mob[activeCompIdx] = extendEval(intQuants.mobility(phaseIdx));
1994 }
1995 if constexpr (has_solvent) {
1996 mob[Indices::contiSolventEqIdx] = extendEval(intQuants.solventMobility());
1997 }
1998 } else {
1999 const auto& paramsCell = materialLawManager->connectionMaterialLawParams(satid, cell_idx);
2000 auto relativePerms = relpermArray();
2001 MaterialLaw::relativePermeabilities(relativePerms, paramsCell, intQuants.fluidState());
2002
2003 // reset the satnumvalue back to original
2004 materialLawManager->connectionMaterialLawParams(satid_elem, cell_idx);
2005
2006 // compute the mobility
2007 for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
2008 if (!FluidSystem::phaseIsActive(phaseIdx)) {
2009 continue;
2010 }
2011
2012 const unsigned activeCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::solventComponentIndex(phaseIdx));
2013 mob[activeCompIdx] = extendEval(relativePerms[phaseIdx] / intQuants.fluidState().viscosity(phaseIdx));
2014 }
2015
2016 // this may not work if viscosity and relperms has been modified?
2017 if constexpr (has_solvent) {
2018 OPM_DEFLOG_THROW(std::runtime_error, "individual mobility for wells does not work in combination with solvent", deferred_logger);
2019 }
2020 }
2021
2022 if (this->isInjector() && !this->inj_fc_multiplier_.empty()) {
2023 const auto perf_ecl_index = this->perforationData()[local_perf_index].ecl_index;
2024 const auto& connections = this->well_ecl_.getConnections();
2025 const auto& connection = connections[perf_ecl_index];
2026 if (connection.filterCakeActive()) {
2027 std::transform(mob.begin(), mob.end(), mob.begin(),
2028 [mult = this->inj_fc_multiplier_[local_perf_index] ](const auto val)
2029 { return val * mult; });
2030 }
2031 }
2032 }
2033
2034
2035 template<typename TypeTag>
2036 bool
2039 WellStateType& well_state,
2040 DeferredLogger& deferred_logger) const
2041 {
2042 OPM_TIMEFUNCTION();
2043 const auto& summary_state = simulator.vanguard().summaryState();
2044
2045 auto bhp_at_thp_limit = computeBhpAtThpLimitProdWithAlq(
2046 simulator, summary_state, this->getALQ(well_state), deferred_logger, /*iterate_if_no_solution */ false);
2047 if (bhp_at_thp_limit) {
2048 std::vector<Scalar> rates(this->number_of_phases_, 0.0);
2049 if (thp_update_iterations) {
2050 computeWellRatesWithBhpIterations(simulator, *bhp_at_thp_limit,
2051 rates, deferred_logger);
2052 } else {
2053 computeWellRatesWithBhp(simulator, *bhp_at_thp_limit,
2054 rates, deferred_logger);
2055 }
2056 auto& ws = well_state.well(this->name());
2057 ws.surface_rates = rates;
2058 ws.bhp = *bhp_at_thp_limit;
2059 ws.thp = this->getTHPConstraint(summary_state);
2060 return true;
2061 } else {
2062 return false;
2063 }
2064 }
2065
2066 template <typename TypeTag>
2067 void
2070 const std::function<Scalar(const Scalar)>& connPICalc,
2071 const std::vector<Scalar>& mobility,
2072 Scalar* connPI) const
2073 {
2074 const int np = this->number_of_phases_;
2075 for (int p = 0; p < np; ++p) {
2076 // Note: E100's notion of PI value phase mobility includes
2077 // the reciprocal FVF.
2078 const int canonical_phase_idx = FluidSystem::activeToCanonicalPhaseIdx(p);
2079 const auto connMob =
2080 mobility[FluidSystem::activePhaseToActiveCompIdx(p)] * fs.invB(canonical_phase_idx).value();
2081
2082 connPI[p] = connPICalc(connMob);
2083 }
2084
2085 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
2086 FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
2087 {
2088 const auto io = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
2089 const auto ig = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
2090
2091 const auto vapoil = connPI[ig] * fs.Rv().value();
2092 const auto disgas = connPI[io] * fs.Rs().value();
2093
2094 connPI[io] += vapoil;
2095 connPI[ig] += disgas;
2096 }
2097 }
2098
2099
2100 template <typename TypeTag>
2101 void
2104 const Phase preferred_phase,
2105 const std::function<Scalar(const Scalar)>& connIICalc,
2106 const std::vector<Scalar>& mobility,
2107 Scalar* connII,
2108 DeferredLogger& deferred_logger) const
2109 {
2110 auto phase_pos = 0;
2111 if (preferred_phase == Phase::GAS) {
2112 phase_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
2113 }
2114 else if (preferred_phase == Phase::OIL) {
2115 phase_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
2116 }
2117 else if (preferred_phase == Phase::WATER) {
2118 phase_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::waterPhaseIdx);
2119 }
2120 else {
2121 OPM_DEFLOG_THROW(NotImplemented,
2122 fmt::format("Unsupported Injector Type ({}) "
2123 "for well {} during connection I.I. calculation",
2124 static_cast<int>(preferred_phase), this->name()),
2125 deferred_logger);
2126 }
2127
2128 const auto mt = std::accumulate(mobility.begin(), mobility.end(), 0.0);
2129 const int canonicalPhaseIdx = FluidSystem::activeToCanonicalPhaseIdx(phase_pos);
2130 connII[phase_pos] = connIICalc(mt * fs.invB(canonicalPhaseIdx).value());
2131 }
2132
2133 template<typename TypeTag>
2134 template<class GasLiftSingleWell>
2135 std::unique_ptr<GasLiftSingleWell>
2137 initializeGliftWellTest_(const Simulator& simulator,
2138 WellStateType& well_state,
2139 const GroupState<Scalar>& group_state,
2140 GLiftEclWells& ecl_well_map,
2141 DeferredLogger& deferred_logger)
2142 {
2143 // Instantiate group info object (without initialization) since it is needed in GasLiftSingleWell
2144 auto& comm = simulator.vanguard().grid().comm();
2145 ecl_well_map.try_emplace(this->name(), &(this->wellEcl()), this->indexOfWell());
2147 ecl_well_map,
2148 simulator.vanguard().schedule(),
2149 simulator.vanguard().summaryState(),
2150 simulator.episodeIndex(),
2151 simulator.model().newtonMethod().numIterations(),
2152 deferred_logger,
2153 well_state,
2154 group_state,
2155 comm,
2156 false
2157 };
2158
2159 // Return GasLiftSingleWell object to use the wellTestALQ() function
2160 std::set<int> sync_groups;
2161 const auto& summary_state = simulator.vanguard().summaryState();
2162 return std::make_unique<GasLiftSingleWell>(*this,
2163 simulator,
2164 summary_state,
2165 deferred_logger,
2166 well_state,
2167 group_state,
2168 group_info,
2169 sync_groups,
2170 comm,
2171 false);
2172
2173 }
2174
2175} // namespace Opm
2176
2177#endif
#define OPM_DEFLOG_THROW(Exception, message, deferred_logger)
Definition: DeferredLoggingErrorHelpers.hpp:45
Definition: DeferredLogger.hpp:57
void info(const std::string &tag, const std::string &message)
void warning(const std::string &tag, const std::string &message)
void debug(const std::string &tag, const std::string &message)
Definition: GasLiftGroupInfo.hpp:46
Definition: GroupState.hpp:41
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:198
Definition: SingleWellState.hpp:43
Scalar temperature
Definition: SingleWellState.hpp:99
PerfData< Scalar > perf_data
Definition: SingleWellState.hpp:120
Class for computing BHP limits.
Definition: WellBhpThpCalculator.hpp:41
Scalar calculateMinimumBhpFromThp(const WellState< Scalar, IndexTraits > &well_state, const Well &well, const SummaryState &summaryState, const Scalar rho) const
Scalar mostStrictBhpFromBhpLimits(const SummaryState &summaryState) const
Obtain the most strict BHP from BHP limits.
bool isStableSolution(const WellState< Scalar, IndexTraits > &well_state, const Well &well, const std::vector< Scalar > &rates, const SummaryState &summaryState) const
EvalWell calculateBhpFromThp(const WellState< Scalar, IndexTraits > &well_state, const std::vector< EvalWell > &rates, const Well &well, const SummaryState &summaryState, const Scalar rho, DeferredLogger &deferred_logger) const
std::optional< Scalar > estimateStableBhp(const WellState< Scalar, IndexTraits > &well_state, const Well &well, const std::vector< Scalar > &rates, const Scalar rho, const SummaryState &summaryState) const
int number_of_local_perforations_
Definition: WellInterfaceGeneric.hpp:340
FluidSystem::Scalar wsolvent_
Definition: WellInterfaceGeneric.hpp:380
const std::vector< FluidSystem::Scalar > & wellIndex() const
Definition: WellInterfaceGeneric.hpp:152
Definition: WellInterfaceIndices.hpp:34
bool solveWellWithBhp(const Simulator &simulator, const double dt, const Scalar bhp, WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:701
void checkWellOperability(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1039
Scalar woxygen() const
Definition: WellInterface_impl.hpp:165
IndividualOrGroup
Definition: WellInterface.hpp:246
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:81
Scalar computeConnectionDFactor(const int perf, const IntensiveQuantities &intQuants, const SingleWellStateType &ws) const
Definition: WellInterface_impl.hpp:1874
void solveWellEquation(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:805
bool updateWellControlAndStatusLocalIteration(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, const Scalar WQTotal, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false)
Definition: WellInterface_impl.hpp:275
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:103
bool iterateWellEquations(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:503
BlackOilFluidState< Eval, FluidSystem, has_temperature, has_energy, Indices::compositionSwitchIdx >=0, has_watVapor, has_brine, has_saltPrecip, has_disgas_in_water, Indices::numPhases > FluidState
Definition: WellInterface.hpp:136
Scalar wfoam() const
Definition: WellInterface_impl.hpp:127
std::optional< Scalar > estimateOperableBhp(const Simulator &ebos_simulator, const double dt, WellStateType &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:668
std::vector< RateVector > connectionRates_
Definition: WellInterface.hpp:369
void computeConnLevelProdInd(const FluidState &fs, const std::function< Scalar(const Scalar)> &connPICalc, const std::vector< Scalar > &mobility, Scalar *connPI) const
Definition: WellInterface_impl.hpp:2069
void addCellRates(RateVector &rates, int cellIdx) const
Definition: WellInterface_impl.hpp:1003
void gliftBeginTimeStepWellTestUpdateALQ(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, GLiftEclWells &ecl_well_map, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1066
void assembleWellEq(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:867
Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const
Definition: WellInterface_impl.hpp:1019
virtual void init(const std::vector< Scalar > &depth_arg, const Scalar gravity_arg, const std::vector< Scalar > &B_avg, const bool changed_to_open_this_step)
Definition: WellInterface_impl.hpp:95
void getMobility(const Simulator &simulator, const int local_perf_index, std::vector< Value > &mob, Callback &extendEval, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1960
void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:883
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:82
std::vector< Scalar > initialWellRateFractions(const Simulator &ebosSimulator, const WellStateType &well_state) const
Definition: WellInterface_impl.hpp:1636
bool updateWellStateWithTHPTargetProd(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:2038
void updateConnectionDFactor(const Simulator &simulator, SingleWellStateType &ws) const
Definition: WellInterface_impl.hpp:1854
Eval getPerfCellPressure(const FluidState &fs) const
Definition: WellInterface_impl.hpp:1945
void initializeProducerWellState(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1686
bool wellUnderZeroGroupRateTarget(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger, std::optional< bool > group_control=std::nullopt) const
Definition: WellInterface_impl.hpp:1604
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:109
bool updateWellControl(const Simulator &simulator, const IndividualOrGroup iog, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:189
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:83
static constexpr bool has_solvent
Definition: WellInterface.hpp:111
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: WellInterface.hpp:89
void updateConnectionTransmissibilityFactor(const Simulator &simulator, SingleWellStateType &ws) const
Definition: WellInterface_impl.hpp:1916
void computeConnLevelInjInd(const FluidState &fs, const Phase preferred_phase, const std::function< Scalar(const Scalar)> &connIICalc, const std::vector< Scalar > &mobility, Scalar *connII, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:2103
typename GasLiftGroupInfo< Scalar, IndexTraits >::GLiftEclWells GLiftEclWells
Definition: WellInterface.hpp:91
std::unique_ptr< GasLiftSingleWell > initializeGliftWellTest_(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, GLiftEclWells &ecl_well_map, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:2137
void wellTesting(const Simulator &simulator, const double simulation_time, WellStateType &well_state, const GroupState< Scalar > &group_state, WellTestState &welltest_state, GLiftEclWells &ecl_well_map, std::map< std::string, double > &open_times, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:383
Scalar wsalt() const
Definition: WellInterface_impl.hpp:141
void prepareWellBeforeAssembling(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:903
typename Base::Eval Eval
Definition: WellInterface.hpp:95
bool wellUnderZeroRateTarget(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1586
WellInterface(const Well &well, const ParallelWellInfo< Scalar > &pw_info, const int time_step, const ModelParameters &param, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_conservation_quantities, const int num_phases, const int index_of_well, const std::vector< PerforationData< Scalar > > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:58
bool solveWellForTesting(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:759
void updateWellOperability(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1128
bool solveWellWithOperabilityCheck(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:566
Scalar wpolymer() const
Definition: WellInterface_impl.hpp:111
bool updateWellOperabilityFromWellEq(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1165
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
Scalar wurea() const
Definition: WellInterface_impl.hpp:177
bool solveWellWithZeroRate(const Simulator &simulator, const double dt, WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:738
Scalar wmicrobes() const
Definition: WellInterface_impl.hpp:153
virtual void scaleSegmentRatesAndPressure(WellStateType &well_state) const
Definition: WellInterface_impl.hpp:1184
virtual void updateWellStateWithTarget(const Simulator &simulator, const GroupState< Scalar > &group_state, WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1192
static constexpr bool has_zFraction
Definition: WellInterface.hpp:112
bool stoppedOrZeroRateTarget(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1623
Definition: WellState.hpp:66
constexpr int numPhases() const
The number of phases present.
Definition: WellState.hpp:249
const SingleWellState< Scalar, IndexTraits > & well(std::size_t well_index) const
Definition: WellState.hpp:290
bool isOpen(const std::string &name) const
Definition: WellState.hpp:197
@ NONE
Definition: DeferredLogger.hpp:46
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
Definition: blackoilbioeffectsmodules.hh:43
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Static data associated with a well perforation.
Definition: PerforationData.hpp:30