My Project
Loading...
Searching...
No Matches
MultisegmentWell.hpp
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
25#include <opm/simulators/wells/WellInterface.hpp>
26#include <opm/simulators/wells/MultisegmentWellEval.hpp>
27
28namespace Opm {
29
30 class DeferredLogger;
31
32 template<typename TypeTag>
33 class MultisegmentWell : public WellInterface<TypeTag>
34 , public MultisegmentWellEval<GetPropType<TypeTag, Properties::FluidSystem>,
35 GetPropType<TypeTag, Properties::Indices>>
36 {
37 public:
41
42 using typename Base::Simulator;
43 using typename Base::IntensiveQuantities;
44 using typename Base::FluidSystem;
45 using typename Base::ModelParameters;
46 using typename Base::MaterialLaw;
47 using typename Base::Indices;
48 using typename Base::RateConverterType;
49 using typename Base::SparseMatrixAdapter;
50 using typename Base::FluidState;
51
52 using Base::has_solvent;
53 using Base::has_polymer;
54 using Base::Water;
55 using Base::Oil;
56 using Base::Gas;
57
58 using typename Base::Scalar;
59
61 using typename Base::BVector;
62 using typename Base::Eval;
63
65 using typename MSWEval::EvalWell;
66 using typename MSWEval::BVectorWell;
67 using MSWEval::SPres;
68 using typename Base::PressureMatrix;
69
70 MultisegmentWell(const Well& well,
72 const int time_step,
73 const ModelParameters& param,
74 const RateConverterType& rate_converter,
75 const int pvtRegionIdx,
76 const int num_components,
77 const int num_phases,
78 const int index_of_well,
79 const std::vector<PerforationData<Scalar>>& perf_data);
80
81 void init(const PhaseUsage* phase_usage_arg,
82 const std::vector<Scalar>& depth_arg,
83 const Scalar gravity_arg,
84 const std::vector<Scalar>& B_avg,
85 const bool changed_to_open_this_step) override;
86
87 void initPrimaryVariablesEvaluation() override;
88
90 void updateWellStateWithTarget(const Simulator& simulator,
91 const GroupState<Scalar>& group_state,
92 WellState<Scalar>& well_state,
93 DeferredLogger& deferred_logger) const override;
94
96 ConvergenceReport getWellConvergence(const Simulator& simulator,
97 const WellState<Scalar>& well_state,
98 const std::vector<Scalar>& B_avg,
100 const bool relax_tolerance) const override;
101
103 void apply(const BVector& x, BVector& Ax) const override;
105 void apply(BVector& r) const override;
106
109 void recoverWellSolutionAndUpdateWellState(const Simulator& simulator,
110 const BVector& x,
111 WellState<Scalar>& well_state,
113
115 void computeWellPotentials(const Simulator& simulator,
116 const WellState<Scalar>& well_state,
117 std::vector<Scalar>& well_potentials,
119
120 void updatePrimaryVariables(const Simulator& simulator,
121 const WellState<Scalar>& well_state,
123
124 void solveEqAndUpdateWellState(const Simulator& simulator,
125 WellState<Scalar>& well_state,
126 DeferredLogger& deferred_logger) override; // const?
127
128 void calculateExplicitQuantities(const Simulator& simulator,
129 const WellState<Scalar>& well_state,
130 DeferredLogger& deferred_logger) override; // should be const?
131
132 void updateIPRImplicit(const Simulator& simulator,
133 WellState<Scalar>& well_state,
135
136 void updateProductivityIndex(const Simulator& simulator,
138 WellState<Scalar>& well_state,
139 DeferredLogger& deferred_logger) const override;
140
141 Scalar connectionDensity(const int globalConnIdx,
142 const int openConnIdx) const override;
143
144 void addWellContributions(SparseMatrixAdapter& jacobian) const override;
145
146 void addWellPressureEquations(PressureMatrix& mat,
147 const BVector& x,
148 const int pressureVarIndex,
149 const bool use_well_weights,
150 const WellState<Scalar>& well_state) const override;
151
152 std::vector<Scalar>
153 computeCurrentWellRates(const Simulator& simulator,
154 DeferredLogger& deferred_logger) const override;
155
156 std::optional<Scalar>
157 computeBhpAtThpLimitProdWithAlq(const Simulator& simulator,
159 const Scalar alq_value,
160 DeferredLogger& deferred_logger) const override;
161
162 std::vector<Scalar> getPrimaryVars() const override;
163
164 int setPrimaryVars(typename std::vector<Scalar>::const_iterator it) override;
165
166 protected:
167 // regularize msw equation
168 bool regularize_;
169
170 // the intial amount of fluids in each segment under surface condition
171 std::vector<std::vector<Scalar> > segment_fluid_initial_;
172
173 mutable int debug_cost_counter_ = 0;
174
175 // updating the well_state based on well solution dwells
176 void updateWellState(const Simulator& simulator,
177 const BVectorWell& dwells,
178 WellState<Scalar>& well_state,
180 const Scalar relaxation_factor = 1.0);
181
182 // computing the accumulation term for later use in well mass equations
183 void computeInitialSegmentFluids(const Simulator& simulator);
184
185 // compute the pressure difference between the perforation and cell center
186 void computePerfCellPressDiffs(const Simulator& simulator);
187
188 template<class Value>
189 void computePerfRate(const IntensiveQuantities& int_quants,
190 const std::vector<Value>& mob_perfcells,
191 const std::vector<Scalar>& Tw,
192 const int seg,
193 const int perf,
194 const Value& segment_pressure,
195 const bool& allow_cf,
196 std::vector<Value>& cq_s,
197 Value& perf_press,
200
201 template<class Value>
202 void computePerfRate(const Value& pressure_cell,
203 const Value& rs,
204 const Value& rv,
205 const std::vector<Value>& b_perfcells,
206 const std::vector<Value>& mob_perfcells,
207 const std::vector<Scalar>& Tw,
208 const int perf,
209 const Value& segment_pressure,
210 const Value& segment_density,
211 const bool& allow_cf,
212 const std::vector<Value>& cmix_s,
213 std::vector<Value>& cq_s,
214 Value& perf_press,
217
218 // compute the fluid properties, such as densities, viscosities, and so on, in the segments
219 // They will be treated implicitly, so they need to be of Evaluation type
220 void computeSegmentFluidProperties(const Simulator& simulator,
222
223 // get the mobility for specific perforation
224 template<class Value>
225 void getMobility(const Simulator& simulator,
226 const int perf,
227 std::vector<Value>& mob,
229
230 void computeWellRatesAtBhpLimit(const Simulator& simulator,
231 std::vector<Scalar>& well_flux,
233
234 void computeWellRatesWithBhp(const Simulator& simulator,
235 const Scalar& bhp,
236 std::vector<Scalar>& well_flux,
237 DeferredLogger& deferred_logger) const override;
238
239 void computeWellRatesWithBhpIterations(const Simulator& simulator,
240 const Scalar& bhp,
241 std::vector<Scalar>& well_flux,
242 DeferredLogger& deferred_logger) const override;
243
244 std::vector<Scalar>
245 computeWellPotentialWithTHP(const WellState<Scalar>& well_state,
246 const Simulator& simulator,
248
249 bool computeWellPotentialsImplicit(const Simulator& simulator,
250 std::vector<Scalar>& well_potentials,
252
253 Scalar getRefDensity() const override;
254
255 bool iterateWellEqWithControl(const Simulator& simulator,
256 const double dt,
257 const Well::InjectionControls& inj_controls,
258 const Well::ProductionControls& prod_controls,
259 WellState<Scalar>& well_state,
260 const GroupState<Scalar>& group_state,
262
263 bool iterateWellEqWithSwitching(const Simulator& simulator,
264 const double dt,
265 const Well::InjectionControls& inj_controls,
266 const Well::ProductionControls& prod_controls,
267 WellState<Scalar>& well_state,
268 const GroupState<Scalar>& group_state,
270 const bool fixed_control = false,
271 const bool fixed_status = false) override;
272
273 void assembleWellEqWithoutIteration(const Simulator& simulator,
274 const double dt,
275 const Well::InjectionControls& inj_controls,
276 const Well::ProductionControls& prod_controls,
277 WellState<Scalar>& well_state,
278 const GroupState<Scalar>& group_state,
280
281 void updateWaterThroughput(const double dt, WellState<Scalar>& well_state) const override;
282
283 EvalWell getSegmentSurfaceVolume(const Simulator& simulator, const int seg_idx) const;
284
285 // turn on crossflow to avoid singular well equations
286 // when the well is banned from cross-flow and the BHP is not properly initialized,
287 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
288 // well rates, it can cause problem for THP calculation
289 // TODO: looking for better alternative to avoid wrong-signed well rates
290 bool openCrossFlowAvoidSingularity(const Simulator& simulator) const;
291
292 // for a well, when all drawdown are in the wrong direction, then this well will not
293 // be able to produce/inject .
294 bool allDrawDownWrongDirection(const Simulator& simulator) const;
295
296 std::optional<Scalar>
297 computeBhpAtThpLimitProd(const WellState<Scalar>& well_state,
298 const Simulator& ebos_simulator,
301
302 std::optional<Scalar>
303 computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
306
307 Scalar maxPerfPress(const Simulator& simulator) const;
308
309 // check whether the well is operable under BHP limit with current reservoir condition
310 void checkOperabilityUnderBHPLimit(const WellState<Scalar>& well_state,
311 const Simulator& ebos_simulator,
313
314 // check whether the well is operable under THP limit with current reservoir condition
315 void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator,
316 const WellState<Scalar>& well_state,
318
319 // updating the inflow based on the current reservoir condition
320 void updateIPR(const Simulator& ebos_simulator,
321 DeferredLogger& deferred_logger) const override;
322 };
323
324}
325
326#ifndef OPM_MULTISEGMENTWELL_IMPL_HEADER_INCLUDED
327#include "MultisegmentWell_impl.hpp"
328#endif
329
330#endif // OPM_MULTISEGMENTWELL_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 GroupState.hpp:38
Definition MultisegmentWellEval.hpp:46
Definition MultisegmentWell.hpp:36
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:283
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:201
std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const override
Compute well rates based on current reservoir conditions and well variables.
Definition MultisegmentWell_impl.hpp:2181
void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition MultisegmentWell_impl.hpp:227
void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
using the solution x to recover the solution xw for wells and applying xw to update Well State
Definition MultisegmentWell_impl.hpp:262
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:180
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:186
Definition WellInterface.hpp:71
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
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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
Static data associated with a well perforation.
Definition PerforationData.hpp:30
Definition PerforationData.hpp:41
Definition BlackoilPhases.hpp:46