13 #include <dune/common/classname.hh>
14 #include <dune/common/parallel/communication.hh>
15 #include <dune/common/exceptions.hh>
16 #include <dune/common/parallel/mpihelper.hh>
22 #if HAVE_UG || DOXYGEN
25 #include <dune/common/parallel/mpicommunication.hh>
53 #include "uggrid/ugincludes.hh"
54 #undef DUNE_UGINCLUDES_HH
59 #include "uggrid/ugwrapper.hh"
64 #include "uggrid/ug_undefs.hh"
84 #include "uggrid/ugincludes.hh"
85 #undef DUNE_UGINCLUDES_HH
90 #include "uggrid/ugwrapper.hh"
94 #include "uggrid/ug_undefs.hh"
100 #include "uggrid/uggridgeometry.hh"
101 #include "uggrid/uggridlocalgeometry.hh"
102 #include "uggrid/uggridentity.hh"
103 #include "uggrid/uggridentityseed.hh"
104 #include "uggrid/uggridintersections.hh"
105 #include "uggrid/uggridintersectioniterators.hh"
106 #include "uggrid/uggridleveliterator.hh"
107 #include "uggrid/uggridleafiterator.hh"
108 #include "uggrid/uggridhieriterator.hh"
109 #include "uggrid/uggridindexsets.hh"
110 #include <dune/grid/uggrid/uggridviews.hh>
112 #include "uggrid/ugmessagebuffer.hh"
113 #include "uggrid/uglbgatherscatter.hh"
120 template <
class DataHandle,
int Gr
idDim,
int codim>
123 template <
class DataHandle,
int Gr
idDim,
int codim>
124 DataHandle *Dune::UGMessageBuffer<DataHandle,GridDim,codim>::duneDataHandle_ =
nullptr;
126 template <
class DataHandle,
int Gr
idDim,
int codim>
127 int Dune::UGMessageBuffer<DataHandle,GridDim,codim>::level = -1;
145 UGGridLeafIntersection,
146 UGGridLevelIntersection,
147 UGGridLeafIntersectionIterator,
148 UGGridLevelIntersectionIterator,
149 UGGridHierarchicIterator,
151 UGGridLevelIndexSet< const UGGrid<dim> >,
152 UGGridLeafIndexSet< const UGGrid<dim> >,
153 UGGridIdSet< const UGGrid<dim> >,
154 typename UG_NS<dim>::UG_ID_TYPE,
155 UGGridIdSet< const UGGrid<dim> >,
156 typename UG_NS<dim>::UG_ID_TYPE,
158 UGGridLevelGridViewTraits,
159 UGGridLeafGridViewTraits,
208 friend class UGGridGeometry<0,dim,const
UGGrid<dim> >;
209 friend class UGGridGeometry<dim,dim,const
UGGrid<dim> >;
210 friend class UGGridGeometry<1,2,const
UGGrid<dim> >;
211 friend class UGGridGeometry<2,3,const
UGGrid<dim> >;
213 friend class UGGridEntity <0,dim,const
UGGrid<dim> >;
214 friend class UGGridEntity <1,dim,const
UGGrid<dim> >;
215 friend class UGGridEntity <2,dim,const
UGGrid<dim> >;
216 friend class UGGridEntity <dim,dim,const
UGGrid<dim> >;
217 friend class UGGridHierarchicIterator<const
UGGrid<dim> >;
218 friend class UGGridLeafIntersection<const
UGGrid<dim> >;
219 friend class UGGridLevelIntersection<const
UGGrid<dim> >;
220 friend class UGGridLeafIntersectionIterator<const
UGGrid<dim> >;
221 friend class UGGridLevelIntersectionIterator<const
UGGrid<dim> >;
223 friend class UGGridLevelIndexSet<const
UGGrid<dim> >;
224 friend class UGGridLeafIndexSet<const
UGGrid<dim> >;
225 friend class UGGridIdSet<const
UGGrid<dim> >;
226 template <
class Gr
idImp_>
228 template <
class Gr
idImp_>
234 friend class UGLBGatherScatter;
237 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
239 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
243 static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
278 template <typename Seed>
282 const int codim = Seed::codimension;
288 int size (
int level,
int codim)
const;
312 return numBoundarySegments_;
331 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
332 return *levelIndexSets_[level];
338 return leafIndexSet_;
414 typename UG_NS<dim>::RefinementRule rule,
438 template<
class DataHandle>
443 if (dataHandle.contains(dim, 0))
444 UGLBGatherScatter::template gather<0>(this->
leafGridView(), dataHandle);
447 if (dataHandle.contains(dim,dim))
448 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
457 if (dataHandle.contains(dim, 0))
458 UGLBGatherScatter::template scatter<0>(this->
leafGridView(), dataHandle);
461 if (dataHandle.contains(dim,dim))
462 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
506 bool loadBalance(
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel);
517 template<
class DataHandle>
518 bool loadBalance (
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel, DataHandle& dataHandle)
522 if (dataHandle.contains(dim, 0))
523 UGLBGatherScatter::template gather<0>(this->
leafGridView(), dataHandle);
526 if (dataHandle.contains(dim,dim))
527 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
536 if (dataHandle.contains(dim, 0))
537 UGLBGatherScatter::template scatter<0>(this->
leafGridView(), dataHandle);
540 if (dataHandle.contains(dim,dim))
541 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
555 template <
int codim,
class Gr
idView,
class DataHandle>
556 void communicateUG_(
const GridView& gv,
int level,
557 DataHandle &dataHandle,
561 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
564 ugIfDir = UG_NS<dim>::IF_FORWARD();
566 ugIfDir = UG_NS<dim>::IF_BACKWARD();
568 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
569 UGMsgBuf::duneDataHandle_ = &dataHandle;
571 UGMsgBuf::level = level;
573 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs = findDDDInterfaces(iftype, codim);
575 unsigned bufSize = UGMsgBuf::ugBufferSize(gv);
578 UGMsgBuf::grid_ =
this;
579 for (
unsigned i=0; i < ugIfs.size(); ++i)
580 UG_NS<dim>::DDD_IFOneway(multigrid_->dddContext(),
584 &UGMsgBuf::ugGather_,
585 &UGMsgBuf::ugScatter_);
589 std::vector<typename UG_NS<dim>::DDD_IF> findDDDInterfaces(
InterfaceType iftype,
608 std::vector<unsigned char>& childElementSides)
const;
628 refinementType_ = type;
640 const FieldVector<double, dim>& pos);
662 typename UG_NS<dim>::MultiGrid* multigrid_;
672 void setIndices(
bool setLevelZero,
673 std::vector<unsigned int>* nodePermutation);
680 std::vector<std::shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
682 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
686 UGGridIdSet<const UGGrid<dim> > idSet_;
701 static int numOfUGGrids;
708 bool someElementHasBeenMarkedForRefinement_;
715 bool someElementHasBeenMarkedForCoarsening_;
718 std::vector<std::shared_ptr<BoundarySegment<dim> > > boundarySegments_;
725 unsigned int numBoundarySegments_;
729 namespace Capabilities
746 template<
int dim,
int codim>
749 static const bool v =
true;
757 template<
int dim,
int codim>
760 static const bool v =
false;
770 static const bool v =
true;
780 static const bool v =
true;
786 template<
int dim,
int codim>
789 static const bool v = (codim>=0 && codim<=dim);
798 static const bool v =
true;
807 static const bool v =
false;
Base class for grid boundary segments of arbitrary geometry.
The specialization of the generic GridFactory for UGGrid.
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
@ ForwardCommunication
communicate as given in InterfaceType
Definition: gridenums.hh:169
Include standard header files.
Definition: agrid.hh:58
CollectiveCommunication< No_Comm > UGCollectiveCommunication
Definition: uggrid.hh:135
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:130
Specialize with 'true' for all codims that a grid implements entities for. (default=false)
Definition: common/capabilities.hh:56
static const bool v
Definition: common/capabilities.hh:57
specialize with 'true' for all codims that a grid provides an iterator for (default=hasEntity<codim>:...
Definition: common/capabilities.hh:72
static const bool v
Definition: common/capabilities.hh:73
specialize with 'true' for all codims that a grid can communicate data on (default=false)
Definition: common/capabilities.hh:95
static const bool v
Definition: common/capabilities.hh:96
Specialize with 'true' if implementation guarantees conforming level grids. (default=false)
Definition: common/capabilities.hh:104
static const bool v
Definition: common/capabilities.hh:105
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false)
Definition: common/capabilities.hh:113
static const bool v
Definition: common/capabilities.hh:114
Wrapper class for entities.
Definition: common/entity.hh:64
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:18
Definition: common/grid.hh:851
Traits::LeafGridView leafGridView() const
View for the leaf grid for All_Partition.
Definition: common/grid.hh:871
bool loadBalance()
default implementation of load balance does nothing and returns false
Definition: common/grid.hh:937
Index Set Interface base class.
Definition: indexidset.hh:76
auto size(GeometryType type) const
Return total number of entities of given geometry type in entity set .
Definition: indexidset.hh:221
Id Set Interface.
Definition: indexidset.hh:450
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:414
GridFamily::Traits::template Codim< cd >::Entity Entity
A type that is a model of a Dune::Entity<cd,dim,...>.
Definition: common/grid.hh:422
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:984
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:312
Grid view abstract base class.
Definition: common/gridview.hh:63
Definition: uggrid.hh:140
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, UGCollectiveCommunication, UGGridLevelGridViewTraits, UGGridLeafGridViewTraits, UGGridEntitySeed, UGGridLocalGeometry > Traits
Definition: uggrid.hh:162
Front-end for the grid manager of the finite element toolbox UG3.
Definition: uggrid.hh:205
void postAdapt()
Clean up refinement markers.
size_t numBoundarySegments() const
Return the number of boundary segments.
Definition: uggrid.hh:309
void setRefinementType(RefinementType type)
Sets the type of grid refinement.
Definition: uggrid.hh:627
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
Query whether element is marked for refinement.
bool mark(int refCount, const typename Traits::template Codim< 0 >::Entity &e)
Mark element for refinement.
friend class UGGridLeafIterator
Definition: uggrid.hh:238
bool loadBalance(int minlevel=0)
Distributes this grid over the available nodes in a distributed machine.
~UGGrid() noexcept(false)
Destructor.
const UGCollectiveCommunication & comm() const
Definition: uggrid.hh:548
UGGridFamily< dim >::Traits Traits
Definition: uggrid.hh:258
friend class UGGridLevelIterator
Definition: uggrid.hh:240
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:518
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.
friend class UGGridLevelGridView
Definition: uggrid.hh:229
UGGridFamily< dim > GridFamily
type of the used GridFamily for this grid
Definition: uggrid.hh:255
bool loadBalance(const std::vector< Rank > &targetProcessors, unsigned int fromLevel)
Distribute this grid over a distributed machine.
bool mark(const typename Traits::template Codim< 0 >::Entity &e, typename UG_NS< dim >::RefinementRule rule, int side=0)
Mark method accepting a UG refinement rule.
UGGrid(UGCollectiveCommunication comm={})
Default constructor.
void setClosureType(ClosureType type)
Sets the type of grid refinement closure.
Definition: uggrid.hh:632
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: uggrid.hh:316
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: uggrid.hh:297
RefinementType
The different forms of grid refinement that UG supports.
Definition: uggrid.hh:611
@ COPY
New level consists of the refined elements and the unrefined ones, too.
Definition: uggrid.hh:615
@ LOCAL
New level consists only of the refined elements and the closure.
Definition: uggrid.hh:613
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: uggrid.hh:303
void setPosition(const typename Traits::template Codim< dim >::Entity &e, const FieldVector< double, dim > &pos)
Sets a vertex to a new position.
int size(int level, int codim) const
Number of grid entities per level and codim.
UG::DOUBLE ctype
The type used to store coordinates.
Definition: uggrid.hh:261
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: uggrid.hh:322
bool adapt()
Triggers the grid refinement process.
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: uggrid.hh:328
friend class UGGridLeafGridView
Definition: uggrid.hh:227
void loadState(const std::string &filename)
Read entire grid hierarchy from disk.
bool loadBalance(DataHandle &dataHandle)
Distributes the grid and some data over the available nodes in a distributed machine.
Definition: uggrid.hh:439
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Create an Entity from an EntitySeed.
Definition: uggrid.hh:280
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: uggrid.hh:336
unsigned int Rank
The type used for process ranks.
Definition: uggrid.hh:264
void saveState(const std::string &filename) const
Save entire grid hierarchy to disk.
bool preAdapt()
returns true, if some elements might be coarsend during grid adaption, here always returns true
void globalRefine(int n)
Does uniform refinement.
ClosureType
Decide whether to add a green closure to locally refined grid sections or not.
Definition: uggrid.hh:619
@ GREEN
Standard red/green refinement.
Definition: uggrid.hh:621
@ NONE
No closure, results in nonconforming meshes.
Definition: uggrid.hh:623
int size(int codim) const
number of leaf entities per codim in this process
Definition: uggrid.hh:291
A set of traits classes to store static information about grid implementation.
Different resources needed by all grid implementations.