My Project
Loading...
Searching...
No Matches
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{
78public:
80
90 using GLiftOptWells = typename BlackoilWellModel<TypeTag>::GLiftOptWells;
91 using GLiftProdWells = typename BlackoilWellModel<TypeTag>::GLiftProdWells;
92 using GLiftWellStateMap =
93 typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
94 using GLiftSyncGroups = typename GasLiftSingleWellGeneric::GLiftSyncGroups;
95
97
98 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
99 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
100 using Eval = typename Base::Eval;
101 using BVector = Dune::BlockVector<VectorBlockType>;
102 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
103
104 using RateConverterType =
106
110
111 static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
112 static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
113 static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
114 static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
115 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
116 // flag for polymer molecular weight related
117 static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
118 static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
119 static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
120 static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableEvaporation>();
121 static constexpr bool has_disgas_in_water = getPropValue<TypeTag, Properties::EnableDisgasInWater>();
122 static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
123 static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
124
125 // For the conversion between the surface volume rate and reservoir voidage rate
126 using FluidState = BlackOilFluidState<Eval,
128 has_temperature,
129 has_energy,
130 Indices::compositionSwitchIdx >= 0,
131 has_watVapor,
132 has_brine,
133 has_saltPrecip,
134 has_disgas_in_water,
135 Indices::numPhases >;
137 WellInterface(const Well& well,
139 const int time_step,
140 const ModelParameters& param,
141 const RateConverterType& rate_converter,
142 const int pvtRegionIdx,
143 const int num_components,
144 const int num_phases,
145 const int index_of_well,
146 const std::vector<PerforationData>& perf_data);
147
149 virtual ~WellInterface() = default;
150
151 virtual void init(const PhaseUsage* phase_usage_arg,
152 const std::vector<double>& depth_arg,
153 const double gravity_arg,
154 const int num_cells,
155 const std::vector< Scalar >& B_avg,
156 const bool changed_to_open_this_step);
157
158 virtual void initPrimaryVariablesEvaluation() = 0;
159
160 virtual ConvergenceReport getWellConvergence(const SummaryState& summary_state,
161 const WellState& well_state,
162 const std::vector<double>& B_avg,
164 const bool relax_tolerance) const = 0;
165
166 virtual void solveEqAndUpdateWellState(const SummaryState& summary_state,
167 WellState& well_state,
169
170 void assembleWellEq(const Simulator& ebosSimulator,
171 const double dt,
172 WellState& well_state,
173 const GroupState& group_state,
175
176 void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
177 const double dt,
178 WellState& well_state,
179 const GroupState& group_state,
181
182 // TODO: better name or further refactoring the function to make it more clear
183 void prepareWellBeforeAssembling(const Simulator& ebosSimulator,
184 const double dt,
185 WellState& well_state,
186 const GroupState& group_state,
188
189
190 virtual void computeWellRatesWithBhp(
191 const Simulator& ebosSimulator,
192 const double& bhp,
193 std::vector<double>& well_flux,
195 ) const = 0;
196
197 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
198 const Simulator& ebos_simulator,
200 const double alq_value,
202 ) const = 0;
203
207 const BVector& x,
208 WellState& well_state,
210
212 virtual void apply(const BVector& x, BVector& Ax) const = 0;
213
215 virtual void apply(BVector& r) const = 0;
216
217 // TODO: before we decide to put more information under mutable, this function is not const
218 virtual void computeWellPotentials(const Simulator& ebosSimulator,
219 const WellState& well_state,
220 std::vector<double>& well_potentials,
222
223 virtual void updateWellStateWithTarget(const Simulator& ebos_simulator,
224 const GroupState& group_state,
225 WellState& well_state,
227
228 virtual void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
229 const Scalar& bhp,
230 std::vector<double>& well_flux,
232
233 bool updateWellStateWithTHPTargetProd(const Simulator& ebos_simulator,
234 WellState& well_state,
236
237 enum class IndividualOrGroup { Individual, Group, Both };
238 bool updateWellControl(const Simulator& ebos_simulator,
239 const IndividualOrGroup iog,
240 WellState& well_state,
241 const GroupState& group_state,
242 DeferredLogger& deferred_logger) /* const */;
243
244 bool updateWellControlAndStatusLocalIteration(const Simulator& ebos_simulator,
245 WellState& well_state,
246 const GroupState& group_state,
247 const Well::InjectionControls& inj_controls,
248 const Well::ProductionControls& prod_controls,
249 const double WQTotal,
251
252 virtual void updatePrimaryVariables(const SummaryState& summary_state,
253 const WellState& well_state,
255
256 virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
257 const WellState& well_state,
258 DeferredLogger& deferred_logger) = 0; // should be const?
259
260 virtual void updateProductivityIndex(const Simulator& ebosSimulator,
262 WellState& well_state,
264
265 virtual double connectionDensity(const int globalConnIdx,
266 const int openConnIdx) const = 0;
267
270 {
271 return false;
272 }
273
274 // Add well contributions to matrix
275 virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
276
277 virtual void addWellPressureEquations(PressureMatrix& mat,
278 const BVector& x,
279 const int pressureVarIndex,
280 const bool use_well_weights,
281 const WellState& well_state) const = 0;
282
283 void addCellRates(RateVector& rates, int cellIdx) const;
284
285 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
286
287 // TODO: theoretically, it should be a const function
288 // Simulator is not const is because that assembleWellEq is non-const Simulator
289 void wellTesting(const Simulator& simulator,
290 const double simulation_time,
291 /* const */ WellState& well_state, const GroupState& group_state, WellTestState& welltest_state,
293
294 void checkWellOperability(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger);
295
296 bool gliftBeginTimeStepWellTestIterateWellEquations(
297 const Simulator& ebos_simulator,
298 const double dt,
299 WellState& well_state,
300 const GroupState &group_state,
302
303 void gliftBeginTimeStepWellTestUpdateALQ(const Simulator& ebos_simulator,
304 WellState& well_state,
306
307 // check whether the well is operable under the current reservoir condition
308 // mostly related to BHP limit and THP limit
309 void updateWellOperability(const Simulator& ebos_simulator,
310 const WellState& well_state,
312
313 bool updateWellOperabilityFromWellEq(const Simulator& ebos_simulator,
314 const WellState& well_state,
316
317 // update perforation water throughput based on solved water rate
318 virtual void updateWaterThroughput(const double dt, WellState& well_state) const = 0;
319
322 virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
324
328 void updateWellStateRates(const Simulator& ebosSimulator,
329 WellState& well_state,
331
332 void solveWellEquation(const Simulator& ebosSimulator,
333 WellState& well_state,
334 const GroupState& group_state,
336
337 const std::vector<RateVector>& connectionRates() const
338 {
339 return connectionRates_;
340 }
341
342 virtual std::vector<double> getPrimaryVars() const
343 {
344 return {};
345 }
346
347 virtual int setPrimaryVars(std::vector<double>::const_iterator)
348 {
349 return 0;
350 }
351
352protected:
353 // simulation parameters
354 const ModelParameters& param_;
355 std::vector<RateVector> connectionRates_;
356 std::vector< Scalar > B_avg_;
357 bool changed_to_stopped_this_step_ = false;
358 bool thp_update_iterations = false;
359
360 double wpolymer() const;
361
362 double wfoam() const;
363
364 double wsalt() const;
365
366 double wmicrobes() const;
367
368 double woxygen() const;
369
370 double wurea() const;
371
372 virtual double getRefDensity() const = 0;
373
374 // Component fractions for each phase for the well
375 const std::vector<double>& compFrac() const;
376
377 std::vector<double> initialWellRateFractions(const Simulator& ebosSimulator, const WellState& well_state) const;
378
379 // check whether the well is operable under BHP limit with current reservoir condition
380 virtual void checkOperabilityUnderBHPLimit(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) =0;
381
382 // check whether the well is operable under THP limit with current reservoir condition
383 virtual void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) =0;
384
385 virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const=0;
386
387 virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
388 const double dt,
389 const WellInjectionControls& inj_controls,
390 const WellProductionControls& prod_controls,
391 WellState& well_state,
392 const GroupState& group_state,
393 DeferredLogger& deferred_logger) = 0;
394
395 // iterate well equations with the specified control until converged
396 virtual bool iterateWellEqWithControl(const Simulator& ebosSimulator,
397 const double dt,
398 const WellInjectionControls& inj_controls,
399 const WellProductionControls& prod_controls,
400 WellState& well_state,
401 const GroupState& group_state,
402 DeferredLogger& deferred_logger) = 0;
403
404 virtual bool iterateWellEqWithSwitching(const Simulator& ebosSimulator,
405 const double dt,
406 const WellInjectionControls& inj_controls,
407 const WellProductionControls& prod_controls,
408 WellState& well_state,
409 const GroupState& group_state,
410 DeferredLogger& deferred_logger) = 0;
411
412 bool iterateWellEquations(const Simulator& ebosSimulator,
413 const double dt,
414 WellState& well_state,
415 const GroupState& group_state,
416 DeferredLogger& deferred_logger);
417
418 bool solveWellForTesting(const Simulator& ebosSimulator, WellState& well_state, const GroupState& group_state,
419 DeferredLogger& deferred_logger);
420
421 Eval getPerfCellPressure(const FluidState& fs) const;
422
423 // get the mobility for specific perforation
424 template<class Value, class Callback>
425 void getMobility(const Simulator& ebosSimulator,
426 const int perf,
427 std::vector<Value>& mob,
428 Callback& extendEval,
429 [[maybe_unused]] DeferredLogger& deferred_logger) const;
430
431 void computeConnLevelProdInd(const FluidState& fs,
432 const std::function<double(const double)>& connPICalc,
433 const std::vector<Scalar>& mobility,
434 double* connPI) const;
435
436 void computeConnLevelInjInd(const FluidState& fs,
437 const Phase preferred_phase,
438 const std::function<double(const double)>& connIICalc,
439 const std::vector<Scalar>& mobility,
440 double* connII,
441 DeferredLogger& deferred_logger) const;
442};
443
444}
445
446#include "WellInterface_impl.hpp"
447
448#endif // OPM_WELLINTERFACE_HEADER_INCLUDED
Definition AquiferInterface.hpp:35
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 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:35
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:1274
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 bool jacobianContainsWellContributions() const
Wether the Jacobian will also have well contributions in it.
Definition WellInterface.hpp:269
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:451
Definition BlackoilPhases.hpp:46