MultisegmentWell.hpp
Go to the documentation of this file.
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21
22#ifndef OPM_MULTISEGMENTWELL_HEADER_INCLUDED
23#define OPM_MULTISEGMENTWELL_HEADER_INCLUDED
24
26
29
30namespace Opm {
31
32 class DeferredLogger;
33
34 template<typename TypeTag>
35 class MultisegmentWell : public WellInterface<TypeTag>
36 , public MultisegmentWellEval<GetPropType<TypeTag, Properties::FluidSystem>,
37 GetPropType<TypeTag, Properties::Indices>>
38 {
39 public:
43
44 using typename Base::Simulator;
45 using typename Base::IntensiveQuantities;
46 using typename Base::FluidSystem;
47 using typename Base::ModelParameters;
48 using typename Base::MaterialLaw;
49 using typename Base::Indices;
50 using typename Base::RateConverterType;
51 using typename Base::SparseMatrixAdapter;
52 using typename Base::FluidState;
53
56 using Base::Water;
57 using Base::Oil;
58 using Base::Gas;
59
60 using typename Base::Scalar;
61
63 using typename Base::BVector;
64 using typename Base::Eval;
65
66 using typename MSWEval::Equations;
67 using typename MSWEval::EvalWell;
68 using typename MSWEval::BVectorWell;
69 using MSWEval::SPres;
70 using typename Base::PressureMatrix;
71 using FSInfo = std::tuple<Scalar,typename std::decay<decltype(std::declval<decltype(std::declval<const Simulator&>().model().intensiveQuantities(0, 0).fluidState())>().saltConcentration())>::type,int>;
72
73 MultisegmentWell(const Well& well,
74 const ParallelWellInfo<Scalar>& pw_info,
75 const int time_step,
76 const ModelParameters& param,
77 const RateConverterType& rate_converter,
78 const int pvtRegionIdx,
79 const int num_components,
80 const int num_phases,
81 const int index_of_well,
82 const std::vector<PerforationData<Scalar>>& perf_data);
83
84 void init(const PhaseUsage* phase_usage_arg,
85 const std::vector<Scalar>& depth_arg,
86 const Scalar gravity_arg,
87 const std::vector<Scalar>& B_avg,
88 const bool changed_to_open_this_step) override;
89
91 void updateWellStateWithTarget(const Simulator& simulator,
92 const GroupState<Scalar>& group_state,
93 WellState<Scalar>& well_state,
94 DeferredLogger& deferred_logger) const override;
95
98 const WellState<Scalar>& well_state,
99 const std::vector<Scalar>& B_avg,
100 DeferredLogger& deferred_logger,
101 const bool relax_tolerance) const override;
102
104 void apply(const BVector& x, BVector& Ax) const override;
106 void apply(BVector& r) const override;
107
111 const BVector& x,
112 WellState<Scalar>& well_state,
113 DeferredLogger& deferred_logger) override;
114
116 void computeWellPotentials(const Simulator& simulator,
117 const WellState<Scalar>& well_state,
118 std::vector<Scalar>& well_potentials,
119 DeferredLogger& deferred_logger) override;
120
121 void updatePrimaryVariables(const Simulator& simulator,
122 const WellState<Scalar>& well_state,
123 DeferredLogger& deferred_logger) override;
124
125 void solveEqAndUpdateWellState(const Simulator& simulator,
126 WellState<Scalar>& well_state,
127 DeferredLogger& deferred_logger) override; // const?
128
129 void calculateExplicitQuantities(const Simulator& simulator,
130 const WellState<Scalar>& well_state,
131 DeferredLogger& deferred_logger) override; // should be const?
132
133 void updateIPRImplicit(const Simulator& simulator,
134 WellState<Scalar>& well_state,
135 DeferredLogger& deferred_logger) override;
136
137 void updateProductivityIndex(const Simulator& simulator,
138 const WellProdIndexCalculator<Scalar>& wellPICalc,
139 WellState<Scalar>& well_state,
140 DeferredLogger& deferred_logger) const override;
141
142 Scalar connectionDensity(const int globalConnIdx,
143 const int openConnIdx) const override;
144
145 void addWellContributions(SparseMatrixAdapter& jacobian) const override;
146
148 const BVector& x,
149 const int pressureVarIndex,
150 const bool use_well_weights,
151 const WellState<Scalar>& well_state) const override;
152
153 std::vector<Scalar>
154 computeCurrentWellRates(const Simulator& simulator,
155 DeferredLogger& deferred_logger) const override;
156
157 std::optional<Scalar>
159 const SummaryState& summary_state,
160 const Scalar alq_value,
161 DeferredLogger& deferred_logger,
162 bool iterate_if_no_solution) const override;
163
164 std::vector<Scalar> getPrimaryVars() const override;
165
166 int setPrimaryVars(typename std::vector<Scalar>::const_iterator it) override;
167
168 protected:
169 // regularize msw equation
171
172 // the intial amount of fluids in each segment under surface condition
173 std::vector<std::vector<Scalar> > segment_fluid_initial_;
174
175 mutable int debug_cost_counter_ = 0;
176
177 // updating the well_state based on well solution dwells
178 void updateWellState(const Simulator& simulator,
179 const BVectorWell& dwells,
180 WellState<Scalar>& well_state,
181 DeferredLogger& deferred_logger,
182 const Scalar relaxation_factor = 1.0);
183
184 // computing the accumulation term for later use in well mass equations
185 void computeInitialSegmentFluids(const Simulator& simulator);
186
187 // compute the pressure difference between the perforation and cell center
188 void computePerfCellPressDiffs(const Simulator& simulator);
189
190 template<class Value>
191 void computePerfRate(const IntensiveQuantities& int_quants,
192 const std::vector<Value>& mob_perfcells,
193 const std::vector<Scalar>& Tw,
194 const int seg,
195 const int perf,
196 const Value& segment_pressure,
197 const bool& allow_cf,
198 std::vector<Value>& cq_s,
199 Value& perf_press,
200 PerforationRates<Scalar>& perf_rates,
201 DeferredLogger& deferred_logger) const;
202
203 template<class Value>
204 void computePerfRate(const Value& pressure_cell,
205 const Value& rs,
206 const Value& rv,
207 const std::vector<Value>& b_perfcells,
208 const std::vector<Value>& mob_perfcells,
209 const std::vector<Scalar>& Tw,
210 const int perf,
211 const Value& segment_pressure,
212 const Value& segment_density,
213 const bool& allow_cf,
214 const std::vector<Value>& cmix_s,
215 std::vector<Value>& cq_s,
216 Value& perf_press,
217 PerforationRates<Scalar>& perf_rates,
218 DeferredLogger& deferred_logger) const;
219
220 // compute the fluid properties, such as densities, viscosities, and so on, in the segments
221 // They will be treated implicitly, so they need to be of Evaluation type
222 void computeSegmentFluidProperties(const Simulator& simulator,
223 DeferredLogger& deferred_logger);
224
225 // get the mobility for specific perforation
226 template<class Value>
227 void getMobility(const Simulator& simulator,
228 const int local_perf_index,
229 std::vector<Value>& mob,
230 DeferredLogger& deferred_logger) const;
231
232 void computeWellRatesAtBhpLimit(const Simulator& simulator,
233 std::vector<Scalar>& well_flux,
234 DeferredLogger& deferred_logger) const;
235
236 void computeWellRatesWithBhp(const Simulator& simulator,
237 const Scalar& bhp,
238 std::vector<Scalar>& well_flux,
239 DeferredLogger& deferred_logger) const override;
240
241 void computeWellRatesWithBhpIterations(const Simulator& simulator,
242 const Scalar& bhp,
243 std::vector<Scalar>& well_flux,
244 DeferredLogger& deferred_logger) const override;
245
246 std::vector<Scalar>
248 const Simulator& simulator,
249 DeferredLogger& deferred_logger) const;
250
251 bool computeWellPotentialsImplicit(const Simulator& simulator,
252 const WellState<Scalar>& well_state,
253 std::vector<Scalar>& well_potentials,
254 DeferredLogger& deferred_logger) const;
255
256 Scalar getRefDensity() const override;
257
258 bool iterateWellEqWithControl(const Simulator& simulator,
259 const double dt,
260 const Well::InjectionControls& inj_controls,
261 const Well::ProductionControls& prod_controls,
262 WellState<Scalar>& well_state,
263 const GroupState<Scalar>& group_state,
264 DeferredLogger& deferred_logger) override;
265
266 bool iterateWellEqWithSwitching(const Simulator& simulator,
267 const double dt,
268 const Well::InjectionControls& inj_controls,
269 const Well::ProductionControls& prod_controls,
270 WellState<Scalar>& well_state,
271 const GroupState<Scalar>& group_state,
272 DeferredLogger& deferred_logger,
273 const bool fixed_control = false,
274 const bool fixed_status = false) override;
275
276 void assembleWellEqWithoutIteration(const Simulator& simulator,
277 const double dt,
278 const Well::InjectionControls& inj_controls,
279 const Well::ProductionControls& prod_controls,
280 WellState<Scalar>& well_state,
281 const GroupState<Scalar>& group_state,
282 DeferredLogger& deferred_logger) override;
283
284 void updateWaterThroughput(const double dt, WellState<Scalar>& well_state) const override;
285
286 EvalWell getSegmentSurfaceVolume(const Simulator& simulator, const int seg_idx) const;
287
288 // turn on crossflow to avoid singular well equations
289 // when the well is banned from cross-flow and the BHP is not properly initialized,
290 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
291 // well rates, it can cause problem for THP calculation
292 // TODO: looking for better alternative to avoid wrong-signed well rates
293 bool openCrossFlowAvoidSingularity(const Simulator& simulator) const;
294
295 // for a well, when all drawdown are in the wrong direction, then this well will not
296 // be able to produce/inject .
297 bool allDrawDownWrongDirection(const Simulator& simulator) const;
298
299 std::optional<Scalar>
301 const Simulator& ebos_simulator,
302 const SummaryState& summary_state,
303 DeferredLogger& deferred_logger) const;
304
305 std::optional<Scalar>
306 computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
307 const SummaryState& summary_state,
308 DeferredLogger& deferred_logger) const;
309
310 Scalar maxPerfPress(const Simulator& simulator) const;
311
312 // check whether the well is operable under BHP limit with current reservoir condition
314 const Simulator& ebos_simulator,
315 DeferredLogger& deferred_logger) override;
316
317 // check whether the well is operable under THP limit with current reservoir condition
318 void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator,
319 const WellState<Scalar>& well_state,
320 DeferredLogger& deferred_logger) override;
321
322 // updating the inflow based on the current reservoir condition
323 void updateIPR(const Simulator& ebos_simulator,
324 DeferredLogger& deferred_logger) const override;
325
326 FSInfo getFirstPerforationFluidStateInfo(const Simulator& simulator) const;
327 };
328
329}
330
332
333#endif // OPM_MULTISEGMENTWELL_HEADER_INCLUDED
Definition: ConvergenceReport.hpp:38
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:43
Definition: MultisegmentWellEval.hpp:47
MultisegmentWellEquations< Scalar, numWellEq, Indices::numEq > Equations
Definition: MultisegmentWellEval.hpp:55
Definition: MultisegmentWell.hpp:38
std::vector< Scalar > computeWellPotentialWithTHP(const WellState< Scalar > &well_state, const Simulator &simulator, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:471
void computeWellPotentials(const Simulator &simulator, const WellState< Scalar > &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: MultisegmentWell_impl.hpp:288
void calculateExplicitQuantities(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:736
std::optional< Scalar > computeBhpAtThpLimitProdWithAlq(const Simulator &simulator, const SummaryState &summary_state, const Scalar alq_value, DeferredLogger &deferred_logger, bool iterate_if_no_solution) const override
Definition: MultisegmentWell_impl.hpp:2058
Scalar connectionDensity(const int globalConnIdx, const int openConnIdx) const override
Definition: MultisegmentWell_impl.hpp:834
void addWellContributions(SparseMatrixAdapter &jacobian) const override
Definition: MultisegmentWell_impl.hpp:856
std::optional< Scalar > computeBhpAtThpLimitProd(const WellState< Scalar > &well_state, const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:2040
Scalar getRefDensity() const override
Definition: MultisegmentWell_impl.hpp:1183
bool iterateWellEqWithSwitching(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false) override
Definition: MultisegmentWell_impl.hpp:1601
std::vector< Scalar > getPrimaryVars() const override
Definition: MultisegmentWell_impl.hpp:2240
void solveEqAndUpdateWellState(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:592
ConvergenceReport getWellConvergence(const Simulator &simulator, const WellState< Scalar > &well_state, const std::vector< Scalar > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance) const override
check whether the well equations get converged for this well
Definition: MultisegmentWell_impl.hpp:195
void updatePrimaryVariables(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:158
void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:752
void computeWellRatesWithBhpIterations(const Simulator &simulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:396
std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:2196
void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:221
void updateWellState(const Simulator &simulator, const BVectorWell &dwells, WellState< Scalar > &well_state, DeferredLogger &deferred_logger, const Scalar relaxation_factor=1.0)
Definition: MultisegmentWell_impl.hpp:697
void getMobility(const Simulator &simulator, const int local_perf_index, std::vector< Value > &mob, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:1143
void computePerfRate(const IntensiveQuantities &int_quants, const std::vector< Value > &mob_perfcells, const std::vector< Scalar > &Tw, const int seg, const int perf, const Value &segment_pressure, const bool &allow_cf, std::vector< Value > &cq_s, Value &perf_press, PerforationRates< Scalar > &perf_rates, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:1038
bool iterateWellEqWithControl(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1477
void checkOperabilityUnderBHPLimit(const WellState< Scalar > &well_state, const Simulator &ebos_simulator, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1191
bool openCrossFlowAvoidSingularity(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:1947
void computeSegmentFluidProperties(const Simulator &simulator, DeferredLogger &deferred_logger)
Definition: MultisegmentWell_impl.hpp:1111
void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:256
int setPrimaryVars(typename std::vector< Scalar >::const_iterator it) override
Definition: MultisegmentWell_impl.hpp:2258
void computeWellRatesWithBhp(const Simulator &simulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:349
EvalWell getSegmentSurfaceVolume(const Simulator &simulator, const int seg_idx) const
Definition: MultisegmentWell_impl.hpp:2020
bool allDrawDownWrongDirection(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:1956
bool computeWellPotentialsImplicit(const Simulator &simulator, const WellState< Scalar > &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:522
void init(const PhaseUsage *phase_usage_arg, const std::vector< Scalar > &depth_arg, const Scalar gravity_arg, const std::vector< Scalar > &B_avg, const bool changed_to_open_this_step) override
Definition: MultisegmentWell_impl.hpp:122
void updateWellStateWithTarget(const Simulator &simulator, const GroupState< Scalar > &group_state, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const override
updating the well state based the current control mode
Definition: MultisegmentWell_impl.hpp:174
int debug_cost_counter_
Definition: MultisegmentWell.hpp:175
void addWellPressureEquations(PressureMatrix &mat, const BVector &x, const int pressureVarIndex, const bool use_well_weights, const WellState< Scalar > &well_state) const override
Definition: MultisegmentWell_impl.hpp:869
void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1766
void checkOperabilityUnderTHPLimit(const Simulator &ebos_simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1425
void computePerfCellPressDiffs(const Simulator &simulator)
Definition: MultisegmentWell_impl.hpp:621
void updateIPR(const Simulator &ebos_simulator, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:1257
void updateWaterThroughput(const double dt, WellState< Scalar > &well_state) const override
Definition: MultisegmentWell_impl.hpp:2009
std::tuple< Scalar, typename std::decay< decltype(std::declval< decltype(std::declval< const Simulator & >().model().intensiveQuantities(0, 0).fluidState())>().saltConcentration())>::type, int > FSInfo
Definition: MultisegmentWell.hpp:71
FSInfo getFirstPerforationFluidStateInfo(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:2273
void computeWellRatesAtBhpLimit(const Simulator &simulator, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:333
bool regularize_
Definition: MultisegmentWell.hpp:170
MultisegmentWell(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)
Definition: MultisegmentWell_impl.hpp:62
Scalar maxPerfPress(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:2168
std::vector< std::vector< Scalar > > segment_fluid_initial_
Definition: MultisegmentWell.hpp:173
void updateIPRImplicit(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1356
std::optional< Scalar > computeBhpAtThpLimitInj(const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:2114
void computeInitialSegmentFluids(const Simulator &simulator)
Definition: MultisegmentWell_impl.hpp:679
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:195
static constexpr int Oil
Definition: WellInterfaceFluidSystem.hpp:63
static constexpr int Water
Definition: WellInterfaceFluidSystem.hpp:62
static constexpr int Gas
Definition: WellInterfaceFluidSystem.hpp:64
int pvtRegionIdx() const
Definition: WellInterfaceGeneric.hpp:124
Definition: WellInterface.hpp:77
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:82
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:100
BlackOilFluidState< Eval, FluidSystem, has_temperature, has_energy, Indices::compositionSwitchIdx >=0, has_watVapor, has_brine, has_saltPrecip, has_disgas_in_water, Indices::numPhases > FluidState
Definition: WellInterface.hpp:132
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: WellInterface.hpp:97
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:83
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: WellInterface.hpp:87
Dune::BlockVector< VectorBlockType > BVector
Definition: WellInterface.hpp:96
static constexpr bool has_polymer
Definition: WellInterface.hpp:110
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:106
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:84
static constexpr bool has_solvent
Definition: WellInterface.hpp:108
typename Base::Eval Eval
Definition: WellInterface.hpp:95
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: WellInterface.hpp:88
Definition: WellProdIndexCalculator.hpp:37
Definition: WellState.hpp:65
Defines the common properties required by the porous medium multi-phase models.
Definition: blackoilboundaryratevector.hh:39
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:233
Static data associated with a well perforation.
Definition: PerforationData.hpp:30
Definition: PerforationData.hpp:41
Definition: BlackoilPhases.hpp:46