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