My Project
GridPartitioning.hpp
1 //===========================================================================
2 //
3 // File: GridPartitioning.hpp
4 //
5 // Created: Mon Sep 7 10:09:13 2009
6 //
7 // Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8 // Bård Skaflestad <bard.skaflestad@sintef.no>
9 //
10 // $Date$
11 //
12 // $Revision$
13 //
14 //===========================================================================
15 
16 /*
17  Copyright 2009, 2010, 2013 SINTEF ICT, Applied Mathematics.
18  Copyright 2009, 2010 Statoil ASA.
19  Copyright 2013 Dr. Markus Blatt - HPC-Simulation-Software & Services
20  Copyright 2020, 2021 OPM-OP AS
21 
22  This file is part of The Open Porous Media project (OPM).
23 
24  OPM is free software: you can redistribute it and/or modify
25  it under the terms of the GNU General Public License as published by
26  the Free Software Foundation, either version 3 of the License, or
27  (at your option) any later version.
28 
29  OPM is distributed in the hope that it will be useful,
30  but WITHOUT ANY WARRANTY; without even the implied warranty of
31  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
32  GNU General Public License for more details.
33 
34  You should have received a copy of the GNU General Public License
35  along with OPM. If not, see <http://www.gnu.org/licenses/>.
36 */
37 
38 #ifndef OPM_GRIDPARTITIONING_HEADER
39 #define OPM_GRIDPARTITIONING_HEADER
40 
41 #include <vector>
42 #include <array>
43 #include <set>
44 #include <tuple>
45 
46 #include <dune/common/parallel/mpihelper.hh>
47 
48 #include <opm/grid/utility/OpmParserIncludes.hpp>
49 namespace Dune
50 {
51 
52  class CpGrid;
53 
54  struct OrderByFirst
55  {
56  bool operator()(const std::pair<int,int>& o, const std::pair<int,int>& v)
57  {
58  return o.first < v.first;
59  }
60  };
61 
70  void partition(const CpGrid& grid,
71  const std::array<int, 3>& initial_split,
72  int& num_part,
73  std::vector<int>& cell_part,
74  bool recursive = false,
75  bool ensureConnectivity = true);
76 
85  void addOverlapLayer(const CpGrid& grid,
86  const std::vector<int>& cell_part,
87  std::vector<std::set<int> >& cell_overlap,
88  int mypart, int overlapLayers, bool all=false);
89 
102  int addOverlapLayer(const CpGrid& grid, const std::vector<int>& cell_part,
103  std::vector<std::tuple<int,int,char>>& exportList,
104  std::vector<std::tuple<int,int,char,int>>& importList,
105  const CollectiveCommunication<Dune::MPIHelper::MPICommunicator>& cc,
106  bool addCornerCells, const double* trans, int layers = 1);
107 
108 namespace cpgrid
109 {
110 #if HAVE_MPI
127  std::tuple<std::vector<int>, std::vector<std::pair<std::string,bool>>,
128  std::vector<std::tuple<int,int,char> >,
129  std::vector<std::tuple<int,int,char,int> > >
130  createZoltanListsFromParts(const CpGrid& grid, const std::vector<cpgrid::OpmWellType> * wells,
131  const double* transmissibilities, const std::vector<int>& parts,
132  bool allowDistributedWells);
133 
149  std::tuple<std::vector<int>, std::vector<std::pair<std::string,bool>>,
150  std::vector<std::tuple<int,int,char> >,
151  std::vector<std::tuple<int,int,char,int> > >
152  vanillaPartitionGridOnRoot(const CpGrid& grid,
153  const std::vector<cpgrid::OpmWellType> * wells,
154  const double* transmissibilities,
155  bool allowDistributedWells);
156 #endif
157 } // namespace cpgrid
158 } // namespace Dune
159 
160 
161 #endif // OPM_GRIDPARTITIONING_HEADER
[ provides Dune::Grid ]
Definition: CpGrid.hpp:207
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
void partition(const CpGrid &grid, const coord_t &initial_split, int &num_part, std::vector< int > &cell_part, bool recursive, bool ensureConnectivity)
Partition a CpGrid based on (ijk) coordinates, with splitting to ensure that each partition is connec...
Definition: GridPartitioning.cpp:195
Definition: GridPartitioning.hpp:55