Go to the documentation of this file.
13 #include <dune/common/classname.hh>
14 #include <dune/common/parallel/collectivecommunication.hh>
15 #include <dune/common/exceptions.hh>
16 #include <dune/common/parallel/mpihelper.hh>
17 #include <dune/common/deprecated.hh>
23 #if HAVE_UG || DOXYGEN
26 #include <dune/common/parallel/mpicollectivecommunication.hh>
52 #ifdef UG_USE_NEW_DIMENSION_DEFINES
58 #include "uggrid/ugincludes.hh"
63 #include "uggrid/ugwrapper.hh"
68 #include "uggrid/ug_undefs.hh"
69 #ifdef UG_USE_NEW_DIMENSION_DEFINES
91 #ifdef UG_USE_NEW_DIMENSION_DEFINES
97 #include "uggrid/ugincludes.hh"
102 #include "uggrid/ugwrapper.hh"
106 #include "uggrid/ug_undefs.hh"
108 #ifdef UG_USE_NEW_DIMENSION_DEFINES
116 #include "uggrid/uggridgeometry.hh"
117 #include "uggrid/uggridlocalgeometry.hh"
118 #include "uggrid/uggridentity.hh"
119 #include "uggrid/uggridentityseed.hh"
120 #include "uggrid/uggridintersections.hh"
121 #include "uggrid/uggridintersectioniterators.hh"
122 #include "uggrid/uggridleveliterator.hh"
123 #include "uggrid/uggridleafiterator.hh"
124 #include "uggrid/uggridhieriterator.hh"
125 #include "uggrid/uggridindexsets.hh"
126 #include <dune/grid/uggrid/uggridviews.hh>
128 #include "uggrid/ugmessagebuffer.hh"
129 #include "uggrid/uglbgatherscatter.hh"
136 template <
class DataHandle,
int Gr
idDim,
int codim>
137 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
139 template <
class DataHandle,
int Gr
idDim,
int codim>
140 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
147 class CollectiveCommunication<
Dune::UGGrid<dim> > :
148 public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
150 typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
152 CollectiveCommunication()
153 : ParentType(MPIHelper::getCommunicator())
165 UGGridLeafIntersection,
166 UGGridLevelIntersection,
167 UGGridLeafIntersectionIterator,
168 UGGridLevelIntersectionIterator,
169 UGGridHierarchicIterator,
171 UGGridLevelIndexSet< const UGGrid<dim> >,
172 UGGridLeafIndexSet< const UGGrid<dim> >,
173 UGGridIdSet< const UGGrid<dim> >,
174 typename UG_NS<dim>::UG_ID_TYPE,
175 UGGridIdSet< const UGGrid<dim> >,
176 typename UG_NS<dim>::UG_ID_TYPE,
177 CollectiveCommunication<Dune::UGGrid<dim> >,
178 UGGridLevelGridViewTraits,
179 UGGridLeafGridViewTraits,
234 friend class UGGridGeometry<0,dim,const
UGGrid<dim> >;
235 friend class UGGridGeometry<dim,dim,const
UGGrid<dim> >;
236 friend class UGGridGeometry<1,2,const
UGGrid<dim> >;
237 friend class UGGridGeometry<2,3,const
UGGrid<dim> >;
239 friend class UGGridEntity <0,dim,const
UGGrid<dim> >;
240 friend class UGGridEntity <dim,dim,const
UGGrid<dim> >;
241 friend class UGGridHierarchicIterator<const
UGGrid<dim> >;
242 friend class UGGridLeafIntersection<const
UGGrid<dim> >;
243 friend class UGGridLevelIntersection<const
UGGrid<dim> >;
244 friend class UGGridLeafIntersectionIterator<const
UGGrid<dim> >;
245 friend class UGGridLevelIntersectionIterator<const
UGGrid<dim> >;
247 friend class UGGridLevelIndexSet<const
UGGrid<dim> >;
248 friend class UGGridLeafIndexSet<const
UGGrid<dim> >;
249 friend class UGGridIdSet<const
UGGrid<dim> >;
250 template <
class Gr
idImp_>
252 template <
class Gr
idImp_>
258 friend class UGLBGatherScatter;
261 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
263 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
267 static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
305 template <typename Seed>
309 const int codim = Seed::codimension;
315 int size (
int level,
int codim)
const;
339 return numBoundarySegments_;
358 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
359 return *levelIndexSets_[level];
365 return leafIndexSet_;
441 typename UG_NS<dim>::RefinementRule rule,
465 return (codim==0) ? 1 : 0;
475 return (codim==0) ? 1 : 0;
485 template<
class DataHandle>
493 if (dataHandle.contains(dim,dim))
494 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
506 if (dataHandle.contains(dim,dim))
507 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
551 bool loadBalance(
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel);
562 template<
class DataHandle>
563 bool loadBalance (
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel, DataHandle& dataHandle)
570 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
582 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
599 template<
class DataHandle>
605 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
606 if (!dataHandle.contains(dim, curCodim))
610 communicateUG_<LevelGridView, DataHandle, 0>(this->
levelGridView(level), level, dataHandle, iftype, dir);
611 else if (curCodim == dim)
612 communicateUG_<LevelGridView, DataHandle, dim>(this->
levelGridView(level), level, dataHandle, iftype, dir);
613 else if (curCodim == dim - 1)
615 else if (curCodim == 1)
616 communicateUG_<LevelGridView, DataHandle, 1>(this->
levelGridView(level), level, dataHandle, iftype, dir);
618 DUNE_THROW(NotImplemented,
619 className(*
this) <<
"::communicate(): Not "
620 "supported for dim " << dim <<
" and codim " << curCodim);
634 template<
class DataHandle>
640 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
641 if (!dataHandle.contains(dim, curCodim))
645 communicateUG_<LeafGridView, DataHandle, 0>(this->
leafGridView(), level, dataHandle, iftype, dir);
646 else if (curCodim == dim)
647 communicateUG_<LeafGridView, DataHandle, dim>(this->
leafGridView(), level, dataHandle, iftype, dir);
648 else if (curCodim == dim - 1)
650 else if (curCodim == 1)
651 communicateUG_<LeafGridView, DataHandle, 1>(this->
leafGridView(), level, dataHandle, iftype, dir);
653 DUNE_THROW(NotImplemented,
654 className(*
this) <<
"::communicate(): Not "
655 "supported for dim " << dim <<
" and codim " << curCodim);
661 const CollectiveCommunication<UGGrid>&
comm ()
const
668 template <
class Gr
idView,
class DataHandle,
int codim>
669 void communicateUG_(
const GridView& gv,
int level,
670 DataHandle &dataHandle,
674 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
677 ugIfDir = UG_NS<dim>::IF_FORWARD();
679 ugIfDir = UG_NS<dim>::IF_BACKWARD();
681 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
682 UGMsgBuf::duneDataHandle_ = &dataHandle;
684 UGMsgBuf::level = level;
686 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
687 findDDDInterfaces_(ugIfs, iftype, codim);
689 unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
692 for (
unsigned i=0; i < ugIfs.size(); ++i)
693 UG_NS<dim>::DDD_IFOneway(ugIfs[i],
696 &UGMsgBuf::ugGather_,
697 &UGMsgBuf::ugScatter_);
700 void findDDDInterfaces_(std::vector<
typename UG_NS<dim>::DDD_IF > &dddIfaces,
716 dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
723 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
726 DUNE_THROW(GridError,
727 "Element communication not supported for "
728 "interfaces of type "
732 else if (codim == dim)
737 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
740 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
741 dddIfaces.push_back(UG_NS<dim>::NodeIF());
744 dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
747 DUNE_THROW(GridError,
748 "Node communication not supported for "
749 "interfaces of type "
753 else if (codim == dim-1)
758 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
761 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
766 dddIfaces.push_back(UG_NS<dim>::EdgeSymmVHIF());
769 DUNE_THROW(GridError,
770 "Edge communication not supported for "
771 "interfaces of type "
781 dddIfaces.push_back(UG_NS<dim>::BorderVectorSymmIF());
784 DUNE_THROW(GridError,
785 "Face communication not supported for "
786 "interfaces of type "
792 DUNE_THROW(GridError,
793 "Communication for codim "
795 <<
" entities is not yet supported "
796 <<
" by the DUNE UGGrid interface!");
817 std::vector<unsigned char>& childElementSides)
const;
837 refinementType_ = type;
859 const FieldVector<double, dim>& pos);
871 void saveState(
const std::string& filename)
const;
877 void loadState(
const std::string& filename);
881 typename UG_NS<dim>::MultiGrid* multigrid_;
884 CollectiveCommunication<UGGrid> ccobj_;
891 void setIndices(
bool setLevelZero,
892 std::vector<unsigned int>* nodePermutation);
899 std::vector<std::shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
901 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
905 UGGridIdSet<const UGGrid<dim> > idSet_;
920 static int numOfUGGrids;
927 bool someElementHasBeenMarkedForRefinement_;
934 bool someElementHasBeenMarkedForCoarsening_;
940 static unsigned int heapSize_;
943 std::vector<std::shared_ptr<BoundarySegment<dim> > > boundarySegments_;
950 unsigned int numBoundarySegments_;
954 namespace Capabilities
974 static const bool v =
true;
983 static const bool v =
true;
992 static const bool v =
true;
1001 static const bool v =
false;
1008 #endif // HAVE_UG || DOXYGEN
1009 #endif // DUNE_UGGRID_HH
bool loadBalance(DataHandle &dataHandle)
Distributes the grid and some data over the available nodes in a distributed machine.
Definition: uggrid.hh:486
Different resources needed by all grid implementations.
UGGrid()
Default constructor.
const CollectiveCommunication< UGGrid > & comm() const
Definition: uggrid.hh:661
int size(int level, int codim) const
Number of grid entities per level and codim.
New level consists of the refined elements and the unrefined ones, too.
Definition: uggrid.hh:824
send/receive interior and border entities
Definition: gridenums.hh:85
RefinementType
The different forms of grid refinement that UG supports.
Definition: uggrid.hh:820
Traits::LevelGridView levelGridView(int level) const
View for a grid level for All_Partition.
Definition: common/grid.hh:936
The specialization of the generic GridFactory for UGGrid.
No closure, results in nonconforming meshes.
Definition: uggrid.hh:832
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:263
bool loadBalance()
default implementation of load balance does nothing and returns false
Definition: common/grid.hh:1034
Include standard header files.
Definition: agrid.hh:59
UG::DOUBLE ctype
The type used to store coordinates.
Definition: uggrid.hh:285
~UGGrid() noexcept(false)
Destructor.
unsigned int overlapSize(int level, int codim) const
Size of the overlap on a given level.
Definition: uggrid.hh:469
unsigned int ghostSize(int level, int codim) const
Size of the ghost cell layer on a given level.
Definition: uggrid.hh:474
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: uggrid.hh:343
GridTraits< dim, dim, Dune::UGGrid< dim >, UGGridGeometry, UGGridEntity, 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:182
int size(int codim) const
number of leaf entities per codim in this process
Definition: uggrid.hh:318
GridFamily::Traits::LevelGridView LevelGridView
type of view for level grid
Definition: common/grid.hh:406
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: uggrid.hh:363
IdSet< const GridImp, LocalIdSetImp, LIDType > LocalIdSet
The type of the local id set.
Definition: common/grid.hh:1230
Front-end for the grid manager of the finite element toolbox UG.
Definition: uggrid.hh:230
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Create an Entity from an EntitySeed.
Definition: uggrid.hh:307
void communicate(DataHandle &dataHandle, InterfaceType iftype, CommunicationDirection dir) const
The communication interface for all codims on the leaf level.
Definition: uggrid.hh:635
IdSet< const GridImp, GlobalIdSetImp, GIDType > GlobalIdSet
The type of the global id set.
Definition: common/grid.hh:1228
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
friend class UGGridLeafGridView
Definition: uggrid.hh:251
void globalRefine(int n)
Does uniform refinement.
size_t numBoundarySegments() const
Return the number of boundary segments.
Definition: uggrid.hh:336
communicate as given in InterfaceType
Definition: gridenums.hh:169
friend class UGGridLevelIterator
Definition: uggrid.hh:264
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
Query whether element is marked for refinement.
Definition: uggrid.hh:159
Wrapper class for entities.
Definition: common/entity.hh:61
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false)
Definition: common/capabilities.hh:86
void setRefinementType(RefinementType type)
Sets the type of grid refinement.
Definition: uggrid.hh:836
UGGridFamily< dim > GridFamily
type of the used GridFamily for this grid
Definition: uggrid.hh:279
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:414
Specialize with 'true' for all codims that a grid implements entities for. (default=false)
Definition: common/capabilities.hh:55
ClosureType
Decide whether to add a green closure to locally refined grid sections or not.
Definition: uggrid.hh:828
static void setDefaultHeapSize(unsigned size)
Sets the default heap size.
Definition: uggrid.hh:851
unsigned int overlapSize(int codim) const
Size of the overlap on the leaf level.
Definition: uggrid.hh:459
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:1156
unsigned int ghostSize(int codim) const
Size of the ghost cell layer on the leaf level.
Definition: uggrid.hh:464
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
void communicate(DataHandle &dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
The communication interface for all codims on a given level.
Definition: uggrid.hh:600
static const bool v
Definition: common/capabilities.hh:88
send interior and border, receive all entities
Definition: gridenums.hh:86
IndexSet< const GridImp, LeafIndexSetImp > LeafIndexSet
The type of the leaf index set.
Definition: common/grid.hh:1226
Specialize with 'true' if implementation guarantees conforming level grids. (default=false)
Definition: common/capabilities.hh:77
IndexType size(GeometryType type) const
Return total number of entities of given geometry type in entity set .
Definition: indexidset.hh:215
IndexSet< const GridImp, LevelIndexSetImp > LevelIndexSet
The type of the level index set.
Definition: common/grid.hh:1224
A set of traits classes to store static information about grid implementation.
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: uggrid.hh:330
send all and receive all entities
Definition: gridenums.hh:89
Standard red/green refinement.
Definition: uggrid.hh:830
void postAdapt()
Clean up refinement markers.
Definition: common/geometry.hh:24
Traits::LeafGridView leafGridView() const
View for the leaf grid for All_Partition.
Definition: common/grid.hh:944
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:1119
static const bool v
Definition: common/capabilities.hh:57
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: uggrid.hh:355
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: uggrid.hh:349
New level consists only of the refined elements and the closure.
Definition: uggrid.hh:822
void setPosition(const typename Traits::template Codim< dim >::Entity &e, const FieldVector< double, dim > &pos)
Sets a vertex to a new position.
GridFamily::Traits::template Codim< cd >::Entity Entity
A type that is a model of a Dune::Entity<cd,dim,...>.
Definition: common/grid.hh:423
UGGridFamily< dim > ::Traits Traits
The traits of this class.
Definition: common/grid.hh:933
static const bool v
Definition: common/capabilities.hh:79
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: uggrid.hh:324
unsigned int Rank
The type used for process ranks.
Definition: uggrid.hh:288
void setClosureType(ClosureType type)
Sets the type of grid refinement closure.
Definition: uggrid.hh:841
friend class UGGridLeafIterator
Definition: uggrid.hh:262
bool mark(int refCount, const typename Traits::template Codim< 0 >::Entity &e)
Mark element for refinement.
void getChildrenOfSubface(const typename Traits::template Codim< 0 >::Entity &e, int elementSide, int maxl, std::vector< typename Traits::template Codim< 0 >::Entity > &childElements, std::vector< unsigned char > &childElementSides) const
Rudimentary substitute for a hierarchic iterator on faces.
bool preAdapt()
returns true, if some elements might be coarsend during grid adaption, here always returns true
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:16
void loadState(const std::string &filename)
Read entire grid hierarchy from disk.
bool adapt()
Triggers the grid refinement process.
void saveState(const std::string &filename) const
Save entire grid hierarchy to disk.
Base class for grid boundary segments of arbitrary geometry.
Grid view abstract base class.
Definition: common/gridview.hh:59
GridFamily::Traits::LeafGridView LeafGridView
type of view for leaf grid
Definition: common/grid.hh:404
UGGridFamily< dim >::Traits Traits
Definition: uggrid.hh:282
friend class UGGridLevelGridView
Definition: uggrid.hh:253
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
bool loadBalance(const std::vector< Rank > &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:563