My Project
Loading...
Searching...
No Matches
ParallelWellInfo.hpp
1/*
2 Copyright 2020 OPM-OP AS
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#ifndef OPM_PARALLELWELLINFO_HEADER_INCLUDED
20#define OPM_PARALLELWELLINFO_HEADER_INCLUDED
21
22#include <dune/common/parallel/communicator.hh>
23#include <dune/common/parallel/interface.hh>
24#include <dune/common/parallel/mpihelper.hh>
25#include <dune/common/parallel/plocalindex.hh>
26#include <dune/common/parallel/remoteindices.hh>
27
28#include <opm/simulators/utils/ParallelCommunication.hpp>
29
30#include <memory>
31#include <unordered_map>
32
33namespace Opm {
34
35class Well;
36
39template<class Scalar>
41{
42public:
43 enum Attribute {
44 owner=1,
45 overlap=2,
46 // there is a bug in older versions of DUNE that will skip
47 // entries with matching attributes in RemoteIndices that are local
48 // therefore we add one more version for above.
49 ownerAbove = 3,
50 overlapAbove = 4
51 };
52 using LocalIndex = Dune::ParallelLocalIndex<Attribute>;
53 using IndexSet = Dune::ParallelIndexSet<int,LocalIndex,50>;
54#if HAVE_MPI
55 using RI = Dune::RemoteIndices<IndexSet>;
56#endif
57
58 explicit CommunicateAboveBelow(const Parallel::Communication& comm);
66 void pushBackEclIndex(int above, int current, bool owner=true);
67
69 void clear();
70
73 void beginReset();
74
80 int endReset();
81
87 std::vector<Scalar> communicateAbove(Scalar first_value,
88 const Scalar* current,
89 std::size_t size);
90
96 std::vector<Scalar> communicateBelow(Scalar first_value,
97 const Scalar* current,
98 std::size_t size);
99
107 template<class RAIterator>
108 void partialSumPerfValues(RAIterator begin, RAIterator end) const;
109
111 const IndexSet& getIndexSet() const;
112
113 int numLocalPerfs() const;
114
115private:
116 Parallel::Communication comm_;
118 IndexSet current_indices_;
119#if HAVE_MPI
121 IndexSet above_indices_;
123 Dune::Interface interface_;
124 Dune::BufferedCommunicator communicator_;
125#endif
126 std::size_t num_local_perfs_{};
127};
128
135template<class Scalar>
137{
138public:
139 using IndexSet = typename CommunicateAboveBelow<Scalar>::IndexSet;
140 using Attribute = typename CommunicateAboveBelow<Scalar>::Attribute;
141 using GlobalIndex = typename IndexSet::IndexPair::GlobalIndex;
142
146 const Parallel::Communication comm,
147 int num_local_perfs);
148
154 std::vector<Scalar> createGlobal(const std::vector<Scalar>& local_perf_container,
155 std::size_t num_components) const;
156
161 void copyGlobalToLocal(const std::vector<Scalar>& global, std::vector<Scalar>& local,
162 std::size_t num_components) const;
163
164 int numGlobalPerfs() const;
165 int globalToLocal(const int globalIndex) const;
166 int localToGlobal(std::size_t localIndex) const;
167
168private:
169 void buildLocalToGlobalMap() const;
170 void buildGlobalToLocalMap() const;
171 mutable std::unordered_map<std::size_t, int> local_to_global_map_; // Cache for L2G mapping
172 mutable std::unordered_map<int, std::size_t> global_to_local_map_; // Cache for G2L mapping
173 mutable bool l2g_map_built_ = false;
174 mutable bool g2l_map_built_ = false;
175 const IndexSet& local_indices_;
176 Parallel::Communication comm_;
177 int num_global_perfs_;
179 std::vector<int> sizes_;
181 std::vector<int> displ_;
183 std::vector<int> map_received_;
187 std::vector<int> perf_ecl_index_;
188};
189
193template<class Scalar>
195{
196public:
197 static constexpr int INVALID_ECL_INDEX = -1;
198
200 explicit ParallelWellInfo(const std::string& name = {""},
201 bool hasLocalCells = true);
202
209 ParallelWellInfo(const std::pair<std::string,bool>& well_info,
210 Parallel::Communication allComm);
211
212 const Parallel::Communication& communication() const
213 {
214 return *comm_;
215 }
216
219
220 void setActiveToLocalMap(const std::unordered_map<int,int> active_to_local_map) const;
221 int activeToLocal(const int activeIndex) const;
222 int localToActive(std::size_t localIndex) const;
223 int globalToLocal(const int globalIndex) const;
224 int localToGlobal(std::size_t localIndex) const;
225
229 template<class T>
230 T broadcastFirstPerforationValue(const T& t) const;
231
237 std::vector<Scalar> communicateAboveValues(Scalar first_value,
238 const Scalar* current,
239 std::size_t size) const;
240
244 std::vector<Scalar> communicateAboveValues(Scalar first_value,
245 const std::vector<Scalar>& current) const;
246
252 std::vector<Scalar> communicateBelowValues(Scalar last_value,
253 const Scalar* current,
254 std::size_t size) const;
255
259 std::vector<Scalar> communicateBelowValues(Scalar last_value,
260 const std::vector<Scalar>& current) const;
261
269 void pushBackEclIndex(int above, int current);
270
272 const std::string& name() const
273 {
274 return name_;
275 }
276
278 bool hasLocalCells() const
279 {
280 return hasLocalCells_;
281 }
282 bool isOwner() const
283 {
284 return isOwner_;
285 }
286
290 void beginReset();
291
293 void endReset();
294
296 template<typename It>
297 typename It::value_type sumPerfValues(It begin, It end) const;
298
306 template<class RAIterator>
308 {
309 commAboveBelow_->partialSumPerfValues(begin, end);
310 }
311
313 void clear();
314
321
322private:
323
325 struct DestroyComm
326 {
327 void operator()(Parallel::Communication* comm);
328 };
329
330
332 std::string name_;
334 bool hasLocalCells_;
336 bool isOwner_;
338 int rankWithFirstPerf_;
342 std::unique_ptr<Parallel::Communication, DestroyComm> comm_;
343
345 std::unique_ptr<CommunicateAboveBelow<Scalar>> commAboveBelow_;
346
347 std::unique_ptr<GlobalPerfContainerFactory<Scalar>> globalPerfCont_;
348
349 mutable std::unordered_map<int,int> active_to_local_map_; // Cache for active perforation index to local perforation index mapping
350 mutable std::unordered_map<int,int> local_to_active_map_; // Cache for local perforation index to active perforation index mapping
351};
352
356template<class Scalar>
358{
359public:
361 const ParallelWellInfo<Scalar>& info);
362
367 void connectionFound(std::size_t index);
368
369 bool checkAllConnectionsFound();
370
371private:
372 std::vector<std::size_t> foundConnections_;
373 const Well& well_;
374 const ParallelWellInfo<Scalar>& pwinfo_;
375};
376
377template<class Scalar>
378bool operator<(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2);
379
380template<class Scalar>
381bool operator==(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2);
382
383template<class Scalar>
384bool operator!=(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2);
385
386template<class Scalar>
387bool operator<(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well);
388
389template<class Scalar>
390bool operator<( const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair);
391
392template<class Scalar>
393bool operator==(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well);
394
395template<class Scalar>
396bool operator==(const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair);
397
398template<class Scalar>
399bool operator!=(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well);
400
401template<class Scalar>
402bool operator!=(const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair);
403
404} // end namespace Opm
405
406#endif // OPM_PARALLELWELLINFO_HEADER_INCLUDED
Class checking that all connections are on active cells.
Definition ParallelWellInfo.hpp:358
void connectionFound(std::size_t index)
Inidicate that the i-th completion was found.
Definition ParallelWellInfo.cpp:765
Class to facilitate getting values associated with the above/below perforation.
Definition ParallelWellInfo.hpp:41
std::vector< Scalar > communicateBelow(Scalar first_value, const Scalar *current, std::size_t size)
Creates an array of values for the perforation below.
Definition ParallelWellInfo.cpp:417
void clear()
Clear all the parallel information.
Definition ParallelWellInfo.cpp:267
void pushBackEclIndex(int above, int current, bool owner=true)
Adds information about original index of the perforations in ECL Schedule.
Definition ParallelWellInfo.cpp:447
const IndexSet & getIndexSet() const
Get index set for the local perforations.
Definition ParallelWellInfo.cpp:488
int endReset()
Indicates that the index information is complete.
Definition ParallelWellInfo.cpp:292
std::vector< Scalar > communicateAbove(Scalar first_value, const Scalar *current, std::size_t size)
Creates an array of values for the perforation above.
Definition ParallelWellInfo.cpp:387
void beginReset()
Indicates that we will add the index information.
Definition ParallelWellInfo.cpp:279
void partialSumPerfValues(RAIterator begin, RAIterator end) const
Do a (in place) partial sum on values attached to all perforations.
Definition ParallelWellInfo.cpp:313
A factory for creating a global data representation for distributed wells.
Definition ParallelWellInfo.hpp:137
std::vector< Scalar > createGlobal(const std::vector< Scalar > &local_perf_container, std::size_t num_components) const
Creates a container that holds values for all perforations.
Definition ParallelWellInfo.cpp:170
void copyGlobalToLocal(const std::vector< Scalar > &global, std::vector< Scalar > &local, std::size_t num_components) const
Copies the values of the global perforation to the local representation.
Definition ParallelWellInfo.cpp:221
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:195
std::vector< Scalar > communicateBelowValues(Scalar last_value, const Scalar *current, std::size_t size) const
Creates an array of values for the perforation below.
Definition ParallelWellInfo.cpp:664
const std::string & name() const
Name of the well.
Definition ParallelWellInfo.hpp:272
T broadcastFirstPerforationValue(const T &t) const
If the well does not have any open connections the member rankWithFirstPerf is not initialized,...
Definition ParallelWellInfo.cpp:625
std::vector< Scalar > communicateAboveValues(Scalar first_value, const Scalar *current, std::size_t size) const
Creates an array of values for the perforation above.
Definition ParallelWellInfo.cpp:645
void beginReset()
Inidicate that we will reset the ecl index information.
Definition ParallelWellInfo.cpp:590
void pushBackEclIndex(int above, int current)
Adds information about the ecl indices of the perforations.
Definition ParallelWellInfo.cpp:584
bool hasLocalCells() const
Whether local cells are perforated somewhen.
Definition ParallelWellInfo.hpp:278
It::value_type sumPerfValues(It begin, It end) const
Sum all the values of the perforations.
Definition ParallelWellInfo.cpp:608
void partialSumPerfValues(RAIterator begin, RAIterator end) const
Do a (in place) partial sum on values attached to all perforations.
Definition ParallelWellInfo.hpp:307
void clear()
Free data of communication data structures.
Definition ParallelWellInfo.cpp:617
const GlobalPerfContainerFactory< Scalar > & getGlobalPerfContainerFactory() const
Get a factor to create a global representation of peforation data.
Definition ParallelWellInfo.cpp:683
void communicateFirstPerforation(bool hasFirst)
Collectively decide which rank has first perforation.
Definition ParallelWellInfo.cpp:572
void endReset()
Inidicate completion of reset of the ecl index information.
Definition ParallelWellInfo.cpp:596
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