My Project
Loading...
Searching...
No Matches
WellGroupHelpers.hpp
1/*
2 Copyright 2019 Norce.
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20
21#ifndef OPM_WELLGROUPHELPERS_HEADER_INCLUDED
22#define OPM_WELLGROUPHELPERS_HEADER_INCLUDED
23
24#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
25#include <opm/input/eclipse/Schedule/Group/GSatProd.hpp>
26#include <opm/input/eclipse/EclipseState/Grid/FieldPropsManager.hpp>
27#include <opm/simulators/utils/ParallelCommunication.hpp>
28#include <opm/input/eclipse/Schedule/ScheduleState.hpp>
29
30#include <map>
31#include <string>
32#include <vector>
33
34namespace Opm {
35
36class DeferredLogger;
37class Group;
38template<class Scalar> class GroupState;
39namespace Network { class ExtNetwork; }
40struct PhaseUsage;
41class Schedule;
42template<class Scalar> class VFPProdProperties;
43template<class Scalar> class WellState;
44class FieldPropsManager;
45
46namespace Network { class ExtNetwork; }
47
48template<class Scalar>
50{
51public:
52
53 static Scalar sumWellPhaseRates(bool res_rates,
54 const Opm::Group& group,
55 const Opm::Schedule& schedule,
56 const Opm::WellState<Scalar>& wellState,
57 const int reportStepIdx,
58 const int phasePos,
59 const bool injector,
60 const bool network = false);
61
62 static Scalar satelliteProduction(const ScheduleState& sched,
63 const std::vector<std::string>& groups,
64 const GSatProd::GSatProdGroup::Rate rateComp);
65
66 static std::optional<GSatProd::GSatProdGroup::Rate>
67 selectRateComponent(const PhaseUsage& pu, const int phasePos);
68
69 static void setCmodeGroup(const Group& group,
70 const Schedule& schedule,
71 const SummaryState& summaryState,
72 const int reportStepIdx,
73 GroupState<Scalar>& group_state);
74
75 static void accumulateGroupEfficiencyFactor(const Group& group,
76 const Schedule& schedule,
77 const int reportStepIdx,
78 Scalar& factor);
79
80 static Scalar sumWellSurfaceRates(const Group& group,
81 const Schedule& schedule,
82 const WellState<Scalar>& wellState,
83 const int reportStepIdx,
84 const int phasePos,
85 const bool injector);
86
88 static std::pair<std::optional<std::string>, Scalar>
89 worstOffendingWell(const Group& group,
90 const Schedule& schedule,
91 const int reportStepIdx,
92 const Group::ProductionCMode& offendedControl,
93 const PhaseUsage& pu,
94 const Parallel::Communication& comm,
95 const WellState<Scalar>& wellState,
97
98 static Scalar sumWellResRates(const Group& group,
99 const Schedule& schedule,
100 const WellState<Scalar>& wellState,
101 const int reportStepIdx,
102 const int phasePos,
103 const bool injector);
104
105 static Scalar sumSolventRates(const Group& group,
106 const Schedule& schedule,
107 const WellState<Scalar>& wellState,
108 const int reportStepIdx,
109 const bool injector);
110
111 static void updateGroupTargetReduction(const Group& group,
112 const Schedule& schedule,
113 const int reportStepIdx,
114 const bool isInjector,
115 const PhaseUsage& pu,
116 const GuideRate& guide_rate,
117 const WellState<Scalar>& wellState,
118 const SummaryState& summaryState,
119 GroupState<Scalar>& group_state,
120 std::vector<Scalar>& groupTargetReduction);
121
122 static void updateGuideRates(const Group& group,
123 const Schedule& schedule,
125 const PhaseUsage& pu,
126 int report_step,
127 double sim_time,
128 WellState<Scalar>& well_state,
129 const GroupState<Scalar>& group_state,
130 const Parallel::Communication& comm,
132 std::vector<Scalar>& pot,
134
135 static void updateGuideRateForProductionGroups(const Group& group,
136 const Schedule& schedule,
137 const PhaseUsage& pu,
138 const int reportStepIdx,
139 const double& simTime,
140 WellState<Scalar>& wellState,
141 const GroupState<Scalar>& group_state,
142 const Parallel::Communication& comm,
143 GuideRate* guideRate,
144 std::vector<Scalar>& pot);
145
146 static void updateGuideRatesForWells(const Schedule& schedule,
147 const PhaseUsage& pu,
148 const int reportStepIdx,
149 const double& simTime,
150 const WellState<Scalar>& wellState,
151 const Parallel::Communication& comm,
152 GuideRate* guideRate);
153
154 static void updateGuideRatesForInjectionGroups(const Group& group,
155 const Schedule& schedule,
156 const SummaryState& summaryState,
157 const PhaseUsage& pu,
158 const int reportStepIdx,
159 const WellState<Scalar>& wellState,
160 const GroupState<Scalar>& group_state,
161 GuideRate* guideRate,
163
164 static void updateVREPForGroups(const Group& group,
165 const Schedule& schedule,
166 const int reportStepIdx,
167 const WellState<Scalar>& wellState,
168 GroupState<Scalar>& group_state);
169
170 template <class RegionalValues>
171 static void updateGpMaintTargetForGroups(const Group& group,
172 const Schedule& schedule,
174 const int reportStepIdx,
175 const double dt,
176 const WellState<Scalar>& well_state,
177 GroupState<Scalar>& group_state);
178
179 static void updateReservoirRatesInjectionGroups(const Group& group,
180 const Schedule& schedule,
181 const int reportStepIdx,
182 const WellState<Scalar>& wellState,
183 GroupState<Scalar>& group_state);
184
185 static void updateSurfaceRatesInjectionGroups(const Group& group,
186 const Schedule& schedule,
187 const int reportStepIdx,
188 const WellState<Scalar>& wellState,
189 GroupState<Scalar>& group_state);
190
191 static void updateWellRates(const Group& group,
192 const Schedule& schedule,
193 const int reportStepIdx,
195 WellState<Scalar>& wellState);
196
197 static void updateGroupProductionRates(const Group& group,
198 const Schedule& schedule,
199 const int reportStepIdx,
200 const WellState<Scalar>& wellState,
201 GroupState<Scalar>& group_state);
202
203 static void updateNetworkLeafNodeProductionRates(const Schedule& schedule,
204 const int reportStepIdx,
205 const WellState<Scalar>& wellState,
206 GroupState<Scalar>& group_state);
207
208
209 static void updateWellRatesFromGroupTargetScale(const Scalar scale,
210 const Group& group,
211 const Schedule& schedule,
212 const int reportStepIdx,
213 bool isInjector,
214 const GroupState<Scalar>& group_state,
215 WellState<Scalar>& wellState);
216
217 static void updateREINForGroups(const Group& group,
218 const Schedule& schedule,
219 const int reportStepIdx,
220 const PhaseUsage& pu,
221 const SummaryState& st,
222 const WellState<Scalar>& wellState,
223 GroupState<Scalar>& group_state,
224 bool sum_rank);
225
226
227 static std::map<std::string, Scalar>
228 computeNetworkPressures(const Network::ExtNetwork& network,
229 const WellState<Scalar>& well_state,
230 const GroupState<Scalar>& group_state,
232 const Schedule& schedule,
233 const int report_time_step);
234
235 static GuideRate::RateVector
236 getWellRateVector(const WellState<Scalar>& well_state,
237 const PhaseUsage& pu,
238 const std::string& name);
239
240 static GuideRate::RateVector
241 getProductionGroupRateVector(const GroupState<Scalar>& group_state,
242 const PhaseUsage& pu,
243 const std::string& group_name);
244
245 static Scalar getGuideRate(const std::string& name,
246 const Schedule& schedule,
247 const WellState<Scalar>& wellState,
248 const GroupState<Scalar>& group_state,
249 const int reportStepIdx,
250 const GuideRate* guideRate,
251 const GuideRateModel::Target target,
252 const PhaseUsage& pu);
253
254 static Scalar getGuideRateInj(const std::string& name,
255 const Schedule& schedule,
256 const WellState<Scalar>& wellState,
257 const GroupState<Scalar>& group_state,
258 const int reportStepIdx,
259 const GuideRate* guideRate,
260 const GuideRateModel::Target target,
261 const Phase& injectionPhase,
262 const PhaseUsage& pu);
263
264 static int groupControlledWells(const Schedule& schedule,
265 const WellState<Scalar>& well_state,
266 const GroupState<Scalar>& group_state,
268 const GuideRate* guideRate,
269 const int report_step,
270 const std::string& group_name,
271 const std::string& always_included_child,
272 const bool is_production_group,
273 const Phase injection_phase);
274
275 static std::pair<bool, Scalar>
276 checkGroupConstraintsInj(const std::string& name,
277 const std::string& parent,
278 const Group& group,
279 const WellState<Scalar>& wellState,
280 const GroupState<Scalar>& group_state,
281 const int reportStepIdx,
282 const GuideRate* guideRate,
283 const Scalar* rates,
284 Phase injectionPhase,
285 const PhaseUsage& pu,
286 const Scalar efficiencyFactor,
287 const Schedule& schedule,
288 const SummaryState& summaryState,
289 const std::vector<Scalar>& resv_coeff,
291
292 static std::vector<std::string>
293 groupChainTopBot(const std::string& bottom,
294 const std::string& top,
295 const Schedule& schedule,
296 const int report_step);
297
298 static std::string
299 control_group(const Group& group,
300 const GroupState<Scalar>& group_state,
301 const int reportStepIdx,
302 const Schedule& schedule);
303
304 static std::pair<bool, Scalar>
305 checkGroupConstraintsProd(const std::string& name,
306 const std::string& parent,
307 const Group& group,
308 const WellState<Scalar>& wellState,
309 const GroupState<Scalar>& group_state,
310 const int reportStepIdx,
311 const GuideRate* guideRate,
312 const Scalar* rates,
313 const PhaseUsage& pu,
314 const Scalar efficiencyFactor,
315 const Schedule& schedule,
316 const SummaryState& summaryState,
317 const std::vector<Scalar>& resv_coeff,
319
320 template <class AverageRegionalPressureType>
321 static void setRegionAveragePressureCalculator(const Group& group,
322 const Schedule& schedule,
323 const int reportStepIdx,
324 const FieldPropsManager& fp,
325 const PhaseUsage& pu,
326 std::map<std::string, std::unique_ptr<AverageRegionalPressureType>>& regionalAveragePressureCalculator);
327};
328
329} // namespace Opm
330
331#endif
Definition DeferredLogger.hpp:57
Definition GroupState.hpp:43
Class which linearly interpolates BHP as a function of rate, tubing head pressure,...
Definition VFPProdProperties.hpp:38
Definition WellGroupHelpers.hpp:50
static std::pair< std::optional< std::string >, Scalar > worstOffendingWell(const Group &group, const Schedule &schedule, const int reportStepIdx, const Group::ProductionCMode &offendedControl, const PhaseUsage &pu, const Parallel::Communication &comm, const WellState< Scalar > &wellState, DeferredLogger &deferred_logger)
Returns the name of the worst offending well and its fraction (i.e. violated_phase / preferred_phase)
Definition WellGroupHelpers.cpp:1706
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:66
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
Definition BlackoilPhases.hpp:46