My Project
WellInterfaceFluidSystem.hpp
1 /*
2  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2017 Statoil ASA.
4  Copyright 2017 IRIS
5  Copyright 2019 Norce
6 
7  This file is part of the Open Porous Media project (OPM).
8 
9  OPM is free software: you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation, either version 3 of the License, or
12  (at your option) any later version.
13 
14  OPM is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  along with OPM. If not, see <http://www.gnu.org/licenses/>.
21 */
22 
23 
24 #ifndef OPM_WELLINTERFACE_FLUID_SYSTEM_HEADER_INCLUDED
25 #define OPM_WELLINTERFACE_FLUID_SYSTEM_HEADER_INCLUDED
26 
27 #include <opm/simulators/wells/WellInterfaceGeneric.hpp>
28 #include <opm/core/props/BlackoilPhases.hpp>
29 
30 #include <limits>
31 
32 namespace Opm
33 {
34 namespace RateConverter
35 {
36  template <class FluidSystem, class Region> class SurfaceToReservoirVoidage;
37 }
38 
39 class Group;
40 class GroupState;
41 class Schedule;
42 class WellState;
43 class SingleWellState;
44 
45 template<class FluidSystem>
47 protected:
48  using RateConverterType = RateConverter::
49  SurfaceToReservoirVoidage<FluidSystem, std::vector<int>>;
50  // to indicate a invalid completion
51  static constexpr int INVALIDCOMPLETION = std::numeric_limits<int>::max();
52 
53 public:
54  void updateWellTestState(const SingleWellState& ws,
55  const double& simulationTime,
56  const bool& writeMessageToOPMLog,
57  WellTestState& wellTestState,
58  DeferredLogger& deferred_logger) const;
59 
60  int flowPhaseToEbosPhaseIdx(const int phaseIdx) const;
61 
62  static constexpr int Water = BlackoilPhases::Aqua;
63  static constexpr int Oil = BlackoilPhases::Liquid;
64  static constexpr int Gas = BlackoilPhases::Vapour;
65 
66  const RateConverterType& rateConverter() const
67  {
68  return rateConverter_;
69  }
70 
71 protected:
72  WellInterfaceFluidSystem(const Well& well,
73  const ParallelWellInfo& parallel_well_info,
74  const int time_step,
75  const RateConverterType& rate_converter,
76  const int pvtRegionIdx,
77  const int num_components,
78  const int num_phases,
79  const int index_of_well,
80  const std::vector<PerforationData>& perf_data);
81 
82  // updating the voidage rates in well_state when requested
83  void calculateReservoirRates(SingleWellState& ws) const;
84 
85  bool checkIndividualConstraints(SingleWellState& ws,
86  const SummaryState& summaryState,
87  DeferredLogger& deferred_logger) const;
88 
89  Well::InjectorCMode activeInjectionConstraint(const SingleWellState& ws,
90  const SummaryState& summaryState,
91  DeferredLogger& deferred_logger) const;
92 
93  Well::ProducerCMode activeProductionConstraint(const SingleWellState& ws,
94  const SummaryState& summaryState,
95  DeferredLogger& deferred_logger) const;
96 
97  std::pair<bool, double> checkGroupConstraintsInj(const Group& group,
98  const WellState& well_state,
99  const GroupState& group_state,
100  const double efficiencyFactor,
101  const Schedule& schedule,
102  const SummaryState& summaryState,
103  DeferredLogger& deferred_logger) const;
104 
105  std::pair<bool, double> checkGroupConstraintsProd(const Group& group,
106  const WellState& well_state,
107  const GroupState& group_state,
108  const double efficiencyFactor,
109  const Schedule& schedule,
110  const SummaryState& summaryState,
111  DeferredLogger& deferred_logger) const;
112 
113  bool checkGroupConstraints(WellState& well_state,
114  const GroupState& group_state,
115  const Schedule& schedule,
116  const SummaryState& summaryState,
117  DeferredLogger& deferred_logger) const;
118 
119  bool checkConstraints(WellState& well_state,
120  const GroupState& group_state,
121  const Schedule& schedule,
122  const SummaryState& summaryState,
123  DeferredLogger& deferred_logger) const;
124 
125  bool checkRateEconLimits(const WellEconProductionLimits& econ_production_limits,
126  const double* rates_or_potentials,
127  Opm::DeferredLogger& deferred_logger) const;
128 
130  bool ratio_limit_violated = false;
131  int worst_offending_completion = INVALIDCOMPLETION;
132  double violation_extent = 0.0;
133  };
134 
135  void checkMaxWaterCutLimit(const WellEconProductionLimits& econ_production_limits,
136  const SingleWellState& ws,
137  RatioLimitCheckReport& report) const;
138 
139  void checkMaxGORLimit(const WellEconProductionLimits& econ_production_limits,
140  const SingleWellState& ws,
141  RatioLimitCheckReport& report) const;
142 
143  void checkMaxWGRLimit(const WellEconProductionLimits& econ_production_limits,
144  const SingleWellState& ws,
145  RatioLimitCheckReport& report) const;
146 
147  void checkRatioEconLimits(const WellEconProductionLimits& econ_production_limits,
148  const SingleWellState& ws,
149  RatioLimitCheckReport& report,
150  DeferredLogger& deferred_logger) const;
151 
152  void updateWellTestStateEconomic(const SingleWellState& ws,
153  const double simulation_time,
154  const bool write_message_to_opmlog,
155  WellTestState& well_test_state,
156  DeferredLogger& deferred_logger) const;
157 
158  std::optional<double>
159  getGroupInjectionTargetRate(const Group& group,
160  const WellState& well_state,
161  const GroupState& group_state,
162  const Schedule& schedule,
163  const SummaryState& summaryState,
164  const InjectorType& injectorType,
165  double efficiencyFactor,
166  DeferredLogger& deferred_logger) const;
167 
168  double
169  getGroupProductionTargetRate(const Group& group,
170  const WellState& well_state,
171  const GroupState& group_state,
172  const Schedule& schedule,
173  const SummaryState& summaryState,
174  double efficiencyFactor) const;
175 
176  // For the conversion between the surface volume rate and reservoir voidage rate
177  const RateConverterType& rateConverter_;
178 
179 private:
180  template <typename RatioFunc>
181  void checkMaxRatioLimitCompletions(const SingleWellState& ws,
182  const double max_ratio_limit,
183  const RatioFunc& ratioFunc,
184  RatioLimitCheckReport& report) const;
185 
186  template<typename RatioFunc>
187  bool checkMaxRatioLimitWell(const SingleWellState& well_state,
188  const double max_ratio_limit,
189  const RatioFunc& ratioFunc) const;
190 };
191 
192 }
193 
194 #endif // OPM_WELLINTERFACE_FLUID_SYSTEM_HEADER_INCLUDED
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:243
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:68
Definition: SingleWellState.hpp:38
Definition: WellInterfaceFluidSystem.hpp:46
Definition: WellInterfaceGeneric.hpp:51
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: WellInterfaceFluidSystem.hpp:129