My Project
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/Schedule.hpp>
42#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
43#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
44#include <opm/input/eclipse/Schedule/Group/Group.hpp>
45
46#include <opm/simulators/timestepping/SimulatorReport.hpp>
47#include <opm/simulators/flow/countGlobalCells.hpp>
48#include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
49#include <opm/simulators/wells/BlackoilWellModelGuideRates.hpp>
50#include <opm/simulators/wells/GasLiftSingleWell.hpp>
51#include <opm/simulators/wells/GasLiftWellState.hpp>
52#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
53#include <opm/simulators/wells/GasLiftStage2.hpp>
54#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
55#include <opm/simulators/wells/PerforationData.hpp>
56#include <opm/simulators/wells/VFPInjProperties.hpp>
57#include <opm/simulators/wells/VFPProdProperties.hpp>
58#include <opm/simulators/wells/WellState.hpp>
59#include <opm/simulators/wells/WGState.hpp>
62#include <opm/simulators/wells/WellInterface.hpp>
63#include <opm/simulators/wells/StandardWell.hpp>
64#include <opm/simulators/wells/MultisegmentWell.hpp>
65#include <opm/simulators/wells/WellGroupHelpers.hpp>
66#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
67#include <opm/simulators/wells/ParallelWellInfo.hpp>
68#include <opm/simulators/timestepping/gatherConvergenceReport.hpp>
69#include <dune/common/fmatrix.hh>
70#include <dune/istl/bcrsmatrix.hh>
71#include <dune/istl/matrixmatrix.hh>
72
73#include <opm/material/densead/Math.hpp>
74
75#include <opm/simulators/utils/DeferredLogger.hpp>
76
77namespace Opm::Properties {
78
79template<class TypeTag, class MyTypeTag>
81 using type = UndefinedProperty;
82};
83
84} // namespace Opm::Properties
85
86namespace Opm {
87
89 template<typename TypeTag>
90 class BlackoilWellModel : public BaseAuxiliaryModule<TypeTag>
92 {
93 public:
94 // --------- Types ---------
96
97 using Grid = GetPropType<TypeTag, Properties::Grid>;
98 using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
99 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
100 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
101 using Indices = GetPropType<TypeTag, Properties::Indices>;
102 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
103 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
104 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
105 using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
106 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
107 using GasLiftSingleWell = typename WellInterface<TypeTag>::GasLiftSingleWell;
108 using GLiftOptWells = typename BlackoilWellModelGeneric::GLiftOptWells;
109 using GLiftProdWells = typename BlackoilWellModelGeneric::GLiftProdWells;
110 using GLiftWellStateMap =
111 typename BlackoilWellModelGeneric::GLiftWellStateMap;
112 using GLiftEclWells = typename GasLiftGroupInfo::GLiftEclWells;
113 using GLiftSyncGroups = typename GasLiftSingleWellGeneric::GLiftSyncGroups;
114 constexpr static std::size_t pressureVarIndex = GetPropType<TypeTag, Properties::Indices>::pressureSwitchIdx;
115 typedef typename BaseAuxiliaryModule<TypeTag>::NeighborSet NeighborSet;
116
117 static const int numEq = Indices::numEq;
118 static const int solventSaturationIdx = Indices::solventSaturationIdx;
119 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
120 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
121 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
122 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
123
124 // TODO: where we should put these types, WellInterface or Well Model?
125 // or there is some other strategy, like TypeTag
126 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
127 typedef Dune::BlockVector<VectorBlockType> BVector;
128
129 typedef BlackOilPolymerModule<TypeTag> PolymerModule;
130 typedef BlackOilMICPModule<TypeTag> MICPModule;
131
132 // For the conversion between the surface volume rate and resrevoir voidage rate
133 using RateConverterType = RateConverter::
134 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
135
136 // For computing average pressured used by gpmaint
137 using AverageRegionalPressureType = RegionAverageCalculator::
138 AverageRegionalPressure<FluidSystem, std::vector<int> >;
139
140 BlackoilWellModel(Simulator& ebosSimulator);
141
142 void init();
143 void initWellContainer(const int reportStepIdx) override;
144
146 // <eWoms auxiliary module stuff>
148 unsigned numDofs() const override
149 // No extra dofs are inserted for wells. (we use a Schur complement.)
150 { return 0; }
151
152 void addNeighbors(std::vector<NeighborSet>& neighbors) const override;
153
154 void applyInitial() override
155 {}
156
157 void linearize(SparseMatrixAdapter& jacobian, GlobalEqVector& res) override;
158
159 void postSolve(GlobalEqVector& deltaX) override
160 {
161 recoverWellSolutionAndUpdateWellState(deltaX);
162 }
163
165 // </ eWoms auxiliary module stuff>
167
168 template <class Restarter>
169 void deserialize(Restarter& /* res */)
170 {
171 // TODO (?)
172 }
173
178 template <class Restarter>
179 void serialize(Restarter& /* res*/)
180 {
181 // TODO (?)
182 }
183
184 void beginEpisode()
185 {
186 OPM_TIMEBLOCK(beginEpsiode);
187 beginReportStep(ebosSimulator_.episodeIndex());
188 }
189
190 void beginTimeStep();
191
192 void beginIteration()
193 {
194 OPM_TIMEBLOCK(beginIteration);
195 assemble(ebosSimulator_.model().newtonMethod().numIterations(),
196 ebosSimulator_.timeStepSize());
197 }
198
199 void endIteration()
200 { }
201
202 void endTimeStep()
203 {
204 OPM_TIMEBLOCK(endTimeStep);
205 timeStepSucceeded(ebosSimulator_.time(), ebosSimulator_.timeStepSize());
206 }
207
208 void endEpisode()
209 {
210 endReportStep();
211 }
212
213 void computeTotalRatesForDof(RateVector& rate,
214 unsigned globalIdx) const;
215
216 template <class Context>
217 void computeTotalRatesForDof(RateVector& rate,
218 const Context& context,
219 unsigned spaceIdx,
220 unsigned timeIdx) const;
221
222
223 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
224
225 using BlackoilWellModelGeneric::initFromRestartFile;
226 void initFromRestartFile(const RestartValue& restartValues)
227 {
228 initFromRestartFile(restartValues,
229 this->ebosSimulator_.vanguard().transferWTestState(),
230 grid().size(0),
232 }
233
234 using BlackoilWellModelGeneric::prepareDeserialize;
235 void prepareDeserialize(const int report_step)
236 {
237 prepareDeserialize(report_step, grid().size(0),
239 }
240
241 data::Wells wellData() const
242 {
243 auto wsrpt = this->wellState()
244 .report(ebosSimulator_.vanguard().globalCell().data(),
245 [this](const int well_index) -> bool
246 {
247 return this->wasDynamicallyShutThisTimeStep(well_index);
248 });
249
250 this->assignWellTracerRates(wsrpt);
251
252 BlackoilWellModelGuideRates(*this).assignWellGuideRates(wsrpt, this->reportStepIndex());
253 this->assignShutConnections(wsrpt, this->reportStepIndex());
254
255 return wsrpt;
256 }
257
258 // subtract Binv(D)rw from r;
259 void apply( BVector& r) const;
260
261 // subtract B*inv(D)*C * x from A*x
262 void apply(const BVector& x, BVector& Ax) const;
263
264 // accumulate the contributions of all Wells in the WellContributions object
265 void getWellContributions(WellContributions& x) const;
266
267 // apply well model with scaling of alpha
268 void applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const;
269
270 // Check if well equations is converged.
271 ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg, const bool checkWellGroupControls = false) const;
272
273 const SimulatorReportSingle& lastReport() const;
274
275 void addWellContributions(SparseMatrixAdapter& jacobian) const;
276
277 // add source from wells to the reservoir matrix
278 void addReservoirSourceTerms(GlobalEqVector& residual,
279 std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress) const;
280
281 // called at the beginning of a report step
282 void beginReportStep(const int time_step);
283
284 void updatePerforationIntensiveQuantities();
285 // it should be able to go to prepareTimeStep(), however, the updateWellControls() and initPrimaryVariablesEvaluation()
286 // makes it a little more difficult. unless we introduce if (iterationIdx != 0) to avoid doing the above functions
287 // twice at the beginning of the time step
290 void calculateExplicitQuantities(DeferredLogger& deferred_logger) const;
291 // some preparation work, mostly related to group control and RESV,
292 // at the beginning of each time step (Not report step)
293 void prepareTimeStep(DeferredLogger& deferred_logger);
294 void initPrimaryVariablesEvaluation() const;
295 std::tuple<bool, bool, double> updateWellControls(DeferredLogger& deferred_logger);
296
297 void updateAndCommunicate(const int reportStepIdx,
298 const int iterationIdx,
299 DeferredLogger& deferred_logger);
300
301 bool updateGroupControls(const Group& group,
302 DeferredLogger& deferred_logger,
303 const int reportStepIdx,
304 const int iterationIdx);
305
306 WellInterfacePtr getWell(const std::string& well_name) const;
307 bool hasWell(const std::string& well_name) const;
308
309 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
310
311 int numLocalWellsEnd() const;
312
313 void addWellPressureEquations(PressureMatrix& jacobian, const BVector& weights,const bool use_well_weights) const;
314
315 std::vector<std::vector<int>> getMaxWellConnections() const;
316
317 void addWellPressureEquationsStruct(PressureMatrix& jacobian) const;
318
319 void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
320
322 const std::vector<WellInterfacePtr>& localNonshutWells() const
323 {
324 return well_container_;
325 }
326
327 int numLocalNonshutWells() const;
328
329 protected:
330 Simulator& ebosSimulator_;
331
332 // a vector of all the wells.
333 std::vector<WellInterfacePtr> well_container_{};
334
335 std::vector<bool> is_cell_perforated_{};
336
337 void initializeWellState(const int timeStepIdx,
338 const SummaryState& summaryState);
339
340 // create the well container
341 void createWellContainer(const int time_step) override;
342
343 WellInterfacePtr
344 createWellPointer(const int wellID,
345 const int time_step) const;
346
347 template <typename WellType>
348 std::unique_ptr<WellType>
349 createTypedWellPointer(const int wellID,
350 const int time_step) const;
351
352 WellInterfacePtr createWellForWellTest(const std::string& well_name, const int report_step, DeferredLogger& deferred_logger) const;
353
354
355 const ModelParameters param_;
356 std::size_t global_num_cells_{};
357 // the number of the cells in the local grid
358 std::size_t local_num_cells_{};
359 double gravity_{};
360 std::vector<double> depth_{};
361 bool alternative_well_rate_init_{};
362
363 std::unique_ptr<RateConverterType> rateConverter_{};
364 std::map<std::string, std::unique_ptr<AverageRegionalPressureType>> regionalAveragePressureCalculator_{};
365
366
367 SimulatorReportSingle last_report_{};
368
369 // used to better efficiency of calcuation
370 mutable BVector scaleAddRes_{};
371
372 std::vector<Scalar> B_avg_{};
373
374 const Grid& grid() const
375 { return ebosSimulator_.vanguard().grid(); }
376
377 const EquilGrid& equilGrid() const
378 { return ebosSimulator_.vanguard().equilGrid(); }
379
380 const EclipseState& eclState() const
381 { return ebosSimulator_.vanguard().eclState(); }
382
383 // compute the well fluxes and assemble them in to the reservoir equations as source terms
384 // and in the well equations.
385 void assemble(const int iterationIdx,
386 const double dt);
387 bool assembleImpl(const int iterationIdx,
388 const double dt,
389 const std::size_t recursion_level,
390 DeferredLogger& local_deferredLogger);
391
392 // called at the end of a time step
393 void timeStepSucceeded(const double& simulationTime, const double dt);
394
395 // called at the end of a report step
396 void endReportStep();
397
398 // using the solution x to recover the solution xw for wells and applying
399 // xw to update Well State
400 void recoverWellSolutionAndUpdateWellState(const BVector& x);
401
402 // setting the well_solutions_ based on well_state.
403 void updatePrimaryVariables(DeferredLogger& deferred_logger);
404
405 void updateAverageFormationFactor();
406
407 void computePotentials(const std::size_t widx,
408 const WellState& well_state_copy,
409 std::string& exc_msg,
410 ExceptionType::ExcEnum& exc_type,
411 DeferredLogger& deferred_logger) override;
412
413 const std::vector<double>& wellPerfEfficiencyFactors() const;
414
415 void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger& deferred_logger) override;
416 void calculateProductivityIndexValues(DeferredLogger& deferred_logger) override;
417 void calculateProductivityIndexValues(const WellInterface<TypeTag>* wellPtr,
418 DeferredLogger& deferred_logger);
419
420 // The number of components in the model.
421 int numComponents() const;
422
423 int reportStepIndex() const;
424
425 void assembleWellEq(const double dt, DeferredLogger& deferred_logger);
426
427 bool maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
428
429 void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
430 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
431 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
432
433 // cannot be const since it accesses the non-const WellState
434 void gasLiftOptimizationStage1SingleWell(WellInterface<TypeTag> *well,
435 DeferredLogger& deferred_logger,
436 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
437 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map,
438 GLiftSyncGroups& groups_to_sync);
439
440 void extractLegacyCellPvtRegionIndex_();
441
442 void extractLegacyDepth_();
443
445 void updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const;
446
447 void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger& deferred_logger);
448
449 void calcRates(const int fipnum,
450 const int pvtreg,
451 const std::vector<double>& production_rates,
452 std::vector<double>& resv_coeff) override;
453
454 void calcInjRates(const int fipnum,
455 const int pvtreg,
456 std::vector<double>& resv_coeff) override;
457
458 void computeWellTemperature();
459
460 void assignWellTracerRates(data::Wells& wsrpt) const;
461
462 int compressedIndexForInterior(int cartesian_cell_idx) const override {
463 return ebosSimulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
464 }
465
466 private:
467 BlackoilWellModel(Simulator& ebosSimulator, const PhaseUsage& pu);
468 };
469
470
471} // namespace Opm
472
473#include "BlackoilWellModel_impl.hpp"
474#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...
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:73
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:92
void serialize(Restarter &)
This method writes the complete state of the well to the harddisk.
Definition: BlackoilWellModel.hpp:179
const std::vector< WellInterfacePtr > & localNonshutWells() const
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:322
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:1438
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1572
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition: BlackoilWellModel.hpp:462
Definition: GasLiftSingleWell.hpp:38
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
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParametersEbos.hpp:327
bool use_multisegment_well_
Whether to use MultisegmentWell to handle multisegment wells it is something temporary before the mul...
Definition: BlackoilModelParametersEbos.hpp:408
Definition: BlackoilPhases.hpp:46
Definition: BlackoilWellModel.hpp:80