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 
79 namespace Opm::Properties {
80 
81 template<class TypeTag, class MyTypeTag>
83  using type = UndefinedProperty;
84 };
85 
86 } // namespace Opm::Properties
87 
88 namespace 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),
230  param_.use_multisegment_well_);
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