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#ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
24#define OPM_WELLINTERFACE_HEADER_INCLUDED
25
26// NOTE: GasLiftSingleWell.hpp includes StandardWell.hpp which includes ourself
27// (WellInterface.hpp), so we need to forward declare GasLiftSingleWell
28// for it to be defined in this file. Similar for BlackoilWellModel
29namespace Opm {
30 template<typename TypeTag> class GasLiftSingleWell;
31 template<typename TypeTag> class BlackoilWellModel;
32}
33
34#include <opm/common/OpmLog/OpmLog.hpp>
35#include <opm/common/ErrorMacros.hpp>
36#include <opm/common/Exceptions.hpp>
37
38#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
39
40#include <opm/core/props/BlackoilPhases.hpp>
41
42#include <opm/simulators/flow/BlackoilModelParameters.hpp>
43
44#include <opm/simulators/wells/BlackoilWellModel.hpp>
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/PerforationData.hpp>
49#include <opm/simulators/wells/WellInterfaceIndices.hpp>
50#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
51#include <opm/simulators/wells/WellState.hpp>
52
53#include <opm/simulators/timestepping/ConvergenceReport.hpp>
54
55#include <opm/simulators/utils/DeferredLogger.hpp>
56
57#include <dune/common/fmatrix.hh>
58#include <dune/istl/bcrsmatrix.hh>
59#include <dune/istl/matrixmatrix.hh>
60
61#include <opm/material/densead/Evaluation.hpp>
62
63#include <cassert>
64#include <vector>
65
66namespace Opm
67{
68
69class WellInjectionProperties;
70class WellProductionProperties;
71
72template<typename TypeTag>
73class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
74 GetPropType<TypeTag, Properties::Indices>>
75{
77 GetPropType<TypeTag, Properties::Indices>>;
78public:
80
81 using Grid = GetPropType<TypeTag, Properties::Grid>;
82 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
83 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
84 using Indices = GetPropType<TypeTag, Properties::Indices>;
85 using IntensiveQuantities = GetPropType<TypeTag, Properties::IntensiveQuantities>;
86 using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
87 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
88 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
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
96 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
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
107 using WellInterfaceFluidSystem<FluidSystem>::Gas;
108 using WellInterfaceFluidSystem<FluidSystem>::Oil;
109 using WellInterfaceFluidSystem<FluidSystem>::Water;
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::EnableVapwat>();
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,
127 FluidSystem,
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,
138 const ParallelWellInfo& pw_info,
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,
163 DeferredLogger& deferred_logger,
164 const bool relax_tolerance) const = 0;
165
166 virtual void solveEqAndUpdateWellState(const SummaryState& summary_state,
167 WellState& well_state,
168 DeferredLogger& deferred_logger) = 0;
169
170 void assembleWellEq(const Simulator& simulator,
171 const double dt,
172 WellState& well_state,
173 const GroupState& group_state,
174 DeferredLogger& deferred_logger);
175
176 void assembleWellEqWithoutIteration(const Simulator& simulator,
177 const double dt,
178 WellState& well_state,
179 const GroupState& group_state,
180 DeferredLogger& deferred_logger);
181
182 // TODO: better name or further refactoring the function to make it more clear
183 void prepareWellBeforeAssembling(const Simulator& simulator,
184 const double dt,
185 WellState& well_state,
186 const GroupState& group_state,
187 DeferredLogger& deferred_logger);
188
189
190 virtual void computeWellRatesWithBhp(
191 const Simulator& simulator,
192 const double& bhp,
193 std::vector<double>& well_flux,
194 DeferredLogger& deferred_logger
195 ) const = 0;
196
197 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
198 const Simulator& simulator,
199 const SummaryState& summary_state,
200 const double alq_value,
201 DeferredLogger& deferred_logger
202 ) const = 0;
203
206 virtual void recoverWellSolutionAndUpdateWellState(const SummaryState& summary_state,
207 const BVector& x,
208 WellState& well_state,
209 DeferredLogger& deferred_logger) = 0;
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& simulator,
219 const WellState& well_state,
220 std::vector<double>& well_potentials,
221 DeferredLogger& deferred_logger) = 0;
222
223 virtual void updateWellStateWithTarget(const Simulator& simulator,
224 const GroupState& group_state,
225 WellState& well_state,
226 DeferredLogger& deferred_logger) const;
227
228 virtual void computeWellRatesWithBhpIterations(const Simulator& simulator,
229 const Scalar& bhp,
230 std::vector<double>& well_flux,
231 DeferredLogger& deferred_logger) const = 0;
232
233 bool updateWellStateWithTHPTargetProd(const Simulator& simulator,
234 WellState& well_state,
235 DeferredLogger& deferred_logger) const;
236
237 enum class IndividualOrGroup { Individual, Group, Both };
238 bool updateWellControl(const Simulator& 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& 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,
250 DeferredLogger& deferred_logger,
251 const bool fixed_control = false,
252 const bool fixed_status = false);
253
254 virtual void updatePrimaryVariables(const SummaryState& summary_state,
255 const WellState& well_state,
256 DeferredLogger& deferred_logger) = 0;
257
258 virtual void calculateExplicitQuantities(const Simulator& simulator,
259 const WellState& well_state,
260 DeferredLogger& deferred_logger) = 0; // should be const?
261
262 virtual void updateProductivityIndex(const Simulator& simulator,
263 const WellProdIndexCalculator& wellPICalc,
264 WellState& well_state,
265 DeferredLogger& deferred_logger) const = 0;
266
267 virtual double connectionDensity(const int globalConnIdx,
268 const int openConnIdx) const = 0;
269
272 {
273 return false;
274 }
275
276 // Add well contributions to matrix
277 virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
278
279 virtual void addWellPressureEquations(PressureMatrix& mat,
280 const BVector& x,
281 const int pressureVarIndex,
282 const bool use_well_weights,
283 const WellState& well_state) const = 0;
284
285 void addCellRates(RateVector& rates, int cellIdx) const;
286
287 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
288
289 // TODO: theoretically, it should be a const function
290 // Simulator is not const is because that assembleWellEq is non-const Simulator
291 void wellTesting(const Simulator& simulator,
292 const double simulation_time,
293 /* const */ WellState& well_state, const GroupState& group_state, WellTestState& welltest_state,
294 DeferredLogger& deferred_logger);
295
296 void checkWellOperability(const Simulator& simulator,
297 const WellState& well_state,
298 DeferredLogger& deferred_logger);
299
300 bool gliftBeginTimeStepWellTestIterateWellEquations(
301 const Simulator& simulator,
302 const double dt,
303 WellState& well_state,
304 const GroupState &group_state,
305 DeferredLogger& deferred_logger);
306
307 void gliftBeginTimeStepWellTestUpdateALQ(const Simulator& simulator,
308 WellState& well_state,
309 DeferredLogger& deferred_logger);
310
311 // check whether the well is operable under the current reservoir condition
312 // mostly related to BHP limit and THP limit
313 void updateWellOperability(const Simulator& simulator,
314 const WellState& well_state,
315 DeferredLogger& deferred_logger);
316
317 bool updateWellOperabilityFromWellEq(const Simulator& simulator,
318 const WellState& well_state,
319 DeferredLogger& deferred_logger);
320
321 // update perforation water throughput based on solved water rate
322 virtual void updateWaterThroughput(const double dt, WellState& well_state) const = 0;
323
326 virtual std::vector<double> computeCurrentWellRates(const Simulator& simulator,
327 DeferredLogger& deferred_logger) const = 0;
328
332 void updateWellStateRates(const Simulator& simulator,
333 WellState& well_state,
334 DeferredLogger& deferred_logger) const;
335
336 void solveWellEquation(const Simulator& simulator,
337 WellState& well_state,
338 const GroupState& group_state,
339 DeferredLogger& deferred_logger);
340
341 const std::vector<RateVector>& connectionRates() const
342 {
343 return connectionRates_;
344 }
345
346 virtual std::vector<double> getPrimaryVars() const
347 {
348 return {};
349 }
350
351 virtual int setPrimaryVars(std::vector<double>::const_iterator)
352 {
353 return 0;
354 }
355
356 std::vector<double> wellIndex(const int perf, const IntensiveQuantities& intQuants, const double trans_mult, const SingleWellState& ws) const;
357
358 void updateConnectionDFactor(const Simulator& simulator, SingleWellState& ws) const;
359
360 void updateConnectionTransmissibilityFactor(const Simulator& simulator, SingleWellState& ws) const;
361
362
363protected:
364 // simulation parameters
365 const ModelParameters& param_;
366 std::vector<RateVector> connectionRates_;
367 std::vector< Scalar > B_avg_;
368 bool changed_to_stopped_this_step_ = false;
369 bool thp_update_iterations = false;
370
371 double wpolymer() const;
372
373 double wfoam() const;
374
375 double wsalt() const;
376
377 double wmicrobes() const;
378
379 double woxygen() const;
380
381 double wurea() const;
382
383 virtual double getRefDensity() const = 0;
384
385 // Component fractions for each phase for the well
386 const std::vector<double>& compFrac() const;
387
388 std::vector<double> initialWellRateFractions(const Simulator& simulator,
389 const WellState& well_state) const;
390
391 // check whether the well is operable under BHP limit with current reservoir condition
392 virtual void checkOperabilityUnderBHPLimit(const WellState& well_state,
393 const Simulator& simulator,
394 DeferredLogger& deferred_logger) = 0;
395
396 // check whether the well is operable under THP limit with current reservoir condition
397 virtual void checkOperabilityUnderTHPLimit(const Simulator& simulator,
398 const WellState& well_state,
399 DeferredLogger& deferred_logger) = 0;
400
401 virtual void updateIPR(const Simulator& simulator,
402 DeferredLogger& deferred_logger) const=0;
403
404 virtual void assembleWellEqWithoutIteration(const Simulator& simulator,
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 // iterate well equations with the specified control until converged
413 virtual bool iterateWellEqWithControl(const Simulator& simulator,
414 const double dt,
415 const WellInjectionControls& inj_controls,
416 const WellProductionControls& prod_controls,
417 WellState& well_state,
418 const GroupState& group_state,
419 DeferredLogger& deferred_logger) = 0;
420
421 virtual bool iterateWellEqWithSwitching(const Simulator& simulator,
422 const double dt,
423 const WellInjectionControls& inj_controls,
424 const WellProductionControls& prod_controls,
425 WellState& well_state,
426 const GroupState& group_state,
427 DeferredLogger& deferred_logger,
428 const bool fixed_control = false,
429 const bool fixed_status = false) = 0;
430
431 virtual void updateIPRImplicit(const Simulator& simulator,
432 WellState& well_state,
433 DeferredLogger& deferred_logger) = 0;
434
435 bool iterateWellEquations(const Simulator& simulator,
436 const double dt,
437 WellState& well_state,
438 const GroupState& group_state,
439 DeferredLogger& deferred_logger);
440
441 bool solveWellWithTHPConstraint(const Simulator& simulator,
442 const double dt,
443 const Well::InjectionControls& inj_controls,
444 const Well::ProductionControls& prod_controls,
445 WellState& well_state,
446 const GroupState& group_state,
447 DeferredLogger& deferred_logger);
448
449 std::optional<double> estimateOperableBhp(const Simulator& simulator,
450 const double dt,
451 WellState& well_state,
452 const SummaryState& summary_state,
453 DeferredLogger& deferred_logger);
454
455 bool solveWellWithBhp(const Simulator& simulator,
456 const double dt,
457 const double bhp,
458 WellState& well_state,
459 DeferredLogger& deferred_logger);
460
461 bool solveWellWithZeroRate(const Simulator& simulator,
462 const double dt,
463 WellState& well_state,
464 DeferredLogger& deferred_logger);
465
466 bool solveWellForTesting(const Simulator& simulator,
467 WellState& well_state,
468 const GroupState& group_state,
469 DeferredLogger& deferred_logger);
470
471 Eval getPerfCellPressure(const FluidState& fs) const;
472
473 // get the mobility for specific perforation
474 template<class Value, class Callback>
475 void getMobility(const Simulator& simulator,
476 const int perf,
477 std::vector<Value>& mob,
478 Callback& extendEval,
479 [[maybe_unused]] DeferredLogger& deferred_logger) const;
480
481 void computeConnLevelProdInd(const FluidState& fs,
482 const std::function<double(const double)>& connPICalc,
483 const std::vector<Scalar>& mobility,
484 double* connPI) const;
485
486 void computeConnLevelInjInd(const FluidState& fs,
487 const Phase preferred_phase,
488 const std::function<double(const double)>& connIICalc,
489 const std::vector<Scalar>& mobility,
490 double* connII,
491 DeferredLogger& deferred_logger) const;
492
493 double computeConnectionDFactor(const int perf,
494 const IntensiveQuantities& intQuants,
495 const SingleWellState& ws) const;
496
497};
498
499} // namespace Opm
500
501#include "WellInterface_impl.hpp"
502
503#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:35
Definition WellInterface.hpp:75
virtual std::vector< double > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
Compute well rates based on current reservoir conditions and well variables.
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
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:271
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
void updateWellStateRates(const Simulator &simulator, WellState &well_state, DeferredLogger &deferred_logger) const
Modify the well_state's rates if there is only one nonzero rate.
Definition WellInterface_impl.hpp:1448
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:61
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Solver parameters for the BlackoilModel.
Definition BlackoilModelParameters.hpp:484
Definition BlackoilPhases.hpp:46