My Project
PropsCentroidsDataHandle.hpp
1 /*
2  Copyright 2020 Equinor AS.
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 3 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 */
26 #ifndef PROPS_CENTROIDS_DATAHANDLE_HPP
27 #define PROPS_CENTROIDS_DATAHANDLE_HPP
28 
29 #if HAVE_MPI
30 #include <opm/simulators/utils/ParallelEclipseState.hpp>
31 #include <opm/simulators/utils/ParallelRestart.hpp>
32 #include <dune/grid/common/datahandleif.hh>
33 #include <dune/grid/common/mcmgmapper.hh>
34 #include <dune/grid/common/partitionset.hh>
35 #include <dune/common/parallel/mpihelper.hh>
36 #include <unordered_map>
37 #include <iostream>
38 
39 namespace Opm
40 {
41 
47 template<class Grid>
48 class PropsCentroidsDataHandle
49  : public Dune::CommDataHandleIF< PropsCentroidsDataHandle<Grid>, double>
50 {
51 public:
53  using DataType = std::pair<double, unsigned char>;
54 
64  PropsCentroidsDataHandle(const Grid& grid, ParallelEclipseState& eclState,
65  const EclipseGrid* eclGridOnRoot,
66  std::vector<double>& centroids,
67  const typename Dune::CartesianIndexMapper<Grid>& cartMapper)
68  : m_grid(grid),
69  m_distributed_fieldProps(eclState.m_fieldProps),
70  m_centroids(centroids)
71  {
72  // Scatter the keys
73  const Parallel::Communication comm = m_grid.comm();
74  if (comm.rank() == 0)
75  {
76  const FieldPropsManager& globalProps = eclState.globalFieldProps();
77  m_intKeys = globalProps.keys<int>();
78  m_doubleKeys = globalProps.keys<double>();
79  m_distributed_fieldProps.copyTran(globalProps);
80  }
81 
82  EclMpiSerializer ser(comm);
83  ser.broadcast(*this);
84 
85  m_no_data = m_intKeys.size() + m_doubleKeys.size() + Grid::dimensionworld;
86 
87  if (comm.rank() == 0) {
88  const FieldPropsManager& globalProps = eclState.globalFieldProps();
89  const auto& idSet = m_grid.localIdSet();
90  const auto& gridView = m_grid.levelGridView(0);
91  using ElementMapper =
92  Dune::MultipleCodimMultipleGeomTypeMapper<typename Grid::LevelGridView>;
93  ElementMapper elemMapper(gridView, Dune::mcmgElementLayout());
94 
95  for (const auto &element : elements(gridView, Dune::Partitions::interiorBorder))
96  {
97  const auto& id = idSet.id(element);
98  auto index = elemMapper.index(element);
99  auto& data = elementData_[id];
100  data.reserve(m_no_data);
101 
102  for (const auto& intKey : m_intKeys)
103  {
104  const auto& fieldData = globalProps.get_int_field_data(intKey);
105  data.emplace_back(fieldData.data[index],
106  static_cast<unsigned char>(fieldData.value_status[index]));
107  }
108 
109  for (const auto& doubleKey : m_doubleKeys)
110  {
111  // We need to allow unsupported keywords to get the data
112  // for TranCalculator, too.
113  const auto& fieldData = globalProps.get_double_field_data(doubleKey,
114  /* allow_unsupported = */ true);
115  data.emplace_back(fieldData.data[index],
116  static_cast<unsigned char>(fieldData.value_status[index]));
117  }
118 
119  auto cartIndex = cartMapper.cartesianIndex(index);
120  const auto& center = eclGridOnRoot->getCellCenter(cartIndex);
121  for (int dim = 0; dim < Grid::dimensionworld; ++dim)
122  data.emplace_back(center[dim], '1'); // write garbage for value_status
123  }
124  }
125  }
126 
127  ~PropsCentroidsDataHandle()
128  {
129  // distributed grid is now correctly set up.
130  for (const auto& intKey : m_intKeys)
131  {
132  m_distributed_fieldProps.m_intProps[intKey].data.resize(m_grid.size(0));
133  m_distributed_fieldProps.m_intProps[intKey].value_status.resize(m_grid.size(0));
134  }
135 
136  for (const auto& doubleKey : m_doubleKeys)
137  {
138  m_distributed_fieldProps.m_doubleProps[doubleKey].data.resize(m_grid.size(0));
139  m_distributed_fieldProps.m_doubleProps[doubleKey].value_status.resize(m_grid.size(0));
140  }
141 
142  m_centroids.resize(m_grid.size(0) * Grid::dimensionworld);
143 
144  // copy data for the persistent mao to the field properties
145  const auto& idSet = m_grid.localIdSet();
146  const auto& gridView = m_grid.levelGridView(0);
147  using ElementMapper =
148  Dune::MultipleCodimMultipleGeomTypeMapper<typename Grid::LevelGridView>;
149  ElementMapper elemMapper(gridView, Dune::mcmgElementLayout());
150 
151  for (const auto &element : elements( gridView, Dune::Partitions::all))
152  {
153  std::size_t counter{};
154  const auto& id = idSet.id(element);
155  auto index = elemMapper.index(element);
156  auto data = elementData_.find(id);
157  assert(data != elementData_.end());
158 
159  for (const auto& intKey : m_intKeys)
160  {
161  const auto& pair = data->second[counter++];
162  m_distributed_fieldProps.m_intProps[intKey].data[index] = static_cast<int>(pair.first);
163  m_distributed_fieldProps.m_intProps[intKey].value_status[index] = static_cast<value::status>(pair.second);
164  }
165 
166  for (const auto& doubleKey : m_doubleKeys)
167  {
168  const auto& pair = data->second[counter++];
169  m_distributed_fieldProps.m_doubleProps[doubleKey].data[index] = pair.first;
170  m_distributed_fieldProps.m_doubleProps[doubleKey].value_status[index] = static_cast<value::status>(pair.second);
171  }
172 
173  auto centroidIter = m_centroids.begin() + Grid::dimensionworld * index;
174  auto centroidIterEnd = centroidIter + Grid::dimensionworld;
175  for ( ; centroidIter != centroidIterEnd; ++centroidIter )
176  *centroidIter = data->second[counter++].first; // value_status discarded
177  }
178  }
179 
180  bool contains(int /* dim */, int codim)
181  {
182  return codim == 0;
183  }
184 
185  bool fixedsize(int /* dim */, int /* codim */)
186  {
187  return true;
188  }
189  bool fixedSize(int /* dim */, int /* codim */)
190  {
191  return true;
192  }
193 
194  template<class EntityType>
195  std::size_t size(const EntityType /* entity */)
196  {
197  return m_no_data;
198  }
199 
200  template<class BufferType, class EntityType>
201  void gather(BufferType& buffer, const EntityType& e) const
202  {
203  auto iter = elementData_.find(m_grid.localIdSet().id(e));
204  assert(iter != elementData_.end());
205  for (const auto& data : iter->second)
206  {
207  buffer.write(data);
208  }
209  }
210 
211  template<class BufferType, class EntityType>
212  void scatter(BufferType& buffer, const EntityType& e, std::size_t n)
213  {
214  assert(n == m_no_data);
215  auto& array = elementData_[m_grid.localIdSet().id(e)];
216  array.resize(n);
217  for (auto& data : array)
218  {
219  buffer.read(data);
220  }
221  }
222 
223  template<class Serializer>
224  void serializeOp(Serializer& serializer)
225  {
226  serializer(m_intKeys);
227  serializer(m_doubleKeys);
228  m_distributed_fieldProps.serializeOp(serializer);
229  }
230 
231 private:
232  using LocalIdSet = typename Grid::LocalIdSet;
233  const Grid& m_grid;
235  ParallelFieldPropsManager& m_distributed_fieldProps;
237  std::vector<std::string> m_intKeys;
239  std::vector<std::string> m_doubleKeys;
243  std::unordered_map<typename LocalIdSet::IdType, std::vector<std::pair<double,unsigned char> > > elementData_;
245  std::vector<double>& m_centroids;
247  std::size_t m_no_data;
248 };
249 
250 } // end namespace Opm
251 #endif // HAVE_MPI
252 #endif // PROPS_CENTROIDS_DATAHANDLE_HPP
253 
Definition: findOverlapRowsAndColumns.hpp:29
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27