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
26
27#include <opm/simulators/wells/WellInterface.hpp>
28#include <opm/simulators/wells/MultisegmentWellEval.hpp>
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
54 using Base::has_solvent;
55 using Base::has_polymer;
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
67 using typename MSWEval::EvalWell;
68 using typename MSWEval::BVectorWell;
69 using MSWEval::SPres;
70 using typename Base::PressureMatrix;
71
72 MultisegmentWell(const Well& well,
74 const int time_step,
75 const ModelParameters& param,
76 const RateConverterType& rate_converter,
77 const int pvtRegionIdx,
78 const int num_components,
79 const int num_phases,
80 const int index_of_well,
81 const std::vector<PerforationData<Scalar>>& perf_data);
82
83 void init(const PhaseUsage* phase_usage_arg,
84 const std::vector<Scalar>& depth_arg,
85 const Scalar gravity_arg,
86 const std::vector<Scalar>& B_avg,
87 const bool changed_to_open_this_step) 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,
161 bool iterate_if_no_solution) const override;
162
163 std::vector<Scalar> getPrimaryVars() const override;
164
165 int setPrimaryVars(typename std::vector<Scalar>::const_iterator it) override;
166
167 protected:
168 // regularize msw equation
169 bool regularize_;
170
171 // the intial amount of fluids in each segment under surface condition
172 std::vector<std::vector<Scalar> > segment_fluid_initial_;
173
174 mutable int debug_cost_counter_ = 0;
175
176 // updating the well_state based on well solution dwells
177 void updateWellState(const Simulator& simulator,
178 const BVectorWell& dwells,
179 WellState<Scalar>& well_state,
181 const Scalar relaxation_factor = 1.0);
182
183 // computing the accumulation term for later use in well mass equations
184 void computeInitialSegmentFluids(const Simulator& simulator);
185
186 // compute the pressure difference between the perforation and cell center
187 void computePerfCellPressDiffs(const Simulator& simulator);
188
189 template<class Value>
190 void computePerfRate(const IntensiveQuantities& int_quants,
191 const std::vector<Value>& mob_perfcells,
192 const std::vector<Scalar>& Tw,
193 const int seg,
194 const int perf,
195 const Value& segment_pressure,
196 const bool& allow_cf,
197 std::vector<Value>& cq_s,
198 Value& perf_press,
201
202 template<class Value>
203 void computePerfRate(const Value& pressure_cell,
204 const Value& rs,
205 const Value& rv,
206 const std::vector<Value>& b_perfcells,
207 const std::vector<Value>& mob_perfcells,
208 const std::vector<Scalar>& Tw,
209 const int perf,
210 const Value& segment_pressure,
211 const Value& segment_density,
212 const bool& allow_cf,
213 const std::vector<Value>& cmix_s,
214 std::vector<Value>& cq_s,
215 Value& perf_press,
218
219 // compute the fluid properties, such as densities, viscosities, and so on, in the segments
220 // They will be treated implicitly, so they need to be of Evaluation type
221 void computeSegmentFluidProperties(const Simulator& simulator,
223
224 // get the mobility for specific perforation
225 template<class Value>
226 void getMobility(const Simulator& simulator,
227 const int local_perf_index,
228 std::vector<Value>& mob,
230
231 void computeWellRatesAtBhpLimit(const Simulator& simulator,
232 std::vector<Scalar>& well_flux,
234
235 void computeWellRatesWithBhp(const Simulator& simulator,
236 const Scalar& bhp,
237 std::vector<Scalar>& well_flux,
238 DeferredLogger& deferred_logger) const override;
239
240 void computeWellRatesWithBhpIterations(const Simulator& simulator,
241 const Scalar& bhp,
242 std::vector<Scalar>& well_flux,
243 DeferredLogger& deferred_logger) const override;
244
245 std::vector<Scalar>
246 computeWellPotentialWithTHP(const WellState<Scalar>& well_state,
247 const Simulator& simulator,
249
250 bool computeWellPotentialsImplicit(const Simulator& simulator,
251 const WellState<Scalar>& well_state,
252 std::vector<Scalar>& well_potentials,
254
255 Scalar getRefDensity() const override;
256
257 bool iterateWellEqWithControl(const Simulator& simulator,
258 const double dt,
259 const Well::InjectionControls& inj_controls,
260 const Well::ProductionControls& prod_controls,
261 WellState<Scalar>& well_state,
262 const GroupState<Scalar>& group_state,
264
265 bool iterateWellEqWithSwitching(const Simulator& simulator,
266 const double dt,
267 const Well::InjectionControls& inj_controls,
268 const Well::ProductionControls& prod_controls,
269 WellState<Scalar>& well_state,
270 const GroupState<Scalar>& group_state,
272 const bool fixed_control = false,
273 const bool fixed_status = false) override;
274
275 void assembleWellEqWithoutIteration(const Simulator& simulator,
276 const double dt,
277 const Well::InjectionControls& inj_controls,
278 const Well::ProductionControls& prod_controls,
279 WellState<Scalar>& well_state,
280 const GroupState<Scalar>& group_state,
282
283 void updateWaterThroughput(const double dt, WellState<Scalar>& well_state) const override;
284
285 EvalWell getSegmentSurfaceVolume(const Simulator& simulator, const int seg_idx) const;
286
287 // turn on crossflow to avoid singular well equations
288 // when the well is banned from cross-flow and the BHP is not properly initialized,
289 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
290 // well rates, it can cause problem for THP calculation
291 // TODO: looking for better alternative to avoid wrong-signed well rates
292 bool openCrossFlowAvoidSingularity(const Simulator& simulator) const;
293
294 // for a well, when all drawdown are in the wrong direction, then this well will not
295 // be able to produce/inject .
296 bool allDrawDownWrongDirection(const Simulator& simulator) const;
297
298 std::optional<Scalar>
299 computeBhpAtThpLimitProd(const WellState<Scalar>& well_state,
300 const Simulator& ebos_simulator,
303
304 std::optional<Scalar>
305 computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
308
309 Scalar maxPerfPress(const Simulator& simulator) const;
310
311 // check whether the well is operable under BHP limit with current reservoir condition
312 void checkOperabilityUnderBHPLimit(const WellState<Scalar>& well_state,
313 const Simulator& ebos_simulator,
315
316 // check whether the well is operable under THP limit with current reservoir condition
317 void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator,
318 const WellState<Scalar>& well_state,
320
321 // updating the inflow based on the current reservoir condition
322 void updateIPR(const Simulator& ebos_simulator,
323 DeferredLogger& deferred_logger) const override;
324 };
325
326}
327
328#include "MultisegmentWell_impl.hpp"
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:43
Definition MultisegmentWellEval.hpp:47
Definition MultisegmentWell.hpp:38
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:276
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:194
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:2203
void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition MultisegmentWell_impl.hpp:220
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:255
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:173
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:195
Definition WellInterface.hpp:77
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:66
Defines the common properties required by the porous medium multi-phase models.
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