19#ifndef OPM_PARALLELWELLINFO_HEADER_INCLUDED
20#define OPM_PARALLELWELLINFO_HEADER_INCLUDED
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>
28#include <opm/simulators/utils/ParallelCommunication.hpp>
31#include <unordered_map>
52 using LocalIndex = Dune::ParallelLocalIndex<Attribute>;
53 using IndexSet = Dune::ParallelIndexSet<int,LocalIndex,50>;
55 using RI = Dune::RemoteIndices<IndexSet>;
107 template<
class RAIterator>
113 int numLocalPerfs()
const;
116 Parallel::Communication comm_;
118 IndexSet current_indices_;
126 std::size_t num_local_perfs_{};
135template<
class Scalar>
139 using IndexSet =
typename CommunicateAboveBelow<Scalar>::IndexSet;
140 using Attribute =
typename CommunicateAboveBelow<Scalar>::Attribute;
141 using GlobalIndex =
typename IndexSet::IndexPair::GlobalIndex;
146 const Parallel::Communication comm,
161 void copyGlobalToLocal(
const std::vector<Scalar>& global, std::vector<Scalar>& local,
164 int numGlobalPerfs()
const;
166 int localToGlobal(std::size_t localIndex)
const;
169 void buildLocalToGlobalMap()
const;
170 void buildGlobalToLocalMap()
const;
171 mutable std::unordered_map<std::size_t, int> local_to_global_map_;
172 mutable std::unordered_map<int, std::size_t> global_to_local_map_;
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_;
193template<
class Scalar>
197 static constexpr int INVALID_ECL_INDEX = -1;
210 Parallel::Communication
allComm);
212 const Parallel::Communication& communication()
const
221 int localToGlobal(std::size_t localIndex)
const;
236 std::size_t size)
const;
242 const std::vector<Scalar>&
current)
const;
251 std::size_t size)
const;
257 const std::vector<Scalar>&
current)
const;
269 const std::string&
name()
const
277 return hasLocalCells_;
293 template<
typename It>
303 template<
class RAIterator>
306 commAboveBelow_->partialSumPerfValues(begin, end);
324 void operator()(Parallel::Communication* comm);
335 int rankWithFirstPerf_;
339 std::unique_ptr<Parallel::Communication, DestroyComm> comm_;
342 std::unique_ptr<CommunicateAboveBelow<Scalar>> commAboveBelow_;
344 std::unique_ptr<GlobalPerfContainerFactory<Scalar>> globalPerfCont_;
350template<
class Scalar>
363 bool checkAllConnectionsFound();
366 std::vector<std::size_t> foundConnections_;
371template<
class Scalar>
374template<
class Scalar>
377template<
class Scalar>
380template<
class Scalar>
383template<
class Scalar>
386template<
class Scalar>
389template<
class Scalar>
392template<
class Scalar>
395template<
class Scalar>
Class checking that all connections are on active cells.
Definition ParallelWellInfo.hpp:352
void connectionFound(std::size_t index)
Inidicate that the i-th completion was found.
Definition ParallelWellInfo.cpp:736
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:635
const std::string & name() const
Name of the well.
Definition ParallelWellInfo.hpp:269
T broadcastFirstPerforationValue(const T &t) const
If the well does not have any open connections the member rankWithFirstPerf is not initialized,...
Definition ParallelWellInfo.cpp:596
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:616
void beginReset()
Inidicate that we will reset the ecl index information.
Definition ParallelWellInfo.cpp:561
void pushBackEclIndex(int above, int current)
Adds information about the ecl indices of the perforations.
Definition ParallelWellInfo.cpp:555
bool hasLocalCells() const
Whether local cells are perforated somewhen.
Definition ParallelWellInfo.hpp:275
It::value_type sumPerfValues(It begin, It end) const
Sum all the values of the perforations.
Definition ParallelWellInfo.cpp:579
void partialSumPerfValues(RAIterator begin, RAIterator end) const
Do a (in place) partial sum on values attached to all perforations.
Definition ParallelWellInfo.hpp:304
void clear()
Free data of communication data structures.
Definition ParallelWellInfo.cpp:588
const GlobalPerfContainerFactory< Scalar > & getGlobalPerfContainerFactory() const
Get a factor to create a global representation of peforation data.
Definition ParallelWellInfo.cpp:654
void communicateFirstPerforation(bool hasFirst)
Collectively decide which rank has first perforation.
Definition ParallelWellInfo.cpp:543
void endReset()
Inidicate completion of reset of the ecl index information.
Definition ParallelWellInfo.cpp:567
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