My Project
Loading...
Searching...
No Matches
BlackoilWellModel.hpp
1/*
2 Copyright 2016 SINTEF ICT, Applied Mathematics.
3 Copyright 2016 - 2017 Statoil ASA.
4 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5 Copyright 2016 - 2018 IRIS AS
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21*/
22
23
24#ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
25#define OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
26
27#include <ebos/eclproblem.hh>
28#include <opm/common/OpmLog/OpmLog.hpp>
29
30#include <cassert>
31#include <cstddef>
32#include <map>
33#include <memory>
34#include <optional>
35#include <set>
36#include <string>
37#include <tuple>
38#include <unordered_map>
39#include <vector>
40
41#include <opm/input/eclipse/Schedule/Group/Group.hpp>
42#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
43#include <opm/input/eclipse/Schedule/Schedule.hpp>
44#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
45
46#include <opm/simulators/flow/countGlobalCells.hpp>
47#include <opm/simulators/flow/SubDomain.hpp>
48
49#include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
50#include <opm/simulators/wells/BlackoilWellModelGuideRates.hpp>
51#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
52#include <opm/simulators/wells/GasLiftSingleWell.hpp>
53#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
54#include <opm/simulators/wells/GasLiftStage2.hpp>
55#include <opm/simulators/wells/GasLiftWellState.hpp>
56#include <opm/simulators/wells/MultisegmentWell.hpp>
57#include <opm/simulators/wells/ParallelWBPCalculation.hpp>
58#include <opm/simulators/wells/ParallelWellInfo.hpp>
59#include <opm/simulators/wells/PerforationData.hpp>
62#include <opm/simulators/wells/StandardWell.hpp>
63#include <opm/simulators/wells/VFPInjProperties.hpp>
64#include <opm/simulators/wells/VFPProdProperties.hpp>
65#include <opm/simulators/wells/WGState.hpp>
66#include <opm/simulators/wells/WellGroupHelpers.hpp>
67#include <opm/simulators/wells/WellInterface.hpp>
68#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
69#include <opm/simulators/wells/WellState.hpp>
70
71#include <opm/simulators/timestepping/SimulatorReport.hpp>
72#include <opm/simulators/timestepping/gatherConvergenceReport.hpp>
73
74#include <dune/common/fmatrix.hh>
75#include <dune/istl/bcrsmatrix.hh>
76#include <dune/istl/matrixmatrix.hh>
77
78#include <opm/material/densead/Math.hpp>
79
80#include <opm/simulators/utils/DeferredLogger.hpp>
81
82namespace Opm::Properties {
83
84template<class TypeTag, class MyTypeTag>
88
89} // namespace Opm::Properties
90
91namespace Opm {
92
94 template<typename TypeTag>
95 class BlackoilWellModel : public BaseAuxiliaryModule<TypeTag>
97 {
98 public:
99 // --------- Types ---------
101
112 using GasLiftSingleWell = typename WellInterface<TypeTag>::GasLiftSingleWell;
113 using GLiftOptWells = typename BlackoilWellModelGeneric::GLiftOptWells;
114 using GLiftProdWells = typename BlackoilWellModelGeneric::GLiftProdWells;
115 using GLiftWellStateMap =
116 typename BlackoilWellModelGeneric::GLiftWellStateMap;
117 using GLiftEclWells = typename GasLiftGroupInfo::GLiftEclWells;
118 using GLiftSyncGroups = typename GasLiftSingleWellGeneric::GLiftSyncGroups;
119 constexpr static std::size_t pressureVarIndex = GetPropType<TypeTag, Properties::Indices>::pressureSwitchIdx;
120 typedef typename BaseAuxiliaryModule<TypeTag>::NeighborSet NeighborSet;
121
122 static const int numEq = Indices::numEq;
123 static const int solventSaturationIdx = Indices::solventSaturationIdx;
124 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
125 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
126 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
127 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
128
129 // TODO: where we should put these types, WellInterface or Well Model?
130 // or there is some other strategy, like TypeTag
131 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
132 typedef Dune::BlockVector<VectorBlockType> BVector;
133
136
137 // For the conversion between the surface volume rate and resrevoir voidage rate
138 using RateConverterType = RateConverter::
139 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
140
141 // For computing average pressured used by gpmaint
142 using AverageRegionalPressureType = RegionAverageCalculator::
143 AverageRegionalPressure<FluidSystem, std::vector<int> >;
144
145 using Domain = SubDomain<Grid>;
146
147 BlackoilWellModel(Simulator& ebosSimulator);
148
149 void init();
150 void initWellContainer(const int reportStepIdx) override;
151
153 // <eWoms auxiliary module stuff>
155 unsigned numDofs() const override
156 // No extra dofs are inserted for wells. (we use a Schur complement.)
157 { return 0; }
158
159 void addNeighbors(std::vector<NeighborSet>& neighbors) const override;
160
161 void applyInitial() override
162 {}
163
164 void linearize(SparseMatrixAdapter& jacobian, GlobalEqVector& res) override;
165 void linearizeDomain(const Domain& domain, SparseMatrixAdapter& jacobian, GlobalEqVector& res);
166
167 void postSolve(GlobalEqVector& deltaX) override
168 {
169 recoverWellSolutionAndUpdateWellState(deltaX);
170 }
171
172 void postSolveDomain(GlobalEqVector& deltaX, const Domain& domain)
173 {
174 recoverWellSolutionAndUpdateWellStateDomain(deltaX, domain);
175 }
176
178 // </ eWoms auxiliary module stuff>
180
181 template <class Restarter>
182 void deserialize(Restarter& /* res */)
183 {
184 // TODO (?)
185 }
186
191 template <class Restarter>
192 void serialize(Restarter& /* res*/)
193 {
194 // TODO (?)
195 }
196
197 void beginEpisode()
198 {
200 beginReportStep(ebosSimulator_.episodeIndex());
201 }
202
203 void beginTimeStep();
204
205 void beginIteration()
206 {
207 OPM_TIMEBLOCK(beginIteration);
208 assemble(ebosSimulator_.model().newtonMethod().numIterations(),
209 ebosSimulator_.timeStepSize());
210 }
211
212 void endIteration()
213 { }
214
215 void endTimeStep()
216 {
217 OPM_TIMEBLOCK(endTimeStep);
218 timeStepSucceeded(ebosSimulator_.time(), ebosSimulator_.timeStepSize());
219 }
220
221 void endEpisode()
222 {
223 endReportStep();
224 }
225
226 void computeTotalRatesForDof(RateVector& rate,
227 unsigned globalIdx) const;
228
229 template <class Context>
230 void computeTotalRatesForDof(RateVector& rate,
231 const Context& context,
232 unsigned spaceIdx,
233 unsigned timeIdx) const;
234
235
236 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
237
238 using BlackoilWellModelGeneric::initFromRestartFile;
239 void initFromRestartFile(const RestartValue& restartValues)
240 {
241 initFromRestartFile(restartValues,
242 this->ebosSimulator_.vanguard().transferWTestState(),
243 grid().size(0),
244 param_.use_multisegment_well_);
245 }
246
247 using BlackoilWellModelGeneric::prepareDeserialize;
248 void prepareDeserialize(const int report_step)
249 {
250 prepareDeserialize(report_step, grid().size(0),
251 param_.use_multisegment_well_);
252 }
253
254 data::Wells wellData() const
255 {
256 auto wsrpt = this->wellState()
257 .report(ebosSimulator_.vanguard().globalCell().data(),
258 [this](const int well_index) -> bool
259 {
260 return this->wasDynamicallyShutThisTimeStep(well_index);
261 });
262
263 const auto& tracerRates = ebosSimulator_.problem().tracerModel().getWellTracerRates();
264 this->assignWellTracerRates(wsrpt, tracerRates);
265
266
267 BlackoilWellModelGuideRates(*this).assignWellGuideRates(wsrpt, this->reportStepIndex());
268 this->assignShutConnections(wsrpt, this->reportStepIndex());
269
270 return wsrpt;
271 }
272
273 data::WellBlockAveragePressures wellBlockAveragePressures() const
274 {
275 return this->computeWellBlockAveragePressures();
276 }
277
278 // subtract Binv(D)rw from r;
279 void apply( BVector& r) const;
280
281 // subtract B*inv(D)*C * x from A*x
282 void apply(const BVector& x, BVector& Ax) const;
283
284 // accumulate the contributions of all Wells in the WellContributions object
285 void getWellContributions(WellContributions& x) const;
286
287 // apply well model with scaling of alpha
288 void applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const;
289
290 // Check if well equations is converged.
291 ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg, const bool checkWellGroupControls = false) const;
292
293 // Check if well equations are converged locally.
294 ConvergenceReport getDomainWellConvergence(const Domain& domain,
295 const std::vector<Scalar>& B_avg) const;
296
297 const SimulatorReportSingle& lastReport() const;
298
299 void addWellContributions(SparseMatrixAdapter& jacobian) const;
300
301 // add source from wells to the reservoir matrix
302 void addReservoirSourceTerms(GlobalEqVector& residual,
303 std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress) const;
304
305 // called at the beginning of a report step
306 void beginReportStep(const int time_step);
307
308 // it should be able to go to prepareTimeStep(), however, the updateWellControls() and initPrimaryVariablesEvaluation()
309 // makes it a little more difficult. unless we introduce if (iterationIdx != 0) to avoid doing the above functions
310 // twice at the beginning of the time step
313 void calculateExplicitQuantities(DeferredLogger& deferred_logger) const;
314 // some preparation work, mostly related to group control and RESV,
315 // at the beginning of each time step (Not report step)
316 void prepareTimeStep(DeferredLogger& deferred_logger);
317 void initPrimaryVariablesEvaluation() const;
318 void initPrimaryVariablesEvaluationDomain(const Domain& domain) const;
319
320 std::pair<bool, bool>
321 updateWellControls(const bool mandatory_network_balance, DeferredLogger& deferred_logger, const bool relax_network_tolerance = false);
322
323 void updateAndCommunicate(const int reportStepIdx,
324 const int iterationIdx,
325 DeferredLogger& deferred_logger);
326
327 bool updateGroupControls(const Group& group,
328 DeferredLogger& deferred_logger,
329 const int reportStepIdx,
330 const int iterationIdx);
331
332 WellInterfacePtr getWell(const std::string& well_name) const;
333 bool hasWell(const std::string& well_name) const;
334
335 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
336
337 void addWellPressureEquations(PressureMatrix& jacobian, const BVector& weights,const bool use_well_weights) const;
338
339 void addWellPressureEquationsStruct(PressureMatrix& jacobian) const;
340
341 void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
342
344 const std::vector<WellInterfacePtr>& localNonshutWells() const
345 {
346 return well_container_;
347 }
348
349 // prototype for assemble function for ASPIN solveLocal()
350 // will try to merge back to assemble() when done prototyping
351 void assembleDomain(const int iterationIdx,
352 const double dt,
353 const Domain& domain);
354 void updateWellControlsDomain(DeferredLogger& deferred_logger, const Domain& domain);
355
356 void logPrimaryVars() const;
357 std::vector<double> getPrimaryVarsDomain(const Domain& domain) const;
358 void setPrimaryVarsDomain(const Domain& domain, const std::vector<double>& vars);
359
360 void setupDomains(const std::vector<Domain>& domains);
361
362 protected:
363 Simulator& ebosSimulator_;
364
365 // a vector of all the wells.
366 std::vector<WellInterfacePtr> well_container_{};
367
368 std::vector<bool> is_cell_perforated_{};
369
370 void initializeWellState(const int timeStepIdx);
371
372 // create the well container
373 void createWellContainer(const int time_step) override;
374
375 WellInterfacePtr
376 createWellPointer(const int wellID,
377 const int time_step) const;
378
379 template <typename WellType>
380 std::unique_ptr<WellType>
381 createTypedWellPointer(const int wellID,
382 const int time_step) const;
383
384 WellInterfacePtr createWellForWellTest(const std::string& well_name, const int report_step, DeferredLogger& deferred_logger) const;
385
386
387 const ModelParameters param_;
388 std::size_t global_num_cells_{};
389 // the number of the cells in the local grid
390 std::size_t local_num_cells_{};
391 double gravity_{};
392 std::vector<double> depth_{};
393 bool alternative_well_rate_init_{};
394
395 std::unique_ptr<RateConverterType> rateConverter_{};
396 std::map<std::string, std::unique_ptr<AverageRegionalPressureType>> regionalAveragePressureCalculator_{};
397
399 {
400 std::optional<typename std::vector<WellInterfacePtr>::size_type> openWellIdx_{};
401 std::size_t wbpCalcIdx_{};
402 };
403
404 std::vector<WBPCalcID> wbpCalcMap_{};
405
406 SimulatorReportSingle last_report_{};
407
408 // solve to get a good network solution, group and well states might be updated during the process.
409 // the reservoir should stay static during this solution procedure.
410 void balanceNetwork(DeferredLogger& deferred_logger);
411
412 // used to better efficiency of calcuation
413 mutable BVector scaleAddRes_{};
414
415 std::vector<Scalar> B_avg_{};
416
417 // Keep track of the domain of each well, if using subdomains.
418 std::map<std::string, int> well_domain_;
419
420 const Grid& grid() const
421 { return ebosSimulator_.vanguard().grid(); }
422
423 const EquilGrid& equilGrid() const
424 { return ebosSimulator_.vanguard().equilGrid(); }
425
426 const EclipseState& eclState() const
427 { return ebosSimulator_.vanguard().eclState(); }
428
429 // compute the well fluxes and assemble them in to the reservoir equations as source terms
430 // and in the well equations.
431 void assemble(const int iterationIdx,
432 const double dt);
433
434 // well controls and network pressures affect each other and are solved in an iterative manner.
435 // the function handles one iteration of updating well controls and network pressures.
436 // it is possible to decouple the update of well controls and network pressures further.
437 // the returned two booleans are {continue_due_to_network, well_group_control_changed}, respectively
438 std::pair<bool, bool> updateWellControlsAndNetworkIteration(const bool mandatory_network_balance,
439 const bool relax_network_tolerance,
440 const double dt,
441 DeferredLogger& local_deferredLogger);
442
443 bool updateWellControlsAndNetwork(const bool mandatory_network_balance,
444 const double dt,
445 DeferredLogger& local_deferredLogger);
446
455 void initializeLocalWellStructure(const int reportStepIdx,
456 const bool enableWellPIScaling);
457
461 void initializeGroupStructure(const int reportStepIdx);
462
463 // called at the end of a time step
464 void timeStepSucceeded(const double simulationTime, const double dt);
465
466 // called at the end of a report step
467 void endReportStep();
468
469 // using the solution x to recover the solution xw for wells and applying
470 // xw to update Well State
471 void recoverWellSolutionAndUpdateWellState(const BVector& x);
472
473 // using the solution x to recover the solution xw for wells and applying
474 // xw to update Well State
475 void recoverWellSolutionAndUpdateWellStateDomain(const BVector& x, const Domain& domain);
476
477 // setting the well_solutions_ based on well_state.
478 void updatePrimaryVariables(DeferredLogger& deferred_logger);
479
480 void initializeWBPCalculationService();
481
482 data::WellBlockAveragePressures
483 computeWellBlockAveragePressures() const;
484
486 makeWellSourceEvaluatorFactory(const std::vector<Well>::size_type wellIdx) const;
487
488 void registerOpenWellsForWBPCalculation();
489
490 void updateAverageFormationFactor();
491
492 void computePotentials(const std::size_t widx,
493 const WellState& well_state_copy,
494 std::string& exc_msg,
495 ExceptionType::ExcEnum& exc_type,
496 DeferredLogger& deferred_logger) override;
497
498 const std::vector<double>& wellPerfEfficiencyFactors() const;
499
500 void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger& deferred_logger) override;
501 void calculateProductivityIndexValues(DeferredLogger& deferred_logger) override;
502 void calculateProductivityIndexValues(const WellInterface<TypeTag>* wellPtr,
503 DeferredLogger& deferred_logger);
504
505 // The number of components in the model.
506 int numComponents() const;
507
508 int reportStepIndex() const;
509
510 void assembleWellEq(const double dt, DeferredLogger& deferred_logger);
511 void assembleWellEqDomain(const double dt, const Domain& domain, DeferredLogger& deferred_logger);
512
513 void prepareWellsBeforeAssembling(const double dt, DeferredLogger& deferred_logger);
514
515 // TODO: finding a better naming
516 void assembleWellEqWithoutIteration(const double dt, DeferredLogger& deferred_logger);
517
518 bool maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
519
520 void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
521 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
522 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
523
524 // cannot be const since it accesses the non-const WellState
525 void gasLiftOptimizationStage1SingleWell(WellInterface<TypeTag> *well,
526 DeferredLogger& deferred_logger,
527 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
528 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map,
529 GLiftSyncGroups& groups_to_sync);
530
531 void extractLegacyCellPvtRegionIndex_();
532
533 void extractLegacyDepth_();
534
536 void updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const;
537
538 void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger& deferred_logger);
539
540 void calcRates(const int fipnum,
541 const int pvtreg,
542 const std::vector<double>& production_rates,
543 std::vector<double>& resv_coeff) override;
544
545 void calcInjRates(const int fipnum,
546 const int pvtreg,
547 std::vector<double>& resv_coeff) override;
548
549 void computeWellTemperature();
550
552 return ebosSimulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
553 }
554
555 private:
556 BlackoilWellModel(Simulator& ebosSimulator, const PhaseUsage& pu);
557 };
558
559
560} // namespace Opm
561
562#include "BlackoilWellModel_impl.hpp"
563#endif
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Definition AquiferInterface.hpp:35
Class for handling the blackoil well model.
Definition BlackoilWellModelGeneric.hpp:82
Class for handling the blackoil well model.
Definition BlackoilWellModel.hpp:97
void initializeGroupStructure(const int reportStepIdx)
Initialize group control modes/constraints and group solution state.
Definition BlackoilWellModel_impl.hpp:329
void serialize(Restarter &)
This method writes the complete state of the well to the harddisk.
Definition BlackoilWellModel.hpp:192
const std::vector< WellInterfacePtr > & localNonshutWells() const
Get list of local nonshut wells.
Definition BlackoilWellModel.hpp:344
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition BlackoilWellModel_impl.hpp:1821
void initializeLocalWellStructure(const int reportStepIdx, const bool enableWellPIScaling)
Update rank's notion of intersecting wells and their associate solution variables.
Definition BlackoilWellModel_impl.hpp:284
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition BlackoilWellModel_impl.hpp:2128
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition BlackoilWellModel.hpp:551
Definition DeferredLogger.hpp:57
std::function< Evaluator()> EvaluatorFactory
Callback for constructing a source term evaluation function on the current MPI rank.
Definition ParallelWBPCalculation.hpp:62
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition RateConverter.hpp:70
Computes hydrocarbon weighed average pressures over regions.
Definition RegionAverageCalculator.hpp:65
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Definition BlackoilWellModel.hpp:399
Definition BlackoilPhases.hpp:46
Definition BlackoilWellModel.hpp:85
A struct for returning timing data from a simulator to its caller.
Definition SimulatorReport.hpp:34
Representing a part of a grid, in a way suitable for performing local solves.
Definition SubDomain.hpp:46