My Project
WellInterface.hpp
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2017 IRIS
5 Copyright 2019 Norce
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_WELLINTERFACE_HEADER_INCLUDED
25#define OPM_WELLINTERFACE_HEADER_INCLUDED
26
27#include <opm/common/OpmLog/OpmLog.hpp>
28#include <opm/common/ErrorMacros.hpp>
29#include <opm/common/Exceptions.hpp>
30
31#include <opm/input/eclipse/Schedule/Well/Well.hpp>
32#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
33
34#include <opm/core/props/BlackoilPhases.hpp>
35
36#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
37#include <opm/simulators/wells/WellState.hpp>
38// NOTE: GasLiftSingleWell.hpp includes StandardWell.hpp which includes ourself
39// (WellInterface.hpp), so we need to forward declare GasLiftSingleWell
40// for it to be defined in this file. Similar for BlackoilWellModel
41namespace Opm {
42 template<typename TypeTag> class GasLiftSingleWell;
43 template<typename TypeTag> class BlackoilWellModel;
44}
45#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
46#include <opm/simulators/wells/GasLiftSingleWell.hpp>
47#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
48#include <opm/simulators/wells/BlackoilWellModel.hpp>
49#include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
50
51#include <opm/simulators/utils/DeferredLogger.hpp>
52
53#include<dune/common/fmatrix.hh>
54#include<dune/istl/bcrsmatrix.hh>
55#include<dune/istl/matrixmatrix.hh>
56
57#include <opm/material/densead/Evaluation.hpp>
58
59#include <opm/simulators/wells/WellInterfaceIndices.hpp>
60#include <opm/simulators/timestepping/ConvergenceReport.hpp>
61
62#include <cassert>
63#include <vector>
64
65namespace Opm
66{
67
68template<typename TypeTag>
69class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
70 GetPropType<TypeTag, Properties::Indices>,
71 GetPropType<TypeTag, Properties::Scalar>>
72{
73public:
74
76
77 using Grid = GetPropType<TypeTag, Properties::Grid>;
78 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
79 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
80 using Indices = GetPropType<TypeTag, Properties::Indices>;
81 using IntensiveQuantities = GetPropType<TypeTag, Properties::IntensiveQuantities>;
82 using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
83 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
84 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
86 using GLiftOptWells = typename BlackoilWellModel<TypeTag>::GLiftOptWells;
87 using GLiftProdWells = typename BlackoilWellModel<TypeTag>::GLiftProdWells;
88 using GLiftWellStateMap =
89 typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
90 using GLiftSyncGroups = typename GasLiftSingleWellGeneric::GLiftSyncGroups;
91
92 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
93
94 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
95 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
96 using BVector = Dune::BlockVector<VectorBlockType>;
97 using Eval = DenseAd::Evaluation<Scalar, /*size=*/Indices::numEq>;
98 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
99
100 using RateConverterType =
102
103 using WellInterfaceFluidSystem<FluidSystem>::Gas;
104 using WellInterfaceFluidSystem<FluidSystem>::Oil;
105 using WellInterfaceFluidSystem<FluidSystem>::Water;
106 using RatioLimitCheckReport = typename WellInterfaceFluidSystem<FluidSystem>::RatioLimitCheckReport;
107
108 static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
109 static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
110 static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
111 static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
112 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
113 // flag for polymer molecular weight related
114 static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
115 static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
116 static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
117 static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableEvaporation>();
118 static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
119 static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
120
121 // For the conversion between the surface volume rate and reservoir voidage rate
122 using FluidState = BlackOilFluidState<Eval,
123 FluidSystem,
124 has_temperature,
125 has_energy,
126 Indices::compositionSwitchIdx >= 0,
127 has_watVapor,
128 has_brine,
129 has_saltPrecip,
130 Indices::numPhases >;
132 WellInterface(const Well& well,
133 const ParallelWellInfo& pw_info,
134 const int time_step,
135 const ModelParameters& param,
136 const RateConverterType& rate_converter,
137 const int pvtRegionIdx,
138 const int num_components,
139 const int num_phases,
140 const int index_of_well,
141 const std::vector<PerforationData>& perf_data);
142
144 virtual ~WellInterface() = default;
145
146 virtual void init(const PhaseUsage* phase_usage_arg,
147 const std::vector<double>& depth_arg,
148 const double gravity_arg,
149 const int num_cells,
150 const std::vector< Scalar >& B_avg,
151 const bool changed_to_open_this_step);
152
153 virtual void initPrimaryVariablesEvaluation() const = 0;
154
155 virtual ConvergenceReport getWellConvergence(const WellState& well_state, const std::vector<double>& B_avg, DeferredLogger& deferred_logger, const bool relax_tolerance) const = 0;
156
157 virtual void solveEqAndUpdateWellState(WellState& well_state, DeferredLogger& deferred_logger) = 0;
158
159 void assembleWellEq(const Simulator& ebosSimulator,
160 const double dt,
161 WellState& well_state,
162 const GroupState& group_state,
163 DeferredLogger& deferred_logger);
164
165 virtual void computeWellRatesWithBhp(
166 const Simulator& ebosSimulator,
167 const double& bhp,
168 std::vector<double>& well_flux,
169 DeferredLogger& deferred_logger
170 ) const = 0;
171
172 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
173 const Simulator& ebos_simulator,
174 const SummaryState& summary_state,
175 DeferredLogger& deferred_logger,
176 double alq_value
177 ) const = 0;
178
181 virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
182 WellState& well_state,
183 DeferredLogger& deferred_logger) const = 0;
184
186 virtual void apply(const BVector& x, BVector& Ax) const = 0;
187
189 virtual void apply(BVector& r) const = 0;
190
191 // TODO: before we decide to put more information under mutable, this function is not const
192 virtual void computeWellPotentials(const Simulator& ebosSimulator,
193 const WellState& well_state,
194 std::vector<double>& well_potentials,
195 DeferredLogger& deferred_logger) = 0;
196
197 virtual void updateWellStateWithTarget(const Simulator& ebos_simulator,
198 const GroupState& group_state,
199 WellState& well_state,
200 DeferredLogger& deferred_logger) const;
201
202 enum class IndividualOrGroup { Individual, Group, Both };
203 bool updateWellControl(const Simulator& ebos_simulator,
204 const IndividualOrGroup iog,
205 WellState& well_state,
206 const GroupState& group_state,
207 DeferredLogger& deferred_logger) /* const */;
208
209 virtual void updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_logger) const = 0;
210
211 virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
212 const WellState& well_state,
213 DeferredLogger& deferred_logger) = 0; // should be const?
214
215 virtual void updateProductivityIndex(const Simulator& ebosSimulator,
216 const WellProdIndexCalculator& wellPICalc,
217 WellState& well_state,
218 DeferredLogger& deferred_logger) const = 0;
219
222 {
223 return false;
224 }
225
226 // Add well contributions to matrix
227 virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
228
229 virtual bool isPressureControlled(const WellState& well_state) const;
230
231 virtual void addWellPressureEquations(PressureMatrix& mat,
232 const BVector& x,
233 const int pressureVarIndex,
234 const bool use_well_weights,
235 const WellState& well_state) const = 0;
236
237 void addCellRates(RateVector& rates, int cellIdx) const;
238
239 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
240
241 template <class EvalWell>
242 Eval restrictEval(const EvalWell& in) const
243 {
244 Eval out = 0.0;
245 out.setValue(in.value());
246 for (int eqIdx = 0; eqIdx < Indices::numEq; ++eqIdx) {
247 out.setDerivative(eqIdx, in.derivative(eqIdx));
248 }
249 return out;
250 }
251
252 // TODO: theoretically, it should be a const function
253 // Simulator is not const is because that assembleWellEq is non-const Simulator
254 void wellTesting(const Simulator& simulator,
255 const double simulation_time,
256 /* const */ WellState& well_state, const GroupState& group_state, WellTestState& welltest_state,
257 DeferredLogger& deferred_logger);
258
259 void checkWellOperability(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger);
260
261 void gliftBeginTimeStepWellTestUpdateALQ(const Simulator& ebos_simulator,
262 WellState& well_state,
263 DeferredLogger& deferred_logger);
264
265 // check whether the well is operable under the current reservoir condition
266 // mostly related to BHP limit and THP limit
267 void updateWellOperability(const Simulator& ebos_simulator,
268 const WellState& well_state,
269 DeferredLogger& deferred_logger);
270
271
272 // update perforation water throughput based on solved water rate
273 virtual void updateWaterThroughput(const double dt, WellState& well_state) const = 0;
274
277 virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
278 DeferredLogger& deferred_logger) const = 0;
279
283 void updateWellStateRates(const Simulator& ebosSimulator,
284 WellState& well_state,
285 DeferredLogger& deferred_logger) const;
286
287 void solveWellEquation(const Simulator& ebosSimulator,
288 WellState& well_state,
289 const GroupState& group_state,
290 DeferredLogger& deferred_logger);
291
292protected:
293
294 // simulation parameters
295 const ModelParameters& param_;
296
297 std::vector<RateVector> connectionRates_;
298
299 std::vector< Scalar > B_avg_;
300
301 bool changed_to_stopped_this_step_ = false;
302
303 double wpolymer() const;
304
305 double wfoam() const;
306
307 double wsalt() const;
308
309 double wmicrobes() const;
310
311 double woxygen() const;
312
313 double wurea() const;
314
315 virtual double getRefDensity() const = 0;
316
317 // Component fractions for each phase for the well
318 const std::vector<double>& compFrac() const;
319
320 std::vector<double> initialWellRateFractions(const Simulator& ebosSimulator, const WellState& well_state) const;
321
322 // check whether the well is operable under BHP limit with current reservoir condition
323 virtual void checkOperabilityUnderBHPLimit(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) =0;
324
325 // check whether the well is operable under THP limit with current reservoir condition
326 virtual void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) =0;
327
328 virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const=0;
329
330 virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
331 const double dt,
332 const Well::InjectionControls& inj_controls,
333 const Well::ProductionControls& prod_controls,
334 WellState& well_state,
335 const GroupState& group_state,
336 DeferredLogger& deferred_logger) = 0;
337
338 // iterate well equations with the specified control until converged
339 virtual bool iterateWellEqWithControl(const Simulator& ebosSimulator,
340 const double dt,
341 const Well::InjectionControls& inj_controls,
342 const Well::ProductionControls& prod_controls,
343 WellState& well_state,
344 const GroupState& group_state,
345 DeferredLogger& deferred_logger) = 0;
346
347 bool iterateWellEquations(const Simulator& ebosSimulator,
348 const double dt,
349 WellState& well_state,
350 const GroupState& group_state,
351 DeferredLogger& deferred_logger);
352
353 bool solveWellForTesting(const Simulator& ebosSimulator, WellState& well_state, const GroupState& group_state,
354 DeferredLogger& deferred_logger);
355
356 Eval getPerfCellPressure(const FluidState& fs) const;
357
358
359};
360
361}
362
363#include "WellInterface_impl.hpp"
364
365#endif // OPM_WELLINTERFACE_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: GasLiftSingleWell.hpp:39
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:243
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:68
Definition: WellInterfaceFluidSystem.hpp:46
Definition: WellInterfaceIndices.hpp:35
Definition: WellInterface.hpp:72
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
void updateWellStateRates(const Simulator &ebosSimulator, WellState &well_state, DeferredLogger &deferred_logger) const
Modify the well_state's rates if there is only one nonzero rate.
Definition: WellInterface_impl.hpp:1137
virtual std::vector< double > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const =0
Compute well rates based on current reservoir conditions and well variables.
virtual void apply(BVector &r) const =0
r = r - C D^-1 Rw
virtual ~WellInterface()=default
Virtual destructor.
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) const =0
using the solution x to recover the solution xw for wells and applying xw to update Well State
WellInterface(const Well &well, const ParallelWellInfo &pw_info, const int time_step, const ModelParameters &param, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_components, const int num_phases, const int index_of_well, const std::vector< PerforationData > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:35
virtual bool jacobianContainsWellContributions() const
Wether the Jacobian will also have well contributions in it.
Definition: WellInterface.hpp:221
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
Definition: WellInterfaceFluidSystem.hpp:129