All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
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 
26 #include <opm/autodiff/WellInterface.hpp>
27 
28 namespace Opm
29 {
30 
31  template<typename TypeTag>
32  class MultisegmentWell: public WellInterface<TypeTag>
33  {
34  public:
36 
37  using typename Base::WellState;
38  using typename Base::Simulator;
39  using typename Base::IntensiveQuantities;
40  using typename Base::FluidSystem;
41  using typename Base::ModelParameters;
42  using typename Base::MaterialLaw;
43  using typename Base::BlackoilIndices;
44 
46  using Base::numEq;
47 
48  using Base::has_solvent;
49  using Base::has_polymer;
50 
51  // TODO: for now, not considering the polymer, solvent and so on to simplify the development process.
52 
53  // TODO: we need to have order for the primary variables and also the order for the well equations.
54  // sometimes, they are similar, while sometimes, they can have very different forms.
55 
56  // TODO: the following system looks not rather flexible. Looking into all kinds of possibilities
57  // TODO: gas is always there? how about oil water case?
58  // Is it gas oil two phase case?
59  static const bool gasoil = numEq == 2 && (BlackoilIndices::compositionSwitchIdx >= 0);
60  static const int GTotal = 0;
61  static const int WFrac = gasoil? -1000: 1;
62  static const int GFrac = gasoil? 1 : 2;
63  static const int SPres = gasoil? 2 : 3;
64 
66  static const int numWellEq = GET_PROP_VALUE(TypeTag, EnablePolymer)? numEq : numEq + 1;
67 
68  using typename Base::Scalar;
69  using typename Base::ConvergenceReport;
70 
72  using typename Base::Mat;
73  using typename Base::BVector;
74  using typename Base::Eval;
75 
76  // sparsity pattern for the matrices
77  // [A C^T [x = [ res
78  // B D ] x_well] res_well]
79 
80  // the vector type for the res_well and x_well
81  typedef Dune::FieldVector<Scalar, numWellEq> VectorBlockWellType;
82  typedef Dune::BlockVector<VectorBlockWellType> BVectorWell;
83 
84  // the matrix type for the diagonal matrix D
85  typedef Dune::FieldMatrix<Scalar, numWellEq, numWellEq > DiagMatrixBlockWellType;
86  typedef Dune::BCRSMatrix <DiagMatrixBlockWellType> DiagMatWell;
87 
88  // the matrix type for the non-diagonal matrix B and C^T
89  typedef Dune::FieldMatrix<Scalar, numWellEq, numEq> OffDiagMatrixBlockWellType;
90  typedef Dune::BCRSMatrix<OffDiagMatrixBlockWellType> OffDiagMatWell;
91 
92  // TODO: for more efficient implementation, we should have EvalReservoir, EvalWell, and EvalRerservoirAndWell
93  // EvalR (Eval), EvalW, EvalRW
94  // TODO: for now, we only use one type to save some implementation efforts, while improve later.
95  typedef DenseAd::Evaluation<double, /*size=*/numEq + numWellEq> EvalWell;
96 
97  MultisegmentWell(const Well* well, const int time_step, const Wells* wells, const ModelParameters& param);
98 
99  virtual void init(const PhaseUsage* phase_usage_arg,
100  const std::vector<bool>* active_arg,
101  const std::vector<double>& depth_arg,
102  const double gravity_arg,
103  const int num_cells);
104 
105 
106  virtual void initPrimaryVariablesEvaluation() const;
107 
108  virtual void assembleWellEq(Simulator& ebosSimulator,
109  const double dt,
110  WellState& well_state,
111  bool only_wells);
112 
114  // TODO: later will check wheter we need current
115  virtual void updateWellStateWithTarget(const int current,
116  WellState& well_state) const;
117 
119  virtual ConvergenceReport getWellConvergence(const std::vector<double>& B_avg) const;
120 
122  virtual void apply(const BVector& x, BVector& Ax) const;
124  virtual void apply(BVector& r) const;
125 
128  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
129  WellState& well_state) const;
130 
132  virtual void computeWellPotentials(const Simulator& ebosSimulator,
133  const WellState& well_state,
134  std::vector<double>& well_potentials);
135 
136  virtual void updatePrimaryVariables(const WellState& well_state) const;
137 
138  virtual void solveEqAndUpdateWellState(WellState& well_state); // const?
139 
140  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
141  const WellState& well_state); // should be const?
142 
145  int numberOfSegments() const;
146 
147  int numberOfPerforations() const;
148 
149  protected:
150  int number_segments_;
151 
152  // components of the pressure drop to be included
153  WellSegment::CompPressureDropEnum compPressureDrop() const;
154  // multi-phase flow model
155  WellSegment::MultiPhaseModelEnum multiphaseModel() const;
156 
157  // get the SegmentSet from the well_ecl_
158  const SegmentSet& segmentSet() const;
159 
160  // protected member variables from the Base class
161  using Base::well_ecl_;
162  using Base::number_of_perforations_; // TODO: can use well_ecl_?
163  using Base::current_step_;
164  using Base::index_of_well_;
165  using Base::number_of_phases_;
166 
167  // TODO: the current implementation really relies on the order of the
168  // perforation does not change from the parser to Wells structure.
169  using Base::well_cells_;
170  using Base::param_;
171  using Base::well_index_;
172  using Base::well_type_;
173  using Base::first_perf_;
174  using Base::saturation_table_number_;
175  using Base::well_efficiency_factor_;
176  using Base::gravity_;
177  using Base::well_controls_;
178  using Base::perf_depth_;
179 
180  // protected functions from the Base class
181  using Base::active;
182  using Base::phaseUsage;
183  using Base::name;
184  using Base::numComponents;
185  using Base::flowPhaseToEbosPhaseIdx;
186  using Base::flowPhaseToEbosCompIdx;
187  using Base::getAllowCrossFlow;
188 
189  // TODO: trying to use the information from the Well opm-parser as much
190  // as possible, it will possibly be re-implemented later for efficiency reason.
191 
192  // the completions that is related to each segment
193  // the completions's ids are their index in the vector well_index_, well_cell_
194  // This is also assuming the order of the completions in Well is the same with
195  // the order of the completions in wells.
196  // it is for convinience reason. we can just calcuate the inforation for segment once then using it for all the perofrations
197  // belonging to this segment
198  std::vector<std::vector<int> > segment_perforations_;
199 
200  // the inlet segments for each segment. It is for convinience and efficiency reason
201  std::vector<std::vector<int> > segment_inlets_;
202 
203  // segment number is an ID of the segment, it is specified in the deck
204  // get the loation of the segment with a segment number in the segmentSet
205  int segmentNumberToIndex(const int segment_number) const;
206 
207  // TODO, the following should go to a class for computing purpose
208  // two off-diagonal matrices
209  mutable OffDiagMatWell duneB_;
210  mutable OffDiagMatWell duneC_;
211  // diagonal matrix for the well
212  mutable DiagMatWell duneD_;
213 
214  // residuals of the well equations
215  mutable BVectorWell resWell_;
216 
217  // the values for the primary varibles
218  // based on different solutioin strategies, the wells can have different primary variables
219  mutable std::vector<std::array<double, numWellEq> > primary_variables_;
220 
221  // the Evaluation for the well primary variables, which contain derivativles and are used in AD calculation
222  mutable std::vector<std::array<EvalWell, numWellEq> > primary_variables_evaluation_;
223 
224  // depth difference between perforations and the perforated grid cells
225  std::vector<double> cell_perforation_depth_diffs_;
226  // pressure correction due to the different depth of the perforation and
227  // center depth of the grid block
228  std::vector<double> cell_perforation_pressure_diffs_;
229 
230  // depth difference between the segment and the peforation
231  // or in another way, the depth difference between the perforation and
232  // the segment the perforation belongs to
233  std::vector<double> perforation_segment_depth_diffs_;
234 
235  // the intial component compistion of segments
236  std::vector<std::vector<double> > segment_comp_initial_;
237 
238  // the densities of segment fluids
239  // we should not have this member variable
240  std::vector<EvalWell> segment_densities_;
241 
242  // the viscosity of the segments
243  std::vector<EvalWell> segment_viscosities_;
244 
245  // the mass rate of the segments
246  std::vector<EvalWell> segment_mass_rates_;
247 
248  std::vector<double> segment_depth_diffs_;
249 
250  void initMatrixAndVectors(const int num_cells) const;
251 
252  // protected functions
253  // EvalWell getBhp(); this one should be something similar to getSegmentPressure();
254  // EvalWell getQs(); this one should be something similar to getSegmentRates()
255  // EValWell wellVolumeFractionScaled, wellVolumeFraction, wellSurfaceVolumeFraction ... these should have different names, and probably will be needed.
256  // bool crossFlowAllowed(const Simulator& ebosSimulator) const; probably will be needed
257  // xw = inv(D)*(rw - C*x)
258  void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
259 
260  // updating the well_state based on well solution dwells
261  void updateWellState(const BVectorWell& dwells,
262  const bool inner_iteration,
263  WellState& well_state) const;
264 
265  // initialize the segment rates with well rates
266  // when there is no more accurate way to initialize the segment rates, we initialize
267  // the segment rates based on well rates with a simple strategy
268  void initSegmentRatesWithWellRates(WellState& well_state) const;
269 
270  // computing the accumulation term for later use in well mass equations
271  void computeInitialComposition();
272 
273  // compute the pressure difference between the perforation and cell center
274  void computePerfCellPressDiffs(const Simulator& ebosSimulator);
275 
276  // fraction value of the primary variables
277  // should we just use member variables to store them instead of calculating them again and again
278  EvalWell volumeFraction(const int seg, const int comp_idx) const;
279 
280  // F_p / g_p, the basic usage of this value is because Q_p = G_t * F_p / G_p
281  EvalWell volumeFractionScaled(const int seg, const int comp_idx) const;
282 
283  // basically Q_p / \sigma_p Q_p
284  EvalWell surfaceVolumeFraction(const int seg, const int comp_idx) const;
285 
286  void computePerfRate(const IntensiveQuantities& int_quants,
287  const std::vector<EvalWell>& mob_perfcells,
288  const int seg,
289  const int perf,
290  const EvalWell& segment_pressure,
291  const bool& allow_cf,
292  std::vector<EvalWell>& cq_s) const;
293 
294  // convert a Eval from reservoir to contain the derivative related to wells
295  EvalWell extendEval(const Eval& in) const;
296 
297  // compute the fluid properties, such as densities, viscosities, and so on, in the segments
298  // They will be treated implicitly, so they need to be of Evaluation type
299  void computeSegmentFluidProperties(const Simulator& ebosSimulator);
300 
301  EvalWell getSegmentPressure(const int seg) const;
302 
303  EvalWell getSegmentRate(const int seg, const int comp_idx) const;
304 
305  EvalWell getSegmentGTotal(const int seg) const;
306 
307  // get the mobility for specific perforation
308  void getMobility(const Simulator& ebosSimulator,
309  const int perf,
310  std::vector<EvalWell>& mob) const;
311 
312  void assembleControlEq() const;
313 
314  void assemblePressureEq(const int seg) const;
315 
316  // hytrostatic pressure loss
317  EvalWell getHydroPressureLoss(const int seg) const;
318 
319  // frictinal pressure loss
320  EvalWell getFrictionPressureLoss(const int seg) const;
321 
322  void handleAccelerationPressureLoss(const int seg) const;
323 
324  // handling the overshooting and undershooting of the fractions
325  void processFractions(const int seg) const;
326 
327  void updateWellStateFromPrimaryVariables(WellState& well_state) const;
328 
329  double scalingFactor(const int comp_idx) const;
330 
331  bool frictionalPressureLossConsidered() const;
332 
333  bool accelerationalPressureLossConsidered() const;
334 
335  // TODO: try to make ebosSimulator const, as it should be
336  void iterateWellEquations(Simulator& ebosSimulator,
337  const double dt,
338  WellState& well_state);
339 
340  void assembleWellEqWithoutIteration(Simulator& ebosSimulator,
341  const double dt,
342  WellState& well_state,
343  bool only_wells);
344  };
345 
346 }
347 
348 #include "MultisegmentWell_impl.hpp"
349 
350 #endif // OPM_MULTISEGMENTWELL_HEADER_INCLUDED
virtual void apply(const BVector &x, BVector &Ax) const
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:490
virtual void updateWellStateWithTarget(const int current, WellState &well_state) const
updating the well state based the control mode specified with current
Definition: MultisegmentWell_impl.hpp:251
a struct to collect information about the convergence checking
Definition: WellInterface.hpp:117
const std::string & name() const
Well name.
Definition: WellInterface_impl.hpp:138
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
Definition: WellInterface.hpp:60
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state) const
using the solution x to recover the solution xw for wells and applying xw to update Well State ...
Definition: MultisegmentWell_impl.hpp:533
static const int numWellEq
the number of well equations // TODO: it should have a more general strategy for it ...
Definition: MultisegmentWell.hpp:66
virtual ConvergenceReport getWellConvergence(const std::vector< double > &B_avg) const
check whether the well equations get converged for this well
Definition: MultisegmentWell_impl.hpp:409
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials)
computing the well potentials for group control
Definition: MultisegmentWell_impl.hpp:548
int numberOfSegments() const
number of segments for this well int number_of_segments_;
Definition: MultisegmentWell_impl.hpp:824
Definition: MultisegmentWell.hpp:32