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