15#include <dune/common/classname.hh>
16#include <dune/common/parallel/communication.hh>
17#include <dune/common/exceptions.hh>
18#include <dune/common/parallel/mpihelper.hh>
24#if HAVE_DUNE_UGGRID || DOXYGEN
27#include <dune/common/parallel/mpicommunication.hh>
55#include "uggrid/ugincludes.hh"
56#undef DUNE_UGINCLUDES_HH
61#include "uggrid/ugwrapper.hh"
66#include "uggrid/ug_undefs.hh"
86#include "uggrid/ugincludes.hh"
87#undef DUNE_UGINCLUDES_HH
92#include "uggrid/ugwrapper.hh"
96#include "uggrid/ug_undefs.hh"
102#include "uggrid/uggridgeometry.hh"
103#include "uggrid/uggridlocalgeometry.hh"
104#include "uggrid/uggridentity.hh"
105#include "uggrid/uggridentityseed.hh"
106#include "uggrid/uggridintersections.hh"
107#include "uggrid/uggridintersectioniterators.hh"
108#include "uggrid/uggridleveliterator.hh"
109#include "uggrid/uggridleafiterator.hh"
110#include "uggrid/uggridhieriterator.hh"
111#include "uggrid/uggridindexsets.hh"
112#include <dune/grid/uggrid/uggridviews.hh>
114#include "uggrid/ugmessagebuffer.hh"
115#include "uggrid/uglbgatherscatter.hh"
122template <
class DataHandle,
int Gr
idDim,
int codim>
125template <
class DataHandle,
int Gr
idDim,
int codim>
126DataHandle *Dune::UGMessageBuffer<DataHandle,GridDim,codim>::duneDataHandle_ =
nullptr;
128template <
class DataHandle,
int Gr
idDim,
int codim>
129int Dune::UGMessageBuffer<DataHandle,GridDim,codim>::level = -1;
147 UGGridLeafIntersection,
148 UGGridLevelIntersection,
149 UGGridLeafIntersectionIterator,
150 UGGridLevelIntersectionIterator,
151 UGGridHierarchicIterator,
153 UGGridLevelIndexSet< const UGGrid<dim> >,
154 UGGridLeafIndexSet< const UGGrid<dim> >,
155 UGGridIdSet< const UGGrid<dim> >,
156 typename UG_NS<dim>::UG_ID_TYPE,
157 UGGridIdSet< const UGGrid<dim> >,
158 typename UG_NS<dim>::UG_ID_TYPE,
160 UGGridLevelGridViewTraits,
161 UGGridLeafGridViewTraits,
210 friend class UGGridGeometry<0,dim,const
UGGrid<dim> >;
211 friend class UGGridGeometry<dim,dim,const
UGGrid<dim> >;
212 friend class UGGridGeometry<1,2,const
UGGrid<dim> >;
213 friend class UGGridGeometry<2,3,const
UGGrid<dim> >;
215 friend class UGGridEntity <0,dim,const
UGGrid<dim> >;
216 friend class UGGridEntity <1,dim,const
UGGrid<dim> >;
217 friend class UGGridEntity <2,dim,const
UGGrid<dim> >;
218 friend class UGGridEntity <dim,dim,const
UGGrid<dim> >;
219 friend class UGGridHierarchicIterator<const
UGGrid<dim> >;
220 friend class UGGridLeafIntersection<const
UGGrid<dim> >;
221 friend class UGGridLevelIntersection<const
UGGrid<dim> >;
222 friend class UGGridLeafIntersectionIterator<const
UGGrid<dim> >;
223 friend class UGGridLevelIntersectionIterator<const
UGGrid<dim> >;
225 friend class UGGridLevelIndexSet<const
UGGrid<dim> >;
226 friend class UGGridLeafIndexSet<const
UGGrid<dim> >;
227 friend class UGGridIdSet<const
UGGrid<dim> >;
228 template <
class Gr
idImp_>
230 template <
class Gr
idImp_>
236 friend class UGLBGatherScatter;
239 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
241 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
245 static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
280 template <typename Seed>
281 typename
Traits::template Codim<Seed::codimension>::
Entity
284 const int codim = Seed::codimension;
285 return typename Traits::template Codim<codim>::Entity(UGGridEntity<codim,dim,
const UGGrid<dim> >(seed.impl().target(),
this));
290 int size (
int level,
int codim)
const;
314 return numBoundarySegments_;
333 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
334 return *levelIndexSets_[level];
340 return leafIndexSet_;
358 bool mark(
int refCount,
const typename Traits::template Codim<0>::Entity & e );
415 bool mark(
const typename Traits::template Codim<0>::Entity & e,
416 typename UG_NS<dim>::RefinementRule rule,
420 int getMark(
const typename Traits::template Codim<0>::Entity& e)
const;
440 template<
class DataHandle>
445 if (dataHandle.contains(dim, 0))
446 UGLBGatherScatter::template gather<0>(this->
leafGridView(), dataHandle);
449 if (dataHandle.contains(dim,dim))
450 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
459 if (dataHandle.contains(dim, 0))
460 UGLBGatherScatter::template scatter<0>(this->
leafGridView(), dataHandle);
463 if (dataHandle.contains(dim,dim))
464 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
508 bool loadBalance(
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel);
519 template<
class DataHandle>
520 bool loadBalance (
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel, DataHandle& dataHandle)
524 if (dataHandle.contains(dim, 0))
525 UGLBGatherScatter::template gather<0>(this->
leafGridView(), dataHandle);
528 if (dataHandle.contains(dim,dim))
529 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
538 if (dataHandle.contains(dim, 0))
539 UGLBGatherScatter::template scatter<0>(this->
leafGridView(), dataHandle);
542 if (dataHandle.contains(dim,dim))
543 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
557 template <
int codim,
class Gr
idView,
class DataHandle>
558 void communicateUG_(
const GridView& gv,
int level,
559 DataHandle &dataHandle,
563 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
566 ugIfDir = UG_NS<dim>::IF_FORWARD();
568 ugIfDir = UG_NS<dim>::IF_BACKWARD();
570 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
571 UGMsgBuf::duneDataHandle_ = &dataHandle;
573 UGMsgBuf::level = level;
575 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs = findDDDInterfaces(iftype, codim);
577 unsigned bufSize = UGMsgBuf::ugBufferSize(gv);
580 UGMsgBuf::grid_ =
this;
581 for (
unsigned i=0; i < ugIfs.size(); ++i)
582 UG_NS<dim>::DDD_IFOneway(multigrid_->dddContext(),
586 &UGMsgBuf::ugGather_,
587 &UGMsgBuf::ugScatter_);
591 std::vector<typename UG_NS<dim>::DDD_IF> findDDDInterfaces(
InterfaceType iftype,
609 std::vector<
typename Traits::template Codim<0>::Entity>& childElements,
610 std::vector<unsigned char>& childElementSides)
const;
630 refinementType_ = type;
641 void setPosition(
const typename Traits::template Codim<dim>::Entity& e,
642 const FieldVector<double, dim>& pos);
664 typename UG_NS<dim>::MultiGrid* multigrid_;
674 void setIndices(
bool setLevelZero,
675 std::vector<unsigned int>* nodePermutation);
682 std::vector<std::shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
684 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
688 UGGridIdSet<const UGGrid<dim> > idSet_;
703 static int numOfUGGrids;
710 bool someElementHasBeenMarkedForRefinement_;
717 bool someElementHasBeenMarkedForCoarsening_;
720 std::vector<std::shared_ptr<BoundarySegment<dim> > > boundarySegments_;
727 unsigned int numBoundarySegments_;
731 namespace Capabilities
748 template<
int dim,
int codim>
751 static const bool v =
true;
759 template<
int dim,
int codim>
762 static const bool v =
false;
772 static const bool v =
true;
782 static const bool v =
true;
788 template<
int dim,
int codim>
791 static const bool v = (codim>=0 && codim<=dim);
800 static const bool v =
true;
809 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:170
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:86
@ ForwardCommunication
communicate as given in InterfaceType
Definition: gridenums.hh:171
Include standard header files.
Definition: agrid.hh:60
Communication< No_Comm > UGCommunication
Definition: uggrid.hh:137
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
Specialize with 'true' for all codims that a grid implements entities for. (default=false)
Definition: common/capabilities.hh:58
static const bool v
Definition: common/capabilities.hh:59
specialize with 'true' for all codims that a grid provides an iterator for (default=hasEntity<codim>:...
Definition: common/capabilities.hh:74
static const bool v
Definition: common/capabilities.hh:75
specialize with 'true' for all codims that a grid can communicate data on (default=false)
Definition: common/capabilities.hh:97
static const bool v
Definition: common/capabilities.hh:98
Specialize with 'true' if implementation guarantees conforming level grids. (default=false)
Definition: common/capabilities.hh:106
static const bool v
Definition: common/capabilities.hh:107
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false)
Definition: common/capabilities.hh:115
static const bool v
Definition: common/capabilities.hh:116
Wrapper class for entities.
Definition: common/entity.hh:66
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:20
Definition: common/grid.hh:862
Traits::LeafGridView leafGridView() const
View for the leaf grid for All_Partition.
Definition: common/grid.hh:882
bool loadBalance()
default implementation of load balance does nothing and returns false
Definition: common/grid.hh:948
Index Set Interface base class.
Definition: indexidset.hh:78
auto size(GeometryType type) const
Return total number of entities of given geometry type in entity set .
Definition: indexidset.hh:223
Id Set Interface.
Definition: indexidset.hh:452
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:995
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:314
Grid view abstract base class.
Definition: common/gridview.hh:66
Definition: uggrid.hh:142
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, UGCommunication, UGGridLevelGridViewTraits, UGGridLeafGridViewTraits, UGGridEntitySeed, UGGridLocalGeometry > Traits
Definition: uggrid.hh:164
Front-end for the grid manager of the finite element toolbox UG3.
Definition: uggrid.hh:207
void postAdapt()
Clean up refinement markers.
size_t numBoundarySegments() const
Return the number of boundary segments.
Definition: uggrid.hh:311
void setRefinementType(RefinementType type)
Sets the type of grid refinement.
Definition: uggrid.hh:629
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:240
bool loadBalance(int minlevel=0)
Distributes this grid over the available nodes in a distributed machine.
~UGGrid() noexcept(false)
Destructor.
UGGridFamily< dim >::Traits Traits
Definition: uggrid.hh:260
friend class UGGridLevelIterator
Definition: uggrid.hh:242
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:520
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:231
UGGridFamily< dim > GridFamily
type of the used GridFamily for this grid
Definition: uggrid.hh:257
bool loadBalance(const std::vector< Rank > &targetProcessors, unsigned int fromLevel)
Distribute this grid over a distributed machine.
const UGCommunication & comm() const
Definition: uggrid.hh:550
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.
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: uggrid.hh:330
void setClosureType(ClosureType type)
Sets the type of grid refinement closure.
Definition: uggrid.hh:634
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: uggrid.hh:324
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: uggrid.hh:299
RefinementType
The different forms of grid refinement that UG supports.
Definition: uggrid.hh:613
@ COPY
New level consists of the refined elements and the unrefined ones, too.
Definition: uggrid.hh:617
@ LOCAL
New level consists only of the refined elements and the closure.
Definition: uggrid.hh:615
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: uggrid.hh:305
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:263
bool adapt()
Triggers the grid refinement process.
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: uggrid.hh:318
friend class UGGridLeafGridView
Definition: uggrid.hh:229
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:441
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Create an Entity from an EntitySeed.
Definition: uggrid.hh:282
unsigned int Rank
The type used for process ranks.
Definition: uggrid.hh:266
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.
UGGrid(UGCommunication comm={})
Default constructor.
ClosureType
Decide whether to add a green closure to locally refined grid sections or not.
Definition: uggrid.hh:621
@ GREEN
Standard red/green refinement.
Definition: uggrid.hh:623
@ NONE
No closure, results in nonconforming meshes.
Definition: uggrid.hh:625
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: uggrid.hh:338
int size(int codim) const
number of leaf entities per codim in this process
Definition: uggrid.hh:293
A set of traits classes to store static information about grid implementation.
Different resources needed by all grid implementations.