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_>
251 friend class UGGridLeafGridView;
252 template <
class Gr
idImp_>
253 friend class UGGridLevelGridView;
258 friend class UGLBGatherScatter;
261 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
262 friend class UGGridLeafIterator;
263 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
264 friend class UGGridLevelIterator;
267 static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
298 ~
UGGrid() noexcept(false);
302 int maxLevel() const;
305 template <typename Seed>
306 typename Traits::template
Codim<Seed::codimension>::
Entity 307 entity(const Seed& seed)
const 309 const int codim = Seed::codimension;
315 int size (
int level,
int codim)
const;
320 return leafIndexSet().size(codim);
326 return this->levelIndexSet(level).size(type);
332 return this->leafIndexSet().size(type);
339 return numBoundarySegments_;
357 if (level<0 || level>maxLevel())
358 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
359 return *levelIndexSets_[level];
365 return leafIndexSet_;
383 bool mark(
int refCount,
const typename Traits::template
Codim<0>::Entity & e );
441 typename UG_NS<dim>::RefinementRule rule,
465 template<
class DataHandle>
470 if (dataHandle.contains(dim, 0))
471 UGLBGatherScatter::template gather<0>(this->
leafGridView(), dataHandle);
474 if (dataHandle.contains(dim,dim))
475 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
484 if (dataHandle.contains(dim, 0))
485 UGLBGatherScatter::template scatter<0>(this->
leafGridView(), dataHandle);
488 if (dataHandle.contains(dim,dim))
489 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
501 bool loadBalance(
int minlevel=0);
533 bool loadBalance(
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel);
544 template<
class DataHandle>
545 bool loadBalance (
const std::vector<Rank>& targetProcessors,
unsigned int fromLevel, DataHandle& dataHandle)
552 UGLBGatherScatter::template gather<dim>(this->
leafGridView(), dataHandle);
557 loadBalance(targetProcessors,fromLevel);
564 UGLBGatherScatter::template scatter<dim>(this->
leafGridView(), dataHandle);
578 template <
class Gr
idView,
class DataHandle,
int codim>
579 void communicateUG_(
const GridView& gv,
int level,
580 DataHandle &dataHandle,
584 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
587 ugIfDir = UG_NS<dim>::IF_FORWARD();
589 ugIfDir = UG_NS<dim>::IF_BACKWARD();
591 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
592 UGMsgBuf::duneDataHandle_ = &dataHandle;
594 UGMsgBuf::level = level;
596 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
597 findDDDInterfaces_(ugIfs, iftype, codim);
599 unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
602 for (
unsigned i=0; i < ugIfs.size(); ++i)
603 UG_NS<dim>::DDD_IFOneway(ugIfs[i],
606 &UGMsgBuf::ugGather_,
607 &UGMsgBuf::ugScatter_);
610 void findDDDInterfaces_(std::vector<
typename UG_NS<dim>::DDD_IF > &dddIfaces,
626 dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
633 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
637 "Element communication not supported for " 638 "interfaces of type " 642 else if (codim == dim)
647 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
650 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
651 dddIfaces.push_back(UG_NS<dim>::NodeIF());
654 dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
658 "Node communication not supported for " 659 "interfaces of type " 663 else if (codim == dim-1)
668 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
671 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
676 dddIfaces.push_back(UG_NS<dim>::EdgeSymmVHIF());
680 "Edge communication not supported for " 681 "interfaces of type " 691 dddIfaces.push_back(UG_NS<dim>::BorderVectorSymmIF());
695 "Face communication not supported for " 696 "interfaces of type " 703 "Communication for codim " 705 <<
" entities is not yet supported " 706 <<
" by the DUNE UGGrid interface!");
723 void getChildrenOfSubface(
const typename Traits::template
Codim<0>::Entity & e,
727 std::vector<unsigned char>& childElementSides)
const;
747 refinementType_ = type;
769 const FieldVector<double, dim>& pos);
775 void globalRefine(
int n);
781 void saveState(
const std::string& filename)
const;
787 void loadState(
const std::string& filename);
791 typename UG_NS<dim>::MultiGrid* multigrid_;
801 void setIndices(
bool setLevelZero,
802 std::vector<unsigned int>* nodePermutation);
809 std::vector<std::shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
811 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
815 UGGridIdSet<const UGGrid<dim> > idSet_;
830 static int numOfUGGrids;
837 bool someElementHasBeenMarkedForRefinement_;
844 bool someElementHasBeenMarkedForCoarsening_;
850 static unsigned int heapSize_;
853 std::vector<std::shared_ptr<BoundarySegment<dim> > > boundarySegments_;
860 unsigned int numBoundarySegments_;
864 namespace Capabilities
884 static const bool v =
true;
893 static const bool v =
true;
903 static const bool v =
true;
913 static const bool v =
true;
922 static const bool v =
true;
931 static const bool v =
false;
938 #endif // HAVE_UG || DOXYGEN 939 #endif // DUNE_UGGRID_HH Include standard header files.
Definition: agrid.hh:58
bool loadBalance(DataHandle &dataHandle)
Distributes the grid and some data over the available nodes in a distributed machine.
Definition: uggrid.hh:466
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: uggrid.hh:363
Base class for grid boundary segments of arbitrary geometry.
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:414
Grid view abstract base class.
Definition: common/gridview.hh:59
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: uggrid.hh:349
UGGridFamily< dim > GridFamily
type of the used GridFamily for this grid
Definition: uggrid.hh:279
Standard red/green refinement.
Definition: uggrid.hh:740
unsigned int Rank
The type used for process ranks.
Definition: uggrid.hh:288
Specialize with 'true' if implementation guarantees conforming level grids. (default=false) ...
Definition: common/capabilities.hh:103
specialize with 'true' for all codims that a grid provides an iterator for (default=false) ...
Definition: common/capabilities.hh:71
void setRefinementType(RefinementType type)
Sets the type of grid refinement.
Definition: uggrid.hh:746
void setClosureType(ClosureType type)
Sets the type of grid refinement closure.
Definition: uggrid.hh:751
Wrapper class for entities.
Definition: common/entity.hh:63
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:1063
Specialize with 'true' for all codims that a grid implements entities for. (default=false) ...
Definition: common/capabilities.hh:55
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:16
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: uggrid.hh:330
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:545
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:263
size_t numBoundarySegments() const
Return the number of boundary segments.
Definition: uggrid.hh:336
Front-end for the grid manager of the finite element toolbox UG.
Definition: uggrid.hh:230
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: uggrid.hh:324
send all and receive all entities
Definition: gridenums.hh:89
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false) ...
Definition: common/capabilities.hh:112
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: uggrid.hh:343
ClosureType
Decide whether to add a green closure to locally refined grid sections or not.
Definition: uggrid.hh:738
static void setDefaultHeapSize(unsigned size)
Sets the default heap size.
Definition: uggrid.hh:761
Definition: uggrid.hh:159
Definition: common/geometry.hh:24
send interior and border, receive all entities
Definition: gridenums.hh:86
RefinementType
The different forms of grid refinement that UG supports.
Definition: uggrid.hh:730
A set of traits classes to store static information about grid implementation.
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
UG::DOUBLE ctype
The type used to store coordinates.
Definition: uggrid.hh:285
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: uggrid.hh:355
const CollectiveCommunication< UGGrid > & comm() const
Definition: uggrid.hh:571
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
int size(int codim) const
number of leaf entities per codim in this process
Definition: uggrid.hh:318
Different resources needed by all grid implementations.
send/receive interior and border entities
Definition: gridenums.hh:85
UGGridFamily< dim >::Traits Traits
Definition: uggrid.hh:282
communicate as given in InterfaceType
Definition: gridenums.hh:169
Grid< dim, dimworld, ct, GridFamily >::LeafGridView leafGridView(const Grid< dim, dimworld, ct, GridFamily > &grid)
leaf grid view for the given grid
Definition: common/grid.hh:809
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:519
New level consists only of the refined elements and the closure.
Definition: uggrid.hh:732
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
The specialization of the generic GridFactory for UGGrid.