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/simulators/wells/BlackoilWellModel.hpp>
41#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
42#include <opm/simulators/wells/GasLiftSingleWell.hpp>
43#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
44#include <opm/simulators/wells/PerforationData.hpp>
45#include <opm/simulators/wells/WellInterfaceIndices.hpp>
46#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
47#include <opm/simulators/wells/WellState.hpp>
48
49#include <opm/simulators/timestepping/ConvergenceReport.hpp>
50
51#include <opm/simulators/utils/BlackoilPhases.hpp>
52#include <opm/simulators/utils/DeferredLogger.hpp>
53
54#include <dune/common/fmatrix.hh>
55#include <dune/istl/bcrsmatrix.hh>
56#include <dune/istl/matrixmatrix.hh>
57
58#include <opm/material/densead/Evaluation.hpp>
59
60#include <vector>
61
62namespace Opm
63{
64
65class WellInjectionProperties;
66class WellProductionProperties;
67
68template<typename TypeTag>
69class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
70 GetPropType<TypeTag, Properties::Indices>>
71{
74public:
85 using GLiftOptWells = typename BlackoilWellModel<TypeTag>::GLiftOptWells;
86 using GLiftProdWells = typename BlackoilWellModel<TypeTag>::GLiftProdWells;
87 using GLiftWellStateMap =
88 typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
89 using GLiftSyncGroups = typename GasLiftSingleWellGeneric<Scalar>::GLiftSyncGroups;
90
91 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
92 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
93 using Eval = typename Base::Eval;
94 using BVector = Dune::BlockVector<VectorBlockType>;
95 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, 1, 1>>;
96
97 using RateConverterType =
98 typename WellInterfaceFluidSystem<FluidSystem>::RateConverterType;
99
100 using WellInterfaceFluidSystem<FluidSystem>::Gas;
101 using WellInterfaceFluidSystem<FluidSystem>::Oil;
102 using WellInterfaceFluidSystem<FluidSystem>::Water;
103
104 using ModelParameters = typename Base::ModelParameters;
105
106 static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
107 static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
108 static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
109 static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
110 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
111 // flag for polymer molecular weight related
112 static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
113 static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
114 static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
115 static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableVapwat>();
116 static constexpr bool has_disgas_in_water = getPropValue<TypeTag, Properties::EnableDisgasInWater>();
117 static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
118 static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
119
120 // For the conversion between the surface volume rate and reservoir voidage rate
121 using FluidState = BlackOilFluidState<Eval,
122 FluidSystem,
123 has_temperature,
124 has_energy,
125 Indices::compositionSwitchIdx >= 0,
126 has_watVapor,
127 has_brine,
128 has_saltPrecip,
129 has_disgas_in_water,
130 Indices::numPhases >;
132 WellInterface(const Well& well,
133 const ParallelWellInfo<Scalar>& 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<Scalar>>& perf_data);
142
144 virtual ~WellInterface() = default;
145
146 virtual void init(const PhaseUsage* phase_usage_arg,
147 const std::vector<Scalar>& depth_arg,
148 const Scalar gravity_arg,
149 const std::vector<Scalar>& B_avg,
150 const bool changed_to_open_this_step);
151
152 virtual void initPrimaryVariablesEvaluation() = 0;
153
154 virtual ConvergenceReport getWellConvergence(const Simulator& simulator,
155 const WellState<Scalar>& well_state,
156 const std::vector<Scalar>& B_avg,
157 DeferredLogger& deferred_logger,
158 const bool relax_tolerance) const = 0;
159
160 virtual void solveEqAndUpdateWellState(const Simulator& simulator,
161 WellState<Scalar>& well_state,
162 DeferredLogger& deferred_logger) = 0;
163
164 void assembleWellEq(const Simulator& simulator,
165 const double dt,
166 WellState<Scalar>& well_state,
167 const GroupState<Scalar>& group_state,
168 DeferredLogger& deferred_logger);
169
170 void assembleWellEqWithoutIteration(const Simulator& simulator,
171 const double dt,
172 WellState<Scalar>& well_state,
173 const GroupState<Scalar>& group_state,
174 DeferredLogger& deferred_logger);
175
176 // TODO: better name or further refactoring the function to make it more clear
177 void prepareWellBeforeAssembling(const Simulator& simulator,
178 const double dt,
179 WellState<Scalar>& well_state,
180 const GroupState<Scalar>& group_state,
181 DeferredLogger& deferred_logger);
182
183
184 virtual void computeWellRatesWithBhp(const Simulator& ebosSimulator,
185 const Scalar& bhp,
186 std::vector<Scalar>& well_flux,
187 DeferredLogger& deferred_logger) const = 0;
188
189 virtual std::optional<Scalar>
190 computeBhpAtThpLimitProdWithAlq(const Simulator& ebos_simulator,
191 const SummaryState& summary_state,
192 const Scalar alq_value,
193 DeferredLogger& deferred_logger) const = 0;
194
197 virtual void recoverWellSolutionAndUpdateWellState(const Simulator& simulator,
198 const BVector& x,
199 WellState<Scalar>& well_state,
200 DeferredLogger& deferred_logger) = 0;
201
203 virtual void apply(const BVector& x, BVector& Ax) const = 0;
204
206 virtual void apply(BVector& r) const = 0;
207
208 // TODO: before we decide to put more information under mutable, this function is not const
209 virtual void computeWellPotentials(const Simulator& simulator,
210 const WellState<Scalar>& well_state,
211 std::vector<Scalar>& well_potentials,
212 DeferredLogger& deferred_logger) = 0;
213
214 virtual void updateWellStateWithTarget(const Simulator& simulator,
215 const GroupState<Scalar>& group_state,
216 WellState<Scalar>& well_state,
217 DeferredLogger& deferred_logger) const;
218
219 virtual void computeWellRatesWithBhpIterations(const Simulator& simulator,
220 const Scalar& bhp,
221 std::vector<Scalar>& well_flux,
222 DeferredLogger& deferred_logger) const = 0;
223
224 bool wellUnderZeroRateTarget(const Simulator& simulator,
225 const WellState<Scalar>& well_state,
226 DeferredLogger& deferred_logger) const;
227
228 bool wellUnderZeroGroupRateTarget(const Simulator& simulator,
229 const WellState<Scalar>& well_state,
230 DeferredLogger& deferred_logger,
231 std::optional<bool> group_control = std::nullopt) const;
232
233 bool stoppedOrZeroRateTarget(const Simulator& simulator,
234 const WellState<Scalar>& well_state,
235 DeferredLogger& deferred_logger) const;
236
237 bool updateWellStateWithTHPTargetProd(const Simulator& simulator,
238 WellState<Scalar>& well_state,
239 DeferredLogger& deferred_logger) const;
240
241 enum class IndividualOrGroup { Individual, Group, Both };
242 bool updateWellControl(const Simulator& simulator,
243 const IndividualOrGroup iog,
244 WellState<Scalar>& well_state,
245 const GroupState<Scalar>& group_state,
246 DeferredLogger& deferred_logger) /* const */;
247
248 bool updateWellControlAndStatusLocalIteration(const Simulator& simulator,
249 WellState<Scalar>& well_state,
250 const GroupState<Scalar>& group_state,
251 const Well::InjectionControls& inj_controls,
252 const Well::ProductionControls& prod_controls,
253 const Scalar WQTotal,
254 DeferredLogger& deferred_logger,
255 const bool fixed_control = false,
256 const bool fixed_status = false);
257
258 virtual void updatePrimaryVariables(const Simulator& simulator,
259 const WellState<Scalar>& well_state,
260 DeferredLogger& deferred_logger) = 0;
261
262 virtual void calculateExplicitQuantities(const Simulator& simulator,
263 const WellState<Scalar>& well_state,
264 DeferredLogger& deferred_logger) = 0; // should be const?
265
266 virtual void updateProductivityIndex(const Simulator& simulator,
267 const WellProdIndexCalculator<Scalar>& wellPICalc,
268 WellState<Scalar>& well_state,
269 DeferredLogger& deferred_logger) const = 0;
270
271 virtual Scalar connectionDensity(const int globalConnIdx,
272 const int openConnIdx) const = 0;
273
276 {
277 return false;
278 }
279
280 // Add well contributions to matrix
281 virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
282
283 virtual void addWellPressureEquations(PressureMatrix& mat,
284 const BVector& x,
285 const int pressureVarIndex,
286 const bool use_well_weights,
287 const WellState<Scalar>& well_state) const = 0;
288
289 void addCellRates(RateVector& rates, int cellIdx) const;
290
291 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
292
293 // TODO: theoretically, it should be a const function
294 // Simulator is not const is because that assembleWellEq is non-const Simulator
295 void wellTesting(const Simulator& simulator,
296 const double simulation_time,
297 /* const */ WellState<Scalar>& well_state,
298 const GroupState<Scalar>& group_state,
299 WellTestState& welltest_state,
300 DeferredLogger& deferred_logger);
301
302 void checkWellOperability(const Simulator& simulator,
303 const WellState<Scalar>& well_state,
304 DeferredLogger& deferred_logger);
305
306 bool gliftBeginTimeStepWellTestIterateWellEquations(const Simulator& ebos_simulator,
307 const double dt,
308 WellState<Scalar>& well_state,
309 const GroupState<Scalar>& group_state,
310 DeferredLogger& deferred_logger);
311
312 void gliftBeginTimeStepWellTestUpdateALQ(const Simulator& simulator,
313 WellState<Scalar>& well_state,
314 DeferredLogger& deferred_logger);
315
316 // check whether the well is operable under the current reservoir condition
317 // mostly related to BHP limit and THP limit
318 void updateWellOperability(const Simulator& simulator,
319 const WellState<Scalar>& well_state,
320 DeferredLogger& deferred_logger);
321
322 bool updateWellOperabilityFromWellEq(const Simulator& simulator,
323 const WellState<Scalar>& well_state,
324 DeferredLogger& deferred_logger);
325
326 // update perforation water throughput based on solved water rate
327 virtual void updateWaterThroughput(const double dt,
328 WellState<Scalar>& well_state) const = 0;
329
332 virtual std::vector<Scalar>
333 computeCurrentWellRates(const Simulator& simulator,
334 DeferredLogger& deferred_logger) const = 0;
335
339 void updateWellStateRates(const Simulator& simulator,
340 WellState<Scalar>& well_state,
341 DeferredLogger& deferred_logger) const;
342
343 void solveWellEquation(const Simulator& simulator,
344 WellState<Scalar>& well_state,
345 const GroupState<Scalar>& group_state,
346 DeferredLogger& deferred_logger);
347
348 const std::vector<RateVector>& connectionRates() const
349 {
350 return connectionRates_;
351 }
352
353 virtual std::vector<Scalar> getPrimaryVars() const
354 {
355 return {};
356 }
357
358 virtual int setPrimaryVars(typename std::vector<Scalar>::const_iterator)
359 {
360 return 0;
361 }
362
363 std::vector<Scalar> wellIndex(const int perf,
364 const IntensiveQuantities& intQuants,
365 const Scalar trans_mult,
366 const SingleWellState<Scalar>& ws) const;
367
368 void updateConnectionDFactor(const Simulator& simulator,
369 SingleWellState<Scalar>& ws) const;
370
371 void updateConnectionTransmissibilityFactor(const Simulator& simulator,
372 SingleWellState<Scalar>& ws) const;
373
374 virtual bool iterateWellEqWithSwitching(const Simulator& simulator,
375 const double dt,
376 const WellInjectionControls& inj_controls,
377 const WellProductionControls& prod_controls,
378 WellState<Scalar>& well_state,
379 const GroupState<Scalar>& group_state,
380 DeferredLogger& deferred_logger,
381 const bool fixed_control = false,
382 const bool fixed_status = false) = 0;
383protected:
384 // simulation parameters
385 std::vector<RateVector> connectionRates_;
386 std::vector<Scalar> B_avg_;
387 bool changed_to_stopped_this_step_ = false;
388 bool thp_update_iterations = false;
389
390 Scalar wpolymer() const;
391 Scalar wfoam() const;
392 Scalar wsalt() const;
393 Scalar wmicrobes() const;
394 Scalar woxygen() const;
395 Scalar wurea() const;
396
397 virtual Scalar getRefDensity() const = 0;
398
399 // Component fractions for each phase for the well
400 const std::vector<Scalar>& compFrac() const;
401
402 std::vector<Scalar>
403 initialWellRateFractions(const Simulator& ebosSimulator,
404 const WellState<Scalar>& well_state) const;
405
406 // check whether the well is operable under BHP limit with current reservoir condition
407 virtual void checkOperabilityUnderBHPLimit(const WellState<Scalar>& well_state,
408 const Simulator& simulator,
409 DeferredLogger& deferred_logger) = 0;
410
411 // check whether the well is operable under THP limit with current reservoir condition
412 virtual void checkOperabilityUnderTHPLimit(const Simulator& simulator,
413 const WellState<Scalar>& well_state,
414 DeferredLogger& deferred_logger) = 0;
415
416 virtual void updateIPR(const Simulator& simulator,
417 DeferredLogger& deferred_logger) const = 0;
418
419 virtual void assembleWellEqWithoutIteration(const Simulator& simulator,
420 const double dt,
421 const WellInjectionControls& inj_controls,
422 const WellProductionControls& prod_controls,
423 WellState<Scalar>& well_state,
424 const GroupState<Scalar>& group_state,
425 DeferredLogger& deferred_logger) = 0;
426
427 // iterate well equations with the specified control until converged
428 virtual bool iterateWellEqWithControl(const Simulator& simulator,
429 const double dt,
430 const WellInjectionControls& inj_controls,
431 const WellProductionControls& prod_controls,
432 WellState<Scalar>& well_state,
433 const GroupState<Scalar>& group_state,
434 DeferredLogger& deferred_logger) = 0;
435
436 virtual void updateIPRImplicit(const Simulator& simulator,
437 WellState<Scalar>& well_state,
438 DeferredLogger& deferred_logger) = 0;
439
440 bool iterateWellEquations(const Simulator& simulator,
441 const double dt,
442 WellState<Scalar>& well_state,
443 const GroupState<Scalar>& group_state,
444 DeferredLogger& deferred_logger);
445
446 bool solveWellWithTHPConstraint(const Simulator& simulator,
447 const double dt,
448 const Well::InjectionControls& inj_controls,
449 const Well::ProductionControls& prod_controls,
450 WellState<Scalar>& well_state,
451 const GroupState<Scalar>& group_state,
452 DeferredLogger& deferred_logger);
453
454 std::optional<Scalar>
455 estimateOperableBhp(const Simulator& ebos_simulator,
456 const double dt,
457 WellState<Scalar>& well_state,
458 const SummaryState& summary_state,
459 DeferredLogger& deferred_logger);
460
461 bool solveWellWithBhp(const Simulator& simulator,
462 const double dt,
463 const Scalar bhp,
464 WellState<Scalar>& well_state,
465 DeferredLogger& deferred_logger);
466
467 bool solveWellWithZeroRate(const Simulator& simulator,
468 const double dt,
469 WellState<Scalar>& well_state,
470 DeferredLogger& deferred_logger);
471
472 bool solveWellForTesting(const Simulator& simulator,
473 WellState<Scalar>& well_state,
474 const GroupState<Scalar>& group_state,
475 DeferredLogger& deferred_logger);
476
477 Eval getPerfCellPressure(const FluidState& fs) const;
478
479 // get the mobility for specific perforation
480 template<class Value, class Callback>
481 void getMobility(const Simulator& simulator,
482 const int perf,
483 std::vector<Value>& mob,
484 Callback& extendEval,
485 [[maybe_unused]] DeferredLogger& deferred_logger) const;
486
487 void computeConnLevelProdInd(const FluidState& fs,
488 const std::function<Scalar(const Scalar)>& connPICalc,
489 const std::vector<Scalar>& mobility,
490 Scalar* connPI) const;
491
492 void computeConnLevelInjInd(const FluidState& fs,
493 const Phase preferred_phase,
494 const std::function<Scalar(const Scalar)>& connIICalc,
495 const std::vector<Scalar>& mobility,
496 Scalar* connII,
497 DeferredLogger& deferred_logger) const;
498
499 Scalar computeConnectionDFactor(const int perf,
500 const IntensiveQuantities& intQuants,
501 const SingleWellState<Scalar>& ws) const;
502};
503
504} // namespace Opm
505
506#ifndef OPM_WELLINTERFACE_IMPL_HEADER_INCLUDED
507#include "WellInterface_impl.hpp"
508#endif
509
510#endif // OPM_WELLINTERFACE_HEADER_INCLUDED
Class for handling the blackoil well model.
Definition WellInterface.hpp:31
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 WellInterface.hpp:30
Definition WellInterfaceFluidSystem.hpp:43
Class encapsulating some information about parallel wells.
Definition WGState.hpp:31
Manages the initializing and running of time dependent problems.
Definition simulator.hh:92
Definition WellInterfaceFluidSystem.hpp:51
Definition WellInterfaceIndices.hpp:34
Definition WellInterface.hpp:71
virtual std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
Compute well rates based on current reservoir conditions and well variables.
WellInterface(const Well &well, const ParallelWellInfo< Scalar > &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< Scalar > > &perf_data)
Constructor.
Definition WellInterface_impl.hpp:56
virtual void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
using the solution x to recover the solution xw for wells and applying xw to update Well State
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
void updateWellStateRates(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Modify the well_state's rates if there is only one nonzero rate.
Definition WellInterface_impl.hpp:1572
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:275
Collect per-connection static information to enable calculating connection-level or well-level produc...
Definition WellProdIndexCalculator.hpp:37
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:62
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Static data associated with a well perforation.
Definition WellState.hpp:54
Definition BlackoilPhases.hpp:46