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/parser/eclipse/EclipseState/Runspec.hpp>
43 
44 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
45 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
46 #include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
47 #include <opm/parser/eclipse/EclipseState/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 FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
101  using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
102  using Indices = GetPropType<TypeTag, Properties::Indices>;
103  using Simulator = GetPropType<TypeTag, Properties::Simulator>;
104  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
105  using RateVector = GetPropType<TypeTag, Properties::RateVector>;
106  using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
107  using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
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 
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 Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
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() 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  template <class Context>
213  void computeTotalRatesForDof(RateVector& rate,
214  const Context& context,
215  unsigned spaceIdx,
216  unsigned timeIdx) const;
217 
218 
219  using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
220 
221  using BlackoilWellModelGeneric::initFromRestartFile;
222  void initFromRestartFile(const RestartValue& restartValues)
223  {
224  initFromRestartFile(restartValues,
225  this->ebosSimulator_.vanguard().transferWTestState(),
226  UgGridHelpers::numCells(grid()),
227  param_.use_multisegment_well_);
228  }
229 
230  data::Wells wellData() const
231  {
232  auto wsrpt = this->wellState()
233  .report(UgGridHelpers::globalCell(grid()),
234  [this](const int well_ndex) -> bool
235  {
236  return this->wasDynamicallyShutThisTimeStep(well_ndex);
237  });
238 
239  this->assignWellTracerRates(wsrpt);
240 
241  this->assignWellGuideRates(wsrpt);
242  this->assignShutConnections(wsrpt, this->reportStepIndex());
243 
244  return wsrpt;
245  }
246 
247  // subtract Binv(D)rw from r;
248  void apply( BVector& r) const;
249 
250  // subtract B*inv(D)*C * x from A*x
251  void apply(const BVector& x, BVector& Ax) const;
252 
253 #if HAVE_CUDA || HAVE_OPENCL
254  // accumulate the contributions of all Wells in the WellContributions object
255  void getWellContributions(WellContributions& x) const;
256 #endif
257 
258  // apply well model with scaling of alpha
259  void applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const;
260 
261  // Check if well equations is converged.
262  ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg, const bool checkGroupConvergence = false) const;
263 
264  const SimulatorReportSingle& lastReport() const;
265 
266  void addWellContributions(SparseMatrixAdapter& jacobian) const
267  {
268  for ( const auto& well: well_container_ ) {
269  well->addWellContributions(jacobian);
270  }
271  }
272 
273  // called at the beginning of a report step
274  void beginReportStep(const int time_step);
275 
276  void updatePerforationIntensiveQuantities();
277  // it should be able to go to prepareTimeStep(), however, the updateWellControls() and initPrimaryVariablesEvaluation()
278  // makes it a little more difficult. unless we introduce if (iterationIdx != 0) to avoid doing the above functions
279  // twice at the beginning of the time step
282  void calculateExplicitQuantities(DeferredLogger& deferred_logger) const;
283  // some preparation work, mostly related to group control and RESV,
284  // at the beginning of each time step (Not report step)
285  void prepareTimeStep(DeferredLogger& deferred_logger);
286  void initPrimaryVariablesEvaluation() const;
287  void updateWellControls(DeferredLogger& deferred_logger, const bool checkGroupControls);
288 
289  void updateAndCommunicate(const int reportStepIdx,
290  const int iterationIdx,
291  DeferredLogger& deferred_logger);
292 
293  WellInterfacePtr getWell(const std::string& well_name) const;
294  bool hasWell(const std::string& well_name) const;
295 
296  void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
297 
299  const std::vector<WellInterfacePtr>& localNonshutWells()
300  {
301  return well_container_;
302  }
303 
304  protected:
305  Simulator& ebosSimulator_;
306 
307  // a vector of all the wells.
308  std::vector<WellInterfacePtr > well_container_{};
309 
310  std::vector<bool> is_cell_perforated_{};
311 
312  void initializeWellState(const int timeStepIdx,
313  const SummaryState& summaryState);
314 
315  // create the well container
316  void createWellContainer(const int time_step) override;
317 
318  WellInterfacePtr
319  createWellPointer(const int wellID,
320  const int time_step) const;
321 
322  template <typename WellType>
323  std::unique_ptr<WellType>
324  createTypedWellPointer(const int wellID,
325  const int time_step) const;
326 
327  WellInterfacePtr createWellForWellTest(const std::string& well_name, const int report_step, DeferredLogger& deferred_logger) const;
328 
329 
330  const ModelParameters param_;
331  size_t global_num_cells_{};
332  // the number of the cells in the local grid
333  size_t local_num_cells_{};
334  double gravity_{};
335  std::vector<double> depth_{};
336  bool alternative_well_rate_init_{};
337 
338  std::unique_ptr<RateConverterType> rateConverter_{};
339  std::unique_ptr<AverageRegionalPressureType> regionalAveragePressureCalculator_{};
340 
341 
342  SimulatorReportSingle last_report_{};
343 
344  // used to better efficiency of calcuation
345  mutable BVector scaleAddRes_{};
346 
347  std::vector<Scalar> B_avg_{};
348 
349  const Grid& grid() const
350  { return ebosSimulator_.vanguard().grid(); }
351 
352  const EclipseState& eclState() const
353  { return ebosSimulator_.vanguard().eclState(); }
354 
355  // compute the well fluxes and assemble them in to the reservoir equations as source terms
356  // and in the well equations.
357  void assemble(const int iterationIdx,
358  const double dt);
359 
360  // called at the end of a time step
361  void timeStepSucceeded(const double& simulationTime, const double dt);
362 
363  // called at the end of a report step
364  void endReportStep();
365 
366  // using the solution x to recover the solution xw for wells and applying
367  // xw to update Well State
368  void recoverWellSolutionAndUpdateWellState(const BVector& x);
369 
370  // setting the well_solutions_ based on well_state.
371  void updatePrimaryVariables(DeferredLogger& deferred_logger);
372 
373  void setupCartesianToCompressed_(const int* global_cell, int local_num__cells);
374 
375  void updateAverageFormationFactor();
376 
377  void computePotentials(const std::size_t widx,
378  const WellState& well_state_copy,
379  std::string& exc_msg,
380  ExceptionType::ExcEnum& exc_type,
381  DeferredLogger& deferred_logger) override;
382 
383  const std::vector<double>& wellPerfEfficiencyFactors() const;
384 
385  void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger& deferred_logger) override;
386  void calculateProductivityIndexValues(DeferredLogger& deferred_logger) override;
387  void calculateProductivityIndexValues(const WellInterface<TypeTag>* wellPtr,
388  DeferredLogger& deferred_logger);
389 
390  // The number of components in the model.
391  int numComponents() const;
392 
393  int reportStepIndex() const;
394 
395  void assembleWellEq(const double dt, DeferredLogger& deferred_logger);
396 
397  void maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
398 
399  void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
400  GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
401  GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
402 
403  void extractLegacyCellPvtRegionIndex_();
404 
405  void extractLegacyDepth_();
406 
408  void updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const;
409 
410  void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger& deferred_logger);
411 
412  void calcRates(const int fipnum,
413  const int pvtreg,
414  std::vector<double>& resv_coeff) override;
415 
416  void calcInjRates(const int fipnum,
417  const int pvtreg,
418  std::vector<double>& resv_coeff) override;
419 
420  void computeWellTemperature();
421 
422  void assignWellTracerRates(data::Wells& wsrpt) const;
423 
424  private:
425  BlackoilWellModel(Simulator& ebosSimulator, const PhaseUsage& pu);
426  };
427 
428 
429 } // namespace Opm
430 
431 #include "BlackoilWellModel_impl.hpp"
432 #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:69
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()
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:299
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:1198
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1301
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:26
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: BlackoilWellModel.hpp:82