dune-grid  2.4
uggrid.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 
4 #ifndef DUNE_UGGRID_HH
5 #define DUNE_UGGRID_HH
6 
11 #include <dune/common/classname.hh>
12 #include <dune/common/parallel/collectivecommunication.hh>
13 #include <dune/common/exceptions.hh>
14 #include <dune/common/parallel/mpihelper.hh>
15 #include <dune/common/deprecated.hh>
16 
19 #include <dune/grid/common/grid.hh>
20 
21 #if HAVE_UG || DOXYGEN
22 
23 #ifdef ModelP
24 #include <dune/common/parallel/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 #include "uggrid/ug_undefs.hh"
57 #undef _2
58 
59 /* Now we're done with 2d, and we can do the whole thing over again for 3d */
60 
61 /* All macros set by UG have been unset. This includes the macros that ensure
62  single inclusion of headers. We can thus include them again. However, we
63  only want to include those headers again that contain dimension-dependent stuff.
64  Therefore, we set a few single-inclusion defines manually before including
65  ugincludes.hh again.
66  */
67 #define UGTYPES_H
68 #define __HEAPS__
69 #define __UGENV__
70 #define __DEVICESH__
71 #ifdef ModelP
72 #define __PPIF__
73 #endif
74 
75 #define _3
76 #include "uggrid/ugincludes.hh"
77 
78 // Wrap a few large UG macros by functions before they get undef'ed away.
79 // This time it's the 3d-versions.
80 #define UG_DIM 3
81 #include "uggrid/ugwrapper.hh"
82 #undef UG_DIM
83 
84 // undef all macros defined by UG
85 #include "uggrid/ug_undefs.hh"
86 
87 #undef _3
88 #undef FOR_DUNE
89 
90 // The components of the UGGrid interface
91 #include "uggrid/uggridgeometry.hh"
92 #include "uggrid/uggridlocalgeometry.hh"
93 #include "uggrid/uggridentity.hh"
94 #include "uggrid/uggridentitypointer.hh"
95 #include "uggrid/uggridentityseed.hh"
96 #include "uggrid/uggridintersections.hh"
97 #include "uggrid/uggridintersectioniterators.hh"
98 #include "uggrid/uggridleveliterator.hh"
99 #include "uggrid/uggridleafiterator.hh"
100 #include "uggrid/uggridhieriterator.hh"
101 #include "uggrid/uggridindexsets.hh"
102 #include <dune/grid/uggrid/uggridviews.hh>
103 #ifdef ModelP
104 #include "uggrid/ugmessagebuffer.hh"
105 #include "uggrid/uglbgatherscatter.hh"
106 #endif
107 
108 // Not needed here, but included for user convenience
109 #include "uggrid/uggridfactory.hh"
110 
111 #ifdef ModelP
112 template <class DataHandle, int GridDim, int codim>
113 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
114 
115 template <class DataHandle, int GridDim, int codim>
116 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
117 #endif // ModelP
118 
119 namespace Dune {
120 
121 #ifdef ModelP
122  template <int dim>
123  class CollectiveCommunication<Dune::UGGrid<dim> > :
124  public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
125  {
126  typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
127  public:
128  CollectiveCommunication()
129  : ParentType(MPIHelper::getCommunicator())
130  {}
131  };
132 #endif // ModelP
133 
134  template<int dim>
136  {
138  UGGridGeometry,
139  UGGridEntity,
140  UGGridEntityPointer,
141  UGGridLevelIterator,
142  UGGridLeafIntersection,
143  UGGridLevelIntersection,
144  UGGridLeafIntersectionIterator,
145  UGGridLevelIntersectionIterator,
146  UGGridHierarchicIterator,
147  UGGridLeafIterator,
148  UGGridLevelIndexSet< const UGGrid<dim> >,
149  UGGridLeafIndexSet< const UGGrid<dim> >,
150  UGGridIdSet< const UGGrid<dim> >,
151  typename UG_NS<dim>::UG_ID_TYPE,
152  UGGridIdSet< const UGGrid<dim> >,
153  typename UG_NS<dim>::UG_ID_TYPE,
154  CollectiveCommunication<Dune::UGGrid<dim> >,
155  UGGridLevelGridViewTraits,
156  UGGridLeafGridViewTraits,
157  UGGridEntitySeed,
158  UGGridLocalGeometry>
160  };
161 
162 
163  //**********************************************************************
164  //
165  // --UGGrid
166  //
167  //**********************************************************************
168 
206  template <int dim>
207  class UGGrid : public GridDefaultImplementation <dim, dim, double, UGGridFamily<dim> >
208  {
210 
211  friend class UGGridGeometry<0,dim,const UGGrid<dim> >;
212  friend class UGGridGeometry<dim,dim,const UGGrid<dim> >;
213  friend class UGGridGeometry<1,2,const UGGrid<dim> >;
214  friend class UGGridGeometry<2,3,const UGGrid<dim> >;
215 
216  friend class UGGridEntity <0,dim,const UGGrid<dim> >;
217  friend class UGGridEntity <dim,dim,const UGGrid<dim> >;
218  friend class UGGridHierarchicIterator<const UGGrid<dim> >;
219  friend class UGGridLeafIntersection<const UGGrid<dim> >;
220  friend class UGGridLevelIntersection<const UGGrid<dim> >;
221  friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >;
222  friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >;
223 
224  friend class UGGridLevelIndexSet<const UGGrid<dim> >;
225  friend class UGGridLeafIndexSet<const UGGrid<dim> >;
226  friend class UGGridIdSet<const UGGrid<dim> >;
227  template <class GridImp_, PartitionIteratorType PiType_>
228  friend class UGGridLeafGridView;
229  template <class GridImp_, PartitionIteratorType PiType_>
230  friend class UGGridLevelGridView;
231 
232  friend class GridFactory<UGGrid<dim> >;
233 
234 #ifdef ModelP
235  friend class UGLBGatherScatter;
236 #endif
237 
238  template <int codim_, PartitionIteratorType PiType_, class GridImp_>
239  friend class UGGridLeafIterator;
240  template <int codim_, PartitionIteratorType PiType_, class GridImp_>
241  friend class UGGridLevelIterator;
242  template <int codim_, class GridImp_>
243  friend class UGGridEntityPointer;
244 
246  static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!");
247 
248  // The different instantiations are mutual friends so they can access
249  // each others numOfUGGrids field
250  friend class UGGrid<2>;
251  friend class UGGrid<3>;
252 
253  //**********************************************************
254  // The Interface Methods
255  //**********************************************************
256  public:
259 
260  // the Traits
262 
264  typedef UG::DOUBLE ctype;
265 
271  UGGrid();
272 
274  ~UGGrid();
275 
278  int maxLevel() const;
279 
281  template<int codim>
282  typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
283 
285  template<int codim>
286  typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
287 
289  template<int codim, PartitionIteratorType PiType>
290  typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
291 
293  template<int codim, PartitionIteratorType PiType>
294  typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
295 
297  template<int codim>
298  typename Traits::template Codim<codim>::LeafIterator leafbegin() const {
299  return typename Traits::template Codim<codim>::template Partition<All_Partition>::LeafIterator(*this);
300  }
301 
303  template<int codim>
304  typename Traits::template Codim<codim>::LeafIterator leafend() const {
305  return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
306  }
307 
309  template<int codim, PartitionIteratorType PiType>
311  return typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator(*this);
312  }
313 
315  template<int codim, PartitionIteratorType PiType>
317  return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
318  }
319 
321  template <typename Seed>
322  DUNE_DEPRECATED_MSG("entityPointer() is deprecated and will be removed after the release of dune-grid 2.4. Use entity() instead to directly obtain an Entity object.")
323  typename Traits::template Codim<Seed::codimension>::EntityPointer
324  entityPointer(const Seed& seed) const
325  {
326  enum {codim = Seed::codimension};
327  return typename Traits::template Codim<codim>::EntityPointer(UGGridEntityPointer<codim,const UGGrid<dim> >(this->getRealImplementation(seed).target(),this));
328  }
329 
331  template <typename Seed>
332  typename Traits::template Codim<Seed::codimension>::Entity
333  entity(const Seed& seed) const
334  {
335  const int codim = Seed::codimension;
336  return typename Traits::template Codim<codim>::Entity(UGGridEntity<codim,dim,const UGGrid<dim> >(this->getRealImplementation(seed).target(),this));
337  }
338 
341  int size (int level, int codim) const;
342 
344  int size (int codim) const
345  {
346  return leafIndexSet().size(codim);
347  }
348 
350  int size (int level, GeometryType type) const
351  {
352  return this->levelIndexSet(level).size(type);
353  }
354 
356  int size (GeometryType type) const
357  {
358  return this->leafIndexSet().size(type);
359  }
360 
362  size_t numBoundarySegments() const {
363  // The number is stored as a member of UGGrid upon grid creation.
364  // The corresponding data structure is not exported by UG. (It is in ug/dom/std/std_internal.h)
365  return numBoundarySegments_;
366  }
367 
369  const typename Traits::GlobalIdSet& globalIdSet() const
370  {
371  return idSet_;
372  }
373 
375  const typename Traits::LocalIdSet& localIdSet() const
376  {
377  return idSet_;
378  }
379 
381  const typename Traits::LevelIndexSet& levelIndexSet(int level) const
382  {
383  if (level<0 || level>maxLevel())
384  DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
385  return *levelIndexSets_[level];
386  }
387 
389  const typename Traits::LeafIndexSet& leafIndexSet() const
390  {
391  return leafIndexSet_;
392  }
393 
396 
409  bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e );
410 
467  bool mark(const typename Traits::template Codim<0>::Entity & e,
468  typename UG_NS<dim>::RefinementRule rule,
469  int side=0);
470 
472  int getMark(const typename Traits::template Codim<0>::Entity& e) const;
473 
476  bool preAdapt();
477 
479  bool adapt();
480 
482  void postAdapt();
486  unsigned int overlapSize(int codim) const {
487  return 0;
488  }
489 
491  unsigned int ghostSize(int codim) const {
492  return (codim==0) ? 1 : 0;
493  }
494 
496  unsigned int overlapSize(int level, int codim) const {
497  return 0;
498  }
499 
501  unsigned int ghostSize(int level, int codim) const {
502  return (codim==0) ? 1 : 0;
503  }
504 
512  template<class DataHandle>
513  bool loadBalance (DataHandle& dataHandle)
514  {
515 #ifdef ModelP
516  // gather element data
517  // UGLBGatherScatter::template gather<0>(this->leafGridView(), dataHandle);
518 
519  // gather node data
520  UGLBGatherScatter::template gather<dim>(this->leafGridView(), dataHandle);
521 #endif
522 
523  // the load balancing step now also attaches
524  // the data to the entities and distributes it
525  loadBalance();
526 
527 #ifdef ModelP
528  // scatter element data
529  // UGLBGatherScatter::template scatter<0>(this->leafGridView(), dataHandle);
530 
531  // scatter node data
532  UGLBGatherScatter::template scatter<dim>(this->leafGridView(), dataHandle);
533 #endif
534 
535  return true;
536  }
537 
544  bool loadBalance(int minlevel=0);
545 
576  bool loadBalance(const std::vector<unsigned int>& targetProcessors, unsigned int fromLevel);
577 
587  template<class DataHandle>
588  bool loadBalance (const std::vector<unsigned int>& targetProcessors, unsigned int fromLevel, DataHandle& dataHandle)
589  {
590 #ifdef ModelP
591  // gather element data
592  // UGLBGatherScatter::template gather<0>(this->leafGridView(), dataHandle);
593 
594  // gather node data
595  UGLBGatherScatter::template gather<dim>(this->leafGridView(), dataHandle);
596 #endif
597 
598  // the load balancing step now also attaches
599  // the data to the entities and distributes it
600  loadBalance(targetProcessors,fromLevel);
601 
602 #ifdef ModelP
603  // scatter element data
604  // UGLBGatherScatter::template scatter<0>(this->leafGridView(), dataHandle);
605 
606  // scatter node data
607  UGLBGatherScatter::template scatter<dim>(this->leafGridView(), dataHandle);
608 #endif
609 
610  return true;
611  }
612 
613 
624  template<class DataHandle>
625  void communicate (DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
626  {
627 #ifdef ModelP
628  typedef typename UGGrid::LevelGridView LevelGridView;
629 
630  for (int curCodim = 0; curCodim <= dim; ++curCodim) {
631  if (!dataHandle.contains(dim, curCodim))
632  continue;
633 
634  if (curCodim == 0)
635  communicateUG_<LevelGridView, DataHandle, 0>(this->levelGridView(level), level, dataHandle, iftype, dir);
636  else if (curCodim == dim)
637  communicateUG_<LevelGridView, DataHandle, dim>(this->levelGridView(level), level, dataHandle, iftype, dir);
638  else if (curCodim == dim - 1)
639  communicateUG_<LevelGridView, DataHandle, dim-1>(this->levelGridView(level), level, dataHandle, iftype, dir);
640  else if (curCodim == 1)
641  communicateUG_<LevelGridView, DataHandle, 1>(this->levelGridView(level), level, dataHandle, iftype, dir);
642  else
643  DUNE_THROW(NotImplemented,
644  className(*this) << "::communicate(): Not "
645  "supported for dim " << dim << " and codim " << curCodim);
646  }
647 #endif // ModelP
648  }
649 
659  template<class DataHandle>
660  void communicate(DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir) const
661  {
662 #ifdef ModelP
663  typedef typename UGGrid::LeafGridView LeafGridView;
664 
665  for (int curCodim = 0; curCodim <= dim; ++curCodim) {
666  if (!dataHandle.contains(dim, curCodim))
667  continue;
668  int level = -1;
669  if (curCodim == 0)
670  communicateUG_<LeafGridView, DataHandle, 0>(this->leafGridView(), level, dataHandle, iftype, dir);
671  else if (curCodim == dim)
672  communicateUG_<LeafGridView, DataHandle, dim>(this->leafGridView(), level, dataHandle, iftype, dir);
673  else if (curCodim == dim - 1)
674  communicateUG_<LeafGridView, DataHandle, dim-1>(this->leafGridView(), level, dataHandle, iftype, dir);
675  else if (curCodim == 1)
676  communicateUG_<LeafGridView, DataHandle, 1>(this->leafGridView(), level, dataHandle, iftype, dir);
677  else
678  DUNE_THROW(NotImplemented,
679  className(*this) << "::communicate(): Not "
680  "supported for dim " << dim << " and codim " << curCodim);
681  }
682 #endif // ModelP
683  }
684 
687  {
688  return ccobj_;
689  }
690 
691  protected:
692 #ifdef ModelP
693  template <class GridView, class DataHandle, int codim>
694  void communicateUG_(const GridView& gv, int level,
695  DataHandle &dataHandle,
696  InterfaceType iftype,
697  CommunicationDirection dir) const
698  {
699  typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
700  // Translate the communication direction from Dune-Speak to UG-Speak
701  if (dir==ForwardCommunication)
702  ugIfDir = UG_NS<dim>::IF_FORWARD();
703  else
704  ugIfDir = UG_NS<dim>::IF_BACKWARD();
705 
706  typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
707  UGMsgBuf::duneDataHandle_ = &dataHandle;
708 
709  UGMsgBuf::level = level;
710 
711  std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
712  findDDDInterfaces_(ugIfs, iftype, codim);
713 
714  unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
715  if (!bufSize)
716  return; // we don't need to communicate if we don't have any data!
717  for (unsigned i=0; i < ugIfs.size(); ++i)
718  UG_NS<dim>::DDD_IFOneway(ugIfs[i],
719  ugIfDir,
720  bufSize,
721  &UGMsgBuf::ugGather_,
722  &UGMsgBuf::ugScatter_);
723  }
724 
725  void findDDDInterfaces_(std::vector<typename UG_NS<dim>::DDD_IF > &dddIfaces,
726  InterfaceType iftype,
727  int codim) const
728  {
729  dddIfaces.clear();
730  if (codim == 0)
731  {
732  switch (iftype) {
734  // do not communicate anything: Elements can not be in
735  // the interior of two processes at the same time
736  return;
738  // The communicated elements are in the sender's
739  // interior and it does not matter what they are for
740  // the receiver
741  dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
742  return;
743  case All_All_Interface :
744  // It does neither matter what the communicated
745  // elements are for sender nor for the receiver. If
746  // they are seen by these two processes, data is send
747  // and received.
748  dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
749  return;
750  default :
751  DUNE_THROW(GridError,
752  "Element communication not supported for "
753  "interfaces of type "
754  << iftype);
755  }
756  }
757  else if (codim == dim)
758  {
759  switch (iftype)
760  {
762  dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
763  return;
765  dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
766  dddIfaces.push_back(UG_NS<dim>::NodeIF());
767  return;
768  case All_All_Interface :
769  dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
770  return;
771  default :
772  DUNE_THROW(GridError,
773  "Node communication not supported for "
774  "interfaces of type "
775  << iftype);
776  }
777  }
778  else if (codim == dim-1)
779  {
780  switch (iftype)
781  {
783  dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
784  return;
786  dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
787  // Is the following line needed or not?
788  // dddIfaces.push_back(UG_NS<dim>::EdgeIF());
789  return;
790  case All_All_Interface :
791  dddIfaces.push_back(UG_NS<dim>::EdgeSymmVHIF());
792  return;
793  default :
794  DUNE_THROW(GridError,
795  "Edge communication not supported for "
796  "interfaces of type "
797  << iftype);
798  }
799  }
800  else if (codim == 1)
801  {
802  switch (iftype)
803  {
806  dddIfaces.push_back(UG_NS<dim>::BorderVectorSymmIF());
807  return;
808  default :
809  DUNE_THROW(GridError,
810  "Face communication not supported for "
811  "interfaces of type "
812  << iftype);
813  }
814  }
815  else
816  {
817  DUNE_THROW(GridError,
818  "Communication for codim "
819  << codim
820  << " entities is not yet supported "
821  << " by the DUNE UGGrid interface!");
822  }
823  };
824 #endif // ModelP
825 
826  public:
827  // **********************************************************
828  // End of Interface Methods
829  // **********************************************************
830 
831 #ifdef DOXYGEN
832 
844  void getChildrenOfSubface(const typename Traits::template Codim<0>::EntityPointer & e,
845  int elementSide,
846  int maxl,
847  std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements,
848  std::vector<unsigned char>& childElementSides) const;
849 
850 #else
851 
852  // This is the actual implementation of the deprecated method. We need this ugly trick of turning
853  // it into a template to avoid triggering a deprecation warning every time UGGrid is used.
854  template< typename T >
855  DUNE_DEPRECATED_MSG("This version of getChildrenOfSubface() uses EntityPointer and is deprecated. It will be removed after the release of Dune 2.4. Please use the new version with entities instead.")
856  typename std::enable_if<
857  std::is_same<
858  T,
859  typename UGGrid<dim>::Traits::template Codim<0>::EntityPointer
860  >::value
861  >::type
862  getChildrenOfSubface(const T & e,
863  int elementSide,
864  int maxl,
865  std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements,
866  std::vector<unsigned char>& childElementSides) const;
867 
868 #endif
869 
877  void getChildrenOfSubface(const typename Traits::template Codim<0>::Entity & e,
878  int elementSide,
879  int maxl,
880  std::vector<typename Traits::template Codim<0>::Entity>& childElements,
881  std::vector<unsigned char>& childElementSides) const;
882 
889  };
890 
892  enum ClosureType {
897  };
898 
901  refinementType_ = type;
902  }
903 
906  closureType_ = type;
907  }
908 
915  static void setDefaultHeapSize(unsigned size) {
916  heapSize_ = size;
917  }
918 
919 #ifdef DOXYGEN
920 
929  void setPosition(const typename Traits::template Codim<dim>::EntityPointer& e,
930  const FieldVector<double, dim>& pos);
931 
932 #else
933 
934  // This is the actual implementation of the deprecated method. We need this ugly trick of turning
935  // it into a template to avoid triggering a deprecation warning every time UGGrid is used.
936  template< typename T >
937  DUNE_DEPRECATED_MSG("This version of setPosition() uses EntityPointer and is deprecated. It will be removed after the release of Dune 2.4. Please use the new version with entities instead.")
938  typename std::enable_if<
939  std::is_same<
940  T,
941  typename UGGrid<dim>::Traits::template Codim<dim>::EntityPointer
942  >::value
943  >::type
944  setPosition(const T& e,
945  const FieldVector<double, dim>& pos);
946 
947 #endif
948 
952  void setPosition(const typename Traits::template Codim<dim>::Entity& e,
953  const FieldVector<double, dim>& pos);
954 
959  void globalRefine(int n);
960 
965  void saveState(const std::string& filename) const;
966 
971  void loadState(const std::string& filename);
972 
973  private:
975  typename UG_NS<dim>::MultiGrid* multigrid_;
976 
978  CollectiveCommunication<UGGrid> ccobj_;
979 
985  void setIndices(bool setLevelZero,
986  std::vector<unsigned int>* nodePermutation);
987 
988  // Each UGGrid object has a unique name to identify it in the
989  // UG environment structure
990  std::string name_;
991 
992  // Our set of level indices
993  std::vector<shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
994 
995  UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
996 
997  // One id set implementation
998  // Used for both the local and the global UGGrid id sets
999  UGGridIdSet<const UGGrid<dim> > idSet_;
1000 
1002  RefinementType refinementType_;
1003 
1005  ClosureType closureType_;
1006 
1014  static int numOfUGGrids;
1015 
1021  bool someElementHasBeenMarkedForRefinement_;
1022 
1028  bool someElementHasBeenMarkedForCoarsening_;
1029 
1034  static unsigned int heapSize_;
1035 
1037  std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
1038 
1044  unsigned int numBoundarySegments_;
1045 
1046  }; // end Class UGGrid
1047 
1048  namespace Capabilities
1049  {
1065  template<int dim>
1066  struct hasEntity< UGGrid<dim>, 0>
1067  {
1068  static const bool v = true;
1069  };
1070 
1074  template<int dim>
1075  struct hasEntity< UGGrid<dim>, dim>
1076  {
1077  static const bool v = true;
1078  };
1079 
1083  template<int dim>
1084  struct DUNE_DEPRECATED_MSG("Capabilities::isParallel will be removed after dune-grid-2.4.") isParallel< UGGrid<dim> >
1085  {
1086 #ifdef ModelP
1087  static const bool DUNE_DEPRECATED_MSG("Capabilities::isParallel will be removed after dune-grid-2.4.") v = true;
1088 #else // !ModelP
1089  static const bool DUNE_DEPRECATED_MSG("Capabilities::isParallel will be removed after dune-grid-2.4.") v = false;
1090 #endif // !ModelP
1091  };
1092 
1096  template<int dim>
1098  {
1099  static const bool v = true;
1100  };
1101 
1105  template<int dim>
1107  {
1108  static const bool v = false;
1109  };
1110 
1111  }
1112 
1113 } // namespace Dune
1114 
1115 #endif // HAVE_UG || DOXYGEN
1116 #endif // DUNE_UGGRID_HH
Traits::template Partition< pitype >::LevelGridView levelGridView(int level) const
View for a grid level.
Definition: common/grid.hh:1105
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: uggrid.hh:316
UGGridFamily< dim >::Traits Traits
Definition: uggrid.hh:261
friend class UGGridLevelGridView
Definition: uggrid.hh:230
Wrapper class for pointers to entities.
Definition: common/entitypointer.hh:112
Definition: common/geometry.hh:24
send interior and border, receive all entities
Definition: gridenums.hh:86
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:437
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
void saveState(const std::string &filename) const
Save entire grid hierarchy to disk.
static void setDefaultHeapSize(unsigned size)
Sets the default heap size.
Definition: uggrid.hh:915
UG::DOUBLE ctype
The type used to store coordinates.
Definition: uggrid.hh:264
Traits::template Codim< codim >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: uggrid.hh:350
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
void setPosition(const typename Traits::template Codim< dim >::EntityPointer &e, const FieldVector< double, dim > &pos)
Sets a vertex to a new position.
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: uggrid.hh:381
GridFamily::Traits::template Codim< cd >::Entity Entity
A type that is a model of a Dune::Entity.
Definition: common/grid.hh:446
Include standard header files.
Definition: agrid.hh:59
int size(int codim) const
number of leaf entities per codim in this process
Definition: uggrid.hh:344
communicate as given in InterfaceType
Definition: gridenums.hh:169
bool preAdapt()
returns true, if some elements might be coarsend during grid adaption, here always returns true ...
UGGridFamily< dim > GridFamily
type of the used GridFamily for this grid
Definition: uggrid.hh:258
void communicate(DataHandle &dataHandle, InterfaceType iftype, CommunicationDirection dir) const
The communication interface for all codims on the leaf level.
Definition: uggrid.hh:660
GridFamily::Traits::template Codim< cd >::EntityPointer EntityPointer
A type that is a model of Dune::EntityPointer.
Definition: common/grid.hh:449
send/receive interior and border entities
Definition: gridenums.hh:85
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: uggrid.hh:389
void loadState(const std::string &filename)
Read entire grid hierarchy from disk.
send all and receive all entities
Definition: gridenums.hh:89
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:263
void setClosureType(ClosureType type)
Sets the type of grid refinement closure.
Definition: uggrid.hh:905
Traits::template Codim< codim >::LevelIterator lend(int level) const
one past the end on this level
bool mark(int refCount, const typename Traits::template Codim< 0 >::Entity &e)
Mark element for refinement.
friend class UGGridLeafIterator
Definition: uggrid.hh:239
unsigned int overlapSize(int level, int codim) const
Size of the overlap on a given level.
Definition: uggrid.hh:496
friend class UGGridEntityPointer
Definition: uggrid.hh:243
bool loadBalance(const std::vector< unsigned int > &targetProcessors, unsigned int fromLevel, DataHandle &dataHandle)
Distributes the grid over the processes of a parallel machine, and sends data along with it...
Definition: uggrid.hh:588
ClosureType
Decide whether to add a green closure to locally refined grid sections or not.
Definition: uggrid.hh:892
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: uggrid.hh:356
Specialize with 'true' for all codims that a grid implements entities for. (default=false) ...
Definition: common/capabilities.hh:57
void postAdapt()
Clean up refinement markers.
void communicate(DataHandle &dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
The communication interface for all codims on a given level.
Definition: uggrid.hh:625
int size(int level, int codim) const
Number of grid entities per level and codim.
Different resources needed by all grid implementations.
Standard red/green refinement.
Definition: uggrid.hh:894
static std::conditional< std::is_reference< InterfaceType >::value, typename std::add_lvalue_reference< typename ReturnImplementationType< typename std::remove_reference< InterfaceType >::type >::ImplementationType >::type, typename std::remove_const< typename ReturnImplementationType< typename std::remove_reference< InterfaceType >::type >::ImplementationType >::type >::type getRealImplementation(InterfaceType &&i)
return real implementation of interface class
Definition: common/grid.hh:1305
unsigned int overlapSize(int codim) const
Size of the overlap on the leaf level.
Definition: uggrid.hh:486
bool loadBalance(DataHandle &dataHandle)
Distributes the grid and some data over the available nodes in a distributed machine.
Definition: uggrid.hh:513
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Create an Entity from an EntitySeed.
Definition: uggrid.hh:333
Traits::template Codim< codim >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: uggrid.hh:298
Traits::template Codim< Seed::codimension >::EntityPointer entityPointer(const Seed &seed) const
Create an EntityPointer from an EntitySeed.
Definition: uggrid.hh:324
int maxLevel() const
unsigned int ghostSize(int level, int codim) const
Size of the ghost cell layer on a given level.
Definition: uggrid.hh:501
const CollectiveCommunication< UGGrid > & comm() const
Definition: uggrid.hh:686
void setRefinementType(RefinementType type)
Sets the type of grid refinement.
Definition: uggrid.hh:900
New level consists of the refined elements and the unrefined ones, too.
Definition: uggrid.hh:888
The specialization of the generic GridFactory for UGGrid.
unsigned int ghostSize(int codim) const
Size of the ghost cell layer on the leaf level.
Definition: uggrid.hh:491
Id Set Interface.
Definition: common/grid.hh:362
friend class UGGridLevelIterator
Definition: uggrid.hh:241
Specialize with 'true' if implementation guarantees conforming level grids. (default=false) ...
Definition: common/capabilities.hh:98
STL namespace.
bool adapt()
Triggers the grid refinement process.
Index Set Interface base class.
Definition: common/grid.hh:361
Grid view abstract base class.
Definition: common/gridview.hh:58
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: uggrid.hh:369
size_t numBoundarySegments() const
Return the number of boundary segments.
Definition: uggrid.hh:362
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:16
GridFamily::Traits::CollectiveCommunication CollectiveCommunication
A type that is a model of Dune::CollectiveCommunication. It provides a portable way for collective co...
Definition: common/grid.hh:545
Definition: uggrid.hh:135
Partition< All_Partition >::LevelGridView LevelGridView
View types for All_Partition.
Definition: common/grid.hh:428
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false) ...
Definition: common/capabilities.hh:107
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
Query whether element is marked for refinement.
friend class UGGridLeafGridView
Definition: uggrid.hh:228
~UGGrid()
Destructor.
UGGrid()
Default constructor.
static const bool v
Definition: common/capabilities.hh:59
Traits::template Partition< pitype >::LeafGridView leafGridView() const
View for the leaf grid.
Definition: common/grid.hh:1115
GridTraits< dim, dim, Dune::UGGrid< dim >, UGGridGeometry, UGGridEntity, UGGridEntityPointer, UGGridLevelIterator, UGGridLeafIntersection, UGGridLevelIntersection, UGGridLeafIntersectionIterator, UGGridLevelIntersectionIterator, UGGridHierarchicIterator, UGGridLeafIterator, UGGridLevelIndexSet< const UGGrid< dim > >, UGGridLeafIndexSet< const UGGrid< dim > >, UGGridIdSet< const UGGrid< dim > >, typename UG_NS< dim >::UG_ID_TYPE, UGGridIdSet< const UGGrid< dim > >, typename UG_NS< dim >::UG_ID_TYPE, CollectiveCommunication< Dune::UGGrid< dim > >, UGGridLevelGridViewTraits, UGGridLeafGridViewTraits, UGGridEntitySeed, UGGridLocalGeometry > Traits
Definition: uggrid.hh:159
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: uggrid.hh:375
Types for GridView.
Definition: common/grid.hh:420
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
Specialize with 'true' if implementation supports parallelism. (default=false)
Definition: common/capabilities.hh:68
void getChildrenOfSubface(const typename Traits::template Codim< 0 >::EntityPointer &e, int elementSide, int maxl, std::vector< typename Traits::template Codim< 0 >::EntityPointer > &childElements, std::vector< unsigned char > &childElementSides) const
Rudimentary substitute for a hierarchic iterator on faces.
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: uggrid.hh:310
Partition< All_Partition >::LeafGridView LeafGridView
Definition: common/grid.hh:429
Traits::template Codim< codim >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: uggrid.hh:304
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:1343
Base class for grid boundary segments of arbitrary geometry.
Front-end for the grid manager of the finite element toolbox UG.
Definition: uggrid.hh:207
A set of traits classes to store static information about grid implementation.
bool loadBalance()
default implementation of load balance does nothing and returns false
Definition: common/grid.hh:1220
IndexType size(GeometryType type) const
Return total number of entities of given geometry type in entity set .
Definition: indexidset.hh:232
RefinementType
The different forms of grid refinement that UG supports.
Definition: uggrid.hh:884
No closure, results in nonconforming meshes.
Definition: uggrid.hh:896
New level consists only of the refined elements and the closure.
Definition: uggrid.hh:886
void globalRefine(int n)
Does uniform refinement.