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