All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
StandardWell.hpp
1 /*
2  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2017 Statoil ASA.
4  Copyright 2016 - 2017 IRIS AS.
5 
6  This file is part of the Open Porous Media project (OPM).
7 
8  OPM is free software: you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation, either version 3 of the License, or
11  (at your option) any later version.
12 
13  OPM is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with OPM. If not, see <http://www.gnu.org/licenses/>.
20 */
21 
22 
23 #ifndef OPM_STANDARDWELL_HEADER_INCLUDED
24 #define OPM_STANDARDWELL_HEADER_INCLUDED
25 
26 
27 #include <opm/autodiff/WellInterface.hpp>
28 
29 namespace Opm
30 {
31 
32  template<typename TypeTag>
33  class StandardWell: public WellInterface<TypeTag>
34  {
35 
36  public:
38  // TODO: some functions working with AD variables handles only with values (double) without
39  // dealing with derivatives. It can be beneficial to make functions can work with either AD or scalar value.
40  // And also, it can also be beneficial to make these functions hanle different types of AD variables.
41  using typename Base::Simulator;
42  using typename Base::WellState;
43  using typename Base::IntensiveQuantities;
44  using typename Base::FluidSystem;
45  using typename Base::MaterialLaw;
46  using typename Base::ModelParameters;
47  using typename Base::BlackoilIndices;
48  using typename Base::PolymerModule;
49 
50  using Base::numEq;
51 
52  // the positions of the primary variables for StandardWell
53  // there are three primary variables, the second and the third ones are F_w and F_g
54  // the first one can be total rate (G_t) or bhp, based on the control
55 
56  static const bool gasoil = numEq == 2 && (BlackoilIndices::compositionSwitchIdx >= 0);
57  static const int XvarWell = 0;
58  static const int WFrac = gasoil? -1000: 1;
59  static const int GFrac = gasoil? 1: 2;
60  static const int SFrac = 3;
61 
62  using typename Base::Scalar;
63  using typename Base::ConvergenceReport;
64 
65 
66  using Base::has_solvent;
67  using Base::has_polymer;
68  using Base::name;
69 
70  // TODO: with flow_ebos,for a 2P deck, // TODO: for the 2p deck, numEq will be 3, a dummy phase is already added from the reservoir side.
71  // it will cause problem here without processing the dummy phase.
72  static const int numWellEq = GET_PROP_VALUE(TypeTag, EnablePolymer)? numEq-1 : numEq; // number of wellEq is only numEq - 1 for polymer
73  using typename Base::Mat;
74  using typename Base::BVector;
75  using typename Base::Eval;
76 
77  // sparsity pattern for the matrices
78  //[A C^T [x = [ res
79  // B D ] x_well] res_well]
80 
81  // the vector type for the res_well and x_well
82  typedef Dune::FieldVector<Scalar, numWellEq> VectorBlockWellType;
83  typedef Dune::BlockVector<VectorBlockWellType> BVectorWell;
84 
85  // the matrix type for the diagonal matrix D
86  typedef Dune::FieldMatrix<Scalar, numWellEq, numWellEq > DiagMatrixBlockWellType;
87  typedef Dune::BCRSMatrix <DiagMatrixBlockWellType> DiagMatWell;
88 
89  // the matrix type for the non-diagonal matrix B and C^T
90  typedef Dune::FieldMatrix<Scalar, numWellEq, numEq> OffDiagMatrixBlockWellType;
91  typedef Dune::BCRSMatrix<OffDiagMatrixBlockWellType> OffDiagMatWell;
92 
93  typedef DenseAd::Evaluation<double, /*size=*/numEq + numWellEq> EvalWell;
94 
95  // TODO: should these go to WellInterface?
96  static const int contiSolventEqIdx = BlackoilIndices::contiSolventEqIdx;
97  static const int contiPolymerEqIdx = BlackoilIndices::contiPolymerEqIdx;
98  static const int solventSaturationIdx = BlackoilIndices::solventSaturationIdx;
99  static const int polymerConcentrationIdx = BlackoilIndices::polymerConcentrationIdx;
100 
101 
102  StandardWell(const Well* well, const int time_step, const Wells* wells, const ModelParameters& param);
103 
104  virtual void init(const PhaseUsage* phase_usage_arg,
105  const std::vector<bool>* active_arg,
106  const std::vector<double>& depth_arg,
107  const double gravity_arg,
108  const int num_cells);
109 
110 
111  virtual void initPrimaryVariablesEvaluation() const;
112 
113  virtual void assembleWellEq(Simulator& ebosSimulator,
114  const double dt,
115  WellState& well_state,
116  bool only_wells);
117 
119  // TODO: later will check wheter we need current
120  virtual void updateWellStateWithTarget(const int current,
121  WellState& xw) const;
122 
124  virtual ConvergenceReport getWellConvergence(const std::vector<double>& B_avg) const;
125 
127  virtual void apply(const BVector& x, BVector& Ax) const;
129  virtual void apply(BVector& r) const;
130 
133  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
134  WellState& well_state) const;
135 
137  virtual void computeWellPotentials(const Simulator& ebosSimulator,
138  const WellState& well_state,
139  std::vector<double>& well_potentials) /* const */;
140 
141  virtual void updatePrimaryVariables(const WellState& well_state) const;
142 
143  virtual void solveEqAndUpdateWellState(WellState& well_state);
144 
145  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
146  const WellState& well_state); // should be const?
147  protected:
148 
149  // protected functions from the Base class
150  using Base::getAllowCrossFlow;
151  using Base::phaseUsage;
152  using Base::active;
153  using Base::flowPhaseToEbosPhaseIdx;
154  using Base::flowPhaseToEbosCompIdx;
155  using Base::numComponents;
156  using Base::wsolvent;
157  using Base::wpolymer;
158  using Base::wellHasTHPConstraints;
159  using Base::mostStrictBhpFromBhpLimits;
160 
161  // protected member variables from the Base class
162  using Base::vfp_properties_;
163  using Base::gravity_;
164  using Base::param_;
165  using Base::well_efficiency_factor_;
166  using Base::first_perf_;
167  using Base::ref_depth_;
168  using Base::perf_depth_;
169  using Base::well_cells_;
170  using Base::number_of_perforations_;
171  using Base::number_of_phases_;
172  using Base::saturation_table_number_;
173  using Base::comp_frac_;
174  using Base::well_index_;
175  using Base::index_of_well_;
176  using Base::well_controls_;
177  using Base::well_type_;
178 
179  using Base::perf_rep_radius_;
180  using Base::perf_length_;
181  using Base::bore_diameters_;
182 
183  // densities of the fluid in each perforation
184  std::vector<double> perf_densities_;
185  // pressure drop between different perforations
186  std::vector<double> perf_pressure_diffs_;
187 
188  // residuals of the well equations
189  BVectorWell resWell_;
190 
191  // two off-diagonal matrices
192  OffDiagMatWell duneB_;
193  OffDiagMatWell duneC_;
194  // diagonal matrix for the well
195  DiagMatWell invDuneD_;
196 
197  // several vector used in the matrix calculation
198  mutable BVectorWell Bx_;
199  mutable BVectorWell invDrw_;
200 
201  // the values for the primary varibles
202  // based on different solutioin strategies, the wells can have different primary variables
203  mutable std::vector<double> primary_variables_;
204 
205  // the Evaluation for the well primary variables, which contain derivativles and are used in AD calculation
206  mutable std::vector<EvalWell> primary_variables_evaluation_;
207 
208  // the saturations in the well bore under surface conditions at the beginning of the time step
209  std::vector<double> F0_;
210 
211  // TODO: this function should be moved to the base class.
212  // while it faces chanllenges for MSWell later, since the calculation of bhp
213  // based on THP is never implemented for MSWell yet.
214  EvalWell getBhp() const;
215 
216  // TODO: it is also possible to be moved to the base class.
217  EvalWell getQs(const int comp_idx) const;
218 
219  EvalWell wellVolumeFractionScaled(const int phase) const;
220 
221  EvalWell wellVolumeFraction(const int phase) const;
222 
223  EvalWell wellSurfaceVolumeFraction(const int phase) const;
224 
225  EvalWell extendEval(const Eval& in) const;
226 
227  bool crossFlowAllowed(const Simulator& ebosSimulator) const;
228 
229  // xw = inv(D)*(rw - C*x)
230  void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
231 
232  // updating the well_state based on well solution dwells
233  void updateWellState(const BVectorWell& dwells,
234  WellState& well_state) const;
235 
236  // calculate the properties for the well connections
237  // to calulate the pressure difference between well connections.
238  void computePropertiesForWellConnectionPressures(const Simulator& ebosSimulator,
239  const WellState& xw,
240  std::vector<double>& b_perf,
241  std::vector<double>& rsmax_perf,
242  std::vector<double>& rvmax_perf,
243  std::vector<double>& surf_dens_perf) const;
244 
245  // TODO: not total sure whether it is a good idea to put this function here
246  // the major reason to put here is to avoid the usage of Wells struct
247  void computeConnectionDensities(const std::vector<double>& perfComponentRates,
248  const std::vector<double>& b_perf,
249  const std::vector<double>& rsmax_perf,
250  const std::vector<double>& rvmax_perf,
251  const std::vector<double>& surf_dens_perf);
252 
253  void computeConnectionPressureDelta();
254 
255  void computeWellConnectionDensitesPressures(const WellState& xw,
256  const std::vector<double>& b_perf,
257  const std::vector<double>& rsmax_perf,
258  const std::vector<double>& rvmax_perf,
259  const std::vector<double>& surf_dens_perf);
260 
261  // computing the accumulation term for later use in well mass equations
262  void computeAccumWell();
263 
264  void computeWellConnectionPressures(const Simulator& ebosSimulator,
265  const WellState& xw);
266 
267  // TODO: to check whether all the paramters are required
268  void computePerfRate(const IntensiveQuantities& intQuants,
269  const std::vector<EvalWell>& mob_perfcells_dense,
270  const double Tw, const EvalWell& bhp, const double& cdp,
271  const bool& allow_cf, std::vector<EvalWell>& cq_s) const;
272 
273  // TODO: maybe we should provide a light version of computePerfRate, which does not include the
274  // calculation of the derivatives
275  void computeWellRatesWithBhp(const Simulator& ebosSimulator,
276  const EvalWell& bhp,
277  std::vector<double>& well_flux) const;
278 
279  std::vector<double> computeWellPotentialWithTHP(const Simulator& ebosSimulator,
280  const double initial_bhp, // bhp from BHP constraints
281  const std::vector<double>& initial_potential) const;
282 
283  template <class ValueType>
284  ValueType calculateBhpFromThp(const std::vector<ValueType>& rates, const int control_index) const;
285 
286  double calculateThpFromBhp(const std::vector<double>& rates, const int control_index, const double bhp) const;
287 
288  // get the mobility for specific perforation
289  void getMobility(const Simulator& ebosSimulator,
290  const int perf,
291  std::vector<EvalWell>& mob) const;
292 
293  double scalingFactor(const int comp_idx) const;
294  };
295 
296 }
297 
298 #include "StandardWell_impl.hpp"
299 
300 #endif // OPM_STANDARDWELL_HEADER_INCLUDED
virtual ConvergenceReport getWellConvergence(const std::vector< double > &B_avg) const
check whether the well equations get converged for this well
Definition: StandardWell_impl.hpp:1366
Definition: StandardWell.hpp:33
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 apply(const BVector &x, BVector &Ax) const
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1538
virtual void updateWellStateWithTarget(const int current, WellState &xw) const
updating the well state based the control mode specified with current
Definition: StandardWell_impl.hpp:1013
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials)
computing the well potentials for group control
Definition: StandardWell_impl.hpp:1744
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
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: StandardWell_impl.hpp:1594