6 template<
typename TypeTag>
7 BlackoilWellModel<TypeTag>::
8 BlackoilWellModel(
const Wells* wells_arg,
9 WellCollection* well_collection,
10 const std::vector< const Well* >& wells_ecl,
11 const ModelParameters& param,
12 const RateConverterType& rate_converter,
13 const bool terminal_output,
14 const int current_timeIdx,
15 const std::vector<int>& pvt_region_idx)
16 : wells_active_(wells_arg!=nullptr)
18 , wells_ecl_(wells_ecl)
19 , number_of_wells_(wells_arg ? (wells_arg->number_of_wells) : 0)
20 , number_of_phases_(wells_arg ? (wells_arg->number_of_phases) : 0)
22 , well_container_(createWellContainer(wells_arg, wells_ecl, param.use_multisegment_well_, current_timeIdx, param) )
23 , well_collection_(well_collection)
24 , terminal_output_(terminal_output)
25 , has_solvent_(GET_PROP_VALUE(TypeTag, EnableSolvent))
26 , has_polymer_(GET_PROP_VALUE(TypeTag, EnablePolymer))
27 , current_timeIdx_(current_timeIdx)
28 , rate_converter_(rate_converter)
29 , pvt_region_idx_(pvt_region_idx)
37 template<
typename TypeTag>
39 BlackoilWellModel<TypeTag>::
40 init(
const PhaseUsage phase_usage_arg,
41 const std::vector<bool>& active_arg,
42 const double gravity_arg,
43 const std::vector<double>& depth_arg,
48 global_nc_ = global_nc;
50 phase_usage_ = phase_usage_arg;
53 if ( ! localWellsActive() ) {
57 calculateEfficiencyFactors();
60 const auto& pu = phase_usage_;
61 const int np = pu.num_phases;
65 assert (np == 3 || (np == 2 && !pu.phase_used[Gas]) );
70 if (PolymerModule::hasPlyshlog()) {
71 computeRepRadiusPerfLength(grid);
75 number_of_cells_ = Opm::UgGridHelpers::numCells(grid);
79 for (
auto& well : well_container_) {
80 well->init(&phase_usage_, &active_, depth_arg, gravity_arg, number_of_cells_);
88 template<
typename TypeTag>
90 BlackoilWellModel<TypeTag>::
91 setVFPProperties(
const VFPProperties* vfp_properties_arg)
93 for (
auto& well : well_container_) {
94 well->setVFPProperties(vfp_properties_arg);
103 template<
typename TypeTag>
105 BlackoilWellModel<TypeTag>::
108 return number_of_wells_;
115 template<
typename TypeTag>
116 std::vector<typename BlackoilWellModel<TypeTag>::WellInterfacePtr >
117 BlackoilWellModel<TypeTag>::
118 createWellContainer(
const Wells* wells,
119 const std::vector< const Well* >& wells_ecl,
120 const bool use_multisegment_well,
122 const ModelParameters& param)
124 std::vector<WellInterfacePtr> well_container;
126 const int nw = wells ? (wells->number_of_wells) : 0;
129 well_container.reserve(nw);
133 for (
int w = 0; w < nw; ++w) {
134 const std::string well_name = std::string(wells->name[w]);
137 const int nw_wells_ecl = wells_ecl.size();
139 for (; index_well < nw_wells_ecl; ++index_well) {
140 if (well_name == wells_ecl[index_well]->name()) {
146 if (index_well == nw_wells_ecl) {
147 OPM_THROW(std::logic_error,
"Could not find well " << well_name <<
" in wells_ecl ");
150 const Well* well_ecl = wells_ecl[index_well];
152 if ( !well_ecl->isMultiSegment(time_step) || !use_multisegment_well) {
153 well_container.emplace_back(
new StandardWell<TypeTag>(well_ecl, time_step, wells, param) );
155 well_container.emplace_back(
new MultisegmentWell<TypeTag>(well_ecl, time_step, wells, param) );
159 return well_container;
166 template<
typename TypeTag>
168 BlackoilWellModel<TypeTag>::
169 assemble(Simulator& ebosSimulator,
170 const int iterationIdx,
172 WellState& well_state)
174 SimulatorReport report;
175 if ( ! wellsActive() ) {
179 if (iterationIdx == 0) {
180 prepareTimeStep(ebosSimulator, well_state);
183 updateWellControls(well_state);
185 initPrimaryVariablesEvaluation();
187 if (iterationIdx == 0) {
188 calculateExplicitQuantities(ebosSimulator, well_state);
191 if (param_.solve_welleq_initially_ && iterationIdx == 0) {
193 report = solveWellEq(ebosSimulator, dt, well_state);
195 assembleWellEq(ebosSimulator, dt, well_state,
false);
197 report.converged =
true;
205 template<
typename TypeTag>
207 BlackoilWellModel<TypeTag>::
208 assembleWellEq(Simulator& ebosSimulator,
210 WellState& well_state,
211 bool only_wells)
const 213 for (
int w = 0; w < number_of_wells_; ++w) {
214 well_container_[w]->assembleWellEq(ebosSimulator, dt, well_state, only_wells);
224 template<
typename TypeTag>
226 BlackoilWellModel<TypeTag>::
227 apply( BVector& r)
const 229 if ( ! localWellsActive() ) {
233 for (
auto& well : well_container_) {
243 template<
typename TypeTag>
245 BlackoilWellModel<TypeTag>::
246 apply(
const BVector& x, BVector& Ax)
const 249 if ( ! localWellsActive() ) {
253 for (
auto& well : well_container_) {
263 template<
typename TypeTag>
265 BlackoilWellModel<TypeTag>::
266 applyScaleAdd(
const Scalar alpha,
const BVector& x, BVector& Ax)
const 268 if ( ! localWellsActive() ) {
272 if( scaleAddRes_.size() != Ax.size() ) {
273 scaleAddRes_.resize( Ax.size() );
278 apply( x, scaleAddRes_ );
280 Ax.axpy( alpha, scaleAddRes_ );
287 template<
typename TypeTag>
289 BlackoilWellModel<TypeTag>::
290 recoverWellSolutionAndUpdateWellState(
const BVector& x, WellState& well_state)
const 292 for (
auto& well : well_container_) {
293 well->recoverWellSolutionAndUpdateWellState(x, well_state);
301 template<
typename TypeTag>
303 BlackoilWellModel<TypeTag>::
304 flowPhaseToEbosPhaseIdx(
const int phaseIdx )
const 306 const auto& pu = phase_usage_;
307 if (active_[Water] && pu.phase_pos[Water] == phaseIdx)
308 return FluidSystem::waterPhaseIdx;
309 if (active_[Oil] && pu.phase_pos[Oil] == phaseIdx)
310 return FluidSystem::oilPhaseIdx;
311 if (active_[Gas] && pu.phase_pos[Gas] == phaseIdx)
312 return FluidSystem::gasPhaseIdx;
314 assert(phaseIdx < 3);
323 template<
typename TypeTag>
325 BlackoilWellModel<TypeTag>::
328 return number_of_phases_;
335 template<
typename TypeTag>
337 BlackoilWellModel<TypeTag>::
338 resetWellControlFromState(
const WellState& xw)
const 340 const int nw = wells_->number_of_wells;
341 for (
int w = 0; w < nw; ++w) {
342 WellControls* wc = wells_->ctrls[w];
343 well_controls_set_current( wc, xw.currentControls()[w]);
351 template<
typename TypeTag>
356 return wells_active_;
363 template<
typename TypeTag>
368 wells_active_ = wells_active;
375 template<
typename TypeTag>
380 return number_of_wells_ > 0;
387 template<
typename TypeTag>
392 for (
auto& well : well_container_) {
393 well->initPrimaryVariablesEvaluation();
401 template<
typename TypeTag>
403 BlackoilWellModel<TypeTag>::
404 solveWellEq(Simulator& ebosSimulator,
406 WellState& well_state)
const 408 const int nw = number_of_wells_;
409 WellState well_state0 = well_state;
411 const int numComp = numComponents();
412 std::vector< Scalar > B_avg( numComp, Scalar() );
413 computeAverageFormationFactor(ebosSimulator, B_avg);
415 const int max_iter = param_.max_welleq_iter_;
420 assembleWellEq(ebosSimulator, dt, well_state,
true);
422 converged = getWellConvergence(ebosSimulator, B_avg);
425 if (wellCollection()->groupControlActive()) {
426 converged = converged && wellCollection()->groupTargetConverged(well_state.wellRates());
434 if( localWellsActive() )
436 for (
auto& well : well_container_) {
437 well->solveEqAndUpdateWellState(well_state);
445 updateWellControls(well_state);
446 initPrimaryVariablesEvaluation();
448 }
while (it < max_iter);
451 if ( terminal_output_ ) {
452 OpmLog::debug(
"Well equation solution gets converged with " + std::to_string(it) +
" iterations");
455 if ( terminal_output_ ) {
456 OpmLog::debug(
"Well equation solution failed in getting converged with " + std::to_string(it) +
" iterations");
459 well_state = well_state0;
460 updatePrimaryVariables(well_state);
462 for (
int w = 0; w < nw; ++w) {
463 WellControls* wc = well_container_[w]->wellControls();
464 well_controls_set_current(wc, well_state.currentControls()[w]);
468 SimulatorReport report;
469 report.converged = converged;
470 report.total_well_iterations = it;
478 template<
typename TypeTag>
480 BlackoilWellModel<TypeTag>::
481 getWellConvergence(
const Simulator& ebosSimulator,
482 const std::vector<Scalar>& B_avg)
const 484 ConvergenceReport report;
486 for (
const auto& well : well_container_) {
487 report += well->getWellConvergence(B_avg);
492 bool nan_residual_found = report.nan_residual_found;
493 const auto& grid = ebosSimulator.gridManager().grid();
494 int value = nan_residual_found ? 1 : 0;
496 nan_residual_found = grid.comm().max(value);
498 if (nan_residual_found) {
499 for (
const auto& well : report.nan_residual_wells) {
500 OpmLog::debug(
"NaN residual found with phase " + well.phase_name +
" for well " + well.well_name);
502 OPM_THROW(Opm::NumericalProblem,
"NaN residual found!");
508 bool too_large_residual_found = report.too_large_residual_found;
509 const auto& grid = ebosSimulator.gridManager().grid();
510 int value = too_large_residual_found ? 1 : 0;
512 too_large_residual_found = grid.comm().max(value);
513 if (too_large_residual_found) {
514 for (
const auto& well : report.too_large_residual_wells) {
515 OpmLog::debug(
"Too large residual found with phase " + well.phase_name +
" fow well " + well.well_name);
517 OPM_THROW(Opm::NumericalProblem,
"Too large residual found!");
522 bool converged_well = report.converged;
524 const auto& grid = ebosSimulator.gridManager().grid();
525 int value = converged_well ? 1 : 0;
527 converged_well = grid.comm().min(value);
530 return converged_well;
537 template<
typename TypeTag>
543 for (
auto& well : well_container_) {
544 well->calculateExplicitQuantities(ebosSimulator, xw);
552 template<
typename TypeTag>
561 if( !wellsActive() ) return ;
565 for (
const auto& well : well_container_) {
566 well->updateWellControl(xw, logger);
569 updateGroupControls(xw);
576 template<
typename TypeTag>
580 const int current_step,
581 const Wells* wells_struct,
583 DynamicListEconLimited& list_econ_limited)
const 585 for (
const auto& well : well_container_) {
586 well->updateListEconLimited(well_state, list_econ_limited);
594 template<
typename TypeTag>
598 const WellState& well_state,
599 std::vector<double>& well_potentials)
const 602 const int nw = number_of_wells_;
603 const int np = number_of_phases_;
604 well_potentials.resize(nw * np, 0.0);
606 for (
int w = 0; w < nw; ++w) {
607 std::vector<double> potentials;
608 well_container_[w]->computeWellPotentials(ebosSimulator, well_state, potentials);
611 for (
int p = 0; p < np; ++p) {
612 well_potentials[w * np + p] = std::abs(potentials[p]);
621 template<
typename TypeTag>
623 BlackoilWellModel<TypeTag>::
624 prepareTimeStep(
const Simulator& ebos_simulator,
625 WellState& well_state)
const 632 resetWellControlFromState(well_state);
635 prepareGroupControl(ebos_simulator, well_state);
638 for (
int w = 0; w < number_of_wells_; ++w) {
639 WellControls* wc = well_container_[w]->wellControls();
640 const int control = well_controls_get_current(wc);
641 well_state.currentControls()[w] = control;
644 well_container_[w]->updateWellStateWithTarget(control, well_state);
648 if (well_state.isNewWell(w) ) {
649 well_state.setNewWell(w,
false);
658 template<
typename TypeTag>
660 BlackoilWellModel<TypeTag>::
661 prepareGroupControl(
const Simulator& ebos_simulator,
662 WellState& well_state)
const 665 if (well_collection_->groupControlActive()) {
666 for (
int w = 0; w < number_of_wells_; ++w) {
667 WellControls* wc = well_container_[w]->wellControls();
668 WellNode& well_node = well_collection_->findWellNode(well_container_[w]->name());
673 int ctrl_index = well_controls_get_current(wc);
675 const int group_control_index = well_node.groupControlIndex();
676 if (group_control_index >= 0 && ctrl_index < 0) {
678 well_controls_set_current(wc, group_control_index);
679 well_state.currentControls()[w] = group_control_index;
684 ctrl_index = well_controls_get_current(wc);
685 if (well_node.groupControlIndex() >= 0 && ctrl_index == well_node.groupControlIndex()) {
687 well_node.setIndividualControl(
false);
690 well_node.setIndividualControl(
true);
694 if (well_collection_->requireWellPotentials()) {
697 std::vector<double> well_potentials;
698 computeWellPotentials(ebos_simulator, well_state, well_potentials);
703 well_collection_->setGuideRatesWithPotentials(wells_, phase_usage_, well_potentials);
706 applyVREPGroupControl(well_state);
708 if (!wellCollection()->groupControlApplied()) {
709 wellCollection()->applyGroupControls();
711 wellCollection()->updateWellTargets(well_state.wellRates());
720 template<
typename TypeTag>
722 BlackoilWellModel<TypeTag>::
723 wellCollection()
const 725 return well_collection_;
732 template<
typename TypeTag>
734 BlackoilWellModel<TypeTag>::
735 calculateEfficiencyFactors()
737 if ( !localWellsActive() ) {
741 const int nw = number_of_wells_;
743 for (
int w = 0; w < nw; ++w) {
744 const std::string well_name = well_container_[w]->name();
745 const WellNode& well_node = wellCollection()->findWellNode(well_name);
747 const double well_efficiency_factor = well_node.getAccumulativeEfficiencyFactor();
749 well_container_[w]->setWellEfficiencyFactor(well_efficiency_factor);
757 template<
typename TypeTag>
759 BlackoilWellModel<TypeTag>::
760 computeWellVoidageRates(
const WellState& well_state,
761 std::vector<double>& well_voidage_rates,
762 std::vector<double>& voidage_conversion_coeffs)
const 764 if ( !localWellsActive() ) {
772 const int nw = numWells();
773 const int np = numPhases();
776 well_voidage_rates.resize(nw, 0);
778 voidage_conversion_coeffs.resize(nw * np, 1.0);
780 std::vector<double> well_rates(np, 0.0);
781 std::vector<double> convert_coeff(np, 1.0);
783 for (
int w = 0; w < nw; ++w) {
784 const bool is_producer = well_container_[w]->wellType() == PRODUCER;
785 const int well_cell_top = well_container_[w]->cells()[0];
786 const int pvtRegionIdx = pvt_region_idx_[well_cell_top];
790 std::transform(well_state.wellRates().begin() + np * w,
791 well_state.wellRates().begin() + np * (w + 1),
792 well_rates.begin(), std::negate<double>());
795 const int fipreg = 0;
797 rate_converter_.calcCoeff(fipreg, pvtRegionIdx, convert_coeff);
798 well_voidage_rates[w] = std::inner_product(well_rates.begin(), well_rates.end(),
799 convert_coeff.begin(), 0.0);
803 std::copy(well_state.wellRates().begin() + np * w,
804 well_state.wellRates().begin() + np * (w + 1),
807 const int fipreg = 0;
808 rate_converter_.calcCoeff(fipreg, pvtRegionIdx, convert_coeff);
809 std::copy(convert_coeff.begin(), convert_coeff.end(),
810 voidage_conversion_coeffs.begin() + np * w);
819 template<
typename TypeTag>
821 BlackoilWellModel<TypeTag>::
822 applyVREPGroupControl(WellState& well_state)
const 824 if ( wellCollection()->havingVREPGroups() ) {
825 std::vector<double> well_voidage_rates;
826 std::vector<double> voidage_conversion_coeffs;
827 computeWellVoidageRates(well_state, well_voidage_rates, voidage_conversion_coeffs);
828 wellCollection()->applyVREPGroupControls(well_voidage_rates, voidage_conversion_coeffs);
831 for (
const WellNode* well_node : wellCollection()->getLeafNodes()) {
832 if (well_node->isInjector() && !well_node->individualControl()) {
833 const int well_index = well_node->selfIndex();
834 well_state.currentControls()[well_index] = well_node->groupControlIndex();
836 WellControls* wc = well_container_[well_index]->wellControls();
837 well_controls_set_current(wc, well_node->groupControlIndex());
847 template<
typename TypeTag>
849 BlackoilWellModel<TypeTag>::
850 updateGroupControls(WellState& well_state)
const 853 if (wellCollection()->groupControlActive()) {
854 for (
int w = 0; w < number_of_wells_; ++w) {
857 WellNode& well_node = well_collection_->findWellNode(well_container_[w]->name());
860 const int current = well_state.currentControls()[w];
861 if (well_node.groupControlIndex() >= 0 && current == well_node.groupControlIndex()) {
863 well_node.setIndividualControl(
false);
866 well_node.setIndividualControl(
true);
870 applyVREPGroupControl(well_state);
873 wellCollection()->updateWellTargets(well_state.wellRates());
875 for (
int w = 0; w < number_of_wells_; ++w) {
880 const int current = well_state.currentControls()[w];
881 well_container_[w]->updateWellStateWithTarget(current, well_state);
890 template<
typename TypeTag>
892 BlackoilWellModel<TypeTag>::
893 setupCompressedToCartesian(
const int* global_cell,
int number_of_cells, std::map<int,int>& cartesian_to_compressed )
const 896 for (
int i = 0; i < number_of_cells; ++i) {
897 cartesian_to_compressed.insert(std::make_pair(global_cell[i], i));
901 for (
int i = 0; i < number_of_cells; ++i) {
902 cartesian_to_compressed.insert(std::make_pair(i, i));
908 template<
typename TypeTag>
910 BlackoilWellModel<TypeTag>::
911 computeRepRadiusPerfLength(
const Grid& grid)
915 const int* global_cell = Opm::UgGridHelpers::globalCell(grid);
917 std::map<int,int> cartesian_to_compressed;
918 setupCompressedToCartesian(global_cell, number_of_cells_,
919 cartesian_to_compressed);
921 for (
const auto& well : well_container_) {
922 well->computeRepRadiusPerfLength(grid, cartesian_to_compressed);
930 template<
typename TypeTag>
932 BlackoilWellModel<TypeTag>::
933 computeAverageFormationFactor(
const Simulator& ebosSimulator,
934 std::vector<double>& B_avg)
const 936 const int np = numPhases();
938 const auto& grid = ebosSimulator.gridManager().grid();
939 const auto& gridView = grid.leafGridView();
940 ElementContext elemCtx(ebosSimulator);
941 const auto& elemEndIt = gridView.template end<0, Dune::Interior_Partition>();
943 for (
auto elemIt = gridView.template begin</*codim=*/0, Dune::Interior_Partition>();
944 elemIt != elemEndIt; ++elemIt)
946 elemCtx.updatePrimaryStencil(*elemIt);
947 elemCtx.updatePrimaryIntensiveQuantities(0);
949 const auto& intQuants = elemCtx.intensiveQuantities(0, 0);
950 const auto& fs = intQuants.fluidState();
952 for (
int phaseIdx = 0; phaseIdx < np; ++phaseIdx )
954 auto& B = B_avg[ phaseIdx ];
955 const int ebosPhaseIdx = flowPhaseToEbosPhaseIdx(phaseIdx);
957 B += 1 / fs.invB(ebosPhaseIdx).value();
960 auto& B = B_avg[solventSaturationIdx];
961 B += 1 / intQuants.solventInverseFormationVolumeFactor().value();
966 grid.comm().sum(B_avg.data(), B_avg.size());
967 for(
auto& bval: B_avg)
977 template<
typename TypeTag>
979 BlackoilWellModel<TypeTag>::
980 updatePrimaryVariables(
const WellState& well_state)
const 982 for (
const auto& well : well_container_) {
983 well->updatePrimaryVariables(well_state);
bool localWellsActive() const
return true if wells are available on this process
Definition: BlackoilWellModel_impl.hpp:378
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:66
Utility class to handle the log messages about well switching.
Definition: WellSwitchingLogger.hpp:44
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilWellModel_impl.hpp:354
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
void calculateExplicitQuantities(const Simulator &ebosSimulator, const WellState &xw) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:540
void updateListEconLimited(const Schedule &schedule, const int current_step, const Wells *wells_struct, const WellState &well_state, DynamicListEconLimited &list_econ_limited) const
upate the dynamic lists related to economic limits
Definition: BlackoilWellModel_impl.hpp:579