My Project
Loading...
Searching...
No Matches
GenericCpGridVanguard.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
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 2 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 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_GENERIC_CPGRID_VANGUARD_HPP
28#define OPM_GENERIC_CPGRID_VANGUARD_HPP
29
30#include <opm/grid/CpGrid.hpp>
31
33
34#include <functional>
35#include <memory>
36#include <optional>
37#include <vector>
38
39#if HAVE_MPI
40#include <filesystem>
41#endif
42
43namespace Opm {
44 class EclipseState;
45 class Schedule;
46 class Well;
47 class ParallelEclipseState;
48 class LgrCollection;
49}
50
51namespace Opm {
52
53#if HAVE_MPI
54namespace details {
55class MPIPartitionFromFile
56{
57public:
58 explicit MPIPartitionFromFile(const std::filesystem::path& partitionFile)
59 : partitionFile_(partitionFile)
60 {}
61
62 std::vector<int> operator()(const Dune::CpGrid& grid) const;
63
64private:
65 std::filesystem::path partitionFile_{};
66};
67
68} // namespace Opm::details
69#endif // HAVE_MPI
70
71
75extern std::optional<std::function<std::vector<int> (const Dune::CpGrid&)>> externalLoadBalancer;
76
77template<class ElementMapper, class GridView, class Scalar>
79protected:
81 using Element = typename GridView::template Codim<0>::Entity;
82
83public:
85
86 virtual ~GenericCpGridVanguard() = default;
87
91 Dune::CpGrid& grid()
92 { return *grid_; }
93
97 const Dune::CpGrid& grid() const
98 { return *grid_; }
99
109 const Dune::CpGrid& equilGrid() const;
110
118 void releaseEquilGrid();
119
123 static void setExternalLoadBalancer(const std::function<std::vector<int> (const Dune::CpGrid&)>& loadBalancer)
124 {
125 externalLoadBalancer = loadBalancer;
126 }
127
132 const CartesianIndexMapper& cartesianIndexMapper() const;
133
137 const CartesianIndexMapper& equilCartesianIndexMapper() const;
138
139 const std::vector<int>& cellPartition() const
140 {
141 return this->cell_part_;
142 }
143
144protected:
150#if HAVE_MPI
151 void doLoadBalance_(const Dune::EdgeWeightMethod edgeWeightsMethod,
152 const bool ownersFirst,
153 const Dune::PartitionMethod partitionMethod,
154 const bool serialPartitioning,
155 const bool enableDistributedWells,
156 const double imbalanceTol,
157 const GridView& gridView,
158 const Schedule& schedule,
159 EclipseState& eclState,
160 FlowGenericVanguard::ParallelWellStruct& parallelWells,
161 const int numJacobiBlocks);
162
163 void distributeFieldProps_(EclipseState& eclState);
164
165private:
166 std::vector<double> extractFaceTrans(const GridView& gridView) const;
167
168 void distributeGrid(const Dune::EdgeWeightMethod edgeWeightsMethod,
169 const bool ownersFirst,
170 const Dune::PartitionMethod partitionMethod,
171 const bool serialPartitioning,
172 const bool enableDistributedWells,
173 const double imbalanceTol,
174 const bool loadBalancerSet,
175 const std::vector<double>& faceTrans,
176 const std::vector<Well>& wells,
177 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
178 EclipseState& eclState,
179 FlowGenericVanguard::ParallelWellStruct& parallelWells);
180
181 void distributeGrid(const Dune::EdgeWeightMethod edgeWeightsMethod,
182 const bool ownersFirst,
183 const Dune::PartitionMethod partitionMethod,
184 const bool serialPartitioning,
185 const bool enableDistributedWells,
186 const double imbalanceTol,
187 const bool loadBalancerSet,
188 const std::vector<double>& faceTrans,
189 const std::vector<Well>& wells,
190 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
191 ParallelEclipseState* eclState,
192 FlowGenericVanguard::ParallelWellStruct& parallelWells);
193
194protected:
195 virtual const std::string& zoltanParams() const = 0;
196 virtual const std::string& metisParams() const = 0;
197
198#endif // HAVE_MPI
199
201
202 void doCreateGrids_(EclipseState& eclState);
203 void addLgrsUpdateLeafView(const LgrCollection& lgrCollection, const int lgrsSize);
204
205 virtual void allocTrans() = 0;
206 virtual double getTransmissibility(unsigned I, unsigned J) const = 0;
207
208 // removing some connection located in inactive grid cells
209 void doFilterConnections_(Schedule& schedule);
210
211 Scalar computeCellThickness(const Element& element) const;
212
213 std::unique_ptr<Dune::CpGrid> grid_;
214 std::unique_ptr<Dune::CpGrid> equilGrid_;
215 std::unique_ptr<CartesianIndexMapper> cartesianIndexMapper_;
216 std::unique_ptr<CartesianIndexMapper> equilCartesianIndexMapper_;
217
218 int mpiRank;
219 std::vector<int> cell_part_{};
220};
221
222} // namespace Opm
223
224#endif // OPM_GENERIC_CPGRID_VANGUARD_HPP
Helper class for grid instantiation of ECL file-format using problems.
Definition findOverlapRowsAndColumns.hpp:31
Definition GenericCpGridVanguard.hpp:78
void releaseEquilGrid()
Indicates that the initial condition has been computed and the memory used by the EQUIL grid can be r...
Definition GenericCpGridVanguard.cpp:138
const Dune::CpGrid & equilGrid() const
Returns a refefence to the grid which should be used by the EQUIL initialization code.
Definition GenericCpGridVanguard.cpp:560
const CartesianIndexMapper & cartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition GenericCpGridVanguard.cpp:568
const CartesianIndexMapper & equilCartesianIndexMapper() const
Returns mapper from compressed to cartesian indices for the EQUIL grid.
Definition GenericCpGridVanguard.cpp:575
const Dune::CpGrid & grid() const
Return a reference to the simulation grid.
Definition GenericCpGridVanguard.hpp:97
static void setExternalLoadBalancer(const std::function< std::vector< int >(const Dune::CpGrid &)> &loadBalancer)
Sets a function that returns external load balancing information when passed the grid.
Definition GenericCpGridVanguard.hpp:123
Dune::CpGrid & grid()
Return a reference to the simulation grid.
Definition GenericCpGridVanguard.hpp:91
void allocCartMapper()
Distribute the simulation grid over multiple processes.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
std::optional< std::function< std::vector< int >(const Dune::CpGrid &)> > externalLoadBalancer
optional functor returning external load balancing information
Definition GenericCpGridVanguard.cpp:125