00001
00002
00003
00004 namespace Opm {
00005
00006 template<typename TypeTag>
00007 BlackoilWellModel<TypeTag>::
00008 BlackoilWellModel(const Wells* wells_arg,
00009 WellCollection* well_collection,
00010 const std::vector< const Well* >& wells_ecl,
00011 const ModelParameters& param,
00012 const RateConverterType& rate_converter,
00013 const bool terminal_output,
00014 const int current_timeIdx,
00015 const std::vector<int>& pvt_region_idx)
00016 : wells_active_(wells_arg!=nullptr)
00017 , wells_(wells_arg)
00018 , wells_ecl_(wells_ecl)
00019 , number_of_wells_(wells_arg ? (wells_arg->number_of_wells) : 0)
00020 , number_of_phases_(wells_arg ? (wells_arg->number_of_phases) : 0)
00021 , param_(param)
00022 , well_container_(createWellContainer(wells_arg, wells_ecl, param.use_multisegment_well_, current_timeIdx, param) )
00023 , well_collection_(well_collection)
00024 , terminal_output_(terminal_output)
00025 , has_solvent_(GET_PROP_VALUE(TypeTag, EnableSolvent))
00026 , has_polymer_(GET_PROP_VALUE(TypeTag, EnablePolymer))
00027 , current_timeIdx_(current_timeIdx)
00028 , rate_converter_(rate_converter)
00029 , pvt_region_idx_(pvt_region_idx)
00030 {
00031 }
00032
00033
00034
00035
00036
00037 template<typename TypeTag>
00038 void
00039 BlackoilWellModel<TypeTag>::
00040 init(const PhaseUsage phase_usage_arg,
00041 const std::vector<bool>& active_arg,
00042 const double gravity_arg,
00043 const std::vector<double>& depth_arg,
00044 long int global_nc,
00045 const Grid& grid)
00046 {
00047
00048 global_nc_ = global_nc;
00049
00050 phase_usage_ = phase_usage_arg;
00051 active_ = active_arg;
00052
00053 if ( ! localWellsActive() ) {
00054 return;
00055 }
00056
00057 calculateEfficiencyFactors();
00058
00059 #ifndef NDEBUG
00060 const auto& pu = phase_usage_;
00061 const int np = pu.num_phases;
00062
00063
00064
00065 assert (np == 3 || (np == 2 && !pu.phase_used[Gas]) );
00066 #endif
00067
00068 if (has_polymer_)
00069 {
00070 if (PolymerModule::hasPlyshlog()) {
00071 computeRepRadiusPerfLength(grid);
00072 }
00073 }
00074
00075 number_of_cells_ = Opm::UgGridHelpers::numCells(grid);
00076
00077
00078
00079 for (auto& well : well_container_) {
00080 well->init(&phase_usage_, &active_, depth_arg, gravity_arg, number_of_cells_);
00081 }
00082 }
00083
00084
00085
00086
00087
00088 template<typename TypeTag>
00089 void
00090 BlackoilWellModel<TypeTag>::
00091 setVFPProperties(const VFPProperties* vfp_properties_arg)
00092 {
00093 for (auto& well : well_container_) {
00094 well->setVFPProperties(vfp_properties_arg);
00095 }
00096 }
00097
00098
00099
00100
00101
00102
00103 template<typename TypeTag>
00104 int
00105 BlackoilWellModel<TypeTag>::
00106 numWells() const
00107 {
00108 return number_of_wells_;
00109 }
00110
00111
00112
00113
00114
00115 template<typename TypeTag>
00116 std::vector<typename BlackoilWellModel<TypeTag>::WellInterfacePtr >
00117 BlackoilWellModel<TypeTag>::
00118 createWellContainer(const Wells* wells,
00119 const std::vector< const Well* >& wells_ecl,
00120 const bool use_multisegment_well,
00121 const int time_step,
00122 const ModelParameters& param)
00123 {
00124 std::vector<WellInterfacePtr> well_container;
00125
00126 const int nw = wells ? (wells->number_of_wells) : 0;
00127
00128 if (nw > 0) {
00129 well_container.reserve(nw);
00130
00131
00132
00133 for (int w = 0; w < nw; ++w) {
00134 const std::string well_name = std::string(wells->name[w]);
00135
00136
00137 const int nw_wells_ecl = wells_ecl.size();
00138 int index_well = 0;
00139 for (; index_well < nw_wells_ecl; ++index_well) {
00140 if (well_name == wells_ecl[index_well]->name()) {
00141 break;
00142 }
00143 }
00144
00145
00146 if (index_well == nw_wells_ecl) {
00147 OPM_THROW(std::logic_error, "Could not find well " << well_name << " in wells_ecl ");
00148 }
00149
00150 const Well* well_ecl = wells_ecl[index_well];
00151
00152 if ( !well_ecl->isMultiSegment(time_step) || !use_multisegment_well) {
00153 well_container.emplace_back(new StandardWell<TypeTag>(well_ecl, time_step, wells, param) );
00154 } else {
00155 well_container.emplace_back(new MultisegmentWell<TypeTag>(well_ecl, time_step, wells, param) );
00156 }
00157 }
00158 }
00159 return well_container;
00160 }
00161
00162
00163
00164
00165
00166 template<typename TypeTag>
00167 SimulatorReport
00168 BlackoilWellModel<TypeTag>::
00169 assemble(Simulator& ebosSimulator,
00170 const int iterationIdx,
00171 const double dt,
00172 WellState& well_state)
00173 {
00174 SimulatorReport report;
00175 if ( ! wellsActive() ) {
00176 return report;
00177 }
00178
00179 if (iterationIdx == 0) {
00180 prepareTimeStep(ebosSimulator, well_state);
00181 }
00182
00183 updateWellControls(well_state);
00184
00185 initPrimaryVariablesEvaluation();
00186
00187 if (iterationIdx == 0) {
00188 calculateExplicitQuantities(ebosSimulator, well_state);
00189 }
00190
00191 if (param_.solve_welleq_initially_ && iterationIdx == 0) {
00192
00193 report = solveWellEq(ebosSimulator, dt, well_state);
00194 }
00195 assembleWellEq(ebosSimulator, dt, well_state, false);
00196
00197 report.converged = true;
00198 return report;
00199 }
00200
00201
00202
00203
00204
00205 template<typename TypeTag>
00206 void
00207 BlackoilWellModel<TypeTag>::
00208 assembleWellEq(Simulator& ebosSimulator,
00209 const double dt,
00210 WellState& well_state,
00211 bool only_wells) const
00212 {
00213 for (int w = 0; w < number_of_wells_; ++w) {
00214 well_container_[w]->assembleWellEq(ebosSimulator, dt, well_state, only_wells);
00215 }
00216 }
00217
00218
00219
00220
00221
00222
00223
00224 template<typename TypeTag>
00225 void
00226 BlackoilWellModel<TypeTag>::
00227 apply( BVector& r) const
00228 {
00229 if ( ! localWellsActive() ) {
00230 return;
00231 }
00232
00233 for (auto& well : well_container_) {
00234 well->apply(r);
00235 }
00236 }
00237
00238
00239
00240
00241
00242
00243 template<typename TypeTag>
00244 void
00245 BlackoilWellModel<TypeTag>::
00246 apply(const BVector& x, BVector& Ax) const
00247 {
00248
00249 if ( ! localWellsActive() ) {
00250 return;
00251 }
00252
00253 for (auto& well : well_container_) {
00254 well->apply(x, Ax);
00255 }
00256 }
00257
00258
00259
00260
00261
00262
00263 template<typename TypeTag>
00264 void
00265 BlackoilWellModel<TypeTag>::
00266 applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const
00267 {
00268 if ( ! localWellsActive() ) {
00269 return;
00270 }
00271
00272 if( scaleAddRes_.size() != Ax.size() ) {
00273 scaleAddRes_.resize( Ax.size() );
00274 }
00275
00276 scaleAddRes_ = 0.0;
00277
00278 apply( x, scaleAddRes_ );
00279
00280 Ax.axpy( alpha, scaleAddRes_ );
00281 }
00282
00283
00284
00285
00286
00287 template<typename TypeTag>
00288 void
00289 BlackoilWellModel<TypeTag>::
00290 recoverWellSolutionAndUpdateWellState(const BVector& x, WellState& well_state) const
00291 {
00292 for (auto& well : well_container_) {
00293 well->recoverWellSolutionAndUpdateWellState(x, well_state);
00294 }
00295 }
00296
00297
00298
00299
00300
00301 template<typename TypeTag>
00302 int
00303 BlackoilWellModel<TypeTag>::
00304 flowPhaseToEbosPhaseIdx( const int phaseIdx ) const
00305 {
00306 const auto& pu = phase_usage_;
00307 if (active_[Water] && pu.phase_pos[Water] == phaseIdx)
00308 return FluidSystem::waterPhaseIdx;
00309 if (active_[Oil] && pu.phase_pos[Oil] == phaseIdx)
00310 return FluidSystem::oilPhaseIdx;
00311 if (active_[Gas] && pu.phase_pos[Gas] == phaseIdx)
00312 return FluidSystem::gasPhaseIdx;
00313
00314 assert(phaseIdx < 3);
00315
00316 return phaseIdx;
00317 }
00318
00319
00320
00321
00322
00323 template<typename TypeTag>
00324 int
00325 BlackoilWellModel<TypeTag>::
00326 numPhases() const
00327 {
00328 return number_of_phases_;
00329 }
00330
00331
00332
00333
00334
00335 template<typename TypeTag>
00336 void
00337 BlackoilWellModel<TypeTag>::
00338 resetWellControlFromState(const WellState& xw) const
00339 {
00340 const int nw = wells_->number_of_wells;
00341 for (int w = 0; w < nw; ++w) {
00342 WellControls* wc = wells_->ctrls[w];
00343 well_controls_set_current( wc, xw.currentControls()[w]);
00344 }
00345 }
00346
00347
00348
00349
00350
00351 template<typename TypeTag>
00352 bool
00353 BlackoilWellModel<TypeTag>::
00354 wellsActive() const
00355 {
00356 return wells_active_;
00357 }
00358
00359
00360
00361
00362
00363 template<typename TypeTag>
00364 void
00365 BlackoilWellModel<TypeTag>::
00366 setWellsActive(const bool wells_active)
00367 {
00368 wells_active_ = wells_active;
00369 }
00370
00371
00372
00373
00374
00375 template<typename TypeTag>
00376 bool
00377 BlackoilWellModel<TypeTag>::
00378 localWellsActive() const
00379 {
00380 return number_of_wells_ > 0;
00381 }
00382
00383
00384
00385
00386
00387 template<typename TypeTag>
00388 void
00389 BlackoilWellModel<TypeTag>::
00390 initPrimaryVariablesEvaluation() const
00391 {
00392 for (auto& well : well_container_) {
00393 well->initPrimaryVariablesEvaluation();
00394 }
00395 }
00396
00397
00398
00399
00400
00401 template<typename TypeTag>
00402 SimulatorReport
00403 BlackoilWellModel<TypeTag>::
00404 solveWellEq(Simulator& ebosSimulator,
00405 const double dt,
00406 WellState& well_state) const
00407 {
00408 const int nw = number_of_wells_;
00409 WellState well_state0 = well_state;
00410
00411 const int numComp = numComponents();
00412 std::vector< Scalar > B_avg( numComp, Scalar() );
00413 computeAverageFormationFactor(ebosSimulator, B_avg);
00414
00415 const int max_iter = param_.max_welleq_iter_;
00416
00417 int it = 0;
00418 bool converged;
00419 do {
00420 assembleWellEq(ebosSimulator, dt, well_state, true);
00421
00422 converged = getWellConvergence(ebosSimulator, B_avg);
00423
00424
00425 if (wellCollection()->groupControlActive()) {
00426 converged = converged && wellCollection()->groupTargetConverged(well_state.wellRates());
00427 }
00428
00429 if (converged) {
00430 break;
00431 }
00432
00433 ++it;
00434 if( localWellsActive() )
00435 {
00436 for (auto& well : well_container_) {
00437 well->solveEqAndUpdateWellState(well_state);
00438 }
00439 }
00440
00441
00442
00443 if( wellsActive() )
00444 {
00445 updateWellControls(well_state);
00446 initPrimaryVariablesEvaluation();
00447 }
00448 } while (it < max_iter);
00449
00450 if (converged) {
00451 if ( terminal_output_ ) {
00452 OpmLog::debug("Well equation solution gets converged with " + std::to_string(it) + " iterations");
00453 }
00454 } else {
00455 if ( terminal_output_ ) {
00456 OpmLog::debug("Well equation solution failed in getting converged with " + std::to_string(it) + " iterations");
00457 }
00458
00459 well_state = well_state0;
00460 updatePrimaryVariables(well_state);
00461
00462 for (int w = 0; w < nw; ++w) {
00463 WellControls* wc = well_container_[w]->wellControls();
00464 well_controls_set_current(wc, well_state.currentControls()[w]);
00465 }
00466 }
00467
00468 SimulatorReport report;
00469 report.converged = converged;
00470 report.total_well_iterations = it;
00471 return report;
00472 }
00473
00474
00475
00476
00477
00478 template<typename TypeTag>
00479 bool
00480 BlackoilWellModel<TypeTag>::
00481 getWellConvergence(const Simulator& ebosSimulator,
00482 const std::vector<Scalar>& B_avg) const
00483 {
00484 ConvergenceReport report;
00485
00486 for (const auto& well : well_container_) {
00487 report += well->getWellConvergence(B_avg);
00488 }
00489
00490
00491 {
00492 bool nan_residual_found = report.nan_residual_found;
00493 const auto& grid = ebosSimulator.gridManager().grid();
00494 int value = nan_residual_found ? 1 : 0;
00495
00496 nan_residual_found = grid.comm().max(value);
00497
00498 if (nan_residual_found) {
00499 for (const auto& well : report.nan_residual_wells) {
00500 OpmLog::debug("NaN residual found with phase " + well.phase_name + " for well " + well.well_name);
00501 }
00502 OPM_THROW(Opm::NumericalProblem, "NaN residual found!");
00503 }
00504 }
00505
00506
00507 {
00508 bool too_large_residual_found = report.too_large_residual_found;
00509 const auto& grid = ebosSimulator.gridManager().grid();
00510 int value = too_large_residual_found ? 1 : 0;
00511
00512 too_large_residual_found = grid.comm().max(value);
00513 if (too_large_residual_found) {
00514 for (const auto& well : report.too_large_residual_wells) {
00515 OpmLog::debug("Too large residual found with phase " + well.phase_name + " fow well " + well.well_name);
00516 }
00517 OPM_THROW(Opm::NumericalProblem, "Too large residual found!");
00518 }
00519 }
00520
00521
00522 bool converged_well = report.converged;
00523 {
00524 const auto& grid = ebosSimulator.gridManager().grid();
00525 int value = converged_well ? 1 : 0;
00526
00527 converged_well = grid.comm().min(value);
00528 }
00529
00530 return converged_well;
00531 }
00532
00533
00534
00535
00536
00537 template<typename TypeTag>
00538 void
00539 BlackoilWellModel<TypeTag>::
00540 calculateExplicitQuantities(const Simulator& ebosSimulator,
00541 const WellState& xw) const
00542 {
00543 for (auto& well : well_container_) {
00544 well->calculateExplicitQuantities(ebosSimulator, xw);
00545 }
00546 }
00547
00548
00549
00550
00551
00552 template<typename TypeTag>
00553 void
00554 BlackoilWellModel<TypeTag>::
00555 updateWellControls(WellState& xw) const
00556 {
00557
00558
00559
00560
00561 if( !wellsActive() ) return ;
00562
00563 wellhelpers::WellSwitchingLogger logger;
00564
00565 for (const auto& well : well_container_) {
00566 well->updateWellControl(xw, logger);
00567 }
00568
00569 updateGroupControls(xw);
00570 }
00571
00572
00573
00574
00575
00576 template<typename TypeTag>
00577 void
00578 BlackoilWellModel<TypeTag>::
00579 updateListEconLimited(const Schedule& schedule,
00580 const int current_step,
00581 const Wells* wells_struct,
00582 const WellState& well_state,
00583 DynamicListEconLimited& list_econ_limited) const
00584 {
00585 for (const auto& well : well_container_) {
00586 well->updateListEconLimited(well_state, list_econ_limited);
00587 }
00588 }
00589
00590
00591
00592
00593
00594 template<typename TypeTag>
00595 void
00596 BlackoilWellModel<TypeTag>::
00597 computeWellPotentials(const Simulator& ebosSimulator,
00598 const WellState& well_state,
00599 std::vector<double>& well_potentials) const
00600 {
00601
00602 const int nw = number_of_wells_;
00603 const int np = number_of_phases_;
00604 well_potentials.resize(nw * np, 0.0);
00605
00606 for (int w = 0; w < nw; ++w) {
00607 std::vector<double> potentials;
00608 well_container_[w]->computeWellPotentials(ebosSimulator, well_state, potentials);
00609
00610
00611 for (int p = 0; p < np; ++p) {
00612 well_potentials[w * np + p] = std::abs(potentials[p]);
00613 }
00614 }
00615 }
00616
00617
00618
00619
00620
00621 template<typename TypeTag>
00622 void
00623 BlackoilWellModel<TypeTag>::
00624 prepareTimeStep(const Simulator& ebos_simulator,
00625 WellState& well_state) const
00626 {
00627
00628
00629
00630
00631
00632 resetWellControlFromState(well_state);
00633
00634
00635 prepareGroupControl(ebos_simulator, well_state);
00636
00637
00638 for (int w = 0; w < number_of_wells_; ++w) {
00639 WellControls* wc = well_container_[w]->wellControls();
00640 const int control = well_controls_get_current(wc);
00641 well_state.currentControls()[w] = control;
00642
00643
00644 well_container_[w]->updateWellStateWithTarget(control, well_state);
00645
00646
00647
00648 if (well_state.isNewWell(w) ) {
00649 well_state.setNewWell(w, false);
00650 }
00651 }
00652 }
00653
00654
00655
00656
00657
00658 template<typename TypeTag>
00659 void
00660 BlackoilWellModel<TypeTag>::
00661 prepareGroupControl(const Simulator& ebos_simulator,
00662 WellState& well_state) const
00663 {
00664
00665 if (well_collection_->groupControlActive()) {
00666 for (int w = 0; w < number_of_wells_; ++w) {
00667 WellControls* wc = well_container_[w]->wellControls();
00668 WellNode& well_node = well_collection_->findWellNode(well_container_[w]->name());
00669
00670
00671
00672
00673 int ctrl_index = well_controls_get_current(wc);
00674
00675 const int group_control_index = well_node.groupControlIndex();
00676 if (group_control_index >= 0 && ctrl_index < 0) {
00677
00678 well_controls_set_current(wc, group_control_index);
00679 well_state.currentControls()[w] = group_control_index;
00680 }
00681
00682
00683
00684 ctrl_index = well_controls_get_current(wc);
00685 if (well_node.groupControlIndex() >= 0 && ctrl_index == well_node.groupControlIndex()) {
00686
00687 well_node.setIndividualControl(false);
00688 } else {
00689
00690 well_node.setIndividualControl(true);
00691 }
00692 }
00693
00694 if (well_collection_->requireWellPotentials()) {
00695
00696
00697 std::vector<double> well_potentials;
00698 computeWellPotentials(ebos_simulator, well_state, well_potentials);
00699
00700
00701
00702
00703 well_collection_->setGuideRatesWithPotentials(wells_, phase_usage_, well_potentials);
00704 }
00705
00706 applyVREPGroupControl(well_state);
00707
00708 if (!wellCollection()->groupControlApplied()) {
00709 wellCollection()->applyGroupControls();
00710 } else {
00711 wellCollection()->updateWellTargets(well_state.wellRates());
00712 }
00713 }
00714 }
00715
00716
00717
00718
00719
00720 template<typename TypeTag>
00721 WellCollection*
00722 BlackoilWellModel<TypeTag>::
00723 wellCollection() const
00724 {
00725 return well_collection_;
00726 }
00727
00728
00729
00730
00731
00732 template<typename TypeTag>
00733 void
00734 BlackoilWellModel<TypeTag>::
00735 calculateEfficiencyFactors()
00736 {
00737 if ( !localWellsActive() ) {
00738 return;
00739 }
00740
00741 const int nw = number_of_wells_;
00742
00743 for (int w = 0; w < nw; ++w) {
00744 const std::string well_name = well_container_[w]->name();
00745 const WellNode& well_node = wellCollection()->findWellNode(well_name);
00746
00747 const double well_efficiency_factor = well_node.getAccumulativeEfficiencyFactor();
00748
00749 well_container_[w]->setWellEfficiencyFactor(well_efficiency_factor);
00750 }
00751 }
00752
00753
00754
00755
00756
00757 template<typename TypeTag>
00758 void
00759 BlackoilWellModel<TypeTag>::
00760 computeWellVoidageRates(const WellState& well_state,
00761 std::vector<double>& well_voidage_rates,
00762 std::vector<double>& voidage_conversion_coeffs) const
00763 {
00764 if ( !localWellsActive() ) {
00765 return;
00766 }
00767
00768
00769
00770
00771
00772 const int nw = numWells();
00773 const int np = numPhases();
00774
00775
00776 well_voidage_rates.resize(nw, 0);
00777
00778 voidage_conversion_coeffs.resize(nw * np, 1.0);
00779
00780 std::vector<double> well_rates(np, 0.0);
00781 std::vector<double> convert_coeff(np, 1.0);
00782
00783 for (int w = 0; w < nw; ++w) {
00784 const bool is_producer = well_container_[w]->wellType() == PRODUCER;
00785 const int well_cell_top = well_container_[w]->cells()[0];
00786 const int pvtRegionIdx = pvt_region_idx_[well_cell_top];
00787
00788
00789 if (is_producer) {
00790 std::transform(well_state.wellRates().begin() + np * w,
00791 well_state.wellRates().begin() + np * (w + 1),
00792 well_rates.begin(), std::negate<double>());
00793
00794
00795 const int fipreg = 0;
00796
00797 rate_converter_.calcCoeff(fipreg, pvtRegionIdx, convert_coeff);
00798 well_voidage_rates[w] = std::inner_product(well_rates.begin(), well_rates.end(),
00799 convert_coeff.begin(), 0.0);
00800 } else {
00801
00802
00803 std::copy(well_state.wellRates().begin() + np * w,
00804 well_state.wellRates().begin() + np * (w + 1),
00805 well_rates.begin());
00806
00807 const int fipreg = 0;
00808 rate_converter_.calcCoeff(fipreg, pvtRegionIdx, convert_coeff);
00809 std::copy(convert_coeff.begin(), convert_coeff.end(),
00810 voidage_conversion_coeffs.begin() + np * w);
00811 }
00812 }
00813 }
00814
00815
00816
00817
00818
00819 template<typename TypeTag>
00820 void
00821 BlackoilWellModel<TypeTag>::
00822 applyVREPGroupControl(WellState& well_state) const
00823 {
00824 if ( wellCollection()->havingVREPGroups() ) {
00825 std::vector<double> well_voidage_rates;
00826 std::vector<double> voidage_conversion_coeffs;
00827 computeWellVoidageRates(well_state, well_voidage_rates, voidage_conversion_coeffs);
00828 wellCollection()->applyVREPGroupControls(well_voidage_rates, voidage_conversion_coeffs);
00829
00830
00831 for (const WellNode* well_node : wellCollection()->getLeafNodes()) {
00832 if (well_node->isInjector() && !well_node->individualControl()) {
00833 const int well_index = well_node->selfIndex();
00834 well_state.currentControls()[well_index] = well_node->groupControlIndex();
00835
00836 WellControls* wc = well_container_[well_index]->wellControls();
00837 well_controls_set_current(wc, well_node->groupControlIndex());
00838 }
00839 }
00840 }
00841 }
00842
00843
00844
00845
00846
00847 template<typename TypeTag>
00848 void
00849 BlackoilWellModel<TypeTag>::
00850 updateGroupControls(WellState& well_state) const
00851 {
00852
00853 if (wellCollection()->groupControlActive()) {
00854 for (int w = 0; w < number_of_wells_; ++w) {
00855
00856
00857 WellNode& well_node = well_collection_->findWellNode(well_container_[w]->name());
00858
00859
00860 const int current = well_state.currentControls()[w];
00861 if (well_node.groupControlIndex() >= 0 && current == well_node.groupControlIndex()) {
00862
00863 well_node.setIndividualControl(false);
00864 } else {
00865
00866 well_node.setIndividualControl(true);
00867 }
00868 }
00869
00870 applyVREPGroupControl(well_state);
00871
00872
00873 wellCollection()->updateWellTargets(well_state.wellRates());
00874
00875 for (int w = 0; w < number_of_wells_; ++w) {
00876
00877
00878
00879
00880 const int current = well_state.currentControls()[w];
00881 well_container_[w]->updateWellStateWithTarget(current, well_state);
00882 }
00883 }
00884 }
00885
00886
00887
00888
00889
00890 template<typename TypeTag>
00891 void
00892 BlackoilWellModel<TypeTag>::
00893 setupCompressedToCartesian(const int* global_cell, int number_of_cells, std::map<int,int>& cartesian_to_compressed ) const
00894 {
00895 if (global_cell) {
00896 for (int i = 0; i < number_of_cells; ++i) {
00897 cartesian_to_compressed.insert(std::make_pair(global_cell[i], i));
00898 }
00899 }
00900 else {
00901 for (int i = 0; i < number_of_cells; ++i) {
00902 cartesian_to_compressed.insert(std::make_pair(i, i));
00903 }
00904 }
00905
00906 }
00907
00908 template<typename TypeTag>
00909 void
00910 BlackoilWellModel<TypeTag>::
00911 computeRepRadiusPerfLength(const Grid& grid)
00912 {
00913
00914
00915 const int* global_cell = Opm::UgGridHelpers::globalCell(grid);
00916
00917 std::map<int,int> cartesian_to_compressed;
00918 setupCompressedToCartesian(global_cell, number_of_cells_,
00919 cartesian_to_compressed);
00920
00921 for (const auto& well : well_container_) {
00922 well->computeRepRadiusPerfLength(grid, cartesian_to_compressed);
00923 }
00924 }
00925
00926
00927
00928
00929
00930 template<typename TypeTag>
00931 void
00932 BlackoilWellModel<TypeTag>::
00933 computeAverageFormationFactor(const Simulator& ebosSimulator,
00934 std::vector<double>& B_avg) const
00935 {
00936 const int np = numPhases();
00937
00938 const auto& grid = ebosSimulator.gridManager().grid();
00939 const auto& gridView = grid.leafGridView();
00940 ElementContext elemCtx(ebosSimulator);
00941 const auto& elemEndIt = gridView.template end<0, Dune::Interior_Partition>();
00942
00943 for (auto elemIt = gridView.template begin</*codim=*/0, Dune::Interior_Partition>();
00944 elemIt != elemEndIt; ++elemIt)
00945 {
00946 elemCtx.updatePrimaryStencil(*elemIt);
00947 elemCtx.updatePrimaryIntensiveQuantities(0);
00948
00949 const auto& intQuants = elemCtx.intensiveQuantities(0, 0);
00950 const auto& fs = intQuants.fluidState();
00951
00952 for ( int phaseIdx = 0; phaseIdx < np; ++phaseIdx )
00953 {
00954 auto& B = B_avg[ phaseIdx ];
00955 const int ebosPhaseIdx = flowPhaseToEbosPhaseIdx(phaseIdx);
00956
00957 B += 1 / fs.invB(ebosPhaseIdx).value();
00958 }
00959 if (has_solvent_) {
00960 auto& B = B_avg[solventSaturationIdx];
00961 B += 1 / intQuants.solventInverseFormationVolumeFactor().value();
00962 }
00963 }
00964
00965
00966 grid.comm().sum(B_avg.data(), B_avg.size());
00967 for(auto& bval: B_avg)
00968 {
00969 bval/=global_nc_;
00970 }
00971 }
00972
00973
00974
00975
00976
00977 template<typename TypeTag>
00978 void
00979 BlackoilWellModel<TypeTag>::
00980 updatePrimaryVariables(const WellState& well_state) const
00981 {
00982 for (const auto& well : well_container_) {
00983 well->updatePrimaryVariables(well_state);
00984 }
00985 }
00986
00987 }