My Project
Loading...
Searching...
No Matches
CollectDataOnIORank.hpp
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*/
23#ifndef OPM_COLLECT_DATA_ON_IO_RANK_HPP
24#define OPM_COLLECT_DATA_ON_IO_RANK_HPP
25
26#include <opm/grid/common/p2pcommunicator.hh>
27
28#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
29
30#include <opm/output/data/Aquifer.hpp>
31#include <opm/output/data/Cells.hpp>
32#include <opm/output/data/Groups.hpp>
33#include <opm/output/data/Solution.hpp>
34#include <opm/output/data/Wells.hpp>
35
36#include <opm/simulators/flow/FlowsData.hpp>
38
39#include <dune/common/version.hh>
40
41#include <array>
42#include <cstddef>
43#include <map>
44#include <set>
45#include <string>
46#include <type_traits>
47#include <utility>
48#include <vector>
49
50namespace Dune {
51template<class Grid> class CartesianIndexMapper;
52}
53
54namespace Opm {
55
56template <class Grid, class EquilGrid, class GridView>
58{
59public:
60#if DUNE_VERSION_GTE(DUNE_GRID, 2, 9)
61 using CollectiveCommunication = typename Grid::Communication;
62#else
63 using CollectiveCommunication = typename Grid::CollectiveCommunication;
64#endif
65 using P2PCommunicatorType = Dune::Point2PointCommunicator<Dune::SimpleMessageBuffer>;
66 using IndexMapType = std::vector<int>;
67 using IndexMapStorageType = std::vector<IndexMapType>;
68
69 static constexpr int dimension = Grid::dimension;
70
71 enum { ioRank = 0 };
72
73 static const bool needsReordering =
74 !std::is_same<Grid, EquilGrid>::value;
75
76 CollectDataOnIORank(const Grid& grid,
77 const EquilGrid* equilGrid,
78 const GridView& gridView,
81 const std::set<std::string>& fipRegionsInterregFlow = {});
82
83 // gather solution to rank 0 for EclipseWriter
84 void collect(const data::Solution& localCellData,
85 const std::map<std::pair<std::string, int>, double>& localBlockData,
86 const data::Wells& localWellData,
87 const data::WellBlockAveragePressures& localWBPData,
88 const data::GroupAndNetworkValues& localGroupAndNetworkData,
89 const data::Aquifers& localAquiferData,
92 const std::array<FlowsData<double>, 3>& localFlowsn,
93 const std::array<FlowsData<double>, 3>& localFloresn);
94
95 const std::map<std::pair<std::string, int>, double>& globalBlockData() const
96 { return globalBlockData_; }
97
98 const data::Solution& globalCellData() const
99 { return globalCellData_; }
100
101 data::Solution& globalCellData()
102 { return globalCellData_; }
103
104 const data::Wells& globalWellData() const
105 { return globalWellData_; }
106
107 const data::WellBlockAveragePressures& globalWBPData() const
108 { return this->globalWBPData_; }
109
110 const data::GroupAndNetworkValues& globalGroupAndNetworkData() const
111 { return globalGroupAndNetworkData_; }
112
113 const data::Aquifers& globalAquiferData() const
114 { return globalAquiferData_; }
115
116 const WellTestState& globalWellTestState() const
117 { return this->globalWellTestState_; }
118
119 InterRegFlowMap& globalInterRegFlows()
120 { return this->globalInterRegFlows_; }
121
122 const InterRegFlowMap& globalInterRegFlows() const
123 { return this->globalInterRegFlows_; }
124
125 const std::array<FlowsData<double>, 3>& globalFlowsn() const
126 { return globalFlowsn_; }
127
128 const std::array<FlowsData<double>, 3>& globalFloresn() const
129 { return globalFloresn_; }
130
131 bool isIORank() const
132 { return toIORankComm_.rank() == ioRank; }
133
134 bool isParallel() const
135 { return toIORankComm_.size() > 1; }
136
137 int localIdxToGlobalIdx(unsigned localIdx) const;
138
139 const std::vector<int>& localIdxToGlobalIdxMapping() const
140 {
141 return localIdxToGlobalIdx_;
142 }
143
144 bool doesNeedReordering() const
145 { return needsReordering;}
146
147 std::size_t numCells () const
148 { return globalCartesianIndex_.size(); }
149
150 const std::vector<int>& globalRanks() const
151 { return globalRanks_; }
152
153 bool isCartIdxOnThisRank(int cartIdx) const;
154
155protected:
156 P2PCommunicatorType toIORankComm_;
157 InterRegFlowMap globalInterRegFlows_;
158 IndexMapType globalCartesianIndex_;
159 IndexMapType localIndexMap_;
160 IndexMapStorageType indexMaps_;
161 std::vector<int> globalRanks_;
162 data::Solution globalCellData_;
163 std::map<std::pair<std::string, int>, double> globalBlockData_;
164 data::Wells globalWellData_;
165 data::WellBlockAveragePressures globalWBPData_;
166 data::GroupAndNetworkValues globalGroupAndNetworkData_;
167 data::Aquifers globalAquiferData_;
168 WellTestState globalWellTestState_;
169 std::vector<int> localIdxToGlobalIdx_;
170 std::array<FlowsData<double>, 3> globalFlowsn_;
171 std::array<FlowsData<double>, 3> globalFloresn_;
175 std::vector<int> sortedCartesianIdx_;
176};
177
178} // end namespace Opm
179
180#endif // OPM_COLLECT_DATA_ON_IO_RANK_HPP
MPI-aware facility for converting collection of tuples of region ID pairs and associate flow rates in...
Definition CollectDataOnIORank.hpp:51
Definition CollectDataOnIORank.hpp:58
std::vector< int > sortedCartesianIdx_
sorted list of cartesian indices present-
Definition CollectDataOnIORank.hpp:175
Inter-region flow accumulation maps for all region definition arrays.
Definition InterRegFlows.hpp:179
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
Simple container for FLOWS data.
Definition FlowsData.hpp:32