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#include <opm/grid/cpgrid/LevelCartesianIndexMapper.hpp>
32
34
35#include <functional>
36#include <memory>
37#include <optional>
38#include <vector>
39
40#if HAVE_MPI
41#include <filesystem>
42#endif
43
44namespace Opm {
45 class EclipseState;
46 class Schedule;
47 class Well;
48 class ParallelEclipseState;
49 class LgrCollection;
50}
51
52namespace Opm {
53
54#if HAVE_MPI
55namespace details {
57{
58public:
59 explicit MPIPartitionFromFile(const std::filesystem::path& partitionFile)
60 : partitionFile_(partitionFile)
61 {}
62
63 std::vector<int> operator()(const Dune::CpGrid& grid) const;
64
65private:
66 std::filesystem::path partitionFile_{};
67};
68
69} // namespace Opm::details
70#endif // HAVE_MPI
71
72
76extern std::optional<std::function<std::vector<int> (const Dune::CpGrid&)>> externalLoadBalancer;
77
78template<class ElementMapper, class GridView, class Scalar>
80protected:
83 using Element = typename GridView::template Codim<0>::Entity;
84
85public:
87
88 virtual ~GenericCpGridVanguard() = default;
89
93 Dune::CpGrid& grid()
94 { return *grid_; }
95
99 const Dune::CpGrid& grid() const
100 { return *grid_; }
101
111 const Dune::CpGrid& equilGrid() const;
112
120 void releaseEquilGrid();
121
125 static void setExternalLoadBalancer(const std::function<std::vector<int> (const Dune::CpGrid&)>& loadBalancer)
126 {
128 }
129
134 const CartesianIndexMapper& cartesianIndexMapper() const;
135
141
145 const CartesianIndexMapper& equilCartesianIndexMapper() const;
146
147 const std::vector<int>& cellPartition() const
148 {
149 return this->cell_part_;
150 }
151
152protected:
158#if HAVE_MPI
159 void doLoadBalance_(const Dune::EdgeWeightMethod edgeWeightsMethod,
160 const bool ownersFirst,
161 const bool addCorners,
162 const int numOverlap,
163 const Dune::PartitionMethod partitionMethod,
164 const bool serialPartitioning,
165 const bool enableDistributedWells,
167 const double imbalanceTol,
168 const GridView& gridView,
169 const Schedule& schedule,
170 EclipseState& eclState,
171 FlowGenericVanguard::ParallelWellStruct& parallelWells,
172 const int numJacobiBlocks,
173 const bool enableEclOutput);
174
175 void distributeFieldProps_(EclipseState& eclState);
176
177private:
178 std::vector<double> extractFaceTrans(const GridView& gridView) const;
179
180 void distributeGrid(const Dune::EdgeWeightMethod edgeWeightsMethod,
181 const bool ownersFirst,
182 const bool addCorners,
183 const int numOverlap,
184 const Dune::PartitionMethod partitionMethod,
185 const bool serialPartitioning,
186 const bool enableDistributedWells,
187 const double imbalanceTol,
188 const bool loadBalancerSet,
189 const std::vector<double>& faceTrans,
190 const std::vector<Well>& wells,
191 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
192 EclipseState& eclState,
193 FlowGenericVanguard::ParallelWellStruct& parallelWells);
194
195 void distributeGrid(const Dune::EdgeWeightMethod edgeWeightsMethod,
196 const bool ownersFirst,
197 const bool addCorners,
198 const int numOverlap,
199 const Dune::PartitionMethod partitionMethod,
200 const bool serialPartitioning,
201 const bool enableDistributedWells,
202 const double imbalanceTol,
203 const bool loadBalancerSet,
204 const std::vector<double>& faceTrans,
205 const std::vector<Well>& wells,
206 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
207 ParallelEclipseState* eclState,
208 FlowGenericVanguard::ParallelWellStruct& parallelWells);
209
210protected:
211 virtual const std::string& zoltanParams() const = 0;
212 virtual double zoltanPhgEdgeSizeThreshold() const = 0;
213 virtual const std::string& metisParams() const = 0;
214
215#endif // HAVE_MPI
216
218
219 void doCreateGrids_(EclipseState& eclState);
220 void addLgrsUpdateLeafView(const LgrCollection& lgrCollection,
221 const int lgrsSize,
222 Dune::CpGrid& grid);
223
224 virtual void allocTrans() = 0;
225 virtual double getTransmissibility(unsigned I, unsigned J) const = 0;
226
227 // removing some connection located in inactive grid cells
228 void doFilterConnections_(Schedule& schedule);
229
230 Scalar computeCellThickness(const Element& element) const;
231
232 std::unique_ptr<Dune::CpGrid> grid_;
233 std::unique_ptr<Dune::CpGrid> equilGrid_;
234 std::unique_ptr<CartesianIndexMapper> cartesianIndexMapper_;
235 std::unique_ptr<CartesianIndexMapper> equilCartesianIndexMapper_;
236
237 int mpiRank;
238 std::vector<int> cell_part_{};
239};
240
241} // namespace Opm
242
243#endif // OPM_GENERIC_CPGRID_VANGUARD_HPP
Helper class for grid instantiation of ECL file-format using problems.
Definition CollectDataOnIORank.hpp:49
Definition GenericCpGridVanguard.hpp:79
void releaseEquilGrid()
Indicates that the initial condition has been computed and the memory used by the EQUIL grid can be r...
Definition GenericCpGridVanguard.cpp:136
const Dune::CpGrid & equilGrid() const
Returns a refefence to the grid which should be used by the EQUIL initialization code.
Definition GenericCpGridVanguard.cpp:615
const CartesianIndexMapper & cartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition GenericCpGridVanguard.cpp:623
const CartesianIndexMapper & equilCartesianIndexMapper() const
Returns mapper from compressed to cartesian indices for the EQUIL grid.
Definition GenericCpGridVanguard.cpp:637
const Dune::CpGrid & grid() const
Return a reference to the simulation grid.
Definition GenericCpGridVanguard.hpp:99
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:125
Dune::CpGrid & grid()
Return a reference to the simulation grid.
Definition GenericCpGridVanguard.hpp:93
void allocCartMapper()
Distribute the simulation grid over multiple processes.
const LevelCartesianIndexMapper levelCartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition GenericCpGridVanguard.cpp:630
Definition RelpermDiagnostics.hpp:31
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:123
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242