My Project
WellConnections.hpp
1 /*
2  Copyright 2016 Dr. Blatt - HPC-Simulation-Software & Services.
3  Copyright 2016 Statoil AS
4 
5  This file is part of The Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 #ifndef DUNE_CPGRID_WELL_CONNECTIONS_HEADER_INCLUDED
21 #define DUNE_CPGRID_WELL_CONNECTIONS_HEADER_INCLUDED
22 
23 #include <set>
24 #include <unordered_set>
25 #include <vector>
26 #include <functional>
27 
28 #ifdef HAVE_MPI
29 #include <opm/grid/utility/platform_dependent/disable_warnings.h>
30 #include "mpi.h"
31 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
32 #endif
33 
34 #include <dune/common/version.hh>
35 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 7)
36 #include <dune/common/parallel/communication.hh>
37 #else
38 #include <dune/common/parallel/mpicollectivecommunication.hh>
39 #endif
40 
41 
42 #include <opm/grid/utility/OpmParserIncludes.hpp>
43 #include <opm/grid/CpGrid.hpp>
44 
45 namespace Dune
46 {
47 namespace cpgrid
48 {
49 
57 {
58 
59 public:
61  typedef std::vector<std::set<int> >::const_iterator const_iterator;
62 
65 
66  WellConnections() = default;
67 
74  WellConnections(const std::vector<OpmWellType>& wells,
75  const std::array<int, 3>& cartesianSize,
76  const std::vector<int>& cartesian_to_compressed);
77 
81  WellConnections(const std::vector<OpmWellType>& wells,
82  const Dune::CpGrid& cpGrid);
83 
90  void init(const std::vector<OpmWellType>& wells,
91  const std::array<int, 3>& cartesianSize,
92  const std::vector<int>& cartesian_to_compressed);
93 
98  const std::set<int>& operator[](std::size_t i) const
99  {
100  return well_indices_[i];
101  }
102 
105  {
106  return well_indices_.begin();
107  }
108 
111  {
112  return well_indices_.end();
113  }
114 
116  std::size_t size() const
117  {
118  return well_indices_.size();
119  }
120 private:
123  std::vector<std::set<int> > well_indices_;
124 };
125 
126 
127 #ifdef HAVE_MPI
139 std::vector<std::vector<int> >
140 perforatingWellIndicesOnProc(const std::vector<int>& parts,
141  const std::vector<Dune::cpgrid::OpmWellType>& wells,
142  const CpGrid& cpgrid);
143 
161 std::vector<std::vector<int> >
162 postProcessPartitioningForWells(std::vector<int>& parts,
163  std::function<int(int)> gid,
164  const std::vector<OpmWellType>& wells,
165  const WellConnections& well_connections,
166  std::vector<std::tuple<int,int,char>>& exportList,
167  std::vector<std::tuple<int,int,char,int>>& importList,
168 #if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
169  const Communication<MPI_Comm>& cc);
170 #else
171  const CollectiveCommunication<MPI_Comm>& cc);
172 #endif
173 
182 std::vector<std::pair<std::string,bool>>
183 computeParallelWells(const std::vector<std::vector<int> >& wells_on_proc,
184  const std::vector<OpmWellType>& wells,
185 #if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
186  const Communication<MPI_Comm>& cc,
187 #else
188  const CollectiveCommunication<MPI_Comm>& cc,
189 #endif
190  int root);
191 #endif
192 } // end namespace cpgrid
193 } // end namespace Dune
194 #endif
[ provides Dune::Grid ]
Definition: CpGrid.hpp:210
A class calculating and representing all connections of wells.
Definition: WellConnections.hpp:57
const_iterator begin() const
Get a begin iterator.
Definition: WellConnections.hpp:104
const std::set< int > & operator[](std::size_t i) const
Access all connections of a well.
Definition: WellConnections.hpp:98
const_iterator end() const
Get the end iterator.
Definition: WellConnections.hpp:110
void init(const std::vector< OpmWellType > &wells, const std::array< int, 3 > &cartesianSize, const std::vector< int > &cartesian_to_compressed)
Initialze the data of the container.
Definition: WellConnections.cpp:81
const_iterator iterator
The iterator type (always const).
Definition: WellConnections.hpp:64
std::vector< std::set< int > >::const_iterator const_iterator
The const iterator type.
Definition: WellConnections.hpp:61
std::size_t size() const
\breif Get the number of wells
Definition: WellConnections.hpp:116
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10