My Project
GasLiftGroupInfo.hpp
1 /*
2  Copyright 2021 Equinor ASA.
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 #ifndef OPM_GASLIFT_GROUP_INFO_HEADER_INCLUDED
21 #define OPM_GASLIFT_GROUP_INFO_HEADER_INCLUDED
22 
23 #include <dune/common/version.hh>
24 #include <dune/common/parallel/mpihelper.hh>
25 
26 #include <opm/core/props/BlackoilPhases.hpp>
27 #include <opm/models/utils/propertysystem.hh>
28 #include <opm/models/utils/parametersystem.hh>
29 #include <opm/input/eclipse/Schedule/Schedule.hpp>
30 #include <opm/input/eclipse/Schedule/Well/Well.hpp>
31 #include <opm/input/eclipse/Schedule/Group/Group.hpp>
32 #include <opm/input/eclipse/Schedule/GasLiftOpt.hpp>
33 #include <opm/input/eclipse/Schedule/SummaryState.hpp>
34 #include <opm/simulators/wells/GasLiftCommon.hpp>
35 #include <opm/simulators/wells/WellState.hpp>
36 #include <opm/simulators/utils/DeferredLogger.hpp>
37 
38 #include <algorithm>
39 #include <map>
40 #include <string>
41 #include <vector>
42 #include <fmt/format.h>
43 
44 namespace Opm
45 {
47 {
48 protected:
49  class GroupRates;
50  // NOTE: In the Well2GroupMap below, in the std::vector value we store
51  // pairs of group names and efficiency factors. The efficiency factors
52  // are the product of the wells efficiency factor and all the efficiency
53  // factors of the child groups of the group all the way down
54  // to the well group.
55  using Well2GroupMap =
56  std::map<std::string, std::vector<std::pair<std::string,double>>>;
57  using GroupRateMap =
58  std::map<std::string, GroupRates>;
59  using GroupIdxMap = std::map<std::string, int>;
60 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 7)
61  using Communication = Dune::Communication<Dune::MPIHelper::MPICommunicator>;
62 #else
63  using Communication = Dune::CollectiveCommunication<Dune::MPIHelper::MPICommunicator>;
64 #endif
65 
66  // TODO: same definition with WellInterface, and
67  // WellState eventually they should go to a common header file.
68  static const int Water = BlackoilPhases::Aqua;
69  static const int Oil = BlackoilPhases::Liquid;
70  static const int Gas = BlackoilPhases::Vapour;
71 public:
72  enum class Rate {oil, gas, water, liquid};
73 
74  using GLiftEclWells = std::map<std::string,std::pair<const Well *,int>>;
76  GLiftEclWells& ecl_wells,
77  const Schedule& schedule,
78  const SummaryState& summary_state,
79  const int report_step_idx,
80  const int iteration_idx,
81  const PhaseUsage& phase_usage,
82  DeferredLogger& deferred_logger,
83  WellState& well_state,
84  const Parallel::Communication& comm,
85  bool glift_debug
86  );
87  std::vector<std::pair<std::string,double>>& getWellGroups(
88  const std::string& well_name);
89 
90  double alqRate(const std::string& group_name);
91  double gasRate(const std::string& group_name) const;
92  int getGroupIdx(const std::string& group_name);
93  double getRate(Rate rate_type, const std::string& group_name) const;
94  std::tuple<double,double,double,double> getRates(const int group_idx) const;
95  std::optional<double> gasTarget(const std::string& group_name) const;
96  std::optional<double> getTarget(
97  Rate rate_type, const std::string& group_name) const;
98  const std::string& groupIdxToName(int group_idx) const;
99  bool hasAnyTarget(const std::string& group_name) const;
100  bool hasWell(const std::string& well_name);
101  void initialize();
102  std::optional<double> liquidTarget(const std::string& group_name) const;
103  std::optional<double> maxAlq(const std::string& group_name);
104  std::optional<double> maxTotalGasRate(const std::string& group_name);
105  double oilRate(const std::string& group_name) const;
106  std::optional<double> oilTarget(const std::string& group_name) const;
107  static const std::string rateToString(Rate rate);
108  double waterRate(const std::string& group_name) const;
109  std::optional<double> waterTarget(const std::string& group_name) const;
110  void update(const std::string& well_name,
111  double delta_oil, double delta_gas, double delta_water, double delta_alq);
112  void updateRate(int idx, double oil_rate, double gas_rate, double water_rate, double alq);
113  const Well2GroupMap& wellGroupMap() { return well_group_map_; }
114 protected:
115  bool checkDoGasLiftOptimization_(const std::string& well_name);
116  bool checkNewtonIterationIdxOk_(const std::string& well_name);
117  void debugDisplayWellContribution_(
118  const std::string& gr_name, const std::string& well_name,
119  double eff_factor,
120  double well_oil_rate, double well_gas_rate, double well_water_rate,
121  double well_alq,
122  double oil_rate, double gas_rate, double water_rate,
123  double alq
124  ) const;
125  void debugDisplayUpdatedGroupRates(const std::string& name,
126  double oil_rate, double gas_rate, double water_rate, double alq) const;
127  void debugEndInitializeGroup(const std::string& name) const;
128  void debugStartInitializeGroup(const std::string& name) const;
129  void displayDebugMessage_(const std::string& msg) const override;
130  void displayDebugMessage_(const std::string& msg, const std::string& well_name);
131  std::tuple<double, double, double> getProducerWellRates_(const int index);
132  std::tuple<double, double, double, double>
133  initializeGroupRatesRecursive_(const Group &group);
134  void initializeWell2GroupMapRecursive_(
135  const Group& group, std::vector<std::string>& group_names,
136  std::vector<double>& group_efficiency, double cur_efficiency);
137  void updateGroupIdxMap_(const std::string& group_name);
138 
139 
140  class GroupRates {
141  public:
142  GroupRates( double oil_rate, double gas_rate, double water_rate, double alq,
143  std::optional<double> oil_target,
144  std::optional<double> gas_target,
145  std::optional<double> water_target,
146  std::optional<double> liquid_target,
147  std::optional<double> total_gas,
148  std::optional<double> max_alq
149  ) :
150  oil_rate_{oil_rate},
151  gas_rate_{gas_rate},
152  water_rate_{water_rate},
153  alq_{alq},
154  oil_target_{oil_target},
155  gas_target_{gas_target},
156  water_target_{water_target},
157  liquid_target_{liquid_target},
158  total_gas_{total_gas},
159  max_alq_{max_alq}
160  {}
161  double alq() const { return alq_; }
162  void assign(double oil_rate, double gas_rate, double water_rate, double alq)
163  {
164  oil_rate_ = oil_rate;
165  gas_rate_ = gas_rate;
166  water_rate_ = water_rate;
167  alq_ = alq;
168  }
169  double gasRate() const { return gas_rate_; }
170  double waterRate() const { return water_rate_; }
171  std::optional<double> gasTarget() const { return gas_target_; }
172  std::optional<double> waterTarget() const { return water_target_; }
173  std::optional<double> maxAlq() const { return max_alq_; }
174  std::optional<double> maxTotalGasRate() const { return total_gas_; }
175  double oilRate() const { return oil_rate_; }
176  std::optional<double> oilTarget() const { return oil_target_; }
177  std::optional<double> liquidTarget() const { return liquid_target_; }
178 
179  void update(double delta_oil, double delta_gas, double delta_water, double delta_alq)
180  {
181  oil_rate_ += delta_oil;
182  gas_rate_ += delta_gas;
183  water_rate_ += delta_water;
184  alq_ += delta_alq;
185  }
186  private:
187  double oil_rate_;
188  double gas_rate_;
189  double water_rate_;
190  double alq_;
191  std::optional<double> oil_target_;
192  std::optional<double> gas_target_;
193  std::optional<double> water_target_;
194  std::optional<double> liquid_target_;
195  std::optional<double> total_gas_;
196  std::optional<double> max_alq_;
197  };
198 
199  GLiftEclWells &ecl_wells_;
200  const Schedule &schedule_;
201  const SummaryState &summary_state_;
202  const int report_step_idx_;
203  const int iteration_idx_;
204  const PhaseUsage &phase_usage_;
205  const GasLiftOpt& glo_;
206  GroupRateMap group_rate_map_;
207  Well2GroupMap well_group_map_;
208  GroupIdxMap group_idx_;
209  int next_group_idx_ = 0;
210  // Optimize only wells under THP control
211  bool optimize_only_thp_wells_ = false;
212 
213 };
214 
215 } // namespace Opm
216 
217 #endif // OPM_GASLIFT_GROUP_INFO_INCLUDED
Definition: DeferredLogger.hpp:57
Definition: GasLiftCommon.hpp:32
Definition: GasLiftGroupInfo.hpp:140
Definition: GasLiftGroupInfo.hpp:47
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