My Project
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/GPMaint.hpp>
26 #include <opm/input/eclipse/Schedule/Schedule.hpp>
27 #include <opm/simulators/wells/GroupState.hpp>
28 #include <opm/simulators/wells/WellState.hpp>
29 #include <opm/input/eclipse/Schedule/Group/Group.hpp>
30 
31 #include <map>
32 #include <string>
33 #include <vector>
34 
35 namespace Opm
36 {
37 
38 class DeferredLogger;
39 class Group;
40 class GroupState;
41 namespace Network { class ExtNetwork; }
42 struct PhaseUsage;
43 class Schedule;
44 class VFPProdProperties;
45 class WellState;
46 
47 template <typename>
48 class WellContainer;
49 
50 namespace Network { class ExtNetwork; }
51 
52 namespace WellGroupHelpers
53 {
54 
55 
56 
57  void setCmodeGroup(const Group& group,
58  const Schedule& schedule,
59  const SummaryState& summaryState,
60  const int reportStepIdx,
61  WellState& wellState,
62  GroupState& group_state);
63 
64  void accumulateGroupEfficiencyFactor(const Group& group,
65  const Schedule& schedule,
66  const int reportStepIdx,
67  double& factor);
68 
69  double sumWellSurfaceRates(const Group& group,
70  const Schedule& schedule,
71  const WellState& wellState,
72  const int reportStepIdx,
73  const int phasePos,
74  const bool injector);
75 
76  double sumWellResRates(const Group& group,
77  const Schedule& schedule,
78  const WellState& wellState,
79  const int reportStepIdx,
80  const int phasePos,
81  const bool injector);
82 
83  double sumSolventRates(const Group& group,
84  const Schedule& schedule,
85  const WellState& wellState,
86  const int reportStepIdx,
87  const bool injector);
88 
89  void updateGroupTargetReduction(const Group& group,
90  const Schedule& schedule,
91  const int reportStepIdx,
92  const bool isInjector,
93  const PhaseUsage& pu,
94  const GuideRate& guide_rate,
95  const WellState& wellState,
96  GroupState& group_state,
97  std::vector<double>& groupTargetReduction);
98 
99  template <class Comm>
100  void updateGuideRates(const Group& group,
101  const Schedule& schedule,
102  const SummaryState& summary_state,
103  const PhaseUsage& pu,
104  int report_step,
105  double sim_time,
106  WellState& well_state,
107  const GroupState& group_state,
108  const Comm& comm,
109  GuideRate* guide_rate,
110  std::vector<double>& pot,
111  Opm::DeferredLogger& deferred_logge);
112 
113  template <class Comm>
114  void updateGuideRateForProductionGroups(const Group& group,
115  const Schedule& schedule,
116  const PhaseUsage& pu,
117  const int reportStepIdx,
118  const double& simTime,
119  WellState& wellState,
120  const GroupState& group_state,
121  const Comm& comm,
122  GuideRate* guideRate,
123  std::vector<double>& pot);
124 
125  template <class Comm>
126  void updateGuideRatesForWells(const Schedule& schedule,
127  const PhaseUsage& pu,
128  const int reportStepIdx,
129  const double& simTime,
130  const WellState& wellState,
131  const Comm& comm,
132  GuideRate* guideRate);
133 
134  void updateGuideRatesForInjectionGroups(const Group& group,
135  const Schedule& schedule,
136  const SummaryState& summaryState,
137  const Opm::PhaseUsage& pu,
138  const int reportStepIdx,
139  const WellState& wellState,
140  const GroupState& group_state,
141  GuideRate* guideRate,
142  Opm::DeferredLogger& deferred_logger);
143 
144  void updateVREPForGroups(const Group& group,
145  const Schedule& schedule,
146  const int reportStepIdx,
147  const WellState& wellState,
148  GroupState& group_state);
149 
150  void updateReservoirRatesInjectionGroups(const Group& group,
151  const Schedule& schedule,
152  const int reportStepIdx,
153  const WellState& wellState,
154  GroupState& group_state);
155 
156  void updateSurfaceRatesInjectionGroups(const Group& group,
157  const Schedule& schedule,
158  const int reportStepIdx,
159  const WellState& wellState,
160  GroupState& group_state);
161 
162  void updateWellRates(const Group& group,
163  const Schedule& schedule,
164  const int reportStepIdx,
165  const WellState& wellStateNupcol,
166  WellState& wellState);
167 
168  void updateGroupProductionRates(const Group& group,
169  const Schedule& schedule,
170  const int reportStepIdx,
171  const WellState& wellState,
172  GroupState& group_state);
173 
174  void updateWellRatesFromGroupTargetScale(const double scale,
175  const Group& group,
176  const Schedule& schedule,
177  const int reportStepIdx,
178  bool isInjector,
179  const GroupState& group_state,
180  WellState& wellState);
181 
182  void updateREINForGroups(const Group& group,
183  const Schedule& schedule,
184  const int reportStepIdx,
185  const PhaseUsage& pu,
186  const SummaryState& st,
187  const WellState& wellState,
188  GroupState& group_state);
189 
190  template <class RegionalValues>
191  void updateGpMaintTargetForGroups(const Group& group,
192  const Schedule& schedule,
193  const RegionalValues& regional_values,
194  const int reportStepIdx,
195  const double dt,
196  const WellState& well_state,
197  GroupState& group_state)
198  {
199  for (const std::string& groupName : group.groups()) {
200  const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx);
201  updateGpMaintTargetForGroups(groupTmp, schedule, regional_values, reportStepIdx, dt, well_state, group_state);
202  }
203  const auto& gpm = group.gpmaint();
204  if(!gpm)
205  return;
206 
207  const auto [name, number] = *gpm->region();
208  const double error = gpm->pressure_target() - regional_values.pressure(number);
209  double current_rate = 0.0;
210  const auto& pu = well_state.phaseUsage();
211  double sign = 1.0;
212  switch (gpm->flow_target()) {
213  case Opm::GPMaint::FlowTarget::RESV_PROD:
214  {
215  current_rate = -group_state.injection_vrep_rate(group.name());
216  sign = -1.0;
217  break;
218  }
219  case Opm::GPMaint::FlowTarget::RESV_OINJ:
220  {
221  if (pu.phase_used[BlackoilPhases::Liquid])
222  current_rate = group_state.injection_reservoir_rates(group.name())[pu.phase_pos[BlackoilPhases::Liquid]];
223 
224  break;
225  }
226  case Opm::GPMaint::FlowTarget::RESV_WINJ:
227  {
228  if (pu.phase_used[BlackoilPhases::Aqua])
229  current_rate = group_state.injection_reservoir_rates(group.name())[pu.phase_pos[BlackoilPhases::Aqua]];
230 
231  break;
232  }
233  case Opm::GPMaint::FlowTarget::RESV_GINJ:
234  {
235  if (pu.phase_used[BlackoilPhases::Vapour])
236  current_rate = group_state.injection_reservoir_rates(group.name())[pu.phase_pos[BlackoilPhases::Vapour]];
237  break;
238  }
239  case Opm::GPMaint::FlowTarget::SURF_OINJ:
240  {
241  if (pu.phase_used[BlackoilPhases::Liquid])
242  current_rate = group_state.injection_surface_rates(group.name())[pu.phase_pos[BlackoilPhases::Liquid]];
243 
244  break;
245  }
246  case Opm::GPMaint::FlowTarget::SURF_WINJ:
247  {
248  if (pu.phase_used[BlackoilPhases::Aqua])
249  current_rate = group_state.injection_surface_rates(group.name())[pu.phase_pos[BlackoilPhases::Aqua]];
250 
251  break;
252  }
253  case Opm::GPMaint::FlowTarget::SURF_GINJ:
254  {
255  if (pu.phase_used[BlackoilPhases::Vapour])
256  current_rate = group_state.injection_surface_rates(group.name())[pu.phase_pos[BlackoilPhases::Vapour]];
257 
258  break;
259  }
260  default:
261  throw std::invalid_argument("Invalid Flow target type in GPMAINT");
262  }
263 
264  auto& gpmaint_state = group_state.gpmaint(group.name());
265  double rate = gpm->rate(gpmaint_state, current_rate, error, dt);
266  group_state.update_gpmaint_target(group.name(), std::max(0.0, sign * rate));
267  }
268  std::map<std::string, double>
269  computeNetworkPressures(const Opm::Network::ExtNetwork& network,
270  const WellState& well_state,
271  const GroupState& group_state,
272  const VFPProdProperties& vfp_prod_props,
273  const Schedule& schedule,
274  const int report_time_step);
275 
276  GuideRate::RateVector
277  getWellRateVector(const WellState& well_state, const PhaseUsage& pu, const std::string& name);
278 
279  GuideRate::RateVector
280  getProductionGroupRateVector(const GroupState& group_state, const PhaseUsage& pu, const std::string& group_name);
281 
282  double getGuideRate(const std::string& name,
283  const Schedule& schedule,
284  const WellState& wellState,
285  const GroupState& group_state,
286  const int reportStepIdx,
287  const GuideRate* guideRate,
288  const GuideRateModel::Target target,
289  const PhaseUsage& pu);
290 
291 
292  double getGuideRateInj(const std::string& name,
293  const Schedule& schedule,
294  const WellState& wellState,
295  const GroupState& group_state,
296  const int reportStepIdx,
297  const GuideRate* guideRate,
298  const GuideRateModel::Target target,
299  const Phase& injectionPhase,
300  const PhaseUsage& pu);
301 
302  int groupControlledWells(const Schedule& schedule,
303  const WellState& well_state,
304  const GroupState& group_state,
305  const int report_step,
306  const std::string& group_name,
307  const std::string& always_included_child,
308  const bool is_production_group,
309  const Phase injection_phase);
310 
311 
313  {
314  public:
315  FractionCalculator(const Schedule& schedule,
316  const WellState& well_state,
317  const GroupState& group_state,
318  const int report_step,
319  const GuideRate* guide_rate,
320  const GuideRateModel::Target target,
321  const PhaseUsage& pu,
322  const bool is_producer,
323  const Phase injection_phase);
324  double fraction(const std::string& name, const std::string& control_group_name, const bool always_include_this);
325  double localFraction(const std::string& name, const std::string& always_included_child);
326 
327  private:
328  std::string parent(const std::string& name);
329  double guideRateSum(const Group& group, const std::string& always_included_child);
330  double guideRate(const std::string& name, const std::string& always_included_child);
331  int groupControlledWells(const std::string& group_name, const std::string& always_included_child);
332  GuideRate::RateVector getGroupRateVector(const std::string& group_name);
333  const Schedule& schedule_;
334  const WellState& well_state_;
335  const GroupState& group_state_;
336  int report_step_;
337  const GuideRate* guide_rate_;
338  GuideRateModel::Target target_;
339  const PhaseUsage& pu_;
340  bool is_producer_;
341  Phase injection_phase_;
342  };
343 
344 
345  std::pair<bool, double> checkGroupConstraintsInj(const std::string& name,
346  const std::string& parent,
347  const Group& group,
348  const WellState& wellState,
349  const GroupState& group_state,
350  const int reportStepIdx,
351  const GuideRate* guideRate,
352  const double* rates,
353  Phase injectionPhase,
354  const PhaseUsage& pu,
355  const double efficiencyFactor,
356  const Schedule& schedule,
357  const SummaryState& summaryState,
358  const std::vector<double>& resv_coeff,
359  DeferredLogger& deferred_logger);
360 
361 
362 
363 
364 
365 
366  std::vector<std::string> groupChainTopBot(const std::string& bottom,
367  const std::string& top,
368  const Schedule& schedule,
369  const int report_step);
370 
371 
372 
373 
374  std::pair<bool, double> checkGroupConstraintsProd(const std::string& name,
375  const std::string& parent,
376  const Group& group,
377  const WellState& wellState,
378  const GroupState& group_state,
379  const int reportStepIdx,
380  const GuideRate* guideRate,
381  const double* rates,
382  const PhaseUsage& pu,
383  const double efficiencyFactor,
384  const Schedule& schedule,
385  const SummaryState& summaryState,
386  const std::vector<double>& resv_coeff,
387  DeferredLogger& deferred_logger);
388 
389 
390 } // namespace WellGroupHelpers
391 
392 } // namespace Opm
393 
394 #endif
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Definition: WellGroupHelpers.hpp:313
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Definition: BlackoilPhases.hpp:46