dune-grid  2.2.1
uggrid.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 
4 #ifndef DUNE_UGGRID_HH
5 #define DUNE_UGGRID_HH
6 
11 #include <dune/common/classname.hh>
12 #include <dune/common/collectivecommunication.hh>
13 #include <dune/common/exceptions.hh>
14 #include <dune/common/mpihelper.hh>
15 #include <dune/common/static_assert.hh>
16 
19 #include <dune/grid/common/grid.hh>
20 
21 #if HAVE_UG
22 
23 #ifdef ModelP
24 #include <dune/common/mpicollectivecommunication.hh>
25 #endif
26 
27 /* The following lines including the necessary UG headers are somewhat
28  tricky. Here's what's happening:
29  UG can support two- and three-dimensional grids. You choose be setting
30  either _2 oder _3 while compiling. This changes all sorts of stuff, in
31  particular data structures in the headers.
32  UG was never supposed to provide 2d and 3d grids at the same time.
33  However, when compiling it as c++, the dimension-dependent parts are
34  wrapped up cleanly in the namespaces UG::D2 and UG::D3, respectively. That
35  way it is possible to link together the UG lib for 2d and the one for 3d.
36  But we also need the headers twice! Once with _2 set and once with _3!
37  So here we go:*/
38 
39 /* The following define tells the UG headers that we want access to a few
40  special fields, for example the extra index fields in the element data structures. */
41 #define FOR_DUNE
42 
43 // Set UG's space-dimension flag to 2d
44 #define _2
45 // And include all necessary UG headers
46 #include "uggrid/ugincludes.hh"
47 
48 // Wrap a few large UG macros by functions before they get undef'ed away.
49 // Here: The 2d-version of the macros
50 #define UG_DIM 2
51 #include "uggrid/ugwrapper.hh"
52 #undef UG_DIM
53 
54 // UG defines a whole load of preprocessor macros. ug_undefs.hh undefines
55 // them all, so we don't get name clashes.
56 #ifdef UG_LGMDOMAIN
57 #include "uggrid/ug_undefs_lgm_seq.hh"
58 #else
59 #include "uggrid/ug_undefs.hh"
60 #endif
61 #undef _2
62 
63 /* Now we're done with 2d, and we can do the whole thing over again for 3d */
64 
65 /* All macros set by UG have been unset. This includes the macros that ensure
66  single inclusion of headers. We can thus include them again. However, we
67  only want to include those headers again that contain dimension-dependent stuff.
68  Therefore, we set a few single-inclusion defines manually before including
69  ugincludes.hh again.
70 */
71 #define UGTYPES_H
72 #define __HEAPS__
73 #define __UGENV__
74 #define __DEVICESH__
75 
76 #define _3
77 #include "uggrid/ugincludes.hh"
78 
79 // Wrap a few large UG macros by functions before they get undef'ed away.
80 // This time it's the 3d-versions.
81 #define UG_DIM 3
82 #include "uggrid/ugwrapper.hh"
83 #undef UG_DIM
84 
85 // undef all macros defined by UG
86 #ifdef UG_LGMDOMAIN
87 #include "uggrid/ug_undefs_lgm_seq.hh"
88 #else
89 #include "uggrid/ug_undefs.hh"
90 #endif
91 
92 #undef _3
93 #undef FOR_DUNE
94 
95 // The components of the UGGrid interface
96 #include "uggrid/uggridgeometry.hh"
97 #include "uggrid/uggridlocalgeometry.hh"
98 #include "uggrid/uggridentity.hh"
99 #include "uggrid/uggridentitypointer.hh"
100 #include "uggrid/uggridentityseed.hh"
101 #include "uggrid/uggridintersections.hh"
102 #include "uggrid/uggridintersectioniterators.hh"
103 #include "uggrid/uggridleveliterator.hh"
104 #include "uggrid/uggridleafiterator.hh"
105 #include "uggrid/uggridhieriterator.hh"
106 #include "uggrid/uggridindexsets.hh"
107 
108 // Not needed here, but included for user convenience
109 #include "uggrid/uggridfactory.hh"
110 
111 #ifdef ModelP
112 namespace Dune {
113 
114 // converts the UG speak message buffers to DUNE speak and vice-versa
115  template <class DataHandle, int GridDim, int codim>
116  class UGMessageBufferBase {
117  protected:
118  typedef UGMessageBufferBase<DataHandle, GridDim, codim> ThisType;
119  typedef UGGrid<GridDim> GridType;
120  typedef typename DataHandle::DataType DataType;
121 
122  enum {
123  dim = GridDim
124  };
125 
126  UGMessageBufferBase(void *ugData)
127  {
128  ugData_ = static_cast<char*>(ugData);
129  };
130 
131  public:
132  void write(const DataType &t)
133  { this->writeRaw_<DataType>(t); }
134 
135  void read(DataType &t)
136  { this->readRaw_<DataType>(t); }
137 
138  protected:
139  friend class Dune::UGGrid<dim>;
140 
141  template <class ValueType>
142  void writeRaw_(const ValueType &v)
143  {
144  *reinterpret_cast<ValueType*>(ugData_) = v;
145  ugData_ += sizeof(ValueType);
146  }
147 
148  template <class ValueType>
149  void readRaw_(ValueType &v)
150  {
151  v = *reinterpret_cast<ValueType*>(ugData_);
152  ugData_ += sizeof(ValueType);
153  }
154 
155  // called by DDD_IFOneway to serialize the data structure to
156  // be send
157  static int ugGather_(typename UG_NS<dim>::DDD_OBJ obj, void* data)
158  {
159  if (codim == 0) {
162  UGMakeableEntity<0, dim, UGGrid<dim> > e(reinterpret_cast<typename UG_NS<dim>::Element*>(obj),nullptr);
163  // safety check to only communicate what is needed
164  if ((level == -1 && UG_NS<dim>::isLeaf(reinterpret_cast<typename UG_NS<dim>::Element*>(obj))) || e.level() == level)
165  {
166  ThisType msgBuf(static_cast<DataType*>(data));
167  if (!duneDataHandle_->fixedsize(dim, codim))
168  msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
169  duneDataHandle_->gather(msgBuf, e);
170  }
171  }
172  else if (codim == dim) {
175  UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > e(reinterpret_cast<typename UG_NS<dim>::Node*>(obj),nullptr);
176  // safety check to only communicate what is needed
177  if ((level == -1 && UG_NS<dim>::isLeaf(reinterpret_cast<typename UG_NS<dim>::Node*>(obj))) || e.level() == level)
178  {
179  ThisType msgBuf(static_cast<DataType*>(data));
180  if (!duneDataHandle_->fixedsize(dim, codim))
181  msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
182  duneDataHandle_->gather(msgBuf, e);
183  }
184  }
185  else if (codim == dim - 1) {
188  UGMakeableEntity<dim-1, dim, Dune::UGGrid<dim> > e(reinterpret_cast<typename UG_NS<dim>::Edge*>(obj),nullptr);
189  // safety check to only communicate what is needed
190  if ((level == -1 && UG_NS<dim>::isLeaf(reinterpret_cast<typename UG_NS<dim>::Edge*>(obj))) || e.level() == level)
191  {
192  ThisType msgBuf(static_cast<DataType*>(data));
193  if (!duneDataHandle_->fixedsize(dim, codim))
194  msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
195  duneDataHandle_->gather(msgBuf, e);
196  }
197  }
198  else {
199  DUNE_THROW(GridError,
200  "Only node and element wise "
201  "communication is currently "
202  "supported by UGGrid");
203  }
204 
205  return 0;
206  }
207 
208  // called by DDD_IFOneway to deserialize the data structure
209  // that has been received
210  static int ugScatter_(typename UG_NS<dim>::DDD_OBJ obj, void* data)
211  {
212  if (codim == 0) {
213  typedef UGMakeableEntity<0, dim, UGGrid<dim> > Entity;
216  Entity e(reinterpret_cast<typename UG_NS<dim>::Element*>(obj),nullptr);
217  // safety check to only communicate what is needed
218  if ((level == -1 && UG_NS<dim>::isLeaf(reinterpret_cast<typename UG_NS<dim>::Element*>(obj))) || e.level() == level)
219  {
220  ThisType msgBuf(static_cast<DataType*>(data));
221  int n;
222  if (!duneDataHandle_->fixedsize(dim, codim))
223  msgBuf.readRaw_(n);
224  else
225  n = duneDataHandle_->template size<Entity>(e);
226  if (n > 0)
227  duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
228  }
229  }
230  else if (codim == dim) {
231  typedef UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > Entity;
234  Entity e(reinterpret_cast<typename UG_NS<dim>::Node*>(obj),nullptr);
235  // safety check to only communicate what is needed
236  if ((level == -1 && UG_NS<dim>::isLeaf(reinterpret_cast<typename UG_NS<dim>::Node*>(obj))) || e.level() == level)
237  {
238  ThisType msgBuf(static_cast<DataType*>(data));
239  int n;
240  if (!duneDataHandle_->fixedsize(dim, codim))
241  msgBuf.readRaw_(n);
242  else
243  n = duneDataHandle_->template size<Entity>(e);
244  if (n > 0)
245  duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
246  }
247  }
248  else if (codim == dim - 1) { // !!!ALEX!!! Is it possible to send codim 1 in UG<2>?
249  typedef UGMakeableEntity<dim-1, dim, Dune::UGGrid<dim> > Entity;
252  Entity e(reinterpret_cast<typename UG_NS<dim>::Edge*>(obj),nullptr);
253  // safety check to only communicate what is needed
254  if ((level == -1 && UG_NS<dim>::isLeaf(reinterpret_cast<typename UG_NS<dim>::Edge*>(obj))) || e.level() == level)
255  {
256  ThisType msgBuf(static_cast<DataType*>(data));
257  int n;
258  if (!duneDataHandle_->fixedsize(dim, codim))
259  msgBuf.readRaw_(n);
260  else
261  n = duneDataHandle_->template size<Entity>(e);
262  if (n > 0)
263  duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
264  }
265  }
266  else {
267  DUNE_THROW(GridError,
268  "Only node and element wise "
269  "communication is currently "
270  "supported by UGGrid");
271  }
272 
273  return 0;
274  }
275  static DataHandle *duneDataHandle_;
276 
277  static int level;
278 
279  char *ugData_;
280  };
281 
282  template <class DataHandle, int GridDim, int codim>
283  class UGMessageBuffer
284  : public UGMessageBufferBase<DataHandle, GridDim, codim>
285  {
286  typedef typename DataHandle::DataType DataType;
287  typedef UGMessageBufferBase<DataHandle, GridDim, codim> Base;
288  enum { dim = GridDim };
289 
290  protected:
291  friend class Dune::UGGrid<dim>;
292 
293  UGMessageBuffer(void *ugData)
294  : Base(ugData)
295  {}
296 
297  // returns number of bytes required for the UG message buffer
298  template <class GridView>
299  static unsigned ugBufferSize_(const GridView &gv)
300  {
301  if (Base::duneDataHandle_->fixedsize(dim, codim)) {
302  return sizeof(DataType)
303  * Base::duneDataHandle_->size(*gv.template begin<codim,InteriorBorder_Partition>());
304  }
305 
306  typedef typename GridView::template Codim<codim>::Entity Entity;
307 
308  // iterate over all entities, find the maximum size for
309  // the current rank
310  int maxSize = 0;
311  typedef typename
312  GridView
313  ::template Codim<codim>
314  ::template Partition<Dune::All_Partition>
315  ::Iterator Iterator;
316  Iterator it = gv.template begin<codim, Dune::All_Partition>();
317  const Iterator endIt = gv.template end<codim, Dune::All_Partition>();
318  for (; it != endIt; ++it) {
319  maxSize = std::max((int) maxSize,
320  (int) Base::duneDataHandle_->size(*it));
321  }
322 
323  // find maximum size for all ranks
324  maxSize = MPIHelper::getCollectiveCommunication().max(maxSize);
325  if (!maxSize)
326  return 0;
327 
328  // add the size of an unsigned integer to the actual
329  // buffer size. (we somewhere have to store the actual
330  // number of objects for each entity.)
331  return sizeof(unsigned) + sizeof(DataType)*maxSize;
332  }
333  };
334 
335  template <class DataHandle, int GridDim>
336  class UGEdgeMessageBuffer
337  : public UGMessageBufferBase<DataHandle, GridDim, GridDim-1>
338  {
339  enum {codim = GridDim-1,
340  dim = GridDim};
341  typedef typename DataHandle::DataType DataType;
342  typedef UGMessageBufferBase<DataHandle, GridDim, codim> Base;
343  protected:
344  friend class Dune::UGGrid<dim>;
345 
346  UGEdgeMessageBuffer(void *ugData)
347  : Base(ugData)
348  {}
349 
350  // returns number of bytes required for the UG message buffer
351  template <class GridView>
352  static unsigned ugBufferSize_(const GridView &gv)
353  {
354  if (Base::duneDataHandle_->fixedsize(dim, codim)) {
355  typedef typename GridView::template Codim<0>::Entity Element;
356  const Element& element = *gv.template begin<0, InteriorBorder_Partition>();
357  return sizeof(DataType)
358  * Base::duneDataHandle_->size(*element.template subEntity<codim>(0));
359  }
360 
361  typedef typename GridView::template Codim<codim>::Entity Entity;
362 
363  // iterate over all entities, find the maximum size for
364  // the current rank
365  int maxSize = 0;
366  typedef typename
367  GridView
368  ::template Codim<0>
369  ::template Partition<Dune::All_Partition>
370  ::Iterator Iterator;
371  Iterator it = gv.template begin<0, Dune::All_Partition>();
372  const Iterator endIt = gv.template end<0, Dune::All_Partition>();
373  for (; it != endIt; ++it) {
374  int numberOfEdges = it->template count<codim>();
375  for (int k = 0; k < numberOfEdges; k++)
376  {
377  typedef typename GridView::template Codim<0>::Entity Element;
378  typedef typename Element::template Codim<codim>::EntityPointer EdgePointer;
379  const EdgePointer edgePointer(it->template subEntity<codim>(k));
380 
381  maxSize = std::max((int) maxSize,
382  (int) Base::duneDataHandle_->size(*edgePointer));
383  }
384  }
385 
386  // find maximum size for all ranks
387  maxSize = MPIHelper::getCollectiveCommunication().max(maxSize);
388  if (!maxSize)
389  return 0;
390 
391  // add the size of an unsigned integer to the actual
392  // buffer size. (we somewhere have to store the actual
393  // number of objects for each entity.)
394  return sizeof(unsigned) + sizeof(DataType)*maxSize;
395  }
396  };
397 
398  template <class DataHandle>
399  class UGMessageBuffer<DataHandle, 2, 1>
400  : public UGEdgeMessageBuffer<DataHandle, 2>
401  {};
402 
403  template <class DataHandle>
404  class UGMessageBuffer<DataHandle, 3, 2>
405  : public UGEdgeMessageBuffer<DataHandle, 3>
406  {};
407 
408 } // end namespace Dune
409 
410 template <class DataHandle, int GridDim, int codim>
411 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
412 
413 template <class DataHandle, int GridDim, int codim>
414 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
415 #endif // ModelP
416 
417 namespace Dune {
418 
419 #ifdef ModelP
420 template <int dim>
421 class CollectiveCommunication<Dune::UGGrid<dim> > :
422  public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
423 {
424  typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
425 public:
426  CollectiveCommunication()
427  : ParentType(MPIHelper::getCommunicator())
428  {}
429 };
430 #endif // ModelP
431 
432 template<int dim, int dimworld>
434 {
436  UGGridGeometry,
437  UGGridEntity,
438  UGGridEntityPointer,
439  UGGridLevelIterator,
440  UGGridLeafIntersection,
441  UGGridLevelIntersection,
442  UGGridLeafIntersectionIterator,
443  UGGridLevelIntersectionIterator,
444  UGGridHierarchicIterator,
445  UGGridLeafIterator,
446  UGGridLevelIndexSet< const UGGrid<dim> >,
447  UGGridLeafIndexSet< const UGGrid<dim> >,
448  UGGridIdSet< const UGGrid<dim>, false >,
449  unsigned int,
450  UGGridIdSet< const UGGrid<dim>, true >,
451  unsigned int,
452  CollectiveCommunication<Dune::UGGrid<dim> >,
455  UGGridEntitySeed,
456  UGGridLocalGeometry>
458 };
459 
460 
461 //**********************************************************************
462 //
463 // --UGGrid
464 //
465 //**********************************************************************
466 
503 template <int dim>
504 class UGGrid : public GridDefaultImplementation <dim, dim, double, UGGridFamily<dim,dim> >
505 {
506  friend class UGGridGeometry<0,dim,const UGGrid<dim> >;
507  friend class UGGridGeometry<dim,dim,const UGGrid<dim> >;
508  friend class UGGridGeometry<1,2,const UGGrid<dim> >;
509  friend class UGGridGeometry<2,3,const UGGrid<dim> >;
510 
511  friend class UGGridEntity <0,dim,const UGGrid<dim> >;
512  friend class UGGridEntity <dim,dim,const UGGrid<dim> >;
513  friend class UGGridHierarchicIterator<const UGGrid<dim> >;
514  friend class UGGridLeafIntersection<const UGGrid<dim> >;
515  friend class UGGridLevelIntersection<const UGGrid<dim> >;
516  friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >;
517  friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >;
518 
519  friend class UGGridLevelIndexSet<const UGGrid<dim> >;
520  friend class UGGridLeafIndexSet<const UGGrid<dim> >;
521  friend class UGGridIdSet<const UGGrid<dim>, false >;
522  friend class UGGridIdSet<const UGGrid<dim>, true >;
523 
524  friend class GridFactory<UGGrid<dim> >;
525 
526  template <int codim_, PartitionIteratorType PiType_, class GridImp_>
527  friend class UGGridLeafIterator;
528  template <int codim_, PartitionIteratorType PiType_, class GridImp_>
529  friend class UGGridLevelIterator;
530  template <int codim_, class GridImp_>
531  friend class UGGridEntityPointer;
532 
534  dune_static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!");
535 
536  // The different instantiations are mutual friends so they can access
537  // each others numOfUGGrids field
538  friend class UGGrid<2>;
539  friend class UGGrid<3>;
540 
541  //**********************************************************
542  // The Interface Methods
543  //**********************************************************
544 public:
547 
548  // the Traits
550 
552  typedef UG::DOUBLE ctype;
553 
559  UGGrid();
560 
562  ~UGGrid();
563 
566  int maxLevel() const;
567 
569  template<int codim>
570  typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
571 
573  template<int codim>
574  typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
575 
577  template<int codim, PartitionIteratorType PiType>
578  typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
579 
581  template<int codim, PartitionIteratorType PiType>
582  typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
583 
585  template<int codim>
586  typename Traits::template Codim<codim>::LeafIterator leafbegin() const {
587  return typename Traits::template Codim<codim>::template Partition<All_Partition>::LeafIterator(*this);
588  }
589 
591  template<int codim>
592  typename Traits::template Codim<codim>::LeafIterator leafend() const {
593  return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
594  }
595 
597  template<int codim, PartitionIteratorType PiType>
599  return typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator(*this);
600  }
601 
603  template<int codim, PartitionIteratorType PiType>
605  return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
606  }
607 
609  template <int codim>
610  typename Traits::template Codim<codim>::EntityPointer
611  entityPointer(const UGGridEntitySeed<codim, const UGGrid<dim> >& seed) const
612  {
613  return typename Traits::template Codim<codim>::EntityPointer(UGGridEntityPointer<codim,const UGGrid<dim> >(seed.target(),this));
614  }
615 
618  int size (int level, int codim) const;
619 
621  int size (int codim) const
622  {
623  return leafIndexSet().size(codim);
624  }
625 
627  int size (int level, GeometryType type) const
628  {
629  return this->levelIndexSet(level).size(type);
630  }
631 
633  int size (GeometryType type) const
634  {
635  return this->leafIndexSet().size(type);
636  }
637 
639  size_t numBoundarySegments() const {
640  // The number is stored as a member of UGGrid upon grid creation.
641  // The corresponding data structure is not exported by UG. (It is in ug/dom/std/std_internal.h)
642  return numBoundarySegments_;
643  }
644 
646  const typename Traits::GlobalIdSet& globalIdSet() const
647  {
648  return globalIdSet_;
649  }
650 
652  const typename Traits::LocalIdSet& localIdSet() const
653  {
654  return localIdSet_;
655  }
656 
658  const typename Traits::LevelIndexSet& levelIndexSet(int level) const
659  {
660  if (level<0 || level>maxLevel())
661  DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
662  return *levelIndexSets_[level];
663  }
664 
666  const typename Traits::LeafIndexSet& leafIndexSet() const
667  {
668  return leafIndexSet_;
669  }
670 
673 
686  bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e );
687 
695  bool mark(const typename Traits::template Codim<0>::Entity & e,
696  typename UG_NS<dim>::RefinementRule rule,
697  int side=0);
698 
700  int getMark(const typename Traits::template Codim<0>::Entity& e) const;
701 
704  bool preAdapt();
705 
707  bool adapt();
708 
710  void postAdapt();
714  unsigned int overlapSize(int codim) const {
715  return 0;
716  }
717 
719  unsigned int ghostSize(int codim) const {
720  return (codim==0) ? 1 : 0;
721  }
722 
724  unsigned int overlapSize(int level, int codim) const {
725  return 0;
726  }
727 
729  unsigned int ghostSize(int level, int codim) const {
730  return (codim==0) ? 1 : 0;
731  }
732 
737  bool loadBalance() {
738  return loadBalance(0,0,2,32,1);
739  }
740 
746  template<class DataHandle>
747  bool loadBalance (DataHandle& data)
748  {
749  DUNE_THROW(NotImplemented, "load balancing with data attached");
750  }
751 
768  bool loadBalance(int strategy, int minlevel, int depth, int maxlevel, int minelement);
769 
780  template<class DataHandle>
781  void communicate (DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
782  {
783 #ifdef ModelP
784  typedef typename UGGrid::LevelGridView LevelGridView;
785 
786  for (int curCodim = 0; curCodim <= dim; ++curCodim) {
787  if (!dataHandle.contains(dim, curCodim))
788  continue;
789 
790  if (curCodim == 0)
791  communicateUG_<LevelGridView, DataHandle, 0>(this->levelView(level), level, dataHandle, iftype, dir);
792  else if (curCodim == dim)
793  communicateUG_<LevelGridView, DataHandle, dim>(this->levelView(level), level, dataHandle, iftype, dir);
794  else if (curCodim == dim - 1)
795  communicateUG_<LevelGridView, DataHandle, dim-1>(this->levelView(level), level, dataHandle, iftype, dir);
796  else
797  DUNE_THROW(NotImplemented,
798  className(*this) << "::communicate(): Only "
799  "supported for codim=0 and "
800  "codim=dim(=" << dim << "), but "
801  "codim=" << curCodim << " was requested");
802  }
803 #endif // ModelP
804  }
805 
815  template<class DataHandle>
816  void communicate(DataHandle& dataHandle,
817  InterfaceType iftype,
818  CommunicationDirection dir) const
819  {
820 #ifdef ModelP
821  typedef typename UGGrid::LeafGridView LeafGridView;
822 
823  for (int curCodim = 0; curCodim <= dim; ++curCodim) {
824  if (!dataHandle.contains(dim, curCodim))
825  continue;
826  int level = -1;
827  if (curCodim == 0)
828  communicateUG_<LeafGridView, DataHandle, 0>(this->leafView(), level, dataHandle, iftype, dir);
829  else if (curCodim == dim)
830  communicateUG_<LeafGridView, DataHandle, dim>(this->leafView(), level, dataHandle, iftype, dir);
831  else if (curCodim == dim - 1) // !!!ALEX!!! Is it possible to send codim 1 in UG<2>?
832  {
833  communicateUG_<LeafGridView, DataHandle, dim-1>(this->leafView(), level, dataHandle, iftype, dir);
834  }
835  else
836  DUNE_THROW(NotImplemented,
837  className(*this) << "::communicate(): Only "
838  "supported for codim=0 and "
839  "codim=dim(=" << dim << "), but "
840  "codim=" << curCodim << " was requested");
841  }
842 #endif // ModelP
843  }
844 
847  {
848  return ccobj_;
849  }
850 
851 protected:
852 #ifdef ModelP
853  template <class GridView, class DataHandle, int codim>
854  void communicateUG_(const GridView& gv, int level,
855  DataHandle &dataHandle,
856  InterfaceType iftype,
857  CommunicationDirection dir) const
858  {
859  typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
860  // Translate the communication direction from Dune-Speak to UG-Speak
861  if (dir==ForwardCommunication)
862  ugIfDir = UG_NS<dim>::IF_FORWARD();
863  else
864  ugIfDir = UG_NS<dim>::IF_BACKWARD();
865 
866  typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
867  UGMsgBuf::duneDataHandle_ = &dataHandle;
868 
869  UGMsgBuf::level = level;
870 
871  std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
872  findDDDInterfaces_(ugIfs, iftype, codim);
873 
874  unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
875  if (!bufSize)
876  return; // we don't need to communicate if we don't have any data!
877  for (unsigned i=0; i < ugIfs.size(); ++i)
878  UG_NS<dim>::DDD_IFOneway(ugIfs[i],
879  ugIfDir,
880  bufSize,
881  &UGMsgBuf::ugGather_,
882  &UGMsgBuf::ugScatter_);
883  }
884 
885  void findDDDInterfaces_(std::vector<typename UG_NS<dim>::DDD_IF > &dddIfaces,
886  InterfaceType iftype,
887  int codim) const
888  {
889  dddIfaces.clear();
890  switch (codim)
891  {
892  case 0:
893  switch (iftype) {
895  // do not communicate anything: Elements can not be in
896  // the interior of two processes at the same time
897  return;
899  // The communicated elements are in the sender's
900  // interior and it does not matter what they are for
901  // the receiver
902  dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
903  return;
904  case All_All_Interface:
905  // It does neither matter what the communicated
906  // elements are for sender nor for the receiver. If
907  // they are seen by these two processes, data is send
908  // and received.
909  dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
910  return;
911  default:
912  DUNE_THROW(GridError,
913  "Element communication not supported for "
914  "interfaces of type "
915  << iftype);
916  }
917 
918  case dim:
919  switch (iftype)
920  {
922  dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
923  return;
925  dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
926  dddIfaces.push_back(UG_NS<dim>::NodeIF());
927  return;
928  case All_All_Interface:
929  dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
930  return;
931  default:
932  DUNE_THROW(GridError,
933  "Node communication not supported for "
934  "interfaces of type "
935  << iftype);
936  }
937 
938  case dim-1:
939  switch (iftype)
940  {
942  dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
943  return;
945  dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
946  // Is the following line needed or not?
947  // dddIfaces.push_back(UG_NS<dim>::EdgeIF());
948  return;
949  default:
950  DUNE_THROW(GridError,
951  "Edge communication not supported for "
952  "interfaces of type "
953  << iftype);
954  }
955 
956  default:
957  DUNE_THROW(GridError,
958  "Communication for codim "
959  << codim
960  << " entities is not yet supported "
961  << " by the DUNE UGGrid interface!");
962  }
963  };
964 #endif // ModelP
965 
966 public:
967  // **********************************************************
968  // End of Interface Methods
969  // **********************************************************
970 
973 
977  void createLGMGrid(const std::string& name);
978 
988  void getChildrenOfSubface(const typename Traits::template Codim<0>::EntityPointer & e,
989  int elementSide,
990  int maxl,
991  std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements,
992  std::vector<unsigned char>& childElementSides) const;
993 
1000 
1007 
1010  refinementType_ = type;
1011  }
1012 
1015  closureType_ = type;
1016  }
1017 
1024  static void setDefaultHeapSize(unsigned size) {
1025  heapSize_ = size;
1026  }
1027 
1031  void setPosition(const typename Traits::template Codim<dim>::EntityPointer& e,
1032  const FieldVector<double, dim>& pos);
1033 
1038  void globalRefine(int n);
1039 
1044  void saveState(const std::string& filename) const;
1045 
1050  void loadState(const std::string& filename);
1051 
1052 private:
1054  typename UG_NS<dim>::MultiGrid* multigrid_;
1055 
1057  CollectiveCommunication<UGGrid> ccobj_;
1058 
1064  void setIndices(bool setLevelZero,
1065  std::vector<unsigned int>* nodePermutation);
1066 
1067  // Each UGGrid object has a unique name to identify it in the
1068  // UG environment structure
1069  std::string name_;
1070 
1071  // Our set of level indices
1072  std::vector<UGGridLevelIndexSet<const UGGrid<dim> >*> levelIndexSets_;
1073 
1074  UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
1075 
1076  UGGridIdSet<const UGGrid<dim>, false > globalIdSet_;
1077 
1078  UGGridIdSet<const UGGrid<dim>, true > localIdSet_;
1079 
1081  RefinementType refinementType_;
1082 
1084  ClosureType closureType_;
1085 
1093  static int numOfUGGrids;
1094 
1100  bool someElementHasBeenMarkedForRefinement_;
1101 
1107  bool someElementHasBeenMarkedForCoarsening_;
1108 
1113  static unsigned int heapSize_;
1114 
1116  std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
1117 
1123  unsigned int numBoundarySegments_;
1124 
1125 }; // end Class UGGrid
1126 
1127 namespace Capabilities
1128 {
1144  template<int dim>
1145  struct hasEntity< UGGrid<dim>, 0>
1146  {
1147  static const bool v = true;
1148  };
1149 
1153  template<int dim>
1154  struct hasEntity< UGGrid<dim>, dim>
1155  {
1156  static const bool v = true;
1157  };
1158 
1162  template<int dim>
1163  struct isParallel< UGGrid<dim> >
1164  {
1165 #ifdef ModelP
1166  static const bool v = true;
1167 #else // !ModelP
1168  static const bool v = false;
1169 #endif // !ModelP
1170  };
1171 
1175  template<int dim>
1177  {
1178  static const bool v = true;
1179  };
1180 
1184  template<int dim>
1186  {
1187  static const bool v = false;
1188  };
1189 
1190 }
1191 
1192 } // namespace Dune
1193 
1194 #endif // HAVE_UG
1195 #endif // DUNE_UGGRID_HH