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>
21 #if HAVE_UG || DOXYGEN
24 #include <dune/common/parallel/mpicollectivecommunication.hh>
46 #include "uggrid/ugincludes.hh"
51 #include "uggrid/ugwrapper.hh"
56 #include "uggrid/ug_undefs.hh"
73 #include "uggrid/ugincludes.hh"
78 #include "uggrid/ugwrapper.hh"
82 #include "uggrid/ug_undefs.hh"
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"
100 #include "uggrid/ugmessagebuffer.hh"
101 #include "uggrid/uglbgatherscatter.hh"
108 template <
class DataHandle,
int Gr
idDim,
int codim>
109 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
111 template <
class DataHandle,
int Gr
idDim,
int codim>
112 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
119 class CollectiveCommunication<Dune::UGGrid<dim> > :
120 public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
122 typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
124 CollectiveCommunication()
125 : ParentType(MPIHelper::getCommunicator())
138 UGGridLeafIntersection,
139 UGGridLevelIntersection,
140 UGGridLeafIntersectionIterator,
141 UGGridLevelIntersectionIterator,
142 UGGridHierarchicIterator,
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> >,
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> >;
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> >;
219 friend class UGGridLevelIndexSet<const
UGGrid<dim> >;
220 friend class UGGridLeafIndexSet<const
UGGrid<dim> >;
221 friend class UGGridIdSet<const
UGGrid<dim> >;
226 friend class UGLBGatherScatter;
229 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
231 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
233 template <
int codim_,
class Gr
idImp_>
237 dune_static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
280 template<
int codim, PartitionIteratorType PiType>
284 template<
int codim, PartitionIteratorType PiType>
296 return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
300 template<
int codim, PartitionIteratorType PiType>
306 template<
int codim, PartitionIteratorType PiType>
308 return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
312 template <
typename Seed>
313 typename Traits::template Codim<Seed::codimension>::EntityPointer
316 enum {codim = Seed::codimension};
322 int size (
int level,
int codim)
const;
346 return numBoundarySegments_;
365 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
366 return *levelIndexSets_[level];
372 return leafIndexSet_;
390 bool mark(
int refCount,
const typename Traits::template Codim<0>::Entity & e );
448 bool mark(
const typename Traits::template Codim<0>::Entity & e,
449 typename UG_NS<dim>::RefinementRule rule,
453 int getMark(
const typename Traits::template Codim<0>::Entity& e)
const;
473 return (codim==0) ? 1 : 0;
483 return (codim==0) ? 1 : 0;
501 template<
class DataHandle>
505 DUNE_THROW(NotImplemented,
"load balancing with data attached");
512 UGLBGatherScatter::template gather<dim>(this->
leafView(), dataHandle);
524 UGLBGatherScatter::template scatter<dim>(this->
leafView(), dataHandle);
528 #endif // HAVE_UG_PATCH10
547 bool loadBalance(
int strategy,
int minlevel,
int depth,
int maxlevel,
int minelement);
559 template<
class DataHandle>
565 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
566 if (!dataHandle.contains(dim, curCodim))
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)
575 else if (curCodim == 1)
576 communicateUG_<LevelGridView, DataHandle, 1>(this->
levelGridView(level), level, dataHandle, iftype, dir);
578 DUNE_THROW(NotImplemented,
579 className(*
this) <<
"::communicate(): Not "
580 "supported for dim " << dim <<
" and codim " << curCodim);
594 template<
class DataHandle>
600 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
601 if (!dataHandle.contains(dim, curCodim))
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);
613 DUNE_THROW(NotImplemented,
614 className(*
this) <<
"::communicate(): Not "
615 "supported for dim " << dim <<
" and codim " << curCodim);
628 template <
class Gr
idView,
class DataHandle,
int codim>
629 void communicateUG_(
const GridView& gv,
int level,
630 DataHandle &dataHandle,
634 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
637 ugIfDir = UG_NS<dim>::IF_FORWARD();
639 ugIfDir = UG_NS<dim>::IF_BACKWARD();
641 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
642 UGMsgBuf::duneDataHandle_ = &dataHandle;
644 UGMsgBuf::level = level;
646 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
647 findDDDInterfaces_(ugIfs, iftype, codim);
649 unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
652 for (
unsigned i=0; i < ugIfs.size(); ++i)
653 UG_NS<dim>::DDD_IFOneway(ugIfs[i],
656 &UGMsgBuf::ugGather_,
657 &UGMsgBuf::ugScatter_);
660 void findDDDInterfaces_(std::vector<
typename UG_NS<dim>::DDD_IF > &dddIfaces,
676 dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
683 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
686 DUNE_THROW(GridError,
687 "Element communication not supported for "
688 "interfaces of type "
692 else if (codim == dim)
697 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
700 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
701 dddIfaces.push_back(UG_NS<dim>::NodeIF());
704 dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
707 DUNE_THROW(GridError,
708 "Node communication not supported for "
709 "interfaces of type "
713 else if (codim == dim-1)
718 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
721 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
726 DUNE_THROW(GridError,
727 "Edge communication not supported for "
728 "interfaces of type "
738 dddIfaces.push_back(UG_NS<dim>::BorderVectorSymmIF());
741 DUNE_THROW(GridError,
742 "Face communication not supported for "
743 "interfaces of type "
749 DUNE_THROW(GridError,
750 "Communication for codim "
752 <<
" entities is not yet supported "
753 <<
" by the DUNE UGGrid interface!");
773 std::vector<
typename Traits::template Codim<0>::EntityPointer>& childElements,
774 std::vector<unsigned char>& childElementSides)
const;
794 refinementType_ = type;
815 void setPosition(
const typename Traits::template Codim<dim>::EntityPointer& e,
816 const FieldVector<double, dim>& pos);
828 void saveState(
const std::string& filename)
const;
834 void loadState(
const std::string& filename);
838 typename UG_NS<dim>::MultiGrid* multigrid_;
841 CollectiveCommunication<UGGrid> ccobj_;
848 void setIndices(
bool setLevelZero,
849 std::vector<unsigned int>* nodePermutation);
856 std::vector<shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
858 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
862 UGGridIdSet<const UGGrid<dim> > idSet_;
877 static int numOfUGGrids;
884 bool someElementHasBeenMarkedForRefinement_;
891 bool someElementHasBeenMarkedForCoarsening_;
897 static unsigned int heapSize_;
900 std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
907 unsigned int numBoundarySegments_;
911 namespace Capabilities
931 static const bool v =
true;
940 static const bool v =
true;
950 static const bool v =
true;
952 static const bool v =
false;
962 static const bool v =
true;
971 static const bool v =
false;
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
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
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