dune-pdelab  2.4-dev
borderdofexchanger.hh
Go to the documentation of this file.
1 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=8 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_GRIDOPERATOR_COMMON_BORDERDOFEXCHANGER_HH
4 #define DUNE_PDELAB_GRIDOPERATOR_COMMON_BORDERDOFEXCHANGER_HH
5 
6 #include <cstddef>
7 #include <vector>
8 #include <algorithm>
9 #include <unordered_map>
10 #include <unordered_set>
11 
12 #include <dune/common/deprecated.hh>
13 #include <dune/common/parallel/mpihelper.hh>
14 #include <dune/common/tuples.hh>
15 
16 #include <dune/geometry/typeindex.hh>
17 #include <dune/grid/common/gridenums.hh>
18 #include <dune/grid/common/datahandleif.hh>
19 
23 
24 namespace Dune {
25  namespace PDELab {
26 
27 
31 
33 
65  template<typename GridOperator>
67  {
68  typedef typename GridOperator::Traits::Jacobian M;
69  typedef M Matrix;
71  typedef typename GridOperatorTraits::JacobianField Scalar;
72  typedef typename GridOperatorTraits::TrialGridFunctionSpace GFSU;
73  typedef typename GridOperatorTraits::TestGridFunctionSpace GFSV;
74  typedef typename GFSV::Traits::GridViewType GridView;
75  static const int dim = GridView::dimension;
76  typedef typename GridView::Traits::Grid Grid;
77  typedef typename Matrix::block_type BlockType;
78  typedef typename GridView::template Codim<dim>::Iterator VertexIterator;
79  typedef typename Grid::Traits::GlobalIdSet IdSet;
80  typedef typename IdSet::IdType IdType;
81 
84  typename GFSV::Ordering::Traits::DOFIndex::value_type,
85  GFSV::Ordering::Traits::DOFIndex::max_depth,
86  typename GFSV::Traits::GridView::Grid::GlobalIdSet::IdType
88 
89  public:
91  typedef std::unordered_map<
92  typename GFSV::Ordering::Traits::DOFIndex,
93  std::unordered_set<GlobalDOFIndex>
95 
96  private:
97  typedef typename GFSV::Ordering::Traits::DOFIndex RowDOFIndex;
98  typedef typename GFSU::Ordering::Traits::DOFIndex ColDOFIndex;
99 
100  typedef std::pair<
101  typename RowDOFIndex::TreeIndex,
102  typename BorderPattern::mapped_type::value_type
103  > PatternMPIData;
104 
105  typedef std::tuple<
106  typename RowDOFIndex::TreeIndex,
107  typename BorderPattern::mapped_type::value_type,
108  typename M::field_type
109  > ValueMPIData;
110 
111  public:
119  : _communication_cache(std::make_shared<CommunicationCache>(grid_operator))
120  , _grid_view(grid_operator.testGridFunctionSpace().gridView())
121  {}
122 
123  void update(const GridOperator& grid_operator)
124  {
125  _communication_cache = std::make_shared<CommunicationCache>(grid_operator);
126  }
127 
129  : public BorderIndexIdCache<GFSV>
130  {
131 
134 
135  public:
136 
138  : BaseT(go.testGridFunctionSpace())
139  , _gfsu(go.trialGridFunctionSpace())
140  , _initialized(false)
141  , _entity_cache(go.testGridFunctionSpace())
142  {}
143 
144  typedef IdType EntityID;
145  typedef typename GFSU::Ordering::Traits::DOFIndex::TreeIndex ColumnTreeIndex;
146  typedef std::size_t size_type;
147 
148  bool initialized() const
149  {
150  return _initialized;
151  }
152 
154  {
155  _initialized = true;
156  }
157 
158  void update()
159  {
160  BaseT::update();
161  _border_pattern.clear();
162  _initialized = false;
163  }
164 
165 
166  const BorderPattern& pattern() const
167  {
168  assert(initialized());
169  return _border_pattern;
170  }
171 
172  template<typename LFSVCache, typename LFSUCache, typename LocalPattern>
173  void addEntries(const LFSVCache& lfsv_cache, const LFSUCache& lfsu_cache, const LocalPattern& pattern)
174  {
175  assert(!initialized());
176 
177  for (typename LocalPattern::const_iterator it = pattern.begin(),
178  end_it = pattern.end();
179  it != end_it;
180  ++it)
181  {
182  // skip constrained entries for now. TODO: Is this correct??
183  if (lfsv_cache.isConstrained(it->i()) || lfsu_cache.isConstrained(it->j()))
184  continue;
185 
186  const typename LFSVCache::DOFIndex& di = lfsv_cache.dofIndex(it->i());
187  const typename LFSUCache::DOFIndex& dj = lfsu_cache.dofIndex(it->j());
188 
189  size_type row_gt_index = GFSV::Ordering::Traits::DOFIndexAccessor::geometryType(di);
190  size_type row_entity_index = GFSV::Ordering::Traits::DOFIndexAccessor::entityIndex(di);
191 
192  size_type col_gt_index = GFSU::Ordering::Traits::DOFIndexAccessor::geometryType(dj);
193  size_type col_entity_index = GFSU::Ordering::Traits::DOFIndexAccessor::entityIndex(dj);
194 
195  // We are only interested in connections between two border entities.
196  if (!this->isBorderEntity(row_gt_index,row_entity_index) ||
197  !this->isBorderEntity(col_gt_index,col_entity_index))
198  continue;
199 
200  _border_pattern[di].insert(GlobalDOFIndex(this->id(col_gt_index,col_entity_index),dj.treeIndex()));
201  }
202  }
203 
204 
205  template<typename Entity>
206  size_type size(const Entity& e) const
207  {
208  _entity_cache.update(e);
209  size_type n = 0;
210  for (size_type i = 0; i < _entity_cache.size(); ++i)
211  {
212  typename BorderPattern::const_iterator it = _border_pattern.find(_entity_cache.dofIndex(i));
213  if (!transfer_dof(i,it))
214  continue;
215  n += it->second.size();
216  }
217 
218  return n;
219  }
220 
221  template<typename Buffer, typename Entity>
222  void gather_pattern(Buffer& buf, const Entity& e) const
223  {
224  _entity_cache.update(e);
225  for (size_type i = 0; i < _entity_cache.size(); ++i)
226  {
227  typename BorderPattern::const_iterator it = _border_pattern.find(_entity_cache.dofIndex(i));
228  if (!transfer_dof(i,it))
229  continue;
230  for (typename BorderPattern::mapped_type::const_iterator col_it = it->second.begin(),
231  col_end = it->second.end();
232  col_it != col_end;
233  ++col_it)
234  buf.write(std::make_pair(_entity_cache.dofIndex(i).treeIndex(),*col_it));
235  }
236  }
237 
238  template<typename Buffer, typename Entity>
239  void gather_data(Buffer& buf, const Entity& e, const M& matrix) const
240  {
241  _entity_cache.update(e);
242  for (size_type i = 0; i < _entity_cache.size(); ++i)
243  {
244  typename BorderPattern::const_iterator it = _border_pattern.find(_entity_cache.dofIndex(i));
245  if (!transfer_dof(i,it))
246  continue;
247  for (typename BorderPattern::mapped_type::const_iterator col_it = it->second.begin(),
248  col_end = it->second.end();
249  col_it != col_end;
250  ++col_it)
251  {
252  typename BaseT::EntityIndex col_entity = this->index(col_it->entityID());
253 
254  ColDOFIndex dj;
255  GFSU::Ordering::Traits::DOFIndexAccessor::store(dj,col_entity.geometryTypeIndex(),col_entity.entityIndex(),col_it->treeIndex());
256  buf.write(make_tuple(_entity_cache.dofIndex(i).treeIndex(),*col_it,matrix(_entity_cache.containerIndex(i),_gfsu.ordering().mapIndex(dj))));
257  }
258  }
259  }
260 
261  private:
262 
263  bool transfer_dof(size_type i, typename BorderPattern::const_iterator it) const
264  {
265  // not a border DOF
266  if (it == _border_pattern.end())
267  return false;
268  else
269  return true;
270 
271  /* Constraints check moved to addEntry()
272  // check for constraint
273  typename GridOperator::Traits::TrialGridFunctionSpaceConstraints::const_iterator cit = _constraints->find(_entity_cache.containerIndex(i));
274 
275  // transfer if DOF is not constrained
276  // TODO: What about non-Dirichlet constraints??
277  return cit == _constraints->end();
278  */
279  }
280 
281  const GFSU& _gfsu;
282  BorderPattern _border_pattern;
283  bool _initialized;
284  mutable EntityIndexCache<GFSV,true> _entity_cache;
285 
286  };
287 
289  template<typename Pattern>
291  : public CommDataHandleIF<PatternExtender<Pattern>, PatternMPIData>
292  {
293 
294  typedef std::size_t size_type;
295 
296  public:
298  typedef PatternMPIData DataType;
299 
300  bool contains (int dim, int codim) const
301  {
302  return
303  codim > 0 &&
304  (_gfsu.dataHandleContains(codim) ||
305  _gfsv.dataHandleContains(codim));
306  }
307 
308  bool fixedsize (int dim, int codim) const
309  {
310  return false;
311  }
312 
315  template<typename Entity>
316  size_type size (Entity& e) const
317  {
318  if (Entity::codimension == 0)
319  return 0;
320 
321  return _communication_cache.size(e);
322  }
323 
326  template<typename MessageBuffer, typename Entity>
327  void gather (MessageBuffer& buff, const Entity& e) const
328  {
329  if (Entity::codimension == 0)
330  return;
331 
332  _communication_cache.gather_pattern(buff,e);
333  }
334 
337  template<typename MessageBuffer, typename Entity>
338  void scatter (MessageBuffer& buff, const Entity& e, size_t n)
339  {
340  if (Entity::codimension == 0)
341  return;
342 
343  for (size_type i = 0; i < n; ++i)
344  {
345  DataType data;
346  buff.read(data);
347 
348  std::pair<bool,typename CommunicationCache::EntityIndex> col_index = _communication_cache.findIndex(data.second.entityID());
349  if (!col_index.first)
350  continue;
351 
352  RowDOFIndex di;
353  GFSV::Ordering::Traits::DOFIndexAccessor::store(di,
354  e.type(),
355  _grid_view.indexSet().index(e),
356  data.first);
357 
358  ColDOFIndex dj;
359  GFSU::Ordering::Traits::DOFIndexAccessor::store(dj,
360  col_index.second.geometryTypeIndex(),
361  col_index.second.entityIndex(),
362  data.second.treeIndex());
363 
364  _pattern.add_link(_gfsv.ordering().mapIndex(di),_gfsu.ordering().mapIndex(dj));
365  }
366  }
367 
369  const GFSU& gfsu,
370  const GFSV& gfsv,
371  Pattern& pattern)
372  : _communication_cache(dof_exchanger.communicationCache())
373  , _grid_view(dof_exchanger.gridView())
374  , _gfsu(gfsu)
375  , _gfsv(gfsv)
376  , _pattern(pattern)
377  {}
378 
379  private:
380 
381  const CommunicationCache& _communication_cache;
382  GridView _grid_view;
383  const GFSU& _gfsu;
384  const GFSV& _gfsv;
385  Pattern& _pattern;
386 
387  };
388 
391  : public CommDataHandleIF<EntryAccumulator,ValueMPIData>
392  {
393 
394  typedef std::size_t size_type;
395 
396  public:
398  typedef ValueMPIData DataType;
399 
400  bool contains(int dim, int codim) const
401  {
402  return
403  codim > 0 &&
404  (_gfsu.dataHandleContains(codim) ||
405  _gfsv.dataHandleContains(codim));
406  }
407 
408  bool fixedsize(int dim, int codim) const
409  {
410  return false;
411  }
412 
413  template<typename Entity>
414  size_type size(Entity& e) const
415  {
416  if (Entity::codimension == 0)
417  return 0;
418 
419  return _communication_cache.size(e);
420  }
421 
422  template<typename MessageBuffer, typename Entity>
423  void gather(MessageBuffer& buff, const Entity& e) const
424  {
425  if (Entity::codimension == 0)
426  return;
427 
428  _communication_cache.gather_data(buff,e,_matrix);
429  }
430 
433  template<typename MessageBuffer, typename Entity>
434  void scatter(MessageBuffer& buff, const Entity& e, size_type n)
435  {
436  if (Entity::codimension == 0)
437  return;
438 
439  for (size_type i = 0; i < n; ++i)
440  {
441  DataType data;
442  buff.read(data);
443 
444  std::pair<bool,typename CommunicationCache::EntityIndex> col_index = _communication_cache.findIndex(get<1>(data).entityID());
445  if (!col_index.first)
446  continue;
447 
448  RowDOFIndex di;
449  GFSV::Ordering::Traits::DOFIndexAccessor::store(di,
450  e.type(),
451  _grid_view.indexSet().index(e),
452  get<0>(data));
453 
454  ColDOFIndex dj;
455  GFSU::Ordering::Traits::DOFIndexAccessor::store(dj,
456  col_index.second.geometryTypeIndex(),
457  col_index.second.entityIndex(),
458  get<1>(data).treeIndex());
459 
460  _matrix(_gfsv.ordering().mapIndex(di),_gfsu.ordering().mapIndex(dj)) += get<2>(data);
461  }
462  }
463 
464 
466  const GFSU& gfsu,
467  const GFSV& gfsv,
468  Matrix& matrix)
469  : _communication_cache(dof_exchanger.communicationCache())
470  , _grid_view(dof_exchanger.gridView())
471  , _gfsu(gfsu)
472  , _gfsv(gfsv)
473  , _matrix(matrix)
474  {}
475 
476  private:
477 
478  const CommunicationCache& _communication_cache;
479  GridView _grid_view;
480  const GFSU& _gfsu;
481  const GFSV& _gfsv;
482  Matrix& _matrix;
483 
484  };
485 
492  void accumulateBorderEntries(const GridOperator& grid_operator, Matrix& matrix)
493  {
494  if (_grid_view.comm().size() > 1)
495  {
496  EntryAccumulator data_handle(*this,
497  grid_operator.testGridFunctionSpace(),
498  grid_operator.trialGridFunctionSpace(),
499  matrix);
500  _grid_view.communicate(data_handle,
501  InteriorBorder_InteriorBorder_Interface,
502  ForwardCommunication);
503  }
504  }
505 
507  {
508  return *_communication_cache;
509  }
510 
512  {
513  return *_communication_cache;
514  }
515 
516  shared_ptr<CommunicationCache> communicationCacheStorage()
517  {
518  return _communication_cache;
519  }
520 
521  const GridView& gridView() const
522  {
523  return _grid_view;
524  }
525 
526  private:
527 
528  shared_ptr<CommunicationCache> _communication_cache;
529  GridView _grid_view;
530 
531  };
532 
533 
534  template<typename GridOperator>
536  {
537 
538  public:
539 
541 
543  typedef Empty BorderPattern;
544 
546  {}
547 
549  {}
550 
551  void accumulateBorderEntries(const GridOperator& grid_operator, typename GridOperator::Traits::Jacobian& matrix)
552  {}
553 
554  CommunicationCache& communicationCache()
555  {
556  return *this;
557  }
558 
559  const CommunicationCache& communicationCache() const
560  {
561  return *this;
562  }
563 
564  void update(const GridOperator& grid_operator)
565  {}
566 
567  };
568 
569 
570  template<typename GridOperator>
572  public NoDataBorderDOFExchanger<GridOperator>
573  {
574 
575  public:
576 
578  {}
579 
581  {}
582 
583  };
584 
585 
586  } // namespace PDELab
587 } // namespace Dune
588 
589 #endif // DUNE_PDELAB_GRIDOPERATOR_COMMON_BORDERDOFEXCHANGER_HH
CommunicationCache & communicationCache()
Definition: borderdofexchanger.hh:506
void update(const GridOperator &grid_operator)
Definition: borderdofexchanger.hh:123
NoDataBorderDOFExchanger CommunicationCache
Definition: borderdofexchanger.hh:540
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:87
size_type size(Entity &e) const
How many objects of type DataType have to be sent for a given entity.
Definition: borderdofexchanger.hh:316
void update()
Definition: borderdofexchanger.hh:158
std::pair< bool, EntityIndex > findIndex(id_type entity_id) const
Definition: borderindexidcache.hh:135
size_type size(const Entity &e) const
Definition: borderdofexchanger.hh:206
EntityIndex index(id_type entity_id) const
Definition: borderindexidcache.hh:125
A DataHandle class to exchange matrix entries.
Definition: borderdofexchanger.hh:390
A DataHandle class to exchange matrix sparsity patterns.
Definition: borderdofexchanger.hh:290
void gather_pattern(Buffer &buf, const Entity &e) const
Definition: borderdofexchanger.hh:222
Definition: borderindexidcache.hh:25
PatternExtender(const NonOverlappingBorderDOFExchanger &dof_exchanger, const GFSU &gfsu, const GFSV &gfsv, Pattern &pattern)
Definition: borderdofexchanger.hh:368
JF JacobianField
The field type of the jacobian.
Definition: gridoperatorutilities.hh:69
void update()
Definition: borderindexidcache.hh:88
const E & e
Definition: interpolate.hh:172
STL namespace.
size_type size(Entity &e) const
Definition: borderdofexchanger.hh:414
bool isBorderEntity(std::size_t gt_index, std::size_t entity_index) const
Definition: borderindexidcache.hh:110
Definition: globaldofindex.hh:14
Definition: borderdofexchanger.hh:571
Standard grid operator implementation.
Definition: gridoperator.hh:34
const DI & dofIndex(size_type i) const
Definition: entityindexcache.hh:58
PatternMPIData DataType
Export type of data for message buffer.
Definition: borderdofexchanger.hh:298
const BorderPattern & pattern() const
Definition: borderdofexchanger.hh:166
void update(const GridOperator &grid_operator)
Definition: borderdofexchanger.hh:564
bool initialized() const
Definition: borderdofexchanger.hh:148
NonOverlappingBorderDOFExchanger(const GridOperator &grid_operator)
Constructor. Sets up the local to global relations.
Definition: borderdofexchanger.hh:118
void accumulateBorderEntries(const GridOperator &grid_operator, typename GridOperator::Traits::Jacobian &matrix)
Definition: borderdofexchanger.hh:551
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:66
const GridView & gridView() const
Definition: borderdofexchanger.hh:521
void gather(MessageBuffer &buff, const Entity &e) const
Pack data from user to message buffer.
Definition: borderdofexchanger.hh:327
NoDataBorderDOFExchanger()
Definition: borderdofexchanger.hh:545
OverlappingBorderDOFExchanger(const GridOperator &grid_operator)
Definition: borderdofexchanger.hh:580
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:93
void finishInitialization()
Definition: borderdofexchanger.hh:153
IdType EntityID
Definition: borderdofexchanger.hh:144
void gather(MessageBuffer &buff, const Entity &e) const
Definition: borderdofexchanger.hh:423
const CI & containerIndex(size_type i) const
Definition: entityindexcache.hh:64
Definition: borderdofexchanger.hh:535
GFSU TrialGridFunctionSpace
The trial grid function space.
Definition: gridoperatorutilities.hh:37
shared_ptr< CommunicationCache > communicationCacheStorage()
Definition: borderdofexchanger.hh:516
void accumulateBorderEntries(const GridOperator &grid_operator, Matrix &matrix)
Sums up the entries corresponding to border vertices.
Definition: borderdofexchanger.hh:492
void gather_data(Buffer &buf, const Entity &e, const M &matrix) const
Definition: borderdofexchanger.hh:239
Definition: adaptivity.hh:27
void update(const Entity &e)
Definition: entityindexcache.hh:49
const CommunicationCache & communicationCache() const
Definition: borderdofexchanger.hh:559
Empty BorderPattern
Data structure for storing border-border matrix pattern entries in a communication-optimized form...
Definition: borderdofexchanger.hh:543
EntryAccumulator(const NonOverlappingBorderDOFExchanger &dof_exchanger, const GFSU &gfsu, const GFSV &gfsv, Matrix &matrix)
Definition: borderdofexchanger.hh:465
std::unordered_map< typename GFSV::Ordering::Traits::DOFIndex, std::unordered_set< GlobalDOFIndex > > BorderPattern
Data structure for storing border-border matrix pattern entries in a communication-optimized form...
Definition: borderdofexchanger.hh:94
bool contains(int dim, int codim) const
Definition: borderdofexchanger.hh:400
CommunicationCache & communicationCache()
Definition: borderdofexchanger.hh:554
void scatter(MessageBuffer &buff, const Entity &e, size_type n)
Unpack data from message buffer to user.
Definition: borderdofexchanger.hh:434
GFSU::Ordering::Traits::DOFIndex::TreeIndex ColumnTreeIndex
Definition: borderdofexchanger.hh:145
NoDataBorderDOFExchanger(const GridOperator &grid_operator)
Definition: borderdofexchanger.hh:548
bool contains(int dim, int codim) const
Definition: borderdofexchanger.hh:300
void addEntries(const LFSVCache &lfsv_cache, const LFSUCache &lfsu_cache, const LocalPattern &pattern)
Definition: borderdofexchanger.hh:173
OverlappingBorderDOFExchanger()
Definition: borderdofexchanger.hh:577
GFSV TestGridFunctionSpace
The test grid function space.
Definition: gridoperatorutilities.hh:40
const CommunicationCache & communicationCache() const
Definition: borderdofexchanger.hh:511
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
void scatter(MessageBuffer &buff, const Entity &e, size_t n)
Unpack data from message buffer to user.
Definition: borderdofexchanger.hh:338
size_type size() const
Definition: entityindexcache.hh:74
CommunicationCache(const GridOperator &go)
Definition: borderdofexchanger.hh:137
ValueMPIData DataType
Export type of data for message buffer.
Definition: borderdofexchanger.hh:398
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
bool fixedsize(int dim, int codim) const
Definition: borderdofexchanger.hh:408
bool fixedsize(int dim, int codim) const
Definition: borderdofexchanger.hh:308
std::size_t size_type
Definition: borderdofexchanger.hh:146