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