dune-pdelab  2.5-dev
genericdatahandle.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_GRIDFUNCTIONSPACE_GENERICDATAHANDLE_HH
4 #define DUNE_PDELAB_GRIDFUNCTIONSPACE_GENERICDATAHANDLE_HH
5 
6 #include <vector>
7 #include <set>
8 #include <limits>
9 
10 #include<dune/common/exceptions.hh>
11 #include <dune/common/parallel/mpihelper.hh>
12 
13 #include <dune/grid/common/datahandleif.hh>
14 #include <dune/grid/common/gridenums.hh>
15 
18 
19 namespace Dune {
20  namespace PDELab {
21 
23  template<typename E>
25  {
26 
27  typedef char DataType;
28 
30  typedef std::size_t size_type;
31 
32  // Wrap the grid's communication buffer to enable sending leaf ordering sizes along with the data
33  static const bool wrap_buffer = true;
34 
35  // export original data type to fix up size information forwarded to standard gather / scatter functors
36  typedef E OriginalDataType;
37 
38  template<typename GFS>
39  bool contains(const GFS& gfs, int dim, int codim) const
40  {
41  return gfs.dataHandleContains(codim);
42  }
43 
44  template<typename GFS>
45  bool fixedSize(const GFS& gfs, int dim, int codim) const
46  {
47  return gfs.dataHandleFixedSize(codim);
48  }
49 
50  template<typename GFS, typename Entity>
51  std::size_t size(const GFS& gfs, const Entity& e) const
52  {
53  // include size of leaf ordering offsets if necessary
54  return gfs.dataHandleSize(e) * sizeof(E) + (gfs.sendLeafSizes() ? TypeTree::TreeInfo<typename GFS::Ordering>::leafCount * sizeof(size_type) : 0);
55  }
56 
57  };
58 
60  template<typename E>
62  {
63 
64  typedef E DataType;
65 
66  // Data is per entity, so we don't need to send leaf ordering size and thus can avoid wrapping the
67  // grid's communication buffer
68  static const bool wrap_buffer = false;
69 
70  template<typename GFS>
71  bool contains(const GFS& gfs, int dim, int codim) const
72  {
73  return gfs.dataHandleContains(codim);
74  }
75 
76  template<typename GFS>
77  bool fixedSize(const GFS& gfs, int dim, int codim) const
78  {
79  return true;
80  }
81 
82  template<typename GFS, typename Entity>
83  std::size_t size(const GFS& gfs, const Entity& e) const
84  {
85  return gfs.dataHandleContains(Entity::codimension) && gfs.entitySet().contains(e) ? _count : 0;
86  }
87 
88  explicit EntityDataCommunicationDescriptor(std::size_t count = 1)
89  : _count(count)
90  {}
91 
92  private:
93 
94  const std::size_t _count;
95 
96  };
97 
99 
105  template<typename GFS, typename V, typename GatherScatter, typename CommunicationDescriptor = DOFDataCommunicationDescriptor<typename V::ElementType> >
107  : public Dune::CommDataHandleIF<GFSDataHandle<GFS,V,GatherScatter,CommunicationDescriptor>,typename CommunicationDescriptor::DataType>
108  {
109 
110  public:
111 
112  typedef typename CommunicationDescriptor::DataType DataType;
113  typedef typename GFS::Traits::SizeType size_type;
114 
115  static const size_type leaf_count = TypeTree::TreeInfo<typename GFS::Ordering>::leafCount;
116 
117  GFSDataHandle(const GFS& gfs, V& v, GatherScatter gather_scatter = GatherScatter(), CommunicationDescriptor communication_descriptor = CommunicationDescriptor())
118  : _gfs(gfs)
119  , _index_cache(gfs)
120  , _local_view(v)
121  , _gather_scatter(gather_scatter)
122  , _communication_descriptor(communication_descriptor)
123  {}
124 
126  bool contains(int dim, int codim) const
127  {
128  return _communication_descriptor.contains(_gfs,dim,codim);
129  }
130 
132  bool fixedsize(int dim, int codim) const
133  {
134  return _communication_descriptor.fixedSize(_gfs,dim,codim);
135  }
136 
141  template<typename Entity>
142  size_type size(const Entity& e) const
143  {
144  return _communication_descriptor.size(_gfs,e);
145  }
146 
148  template<typename MessageBuffer, typename Entity>
149  typename std::enable_if<
150  CommunicationDescriptor::wrap_buffer && AlwaysTrue<Entity>::value // we can only support this if the buffer is wrapped
151  >::type
152  gather(MessageBuffer& buff, const Entity& e) const
153  {
154  PolymorphicBufferWrapper<MessageBuffer> buf_wrapper(buff);
155  _index_cache.update(e);
156  _local_view.bind(_index_cache);
157  if (_gfs.sendLeafSizes())
158  {
159  // send the leaf ordering offsets as exported by the EntityIndexCache
160  for (auto it = _index_cache.offsets().begin() + 1,
161  end_it = _index_cache.offsets().end();
162  it != end_it;
163  ++it)
164  {
165  buf_wrapper.write(static_cast<typename CommunicationDescriptor::size_type>(*it));
166  }
167  }
168  // send the normal data
169  if (_gather_scatter.gather(buf_wrapper,e,_local_view))
170  _local_view.commit();
171  _local_view.unbind();
172  }
173 
175  template<typename MessageBuffer, typename Entity>
176  typename std::enable_if<
177  !CommunicationDescriptor::wrap_buffer && AlwaysTrue<Entity>::value
178  >::type
179  gather(MessageBuffer& buff, const Entity& e) const
180  {
181  _index_cache.update(e);
182  _local_view.bind(_index_cache);
183  if (_gather_scatter.gather(buff,e,_local_view))
184  _local_view.commit();
185  _local_view.unbind();
186  }
187 
194  template<typename MessageBuffer, typename Entity>
195  typename std::enable_if<
196  CommunicationDescriptor::wrap_buffer && AlwaysTrue<Entity>::value // we require the buffer to be wrapped
197  >::type
198  scatter(MessageBuffer& buff, const Entity& e, size_type n)
199  {
200  PolymorphicBufferWrapper<MessageBuffer> buf_wrapper(buff);
201  _index_cache.update(e);
202  _local_view.bind(_index_cache);
203  bool needs_commit = false;
204  if (_gfs.sendLeafSizes())
205  {
206  // receive leaf ordering offsets and store in local array
207  typename IndexCache::Offsets remote_offsets = {{0}};
208  for (auto it = remote_offsets.begin() + 1,
209  end_it = remote_offsets.end();
210  it != end_it;
211  ++it)
212  {
213  typename CommunicationDescriptor::size_type data = 0;
214  buf_wrapper.read(data);
215  *it = data;
216  }
217  // call special version of scatter() that can handle empty leafs in the ordering tree
218  needs_commit = _gather_scatter.scatter(buf_wrapper,remote_offsets,_index_cache.offsets(),e,_local_view);
219  }
220  else
221  {
222  // call standard version of scatter - make sure to fix the reported communication size
223  needs_commit = _gather_scatter.scatter(buf_wrapper,n / sizeof(typename CommunicationDescriptor::OriginalDataType),e,_local_view);
224  }
225 
226  if (needs_commit)
227  _local_view.commit();
228 
229  _local_view.unbind();
230  }
231 
238  template<typename MessageBuffer, typename Entity>
239  typename std::enable_if<
240  !CommunicationDescriptor::wrap_buffer && AlwaysTrue<Entity>::value
241  >::type
242  scatter(MessageBuffer& buff, const Entity& e, size_type n)
243  {
244  _index_cache.update(e);
245  _local_view.bind(_index_cache);
246 
247  if (_gather_scatter.scatter(buff,n,e,_local_view))
248  _local_view.commit();
249 
250  _local_view.unbind();
251  }
252 
253 
254  private:
255 
257  typedef typename V::template LocalView<IndexCache> LocalView;
258 
259  const GFS& _gfs;
260  mutable IndexCache _index_cache;
261  mutable LocalView _local_view;
262  mutable GatherScatter _gather_scatter;
263  CommunicationDescriptor _communication_descriptor;
264 
265  };
266 
267 
268  template<typename GatherScatter>
270  {
271 
272  public:
273 
274  typedef std::size_t size_type;
275 
276  template<typename MessageBuffer, typename Entity, typename LocalView>
277  bool gather(MessageBuffer& buff, const Entity& e, const LocalView& local_view) const
278  {
279  for (std::size_t i = 0; i < local_view.size(); ++i)
280  _gather_scatter.gather(buff,local_view[i]);
281  return false;
282  }
283 
284  // default scatter - requires function space structure to be identical on sender and receiver side
285  template<typename MessageBuffer, typename Entity, typename LocalView>
286  bool scatter(MessageBuffer& buff, size_type n, const Entity& e, LocalView& local_view) const
287  {
288  if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
289  {
290  if (local_view.size() != n)
291  DUNE_THROW(Exception,"size mismatch in GridFunctionSpace data handle, have " << local_view.size() << "DOFs, but received " << n);
292 
293  for (std::size_t i = 0; i < local_view.size(); ++i)
294  _gather_scatter.scatter(buff,local_view[i]);
295  return true;
296  }
297  else
298  {
299  if (local_view.size() != 0)
300  DUNE_THROW(Exception,"expected no DOFs in partition '" << e.partitionType() << "', but have " << local_view.size());
301 
302  for (std::size_t i = 0; i < local_view.size(); ++i)
303  {
304  typename LocalView::ElementType dummy;
305  buff.read(dummy);
306  }
307  return false;
308  }
309  }
310 
311  // enhanced scatter with support for function spaces with different structure on sender and receiver side
312  template<typename MessageBuffer, typename Offsets, typename Entity, typename LocalView>
313  bool scatter(MessageBuffer& buff, const Offsets& remote_offsets, const Offsets& local_offsets, const Entity& e, LocalView& local_view) const
314  {
315  if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
316  {
317  // the idea here is this:
318  // the compile time structure of the overall function space (and its ordering) will be identical on both sides
319  // of the communication, but one side may be missing information on some leaf orderings, e.g. because the DOF
320  // belongs to a MultiDomain subdomain that only has an active grid cell on one side of the communication.
321  // So we step through the leaves and simply ignore any block where one of the two sides is of size 0.
322  // Otherwise, it's business as usual: we make sure that the sizes from both sides match and then process all
323  // data with the DOF-local gather / scatter functor.
324  size_type remote_i = 0;
325  size_type local_i = 0;
326  bool needs_commit = false;
327  for (size_type block = 1; block < local_offsets.size(); ++block)
328  {
329  // we didn't get any data - just ignore
330  if (remote_offsets[block] == remote_i)
331  {
332  local_i = local_offsets[block];
333  continue;
334  }
335 
336  // we got data for DOFs we don't have - drain buffer
337  if (local_offsets[block] == local_i)
338  {
339  for (; remote_i < remote_offsets[block]; ++remote_i)
340  {
341  typename LocalView::ElementType dummy;
342  buff.read(dummy);
343  }
344  continue;
345  }
346 
347  if (remote_offsets[block] - remote_i != local_offsets[block] - local_i)
348  DUNE_THROW(Exception,"size mismatch in GridFunctionSpace data handle block " << block << ", have " << local_offsets[block] - local_i << "DOFs, but received " << remote_offsets[block] - remote_i);
349 
350  for (; local_i < local_offsets[block]; ++local_i)
351  _gather_scatter.scatter(buff,local_view[local_i]);
352 
353  remote_i = remote_offsets[block];
354  needs_commit = true;
355  }
356  return needs_commit;
357  }
358  else
359  {
360  if (local_view.size() != 0)
361  DUNE_THROW(Exception,"expected no DOFs in partition '" << e.partitionType() << "', but have " << local_view.size());
362 
363  for (std::size_t i = 0; i < remote_offsets.back(); ++i)
364  {
365  typename LocalView::ElementType dummy;
366  buff.read(dummy);
367  }
368  return false;
369  }
370  }
371 
372  DataGatherScatter(GatherScatter gather_scatter = GatherScatter())
373  : _gather_scatter(gather_scatter)
374  {}
375 
376  private:
377 
378  GatherScatter _gather_scatter;
379 
380  };
381 
382 
383  template<typename GatherScatter>
385  {
386 
387  public:
388 
389  typedef std::size_t size_type;
390 
391  template<typename MessageBuffer, typename Entity, typename LocalView>
392  bool gather(MessageBuffer& buff, const Entity& e, const LocalView& local_view) const
393  {
394  for (std::size_t i = 0; i < local_view.size(); ++i)
395  _gather_scatter.gather(buff,e,local_view[i]);
396  return false;
397  }
398 
399  // see documentation in DataGatherScatter for further info on the scatter() implementations
400  template<typename MessageBuffer, typename Entity, typename LocalView>
401  bool scatter(MessageBuffer& buff, size_type n, const Entity& e, LocalView& local_view) const
402  {
403  if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
404  {
405  if (local_view.size() != n)
406  DUNE_THROW(Exception,"size mismatch in GridFunctionSpace data handle, have " << local_view.size() << "DOFs, but received " << n);
407 
408  for (std::size_t i = 0; i < local_view.size(); ++i)
409  _gather_scatter.scatter(buff,e,local_view[i]);
410  return true;
411  }
412  else
413  {
414  if (local_view.size() != 0)
415  DUNE_THROW(Exception,"expected no DOFs in partition '" << e.partitionType() << "', but have " << local_view.size());
416 
417  for (std::size_t i = 0; i < local_view.size(); ++i)
418  {
419  typename LocalView::ElementType dummy;
420  buff.read(dummy);
421  }
422  return false;
423  }
424  }
425 
426  // see documentation in DataGatherScatter for further info on the scatter() implementations
427  template<typename MessageBuffer, typename Offsets, typename Entity, typename LocalView>
428  bool scatter(MessageBuffer& buff, const Offsets& remote_offsets, const Offsets& local_offsets, const Entity& e, LocalView& local_view) const
429  {
430  if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
431  {
432  size_type remote_i = 0;
433  size_type local_i = 0;
434  bool needs_commit = false;
435  for (size_type block = 1; block < local_offsets.size(); ++block)
436  {
437 
438  // we didn't get any data - just ignore
439  if (remote_offsets[block] == remote_i)
440  {
441  local_i = local_offsets[block];
442  continue;
443  }
444 
445  // we got data for DOFs we don't have - drain buffer
446  if (local_offsets[block] == local_i)
447  {
448  for (; remote_i < remote_offsets[block]; ++remote_i)
449  {
450  typename LocalView::ElementType dummy;
451  buff.read(dummy);
452  }
453  continue;
454  }
455 
456  if (remote_offsets[block] - remote_i != local_offsets[block] - local_i)
457  DUNE_THROW(Exception,"size mismatch in GridFunctionSpace data handle block " << block << ", have " << local_offsets[block] - local_i << "DOFs, but received " << remote_offsets[block] - remote_i);
458 
459  for (; local_i < local_offsets[block]; ++local_i)
460  _gather_scatter.scatter(buff,e,local_view[local_i]);
461 
462  remote_i = remote_offsets[block];
463  needs_commit = true;
464  }
465  return needs_commit;
466  }
467  else
468  {
469  if (local_view.size() != 0)
470  DUNE_THROW(Exception,"expected no DOFs in partition '" << e.partitionType() << "', but have " << local_view.size());
471 
472  for (std::size_t i = 0; i < remote_offsets.back(); ++i)
473  {
474  typename LocalView::ElementType dummy;
475  buff.read(dummy);
476  }
477  return false;
478  }
479  }
480 
481  DataEntityGatherScatter(GatherScatter gather_scatter = GatherScatter())
482  : _gather_scatter(gather_scatter)
483  {}
484 
485  private:
486 
487  GatherScatter _gather_scatter;
488 
489  };
490 
491 
492  template<typename GatherScatter>
494  {
495 
496  public:
497 
498  typedef std::size_t size_type;
499 
500  template<typename MessageBuffer, typename Entity, typename LocalView>
501  bool gather(MessageBuffer& buff, const Entity& e, const LocalView& local_view) const
502  {
503  for (std::size_t i = 0; i < local_view.size(); ++i)
504  _gather_scatter.gather(buff,local_view.cache().containerIndex(i),local_view[i]);
505  return false;
506  }
507 
508  // see documentation in DataGatherScatter for further info on the scatter() implementations
509  template<typename MessageBuffer, typename Entity, typename LocalView>
510  bool scatter(MessageBuffer& buff, size_type n, const Entity& e, LocalView& local_view) const
511  {
512  if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
513  {
514  if (local_view.size() != n)
515  DUNE_THROW(Exception,"size mismatch in GridFunctionSpace data handle, have " << local_view.size() << "DOFs, but received " << n);
516 
517  for (std::size_t i = 0; i < local_view.size(); ++i)
518  _gather_scatter.scatter(buff,local_view.cache().containerIndex(i),local_view[i]);
519 
520  return true;
521  }
522  else
523  {
524  if (local_view.size() != 0)
525  DUNE_THROW(Exception,"expected no DOFs in partition '" << e.partitionType() << "', but have " << local_view.size());
526 
527  for (std::size_t i = 0; i < local_view.size(); ++i)
528  {
529  typename LocalView::ElementType dummy;
530  buff.read(dummy);
531  }
532  return false;
533  }
534  }
535 
536  // see documentation in DataGatherScatter for further info on the scatter() implementations
537  template<typename MessageBuffer, typename Offsets, typename Entity, typename LocalView>
538  bool scatter(MessageBuffer& buff, const Offsets& remote_offsets, const Offsets& local_offsets, const Entity& e, LocalView& local_view) const
539  {
540  if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
541  {
542  size_type remote_i = 0;
543  size_type local_i = 0;
544  bool needs_commit = false;
545  for (size_type block = 1; block < local_offsets.size(); ++block)
546  {
547 
548  // we didn't get any data - just ignore
549  if (remote_offsets[block] == remote_i)
550  {
551  local_i = local_offsets[block];
552  continue;
553  }
554 
555  // we got data for DOFs we don't have - drain buffer
556  if (local_offsets[block] == local_i)
557  {
558  for (; remote_i < remote_offsets[block]; ++remote_i)
559  {
560  typename LocalView::ElementType dummy;
561  buff.read(dummy);
562  }
563  continue;
564  }
565 
566  if (remote_offsets[block] - remote_i != local_offsets[block] - local_i)
567  DUNE_THROW(Exception,"size mismatch in GridFunctionSpace data handle block " << block << ", have " << local_offsets[block] - local_i << "DOFs, but received " << remote_offsets[block] - remote_i);
568 
569  for (; local_i < local_offsets[block]; ++local_i)
570  _gather_scatter.scatter(buff,local_view.cache().containerIndex(local_i),local_view[local_i]);
571 
572  remote_i = remote_offsets[block];
573  needs_commit = true;
574  }
575  return needs_commit;
576  }
577  else
578  {
579  if (local_view.size() != 0)
580  DUNE_THROW(Exception,"expected no DOFs in partition '" << e.partitionType() << "', but have " << local_view.size());
581 
582  for (std::size_t i = 0; i < remote_offsets.back(); ++i)
583  {
584  typename LocalView::ElementType dummy;
585  buff.read(dummy);
586  }
587  return false;
588  }
589  }
590 
591 
592  DataContainerIndexGatherScatter(GatherScatter gather_scatter = GatherScatter())
593  : _gather_scatter(gather_scatter)
594  {}
595 
596  private:
597 
598  GatherScatter _gather_scatter;
599 
600  };
601 
602 
604  {
605  public:
606  template<class MessageBuffer, class DataType>
607  void gather (MessageBuffer& buff, DataType& data) const
608  {
609  buff.write(data);
610  }
611 
612  template<class MessageBuffer, class DataType>
613  void scatter (MessageBuffer& buff, DataType& data) const
614  {
615  DataType x;
616  buff.read(x);
617  data += x;
618  }
619  };
620 
621  template<class GFS, class V>
623  : public GFSDataHandle<GFS,V,DataGatherScatter<AddGatherScatter> >
624  {
626 
627  public:
628 
629  AddDataHandle (const GFS& gfs_, V& v_)
630  : BaseT(gfs_,v_)
631  {}
632  };
633 
635  {
636  public:
637  template<class MessageBuffer, class DataType>
638  void gather (MessageBuffer& buff, DataType& data) const
639  {
640  buff.write(data);
641  data = (DataType) 0;
642  }
643 
644  template<class MessageBuffer, class DataType>
645  void scatter (MessageBuffer& buff, DataType& data) const
646  {
647  DataType x;
648  buff.read(x);
649  data += x;
650  }
651  };
652 
653  template<class GFS, class V>
655  : public GFSDataHandle<GFS,V,DataGatherScatter<AddClearGatherScatter> >
656  {
658 
659  public:
660 
661  AddClearDataHandle (const GFS& gfs_, V& v_)
662  : BaseT(gfs_,v_)
663  {}
664  };
665 
667  {
668  public:
669  template<class MessageBuffer, class DataType>
670  void gather (MessageBuffer& buff, DataType& data) const
671  {
672  buff.write(data);
673  }
674 
675  template<class MessageBuffer, class DataType>
676  void scatter (MessageBuffer& buff, DataType& data) const
677  {
678  DataType x;
679  buff.read(x);
680  data = x;
681  }
682  };
683 
684  template<class GFS, class V>
686  : public GFSDataHandle<GFS,V,DataGatherScatter<CopyGatherScatter> >
687  {
689 
690  public:
691 
692  CopyDataHandle (const GFS& gfs_, V& v_)
693  : BaseT(gfs_,v_)
694  {}
695  };
696 
698  {
699  public:
700  template<class MessageBuffer, class DataType>
701  void gather (MessageBuffer& buff, DataType& data) const
702  {
703  buff.write(data);
704  }
705 
706  template<class MessageBuffer, class DataType>
707  void scatter (MessageBuffer& buff, DataType& data) const
708  {
709  DataType x;
710  buff.read(x);
711  data = std::min(data,x);
712  }
713  };
714 
715  template<class GFS, class V>
717  : public GFSDataHandle<GFS,V,DataGatherScatter<MinGatherScatter> >
718  {
720 
721  public:
722 
723  MinDataHandle (const GFS& gfs_, V& v_)
724  : BaseT(gfs_,v_)
725  {}
726  };
727 
729  {
730  public:
731  template<class MessageBuffer, class DataType>
732  void gather (MessageBuffer& buff, DataType& data) const
733  {
734  buff.write(data);
735  }
736 
737  template<class MessageBuffer, class DataType>
738  void scatter (MessageBuffer& buff, DataType& data) const
739  {
740  DataType x;
741  buff.read(x);
742  data = std::max(data,x);
743  }
744  };
745 
746  template<class GFS, class V>
748  : public GFSDataHandle<GFS,V,DataGatherScatter<MaxGatherScatter> >
749  {
751 
752  public:
753 
754  MaxDataHandle (const GFS& gfs_, V& v_)
755  : BaseT(gfs_,v_)
756  {}
757  };
758 
759 
761 
769  {
770  public:
771 
772  template<typename MessageBuffer, typename Entity, typename LocalView>
773  bool gather(MessageBuffer& buff, const Entity& e, LocalView& local_view) const
774  {
775  // Figure out where we are...
776  const bool ghost = e.partitionType()!=Dune::InteriorEntity && e.partitionType()!=Dune::BorderEntity;
777 
778  // ... and send something (doesn't really matter what, we'll throw it away on the receiving side).
779  buff.write(ghost);
780 
781  return false;
782  }
783 
784  template<typename MessageBuffer, typename Entity, typename LocalView>
785  bool scatter(MessageBuffer& buff, std::size_t n, const Entity& e, LocalView& local_view) const
786  {
787  // Figure out where we are - we have to do this again on the receiving side due to the asymmetric
788  // communication interface!
789  const bool ghost = e.partitionType()!=Dune::InteriorEntity && e.partitionType()!=Dune::BorderEntity;
790 
791  // drain buffer
792  bool dummy;
793  buff.read(dummy);
794 
795  for (std::size_t i = 0; i < local_view.size(); ++i)
796  local_view[i] = ghost;
797 
798  return true;
799  }
800 
801  };
802 
804 
811  template<class GFS, class V>
813  : public Dune::PDELab::GFSDataHandle<GFS,
814  V,
815  GhostGatherScatter,
816  EntityDataCommunicationDescriptor<bool> >
817  {
819  GFS,
820  V,
823  > BaseT;
824 
826  "GhostDataHandle expects a vector of bool values");
827 
828  public:
829 
831 
840  GhostDataHandle(const GFS& gfs_, V& v_, bool init_vector = true)
841  : BaseT(gfs_,v_)
842  {
843  if (init_vector)
844  v_ = false;
845  }
846  };
847 
848 
850 
858  template<typename RankIndex>
860  {
861 
862  public:
863 
864  template<typename MessageBuffer, typename Entity, typename LocalView>
865  bool gather(MessageBuffer& buff, const Entity& e, LocalView& local_view) const
866  {
867  // We only gather from interior and border entities, so we can throw in our ownership
868  // claim without any further checks.
869  buff.write(_rank);
870 
871  return true;
872  }
873 
874  template<typename MessageBuffer, typename Entity, typename LocalView>
875  bool scatter(MessageBuffer& buff, std::size_t n, const Entity& e, LocalView& local_view) const
876  {
877  // Value used for DOFs with currently unknown rank.
878  const RankIndex unknown_rank = std::numeric_limits<RankIndex>::max();
879 
880  // We can only own this DOF if it is either on the interior or border partition.
881  const bool is_interior_or_border = (e.partitionType()==Dune::InteriorEntity || e.partitionType()==Dune::BorderEntity);
882 
883  // Receive data.
884  RankIndex received_rank;
885  buff.read(received_rank);
886 
887  for (std::size_t i = 0; i < local_view.size(); ++i)
888  {
889  // Get the currently stored owner rank for this DOF.
890  RankIndex current_rank = local_view[i];
891 
892  // We only gather from interior and border entities, so we need to make sure
893  // we relinquish any ownership claims on overlap and ghost entities on the
894  // receiving side. We also need to make sure not to overwrite any data already
895  // received, so we only blank the rank value if the currently stored value is
896  // equal to our own rank.
897  if (!is_interior_or_border && current_rank == _rank)
898  current_rank = unknown_rank;
899 
900  // Assign DOFs to minimum rank value.
901  local_view[i] = std::min(current_rank,received_rank);
902  }
903  return true;
904  }
905 
907 
911  : _rank(rank)
912  {}
913 
914  private:
915 
916  const RankIndex _rank;
917 
918  };
919 
921 
929  template<class GFS, class V>
931  : public Dune::PDELab::GFSDataHandle<GFS,
932  V,
933  DisjointPartitioningGatherScatter<
934  typename V::ElementType
935  >,
936  EntityDataCommunicationDescriptor<
937  typename V::ElementType
938  >
939  >
940  {
942  GFS,
943  V,
945  typename V::ElementType
946  >,
948  typename V::ElementType
949  >
950  > BaseT;
951 
952  public:
953 
955 
964  DisjointPartitioningDataHandle(const GFS& gfs_, V& v_, bool init_vector = true)
965  : BaseT(gfs_,v_,DisjointPartitioningGatherScatter<typename V::ElementType>(gfs_.gridView().comm().rank()))
966  {
967  if (init_vector)
968  v_ = gfs_.gridView().comm().rank();
969  }
970  };
971 
972 
974 
981  {
982 
983  template<typename MessageBuffer, typename Entity, typename LocalView>
984  bool gather(MessageBuffer& buff, const Entity& e, LocalView& local_view) const
985  {
986  buff.write(local_view.size() > 0);
987  return false;
988  }
989 
990  template<typename MessageBuffer, typename Entity, typename LocalView>
991  bool scatter(MessageBuffer& buff, std::size_t n, const Entity& e, LocalView& local_view) const
992  {
993  bool remote_entity_has_dofs;
994  buff.read(remote_entity_has_dofs);
995 
996  for (std::size_t i = 0; i < local_view.size(); ++i)
997  {
998  local_view[i] |= remote_entity_has_dofs;
999  }
1000  return true;
1001  }
1002 
1003  };
1004 
1005 
1007 
1013  template<class GFS, class V>
1015  : public Dune::PDELab::GFSDataHandle<GFS,
1016  V,
1017  SharedDOFGatherScatter,
1018  EntityDataCommunicationDescriptor<bool> >
1019  {
1021  GFS,
1022  V,
1025  > BaseT;
1026 
1028  "SharedDOFDataHandle expects a vector of bool values");
1029 
1030  public:
1031 
1033 
1042  SharedDOFDataHandle(const GFS& gfs_, V& v_, bool init_vector = true)
1043  : BaseT(gfs_,v_)
1044  {
1045  if (init_vector)
1046  v_ = false;
1047  }
1048  };
1049 
1050 
1052 
1059  template<typename GFS, typename RankIndex>
1061  : public Dune::CommDataHandleIF<GFSNeighborDataHandle<GFS,RankIndex>,RankIndex>
1062  {
1063 
1064  // We deliberately avoid using the GFSDataHandle here, as we don't want to incur the
1065  // overhead of invoking the whole GFS infrastructure.
1066 
1067  public:
1068 
1069  typedef RankIndex DataType;
1070  typedef typename GFS::Traits::SizeType size_type;
1071 
1072  GFSNeighborDataHandle(const GFS& gfs, RankIndex rank, std::set<RankIndex>& neighbors)
1073  : _gfs(gfs)
1074  , _rank(rank)
1075  , _neighbors(neighbors)
1076  {}
1077 
1078  bool contains(int dim, int codim) const
1079  {
1080  // Only create neighbor relations for codims used by the GFS.
1081  return _gfs.dataHandleContains(codim);
1082  }
1083 
1084  bool fixedsize(int dim, int codim) const
1085  {
1086  // We always send a single value, the MPI rank.
1087  return true;
1088  }
1089 
1090  template<typename Entity>
1091  size_type size(Entity& e) const
1092  {
1093  return 1;
1094  }
1095 
1096  template<typename MessageBuffer, typename Entity>
1097  void gather(MessageBuffer& buff, const Entity& e) const
1098  {
1099  buff.write(_rank);
1100  }
1101 
1102  template<typename MessageBuffer, typename Entity>
1103  void scatter(MessageBuffer& buff, const Entity& e, size_type n)
1104  {
1105  RankIndex rank;
1106  buff.read(rank);
1107  _neighbors.insert(rank);
1108  }
1109 
1110  private:
1111 
1112  const GFS& _gfs;
1113  const RankIndex _rank;
1114  std::set<RankIndex>& _neighbors;
1115 
1116  };
1117 
1118 
1119  } // namespace PDELab
1120 } // namespace Dune
1121 
1122 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_GENERICDATAHANDLE_HH
CopyDataHandle(const GFS &gfs_, V &v_)
Definition: genericdatahandle.hh:692
std::size_t size_type
Definition: genericdatahandle.hh:389
std::size_t size(const GFS &gfs, const Entity &e) const
Definition: genericdatahandle.hh:51
static const int dim
Definition: adaptivity.hh:83
E DataType
Definition: genericdatahandle.hh:64
Implement a data handle with a grid function space.
Definition: genericdatahandle.hh:106
GhostDataHandle(const GFS &gfs_, V &v_, bool init_vector=true)
Creates a new GhostDataHandle.
Definition: genericdatahandle.hh:840
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
Data handle for marking shared DOFs.
Definition: genericdatahandle.hh:1014
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
GFSNeighborDataHandle(const GFS &gfs, RankIndex rank, std::set< RankIndex > &neighbors)
Definition: genericdatahandle.hh:1072
bool scatter(MessageBuffer &buff, const Offsets &remote_offsets, const Offsets &local_offsets, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:313
void gather(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:732
E OriginalDataType
Definition: genericdatahandle.hh:36
bool scatter(MessageBuffer &buff, size_type n, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:401
Definition: genericdatahandle.hh:622
void scatter(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:613
AddClearDataHandle(const GFS &gfs_, V &v_)
Definition: genericdatahandle.hh:661
Definition: genericdatahandle.hh:384
std::enable_if< CommunicationDescriptor::wrap_buffer &&AlwaysTrue< Entity >::value >::type gather(MessageBuffer &buff, const Entity &e) const
pack data from user to message buffer - version with support for sending leaf ordering sizes ...
Definition: genericdatahandle.hh:152
GFS::Traits::SizeType size_type
Definition: genericdatahandle.hh:113
Communication descriptor for sending count items of type E per entity with attached DOFs...
Definition: genericdatahandle.hh:61
std::enable_if< !CommunicationDescriptor::wrap_buffer &&AlwaysTrue< Entity >::value >::type scatter(MessageBuffer &buff, const Entity &e, size_type n)
unpack data from message buffer to user
Definition: genericdatahandle.hh:242
bool contains(const GFS &gfs, int dim, int codim) const
Definition: genericdatahandle.hh:39
void scatter(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:707
RankIndex DataType
Definition: genericdatahandle.hh:1069
bool scatter(MessageBuffer &buff, const Offsets &remote_offsets, const Offsets &local_offsets, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:538
bool gather(MessageBuffer &buff, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:773
Definition: genericdatahandle.hh:697
DataGatherScatter(GatherScatter gather_scatter=GatherScatter())
Definition: genericdatahandle.hh:372
Definition: genericdatahandle.hh:654
std::size_t size(const GFS &gfs, const Entity &e) const
Definition: genericdatahandle.hh:83
bool gather(MessageBuffer &buff, const Entity &e, const LocalView &local_view) const
Definition: genericdatahandle.hh:277
void gather(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:607
Definition: genericdatahandle.hh:603
bool contains(int dim, int codim) const
returns true if data for this codim should be communicated
Definition: genericdatahandle.hh:126
void scatter(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:676
bool scatter(MessageBuffer &buff, size_type n, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:286
MinDataHandle(const GFS &gfs_, V &v_)
Definition: genericdatahandle.hh:723
std::enable_if< !CommunicationDescriptor::wrap_buffer &&AlwaysTrue< Entity >::value >::type gather(MessageBuffer &buff, const Entity &e) const
pack data from user to message buffer - version without support for sending leaf ordering sizes ...
Definition: genericdatahandle.hh:179
DataEntityGatherScatter(GatherScatter gather_scatter=GatherScatter())
Definition: genericdatahandle.hh:481
Wrapper for message buffers of grid DataHandles that allows for sending different types of data...
Definition: polymorphicbufferwrapper.hh:26
Definition: genericdatahandle.hh:728
bool gather(MessageBuffer &buff, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:984
void scatter(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:645
Definition: genericdatahandle.hh:685
bool fixedsize(int dim, int codim) const
returns true if size per entity of given dim and codim is a constant
Definition: genericdatahandle.hh:132
char DataType
Definition: genericdatahandle.hh:27
bool scatter(MessageBuffer &buff, std::size_t n, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:785
bool scatter(MessageBuffer &buff, const Offsets &remote_offsets, const Offsets &local_offsets, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:428
bool fixedSize(const GFS &gfs, int dim, int codim) const
Definition: genericdatahandle.hh:77
SharedDOFDataHandle(const GFS &gfs_, V &v_, bool init_vector=true)
Creates a new SharedDOFDataHandle.
Definition: genericdatahandle.hh:1042
DisjointPartitioningDataHandle(const GFS &gfs_, V &v_, bool init_vector=true)
Creates a new DisjointPartitioningDataHandle.
Definition: genericdatahandle.hh:964
void gather(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:701
size_type size(const Entity &e) const
how many objects of type DataType have to be sent for a given entity
Definition: genericdatahandle.hh:142
void scatter(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:738
void read(T &data)
Definition: polymorphicbufferwrapper.hh:40
GatherScatter functor for creating a disjoint DOF partitioning.
Definition: genericdatahandle.hh:859
CommunicationDescriptor::DataType DataType
Definition: genericdatahandle.hh:112
GFS::Traits::SizeType size_type
Definition: genericdatahandle.hh:1070
bool gather(MessageBuffer &buff, const Entity &e, const LocalView &local_view) const
Definition: genericdatahandle.hh:501
Definition: genericdatahandle.hh:747
void scatter(MessageBuffer &buff, const Entity &e, size_type n)
Definition: genericdatahandle.hh:1103
DisjointPartitioningGatherScatter(RankIndex rank)
Create a DisjointPartitioningGatherScatter object.
Definition: genericdatahandle.hh:910
GatherScatter functor for marking shared DOFs.
Definition: genericdatahandle.hh:980
DataContainerIndexGatherScatter(GatherScatter gather_scatter=GatherScatter())
Definition: genericdatahandle.hh:592
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const Entity & e
Definition: localfunctionspace.hh:111
GFSDataHandle(const GFS &gfs, V &v, GatherScatter gather_scatter=GatherScatter(), CommunicationDescriptor communication_descriptor=CommunicationDescriptor())
Definition: genericdatahandle.hh:117
AddDataHandle(const GFS &gfs_, V &v_)
Definition: genericdatahandle.hh:629
Definition: genericdatahandle.hh:269
void write(const T &data)
Definition: polymorphicbufferwrapper.hh:32
size_type size(Entity &e) const
Definition: genericdatahandle.hh:1091
bool scatter(MessageBuffer &buff, std::size_t n, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:875
static const bool wrap_buffer
Definition: genericdatahandle.hh:33
Communication descriptor for sending one item of type E per DOF.
Definition: genericdatahandle.hh:24
void gather(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:670
std::size_t size_type
size type to use if communicating leaf ordering sizes
Definition: genericdatahandle.hh:30
bool contains(const GFS &gfs, int dim, int codim) const
Definition: genericdatahandle.hh:71
MaxDataHandle(const GFS &gfs_, V &v_)
Definition: genericdatahandle.hh:754
Definition: genericdatahandle.hh:634
std::enable_if< CommunicationDescriptor::wrap_buffer &&AlwaysTrue< Entity >::value >::type scatter(MessageBuffer &buff, const Entity &e, size_type n)
unpack data from message buffer to user
Definition: genericdatahandle.hh:198
Data handle for marking ghost DOFs.
Definition: genericdatahandle.hh:812
bool gather(MessageBuffer &buff, const Entity &e, const LocalView &local_view) const
Definition: genericdatahandle.hh:392
std::array< size_type, leaf_count+1 > Offsets
Definition: entityindexcache.hh:38
Definition: genericdatahandle.hh:493
bool fixedsize(int dim, int codim) const
Definition: genericdatahandle.hh:1084
std::size_t size_type
Definition: genericdatahandle.hh:498
bool contains(int dim, int codim) const
Definition: genericdatahandle.hh:1078
bool scatter(MessageBuffer &buff, size_type n, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:510
bool fixedSize(const GFS &gfs, int dim, int codim) const
Definition: genericdatahandle.hh:45
Definition: genericdatahandle.hh:716
Data handle for collecting set of neighboring MPI ranks.
Definition: genericdatahandle.hh:1060
Definition: genericdatahandle.hh:666
bool scatter(MessageBuffer &buff, std::size_t n, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:991
GatherScatter data handle for creating a disjoint DOF partitioning.
Definition: genericdatahandle.hh:930
bool gather(MessageBuffer &buff, const Entity &e, LocalView &local_view) const
Definition: genericdatahandle.hh:865
EntityDataCommunicationDescriptor(std::size_t count=1)
Definition: genericdatahandle.hh:88
void gather(MessageBuffer &buff, const Entity &e) const
Definition: genericdatahandle.hh:1097
GatherScatter functor for marking ghost DOFs.
Definition: genericdatahandle.hh:768
std::size_t size_type
Definition: genericdatahandle.hh:274
void gather(MessageBuffer &buff, DataType &data) const
Definition: genericdatahandle.hh:638