dune-grid  2.3.0
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/static_assert.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 
72 #define _3
73 #include "uggrid/ugincludes.hh"
74 
75 // Wrap a few large UG macros by functions before they get undef'ed away.
76 // This time it's the 3d-versions.
77 #define UG_DIM 3
78 #include "uggrid/ugwrapper.hh"
79 #undef UG_DIM
80 
81 // undef all macros defined by UG
82 #include "uggrid/ug_undefs.hh"
83 
84 #undef _3
85 #undef FOR_DUNE
86 
87 // The components of the UGGrid interface
88 #include "uggrid/uggridgeometry.hh"
89 #include "uggrid/uggridlocalgeometry.hh"
90 #include "uggrid/uggridentity.hh"
91 #include "uggrid/uggridentitypointer.hh"
92 #include "uggrid/uggridentityseed.hh"
93 #include "uggrid/uggridintersections.hh"
94 #include "uggrid/uggridintersectioniterators.hh"
95 #include "uggrid/uggridleveliterator.hh"
96 #include "uggrid/uggridleafiterator.hh"
97 #include "uggrid/uggridhieriterator.hh"
98 #include "uggrid/uggridindexsets.hh"
99 #ifdef ModelP
100 #include "uggrid/ugmessagebuffer.hh"
101 #include "uggrid/uglbgatherscatter.hh"
102 #endif
103 
104 // Not needed here, but included for user convenience
105 #include "uggrid/uggridfactory.hh"
106 
107 #ifdef ModelP
108 template <class DataHandle, int GridDim, int codim>
109 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
110 
111 template <class DataHandle, int GridDim, int codim>
112 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
113 #endif // ModelP
114 
115 namespace Dune {
116 
117 #ifdef ModelP
118  template <int dim>
119  class CollectiveCommunication<Dune::UGGrid<dim> > :
120  public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
121  {
122  typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
123  public:
124  CollectiveCommunication()
125  : ParentType(MPIHelper::getCommunicator())
126  {}
127  };
128 #endif // ModelP
129 
130  template<int dim>
132  {
134  UGGridGeometry,
135  UGGridEntity,
136  UGGridEntityPointer,
137  UGGridLevelIterator,
138  UGGridLeafIntersection,
139  UGGridLevelIntersection,
140  UGGridLeafIntersectionIterator,
141  UGGridLevelIntersectionIterator,
142  UGGridHierarchicIterator,
143  UGGridLeafIterator,
144  UGGridLevelIndexSet< const UGGrid<dim> >,
145  UGGridLeafIndexSet< const UGGrid<dim> >,
146  UGGridIdSet< const UGGrid<dim> >,
147  typename UG_NS<dim>::UG_ID_TYPE,
148  UGGridIdSet< const UGGrid<dim> >,
149  typename UG_NS<dim>::UG_ID_TYPE,
150  CollectiveCommunication<Dune::UGGrid<dim> >,
153  UGGridEntitySeed,
154  UGGridLocalGeometry>
156  };
157 
158 
159  //**********************************************************************
160  //
161  // --UGGrid
162  //
163  //**********************************************************************
164 
201  template <int dim>
202  class UGGrid : public GridDefaultImplementation <dim, dim, double, UGGridFamily<dim> >
203  {
205 
206  friend class UGGridGeometry<0,dim,const UGGrid<dim> >;
207  friend class UGGridGeometry<dim,dim,const UGGrid<dim> >;
208  friend class UGGridGeometry<1,2,const UGGrid<dim> >;
209  friend class UGGridGeometry<2,3,const UGGrid<dim> >;
210 
211  friend class UGGridEntity <0,dim,const UGGrid<dim> >;
212  friend class UGGridEntity <dim,dim,const UGGrid<dim> >;
213  friend class UGGridHierarchicIterator<const UGGrid<dim> >;
214  friend class UGGridLeafIntersection<const UGGrid<dim> >;
215  friend class UGGridLevelIntersection<const UGGrid<dim> >;
216  friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >;
217  friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >;
218 
219  friend class UGGridLevelIndexSet<const UGGrid<dim> >;
220  friend class UGGridLeafIndexSet<const UGGrid<dim> >;
221  friend class UGGridIdSet<const UGGrid<dim> >;
222 
223  friend class GridFactory<UGGrid<dim> >;
224 
225 #ifdef ModelP
226  friend class UGLBGatherScatter;
227 #endif
228 
229  template <int codim_, PartitionIteratorType PiType_, class GridImp_>
230  friend class UGGridLeafIterator;
231  template <int codim_, PartitionIteratorType PiType_, class GridImp_>
232  friend class UGGridLevelIterator;
233  template <int codim_, class GridImp_>
234  friend class UGGridEntityPointer;
235 
237  dune_static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!");
238 
239  // The different instantiations are mutual friends so they can access
240  // each others numOfUGGrids field
241  friend class UGGrid<2>;
242  friend class UGGrid<3>;
243 
244  //**********************************************************
245  // The Interface Methods
246  //**********************************************************
247  public:
250 
251  // the Traits
253 
255  typedef UG::DOUBLE ctype;
256 
262  UGGrid();
263 
265  ~UGGrid();
266 
269  int maxLevel() const;
270 
272  template<int codim>
273  typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
274 
276  template<int codim>
277  typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
278 
280  template<int codim, PartitionIteratorType PiType>
281  typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
282 
284  template<int codim, PartitionIteratorType PiType>
285  typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
286 
288  template<int codim>
289  typename Traits::template Codim<codim>::LeafIterator leafbegin() const {
290  return typename Traits::template Codim<codim>::template Partition<All_Partition>::LeafIterator(*this);
291  }
292 
294  template<int codim>
295  typename Traits::template Codim<codim>::LeafIterator leafend() const {
296  return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
297  }
298 
300  template<int codim, PartitionIteratorType PiType>
302  return typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator(*this);
303  }
304 
306  template<int codim, PartitionIteratorType PiType>
308  return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
309  }
310 
312  template <typename Seed>
313  typename Traits::template Codim<Seed::codimension>::EntityPointer
314  entityPointer(const Seed& seed) const
315  {
316  enum {codim = Seed::codimension};
317  return typename Traits::template Codim<codim>::EntityPointer(UGGridEntityPointer<codim,const UGGrid<dim> >(this->getRealImplementation(seed).target(),this));
318  }
319 
322  int size (int level, int codim) const;
323 
325  int size (int codim) const
326  {
327  return leafIndexSet().size(codim);
328  }
329 
331  int size (int level, GeometryType type) const
332  {
333  return this->levelIndexSet(level).size(type);
334  }
335 
337  int size (GeometryType type) const
338  {
339  return this->leafIndexSet().size(type);
340  }
341 
343  size_t numBoundarySegments() const {
344  // The number is stored as a member of UGGrid upon grid creation.
345  // The corresponding data structure is not exported by UG. (It is in ug/dom/std/std_internal.h)
346  return numBoundarySegments_;
347  }
348 
350  const typename Traits::GlobalIdSet& globalIdSet() const
351  {
352  return idSet_;
353  }
354 
356  const typename Traits::LocalIdSet& localIdSet() const
357  {
358  return idSet_;
359  }
360 
362  const typename Traits::LevelIndexSet& levelIndexSet(int level) const
363  {
364  if (level<0 || level>maxLevel())
365  DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
366  return *levelIndexSets_[level];
367  }
368 
370  const typename Traits::LeafIndexSet& leafIndexSet() const
371  {
372  return leafIndexSet_;
373  }
374 
377 
390  bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e );
391 
448  bool mark(const typename Traits::template Codim<0>::Entity & e,
449  typename UG_NS<dim>::RefinementRule rule,
450  int side=0);
451 
453  int getMark(const typename Traits::template Codim<0>::Entity& e) const;
454 
457  bool preAdapt();
458 
460  bool adapt();
461 
463  void postAdapt();
467  unsigned int overlapSize(int codim) const {
468  return 0;
469  }
470 
472  unsigned int ghostSize(int codim) const {
473  return (codim==0) ? 1 : 0;
474  }
475 
477  unsigned int overlapSize(int level, int codim) const {
478  return 0;
479  }
480 
482  unsigned int ghostSize(int level, int codim) const {
483  return (codim==0) ? 1 : 0;
484  }
485 
490  bool loadBalance() {
491  return loadBalance(0,0,2,32,1);
492  }
493 
501  template<class DataHandle>
502  bool loadBalance (DataHandle& dataHandle)
503  {
504 #if !HAVE_UG_PATCH10
505  DUNE_THROW(NotImplemented, "load balancing with data attached");
506 #else
507 #ifdef ModelP
508  // gather element data
509  // UGLBGatherScatter::template gather<0>(this->leafView(), dataHandle);
510 
511  // gather node data
512  UGLBGatherScatter::template gather<dim>(this->leafView(), dataHandle);
513 #endif
514 
515  // the load balancing step now also attaches
516  // the data to the entities and distributes it
517  loadBalance();
518 
519 #ifdef ModelP
520  // scatter element data
521  // UGLBGatherScatter::template scatter<0>(this->leafView(), dataHandle);
522 
523  // scatter node data
524  UGLBGatherScatter::template scatter<dim>(this->leafView(), dataHandle);
525 #endif
526 
527  return true;
528 #endif // HAVE_UG_PATCH10
529  }
530 
547  bool loadBalance(int strategy, int minlevel, int depth, int maxlevel, int minelement);
548 
559  template<class DataHandle>
560  void communicate (DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
561  {
562 #ifdef ModelP
563  typedef typename UGGrid::LevelGridView LevelGridView;
564 
565  for (int curCodim = 0; curCodim <= dim; ++curCodim) {
566  if (!dataHandle.contains(dim, curCodim))
567  continue;
568 
569  if (curCodim == 0)
570  communicateUG_<LevelGridView, DataHandle, 0>(this->levelGridView(level), level, dataHandle, iftype, dir);
571  else if (curCodim == dim)
572  communicateUG_<LevelGridView, DataHandle, dim>(this->levelGridView(level), level, dataHandle, iftype, dir);
573  else if (curCodim == dim - 1)
574  communicateUG_<LevelGridView, DataHandle, dim-1>(this->levelGridView(level), level, dataHandle, iftype, dir);
575  else if (curCodim == 1)
576  communicateUG_<LevelGridView, DataHandle, 1>(this->levelGridView(level), level, dataHandle, iftype, dir);
577  else
578  DUNE_THROW(NotImplemented,
579  className(*this) << "::communicate(): Not "
580  "supported for dim " << dim << " and codim " << curCodim);
581  }
582 #endif // ModelP
583  }
584 
594  template<class DataHandle>
595  void communicate(DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir) const
596  {
597 #ifdef ModelP
598  typedef typename UGGrid::LeafGridView LeafGridView;
599 
600  for (int curCodim = 0; curCodim <= dim; ++curCodim) {
601  if (!dataHandle.contains(dim, curCodim))
602  continue;
603  int level = -1;
604  if (curCodim == 0)
605  communicateUG_<LeafGridView, DataHandle, 0>(this->leafView(), level, dataHandle, iftype, dir);
606  else if (curCodim == dim)
607  communicateUG_<LeafGridView, DataHandle, dim>(this->leafView(), level, dataHandle, iftype, dir);
608  else if (curCodim == dim - 1)
609  communicateUG_<LeafGridView, DataHandle, dim-1>(this->leafView(), level, dataHandle, iftype, dir);
610  else if (curCodim == 1)
611  communicateUG_<LeafGridView, DataHandle, 1>(this->leafView(), level, dataHandle, iftype, dir);
612  else
613  DUNE_THROW(NotImplemented,
614  className(*this) << "::communicate(): Not "
615  "supported for dim " << dim << " and codim " << curCodim);
616  }
617 #endif // ModelP
618  }
619 
622  {
623  return ccobj_;
624  }
625 
626  protected:
627 #ifdef ModelP
628  template <class GridView, class DataHandle, int codim>
629  void communicateUG_(const GridView& gv, int level,
630  DataHandle &dataHandle,
631  InterfaceType iftype,
632  CommunicationDirection dir) const
633  {
634  typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
635  // Translate the communication direction from Dune-Speak to UG-Speak
636  if (dir==ForwardCommunication)
637  ugIfDir = UG_NS<dim>::IF_FORWARD();
638  else
639  ugIfDir = UG_NS<dim>::IF_BACKWARD();
640 
641  typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
642  UGMsgBuf::duneDataHandle_ = &dataHandle;
643 
644  UGMsgBuf::level = level;
645 
646  std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
647  findDDDInterfaces_(ugIfs, iftype, codim);
648 
649  unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
650  if (!bufSize)
651  return; // we don't need to communicate if we don't have any data!
652  for (unsigned i=0; i < ugIfs.size(); ++i)
653  UG_NS<dim>::DDD_IFOneway(ugIfs[i],
654  ugIfDir,
655  bufSize,
656  &UGMsgBuf::ugGather_,
657  &UGMsgBuf::ugScatter_);
658  }
659 
660  void findDDDInterfaces_(std::vector<typename UG_NS<dim>::DDD_IF > &dddIfaces,
661  InterfaceType iftype,
662  int codim) const
663  {
664  dddIfaces.clear();
665  if (codim == 0)
666  {
667  switch (iftype) {
669  // do not communicate anything: Elements can not be in
670  // the interior of two processes at the same time
671  return;
673  // The communicated elements are in the sender's
674  // interior and it does not matter what they are for
675  // the receiver
676  dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
677  return;
678  case All_All_Interface :
679  // It does neither matter what the communicated
680  // elements are for sender nor for the receiver. If
681  // they are seen by these two processes, data is send
682  // and received.
683  dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
684  return;
685  default :
686  DUNE_THROW(GridError,
687  "Element communication not supported for "
688  "interfaces of type "
689  << iftype);
690  }
691  }
692  else if (codim == dim)
693  {
694  switch (iftype)
695  {
697  dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
698  return;
700  dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
701  dddIfaces.push_back(UG_NS<dim>::NodeIF());
702  return;
703  case All_All_Interface :
704  dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
705  return;
706  default :
707  DUNE_THROW(GridError,
708  "Node communication not supported for "
709  "interfaces of type "
710  << iftype);
711  }
712  }
713  else if (codim == dim-1)
714  {
715  switch (iftype)
716  {
718  dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
719  return;
721  dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
722  // Is the following line needed or not?
723  // dddIfaces.push_back(UG_NS<dim>::EdgeIF());
724  return;
725  default :
726  DUNE_THROW(GridError,
727  "Edge communication not supported for "
728  "interfaces of type "
729  << iftype);
730  }
731  }
732  else if (codim == 1)
733  {
734  switch (iftype)
735  {
738  dddIfaces.push_back(UG_NS<dim>::BorderVectorSymmIF());
739  return;
740  default :
741  DUNE_THROW(GridError,
742  "Face communication not supported for "
743  "interfaces of type "
744  << iftype);
745  }
746  }
747  else
748  {
749  DUNE_THROW(GridError,
750  "Communication for codim "
751  << codim
752  << " entities is not yet supported "
753  << " by the DUNE UGGrid interface!");
754  }
755  };
756 #endif // ModelP
757 
758  public:
759  // **********************************************************
760  // End of Interface Methods
761  // **********************************************************
762 
770  void getChildrenOfSubface(const typename Traits::template Codim<0>::EntityPointer & e,
771  int elementSide,
772  int maxl,
773  std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements,
774  std::vector<unsigned char>& childElementSides) const;
775 
782  };
783 
785  enum ClosureType {
790  };
791 
794  refinementType_ = type;
795  }
796 
799  closureType_ = type;
800  }
801 
808  static void setDefaultHeapSize(unsigned size) {
809  heapSize_ = size;
810  }
811 
815  void setPosition(const typename Traits::template Codim<dim>::EntityPointer& e,
816  const FieldVector<double, dim>& pos);
817 
822  void globalRefine(int n);
823 
828  void saveState(const std::string& filename) const;
829 
834  void loadState(const std::string& filename);
835 
836  private:
838  typename UG_NS<dim>::MultiGrid* multigrid_;
839 
841  CollectiveCommunication<UGGrid> ccobj_;
842 
848  void setIndices(bool setLevelZero,
849  std::vector<unsigned int>* nodePermutation);
850 
851  // Each UGGrid object has a unique name to identify it in the
852  // UG environment structure
853  std::string name_;
854 
855  // Our set of level indices
856  std::vector<shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
857 
858  UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
859 
860  // One id set implementation
861  // Used for both the local and the global UGGrid id sets
862  UGGridIdSet<const UGGrid<dim> > idSet_;
863 
865  RefinementType refinementType_;
866 
868  ClosureType closureType_;
869 
877  static int numOfUGGrids;
878 
884  bool someElementHasBeenMarkedForRefinement_;
885 
891  bool someElementHasBeenMarkedForCoarsening_;
892 
897  static unsigned int heapSize_;
898 
900  std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
901 
907  unsigned int numBoundarySegments_;
908 
909  }; // end Class UGGrid
910 
911  namespace Capabilities
912  {
928  template<int dim>
929  struct hasEntity< UGGrid<dim>, 0>
930  {
931  static const bool v = true;
932  };
933 
937  template<int dim>
938  struct hasEntity< UGGrid<dim>, dim>
939  {
940  static const bool v = true;
941  };
942 
946  template<int dim>
947  struct isParallel< UGGrid<dim> >
948  {
949 #ifdef ModelP
950  static const bool v = true;
951 #else // !ModelP
952  static const bool v = false;
953 #endif // !ModelP
954  };
955 
959  template<int dim>
961  {
962  static const bool v = true;
963  };
964 
968  template<int dim>
970  {
971  static const bool v = false;
972  };
973 
974  }
975 
976 } // namespace Dune
977 
978 #endif // HAVE_UG || DOXYGEN
979 #endif // DUNE_UGGRID_HH
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:435
void loadState(const std::string &filename)
Read entire grid hierarchy from disk.
friend class UGGridEntityPointer
Definition: uggrid.hh:234
size_t numBoundarySegments() const
Return the number of boundary segments.
Definition: uggrid.hh:343
unsigned int overlapSize(int codim) const
Size of the overlap on the leaf level.
Definition: uggrid.hh:467
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: uggrid.hh:337
Partition< All_Partition >::LeafGridView LeafGridView
Definition: common/grid.hh:427
No closure, results in nonconforming meshes.
Definition: uggrid.hh:789
static ReturnImplementationType< InterfaceType >::ImplementationType & getRealImplementation(InterfaceType &i)
return real implementation of interface class
Definition: common/grid.hh:1223
bool loadBalance(DataHandle &dataHandle)
Distributes the grid and some data over the available nodes in a distributed machine.
Definition: uggrid.hh:502
New level consists only of the refined elements and the closure.
Definition: uggrid.hh:779
static const bool v
Definition: common/capabilities.hh:57
Definition: defaultgridview.hh:223
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: uggrid.hh:307
const CollectiveCommunication< UGGrid > & comm() const
Definition: uggrid.hh:621
bool loadBalance()
Default load balancing.
Definition: uggrid.hh:490
send interior and border, receive all entities
Definition: gridenums.hh:82
unsigned int overlapSize(int level, int codim) const
Size of the overlap on a given level.
Definition: uggrid.hh:477
communicate as given in InterfaceType
Definition: gridenums.hh:165
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: uggrid.hh:301
static const bool v
Definition: common/capabilities.hh:66
Index Set Interface base class.
Definition: common/grid.hh:359
void communicate(DataHandle &dataHandle, InterfaceType iftype, CommunicationDirection dir) const
The communication interface for all codims on the leaf level.
Definition: uggrid.hh:595
Definition: uggrid.hh:131
Base class for grid boundary segments of arbitrary geometry.
unsigned int ghostSize(int level, int codim) const
Size of the ghost cell layer on a given level.
Definition: uggrid.hh:482
Traits::template Partition< pitype >::LevelGridView levelGridView(int level) const
View for a grid level.
Definition: common/grid.hh:1064
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:263
void setRefinementType(RefinementType type)
Sets the type of grid refinement.
Definition: uggrid.hh:793
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: uggrid.hh:370
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: uggrid.hh:350
static void setDefaultHeapSize(unsigned size)
Sets the default heap size.
Definition: uggrid.hh:808
Traits::template Codim< codim >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: uggrid.hh:289
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:543
int maxLevel() const
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:164
void globalRefine(int n)
Does uniform refinement.
bool mark(int refCount, const typename Traits::template Codim< 0 >::Entity &e)
Mark element for refinement.
send/receive interior and border entities
Definition: gridenums.hh:81
Types for GridView.
Definition: common/grid.hh:418
~UGGrid()
Destructor.
Specialize with 'true' if implementation guarantees conforming level grids. (default=false) ...
Definition: common/capabilities.hh:86
New level consists of the refined elements and the unrefined ones, too.
Definition: uggrid.hh:781
int size(int level, int codim) const
Number of grid entities per level and codim.
UGGridFamily< dim >::Traits Traits
Definition: uggrid.hh:252
The specialization of the generic GridFactory for UGGrid.
[ provides Dune::Grid ]
Definition: uggrid.hh:202
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.
void setClosureType(ClosureType type)
Sets the type of grid refinement closure.
Definition: uggrid.hh:798
Traits::template Codim< codim >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: uggrid.hh:295
void saveState(const std::string &filename) const
Save entire grid hierarchy to disk.
Standard red/green refinement.
Definition: uggrid.hh:787
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: uggrid.hh:362
void communicate(DataHandle &dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
The communication interface for all codims on a given level.
Definition: uggrid.hh:560
int size(int codim) const
number of leaf entities per codim in this process
Definition: uggrid.hh:325
Grid view abstract base classInterface class for a view on grids. Grids return two types of view...
Definition: common/gridview.hh:56
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: uggrid.hh:331
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:16
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false) ...
Definition: common/capabilities.hh:95
GridFamily::Traits::template Codim< cd >::EntityPointer EntityPointer
A type that is a model of Dune::EntityPointer<cd,dim,...>.
Definition: common/grid.hh:447
ClosureType
Decide whether to add a green closure to locally refined grid sections or not.
Definition: uggrid.hh:785
Partition< All_Partition >::LevelGridView LevelGridView
View types for All_Partition.
Definition: common/grid.hh:426
A set of traits classes to store static information about grid implementation.
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
Query whether element is marked for refinement.
Specialize with 'true' if implementation supports parallelism. (default=false)
Definition: common/capabilities.hh:64
static const bool v
Definition: common/capabilities.hh:88
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 > >, DefaultLevelGridViewTraits, DefaultLeafGridViewTraits, UGGridEntitySeed, UGGridLocalGeometry > Traits
Definition: uggrid.hh:155
void postAdapt()
Clean up refinement markers.
Specialize with 'true' for all codims that a grid implements entities for. (default=false) ...
Definition: common/capabilities.hh:55
bool adapt()
Triggers the grid refinement process.
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: uggrid.hh:356
friend class UGGridLeafIterator
Definition: uggrid.hh:230
Different resources needed by all grid implementations.
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:1260
static const bool v
Definition: common/capabilities.hh:97
RefinementType
The different forms of grid refinement that UG supports.
Definition: uggrid.hh:777
void setPosition(const typename Traits::template Codim< dim >::EntityPointer &e, const FieldVector< double, dim > &pos)
Sets a vertex to a new position.
Definition: common/geometry.hh:24
Traits::template Partition< pitype >::LeafGridView leafView() const
View for the leaf grid.
Definition: common/grid.hh:1039
UGGridFamily< dim > GridFamily
type of the used GridFamily for this grid
Definition: uggrid.hh:249
Traits::template Codim< codim >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:80
Definition: defaultgridview.hh:23
UGGrid()
Default constructor.
UG::DOUBLE ctype
The type used to store coordinates.
Definition: uggrid.hh:255
friend class UGGridLevelIterator
Definition: uggrid.hh:232
Id Set Interface.
Definition: common/grid.hh:360
Traits::template Codim< codim >::LevelIterator lend(int level) const
one past the end on this level
bool preAdapt()
returns true, if some elements might be coarsend during grid adaption, here always returns true ...
unsigned int ghostSize(int codim) const
Size of the ghost cell layer on the leaf level.
Definition: uggrid.hh:472
Traits::template Codim< Seed::codimension >::EntityPointer entityPointer(const Seed &seed) const
Create an EntityPointer from an EntitySeed.
Definition: uggrid.hh:314
IndexType size(GeometryType type) const
Return total number of entities of given geometry type in entity set .
Definition: indexidset.hh:204
send all and receive all entities
Definition: gridenums.hh:85