My Project
Loading...
Searching...
No Matches
BlackoilWellModelGasLift_impl.hpp
1/*
2 Copyright 2016 - 2019 SINTEF Digital, Mathematics & Cybernetics.
3 Copyright 2016 - 2018 Equinor ASA.
4 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5 Copyright 2016 - 2018 Norce AS
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#define OPM_BLACKOILWELLMODEL_GASLIFT_IMPL_HEADER_INCLUDED
23#define OPM_BLACKOILWELLMODEL_GASLIFT_IMPL_HEADER_INCLUDED
24
25// Improve IDE experience
26#ifndef OPM_BLACKOILWELLMODEL_GASLIFT_HEADER_INCLUDED
27#include <config.h>
28#include <opm/simulators/wells/BlackoilWellModelGasLift.hpp>
29#endif
30#include <opm/common/TimingMacros.hpp>
31#include <opm/simulators/wells/GasLiftSingleWell.hpp>
32
33#if HAVE_MPI
34#include <opm/simulators/utils/MPISerializer.hpp>
35#endif
36
37namespace Opm {
38
39template<typename TypeTag>
40bool
41BlackoilWellModelGasLift<TypeTag>::
42maybeDoGasLiftOptimize(const Simulator& simulator,
43 const std::vector<WellInterfacePtr>& well_container,
44 WellState<Scalar>& wellState,
45 GroupState<Scalar>& groupState,
46 DeferredLogger& deferred_logger)
47{
49 const auto& glo = simulator.vanguard().schedule().glo(simulator.episodeIndex());
50 if(!glo.active()) {
51 return false;
52 }
53 bool do_glift_optimization = false;
54 int num_wells_changed = 0;
55 const double simulation_time = simulator.time();
56 const Scalar min_wait = glo.min_wait();
57 // We only optimize if a min_wait time has past.
58 // If all_newton is true we still want to optimize several times pr timestep
59 // i.e. we also optimize if check simulation_time == last_glift_opt_time_
60 // that is when the last_glift_opt_time is already updated with the current time step
61 if (simulation_time == this->last_glift_opt_time_ ||
62 simulation_time >= (this->last_glift_opt_time_ + min_wait))
63 {
65 this->last_glift_opt_time_ = simulation_time;
66 }
67
69 GLiftOptWells glift_wells;
70 GLiftProdWells prod_wells;
71 GLiftWellStateMap state_map;
72 // NOTE: To make GasLiftGroupInfo (see below) independent of the TypeTag
73 // associated with *this (i.e. BlackoilWellModel<TypeTag>) we observe
74 // that GasLiftGroupInfo's only dependence on *this is that it needs to
75 // access the eclipse Wells in the well container (the eclipse Wells
76 // themselves are independent of the TypeTag).
77 // Hence, we extract them from the well container such that we can pass
78 // them to the GasLiftGroupInfo constructor.
79 GLiftEclWells ecl_well_map;
80 initGliftEclWellMap(well_container, ecl_well_map);
81 GasLiftGroupInfo group_info {
83 simulator.vanguard().schedule(),
84 simulator.vanguard().summaryState(),
85 simulator.episodeIndex(),
86 simulator.model().newtonMethod().numIterations(),
87 phase_usage_,
89 wellState,
90 groupState,
91 simulator.vanguard().grid().comm(),
92 this->glift_debug
93 };
94 group_info.initialize();
95
96 gasLiftOptimizationStage1(simulator,
98 wellState,
99 groupState,
103 state_map,
105
106 this->gasLiftOptimizationStage2(simulator.vanguard().gridView().comm(),
107 simulator.vanguard().schedule(),
108 simulator.vanguard().summaryState(),
109 wellState,
110 groupState,
114 state_map,
115 simulator.episodeIndex(),
117
118 if constexpr (glift_debug) {
119 std::vector<WellInterfaceGeneric<Scalar>*> wc;
120 wc.reserve(well_container.size());
121 for (const auto& w : well_container) {
122 wc.push_back(static_cast<WellInterfaceGeneric<Scalar>*>(w.get()));
123 }
124 this->gliftDebugShowALQ(wc,
125 wellState,
127 }
129 }
130 num_wells_changed = simulator.vanguard().gridView().comm().sum(num_wells_changed);
131 return num_wells_changed > 0;
132}
133
134template<typename TypeTag>
135void
136BlackoilWellModelGasLift<TypeTag>::
137gasLiftOptimizationStage1(const Simulator& simulator,
138 const std::vector<WellInterfacePtr>& well_container,
139 WellState<Scalar>& wellState,
140 GroupState<Scalar>& groupState,
141 GLiftProdWells& prod_wells,
142 GLiftOptWells &glift_wells,
144 GLiftWellStateMap& state_map,
145 DeferredLogger& deferred_logger)
146{
148 auto comm = simulator.vanguard().grid().comm();
149 int num_procs = comm.size();
150 // NOTE: Gas lift optimization stage 1 seems to be difficult
151 // to do in parallel since the wells are optimized on different
152 // processes and each process needs to know the current ALQ allocated
153 // to each group it is a memeber of in order to check group limits and avoid
154 // allocating more ALQ than necessary. (Surplus ALQ is removed in
155 // stage 2). In stage1, as each well is adding ALQ, the current group ALQ needs
156 // to be communicated to the other processes. But there is no common
157 // synchronization point that all process will reach in the
158 // runOptimizeLoop_() in GasLiftSingleWell.cpp.
159 //
160 // TODO: Maybe a better solution could be invented by distributing
161 // wells according to certain parent groups. Then updated group rates
162 // might not have to be communicated to the other processors.
163
164 // Currently, the best option seems to be to run this part sequentially
165 // (not in parallel).
166 //
167 // TODO: The simplest approach seems to be if a) one process could take
168 // ownership of all the wells (the union of all the wells in the
169 // well_container_ of each process) then this process could do the
170 // optimization, while the other processes could wait for it to
171 // finish (e.g. comm.barrier()), or alternatively, b) if all
172 // processes could take ownership of all the wells. Then there
173 // would be no need for synchronization here..
174 //
175 for (int i = 0; i< num_procs; i++) {
176 int num_rates_to_sync = 0; // communication variable
177 GLiftSyncGroups groups_to_sync;
178 if (comm.rank() == i) {
179 // Run stage1: Optimize single wells while also checking group limits
180 for (const auto& well : well_container) {
181 // NOTE: Only the wells in "group_info" needs to be optimized
182 if (group_info.hasWell(well->name())) {
183 gasLiftOptimizationStage1SingleWell(well.get(),
184 simulator,
185 wellState,
186 groupState,
190 state_map,
193 }
194 }
196 }
197 {
200 }
201 if (num_rates_to_sync > 0) {
203 std::vector<int> group_indexes;
205 std::vector<Scalar> group_alq_rates;
207 std::vector<Scalar> group_oil_rates;
209 std::vector<Scalar> group_gas_rates;
211 std::vector<Scalar> group_water_rates;
213 if (comm.rank() == i) {
214 for (auto idx : groups_to_sync) {
215 auto [oil_rate, gas_rate, water_rate, alq] = group_info.getRates(idx);
216 group_indexes.push_back(idx);
217 group_oil_rates.push_back(oil_rate);
218 group_gas_rates.push_back(gas_rate);
219 group_water_rates.push_back(water_rate);
220 group_alq_rates.push_back(alq);
221 }
222 } else {
228 }
229#if HAVE_MPI
230 Parallel::MpiSerializer ser(comm);
231 ser.broadcast(i, group_indexes, group_oil_rates,
233#endif
234 if (comm.rank() != i) {
235 for (int j = 0; j < num_rates_to_sync; ++j) {
236 group_info.updateRate(group_indexes[j],
240 group_alq_rates[j]);
241 }
242 }
243 if constexpr (glift_debug) {
244 int counter = 0;
245 if (comm.rank() == i) {
246 counter = wellState.gliftGetDebugCounter();
247 }
248 counter = comm.sum(counter);
249 if (comm.rank() != i) {
250 wellState.gliftSetDebugCounter(counter);
251 }
252 }
253 }
254 }
255}
256
257// NOTE: this method cannot be const since it passes this->wellState()
258// (see below) to the GasLiftSingleWell constructor which accepts WellState
259// as a non-const reference..
260template<typename TypeTag>
261void
262BlackoilWellModelGasLift<TypeTag>::
263gasLiftOptimizationStage1SingleWell(WellInterface<TypeTag>* well,
264 const Simulator& simulator,
265 WellState<Scalar>& wellState,
266 GroupState<Scalar>& groupState,
267 GLiftProdWells& prod_wells,
268 GLiftOptWells& glift_wells,
270 GLiftWellStateMap& state_map,
271 GLiftSyncGroups& sync_groups,
272 DeferredLogger& deferred_logger)
273{
275 const auto& summary_state = simulator.vanguard().summaryState();
276 auto glift = std::make_unique<GasLiftSingleWell<TypeTag>>(*well,
277 simulator,
280 wellState,
281 groupState,
284 simulator.vanguard().gridView().comm(),
285 this->glift_debug);
286 auto state = glift->runOptimize(simulator.model().newtonMethod().numIterations());
287 if (state) {
288 state_map.emplace(well->name(), std::move(state));
289 glift_wells.emplace(well->name(), std::move(glift));
290 return;
291 }
292 prod_wells.insert({well->name(), well});
293}
294
295template<typename TypeTag>
296void
297BlackoilWellModelGasLift<TypeTag>::
298initGliftEclWellMap(const std::vector<WellInterfacePtr>& well_container,
299 GLiftEclWells& ecl_well_map)
300{
301 for (const auto& well : well_container) {
302 ecl_well_map.try_emplace(well->name(), &well->wellEcl(), well->indexOfWell());
303 }
304}
305
306} // namespace Opm
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