My Project
MultisegmentWell.hpp
1 /*
2  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2017 Statoil ASA.
4 
5  This file is part of the Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #ifndef OPM_MULTISEGMENTWELL_HEADER_INCLUDED
23 #define OPM_MULTISEGMENTWELL_HEADER_INCLUDED
24 
25 #include <opm/simulators/wells/WellInterface.hpp>
26 #include <opm/simulators/wells/MultisegmentWellEval.hpp>
27 
28 #include <opm/input/eclipse/EclipseState/Runspec.hpp>
29 
30 namespace Opm
31 {
32  class DeferredLogger;
33 
34  template<typename TypeTag>
35  class MultisegmentWell : public WellInterface<TypeTag>
36  , public MultisegmentWellEval<GetPropType<TypeTag, Properties::FluidSystem>,
37  GetPropType<TypeTag, Properties::Indices>,
38  GetPropType<TypeTag, Properties::Scalar>>
39  {
40  public:
43  GetPropType<TypeTag, Properties::Indices>,
44  GetPropType<TypeTag, Properties::Scalar>>;
45 
46  using typename Base::Simulator;
47  using typename Base::IntensiveQuantities;
48  using typename Base::FluidSystem;
49  using typename Base::ModelParameters;
50  using typename Base::MaterialLaw;
51  using typename Base::Indices;
52  using typename Base::RateConverterType;
53  using typename Base::SparseMatrixAdapter;
54  using typename Base::FluidState;
55 
56  using Base::has_solvent;
57  using Base::has_polymer;
58  using Base::Water;
59  using Base::Oil;
60  using Base::Gas;
61 
62  using typename Base::Scalar;
63 
65  using typename Base::BVector;
66  using typename Base::Eval;
67 
68  using typename MSWEval::EvalWell;
69  using typename MSWEval::BVectorWell;
70  using typename MSWEval::DiagMatWell;
71  using typename MSWEval::OffDiagMatrixBlockWellType;
72  using MSWEval::GFrac;
73  using MSWEval::WFrac;
74  using MSWEval::WQTotal;
75  using MSWEval::SPres;
76  using MSWEval::numWellEq;
77  using typename Base::PressureMatrix;
78 
79  MultisegmentWell(const Well& well,
80  const ParallelWellInfo& pw_info,
81  const int time_step,
82  const ModelParameters& param,
83  const RateConverterType& rate_converter,
84  const int pvtRegionIdx,
85  const int num_components,
86  const int num_phases,
87  const int index_of_well,
88  const std::vector<PerforationData>& perf_data);
89 
90  virtual void init(const PhaseUsage* phase_usage_arg,
91  const std::vector<double>& depth_arg,
92  const double gravity_arg,
93  const int num_cells,
94  const std::vector< Scalar >& B_avg,
95  const bool changed_to_open_this_step) override;
96 
97  virtual void initPrimaryVariablesEvaluation() const override;
98 
100  virtual void updateWellStateWithTarget(const Simulator& ebos_simulator,
101  const GroupState& group_state,
102  WellState& well_state,
103  DeferredLogger& deferred_logger) const override;
104 
106  virtual ConvergenceReport getWellConvergence(const WellState& well_state,
107  const std::vector<double>& B_avg,
108  DeferredLogger& deferred_logger,
109  const bool relax_tolerance = false) const override;
110 
112  virtual void apply(const BVector& x, BVector& Ax) const override;
114  virtual void apply(BVector& r) const override;
115 
118  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
119  WellState& well_state,
120  DeferredLogger& deferred_logger) const override;
121 
123  virtual void computeWellPotentials(const Simulator& ebosSimulator,
124  const WellState& well_state,
125  std::vector<double>& well_potentials,
126  DeferredLogger& deferred_logger) override;
127 
128  virtual void updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_logger) const override;
129 
130  virtual void solveEqAndUpdateWellState(WellState& well_state, DeferredLogger& deferred_logger) override; // const?
131 
132  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
133  const WellState& well_state,
134  DeferredLogger& deferred_logger) override; // should be const?
135 
136  virtual void updateProductivityIndex(const Simulator& ebosSimulator,
137  const WellProdIndexCalculator& wellPICalc,
138  WellState& well_state,
139  DeferredLogger& deferred_logger) const override;
140 
141  virtual void addWellContributions(SparseMatrixAdapter& jacobian) const override;
142 
143  virtual void addWellPressureEquations(PressureMatrix& mat,
144  const BVector& x,
145  const int pressureVarIndex,
146  const bool use_well_weights,
147  const WellState& well_state) const override;
148 
149  virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
150  DeferredLogger& deferred_logger) const override;
151 
152  void computeConnLevelProdInd(const FluidState& fs,
153  const std::function<double(const double)>& connPICalc,
154  const std::vector<Scalar>& mobility,
155  double* connPI) const;
156 
157  void computeConnLevelInjInd(const FluidState& fs,
158  const Phase preferred_phase,
159  const std::function<double(const double)>& connIICalc,
160  const std::vector<Scalar>& mobility,
161  double* connII,
162  DeferredLogger& deferred_logger) const;
163 
164  virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
165  const Simulator& ebos_simulator,
166  const SummaryState& summary_state,
167  DeferredLogger& deferred_logger,
168  double alq_value) const override;
169 
170  protected:
171  int number_segments_;
172 
173  // regularize msw equation
174  bool regularize_;
175 
176  // components of the pressure drop to be included
177  WellSegments::CompPressureDrop compPressureDrop() const;
178  // multi-phase flow model
179  WellSegments::MultiPhaseModel multiphaseModel() const;
180 
181  // the intial amount of fluids in each segment under surface condition
182  std::vector<std::vector<double> > segment_fluid_initial_;
183 
184  mutable int debug_cost_counter_ = 0;
185 
186  // updating the well_state based on well solution dwells
187  void updateWellState(const BVectorWell& dwells,
188  WellState& well_state,
189  DeferredLogger& deferred_logger,
190  const double relaxation_factor=1.0) const;
191 
192 
193  // computing the accumulation term for later use in well mass equations
194  void computeInitialSegmentFluids(const Simulator& ebos_simulator);
195 
196  // compute the pressure difference between the perforation and cell center
197  void computePerfCellPressDiffs(const Simulator& ebosSimulator);
198 
199  void computePerfRateScalar(const IntensiveQuantities& int_quants,
200  const std::vector<Scalar>& mob_perfcells,
201  const double Tw,
202  const int seg,
203  const int perf,
204  const Scalar& segment_pressure,
205  const bool& allow_cf,
206  std::vector<Scalar>& cq_s,
207  DeferredLogger& deferred_logger) const;
208 
209  void computePerfRateEval(const IntensiveQuantities& int_quants,
210  const std::vector<EvalWell>& mob_perfcells,
211  const double Tw,
212  const int seg,
213  const int perf,
214  const EvalWell& segment_pressure,
215  const bool& allow_cf,
216  std::vector<EvalWell>& cq_s,
217  EvalWell& perf_press,
218  double& perf_dis_gas_rate,
219  double& perf_vap_oil_rate,
220  DeferredLogger& deferred_logger) const;
221 
222  template<class Value>
223  void computePerfRate(const Value& pressure_cell,
224  const Value& rs,
225  const Value& rv,
226  const std::vector<Value>& b_perfcells,
227  const std::vector<Value>& mob_perfcells,
228  const double Tw,
229  const int perf,
230  const Value& segment_pressure,
231  const Value& segment_density,
232  const bool& allow_cf,
233  const std::vector<Value>& cmix_s,
234  std::vector<Value>& cq_s,
235  Value& perf_press,
236  double& perf_dis_gas_rate,
237  double& perf_vap_oil_rate,
238  DeferredLogger& deferred_logger) const;
239 
240  // compute the fluid properties, such as densities, viscosities, and so on, in the segments
241  // They will be treated implicitly, so they need to be of Evaluation type
242  void computeSegmentFluidProperties(const Simulator& ebosSimulator,
243  DeferredLogger& deferred_logger);
244 
245  // get the mobility for specific perforation
246  void getMobilityEval(const Simulator& ebosSimulator,
247  const int perf,
248  std::vector<EvalWell>& mob) const;
249 
250  // get the mobility for specific perforation
251  void getMobilityScalar(const Simulator& ebosSimulator,
252  const int perf,
253  std::vector<Scalar>& mob) const;
254 
255  void computeWellRatesAtBhpLimit(const Simulator& ebosSimulator,
256  std::vector<double>& well_flux,
257  DeferredLogger& deferred_logger) const;
258 
259  virtual void computeWellRatesWithBhp(const Simulator& ebosSimulator,
260  const double& bhp,
261  std::vector<double>& well_flux,
262  DeferredLogger& deferred_logger) const override;
263 
264  void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
265  const Scalar& bhp,
266  std::vector<double>& well_flux,
267  DeferredLogger& deferred_logger) const;
268 
269  std::vector<double> computeWellPotentialWithTHP(
270  const WellState& well_state,
271  const Simulator& ebos_simulator,
272  DeferredLogger& deferred_logger) const;
273 
274  virtual double getRefDensity() const override;
275 
276  virtual bool iterateWellEqWithControl(const Simulator& ebosSimulator,
277  const double dt,
278  const Well::InjectionControls& inj_controls,
279  const Well::ProductionControls& prod_controls,
280  WellState& well_state,
281  const GroupState& group_state,
282  DeferredLogger& deferred_logger) override;
283 
284  virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
285  const double dt,
286  const Well::InjectionControls& inj_controls,
287  const Well::ProductionControls& prod_controls,
288  WellState& well_state,
289  const GroupState& group_state,
290  DeferredLogger& deferred_logger) override;
291 
292  virtual void updateWaterThroughput(const double dt, WellState& well_state) const override;
293 
294  EvalWell getSegmentSurfaceVolume(const Simulator& ebos_simulator, const int seg_idx) const;
295 
296  // turn on crossflow to avoid singular well equations
297  // when the well is banned from cross-flow and the BHP is not properly initialized,
298  // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
299  // well rates, it can cause problem for THP calculation
300  // TODO: looking for better alternative to avoid wrong-signed well rates
301  bool openCrossFlowAvoidSingularity(const Simulator& ebos_simulator) const;
302 
303  // for a well, when all drawdown are in the wrong direction, then this well will not
304  // be able to produce/inject .
305  bool allDrawDownWrongDirection(const Simulator& ebos_simulator) const;
306 
307 
308 
309  std::optional<double> computeBhpAtThpLimitProd(
310  const WellState& well_state,
311  const Simulator& ebos_simulator,
312  const SummaryState& summary_state,
313  DeferredLogger& deferred_logger) const;
314 
315  std::optional<double> computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
316  const SummaryState& summary_state,
317  DeferredLogger& deferred_logger) const;
318 
319  double maxPerfPress(const Simulator& ebos_simulator) const;
320 
321  // check whether the well is operable under BHP limit with current reservoir condition
322  virtual void checkOperabilityUnderBHPLimit(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) override;
323 
324  // check whether the well is operable under THP limit with current reservoir condition
325  virtual void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) override;
326 
327  // updating the inflow based on the current reservoir condition
328  virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const override;
329  };
330 
331 }
332 
333 #include "MultisegmentWell_impl.hpp"
334 
335 #endif // OPM_MULTISEGMENTWELL_HEADER_INCLUDED
Represents the convergence status of the whole simulator, to make it possible to query and store the ...
Definition: ConvergenceReport.hpp:36
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Definition: MultisegmentWellEval.hpp:55
Definition: MultisegmentWell.hpp:39
virtual void updateWellStateWithTarget(const Simulator &ebos_simulator, const GroupState &group_state, WellState &well_state, DeferredLogger &deferred_logger) const override
updating the well state based the current control mode
Definition: MultisegmentWell_impl.hpp:153
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) const override
using the solution x to recover the solution xw for wells and applying xw to update Well State
Definition: MultisegmentWell_impl.hpp:237
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:195
virtual std::vector< double > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const override
Compute well rates based on current reservoir conditions and well variables.
Definition: MultisegmentWell_impl.hpp:2003
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: MultisegmentWell_impl.hpp:255
virtual ConvergenceReport getWellConvergence(const WellState &well_state, const std::vector< double > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance=false) const override
check whether the well equations get converged for this well
Definition: MultisegmentWell_impl.hpp:172
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:243
Definition: WellInterface.hpp:72
Collect per-connection static information to enable calculating connection-level or well-level produc...
Definition: WellProdIndexCalculator.hpp:36
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
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
Definition: BlackoilPhases.hpp:46