11 #include <dune/common/classname.hh>
12 #include <dune/common/collectivecommunication.hh>
13 #include <dune/common/exceptions.hh>
14 #include <dune/common/mpihelper.hh>
15 #include <dune/common/static_assert.hh>
24 #include <dune/common/mpicollectivecommunication.hh>
46 #include "uggrid/ugincludes.hh"
51 #include "uggrid/ugwrapper.hh"
57 #include "uggrid/ug_undefs_lgm_seq.hh"
59 #include "uggrid/ug_undefs.hh"
77 #include "uggrid/ugincludes.hh"
82 #include "uggrid/ugwrapper.hh"
87 #include "uggrid/ug_undefs_lgm_seq.hh"
89 #include "uggrid/ug_undefs.hh"
96 #include "uggrid/uggridgeometry.hh"
97 #include "uggrid/uggridlocalgeometry.hh"
98 #include "uggrid/uggridentity.hh"
99 #include "uggrid/uggridentitypointer.hh"
100 #include "uggrid/uggridentityseed.hh"
101 #include "uggrid/uggridintersections.hh"
102 #include "uggrid/uggridintersectioniterators.hh"
103 #include "uggrid/uggridleveliterator.hh"
104 #include "uggrid/uggridleafiterator.hh"
105 #include "uggrid/uggridhieriterator.hh"
106 #include "uggrid/uggridindexsets.hh"
115 template <
class DataHandle,
int Gr
idDim,
int codim>
116 class UGMessageBufferBase {
118 typedef UGMessageBufferBase<DataHandle, GridDim, codim> ThisType;
119 typedef UGGrid<GridDim> GridType;
120 typedef typename DataHandle::DataType DataType;
126 UGMessageBufferBase(
void *ugData)
128 ugData_ =
static_cast<char*
>(ugData);
132 void write(
const DataType &t)
133 { this->writeRaw_<DataType>(t); }
135 void read(DataType &t)
136 { this->readRaw_<DataType>(t); }
141 template <
class ValueType>
142 void writeRaw_(
const ValueType &v)
144 *
reinterpret_cast<ValueType*
>(ugData_) = v;
145 ugData_ +=
sizeof(ValueType);
148 template <
class ValueType>
149 void readRaw_(ValueType &v)
151 v = *
reinterpret_cast<ValueType*
>(ugData_);
152 ugData_ +=
sizeof(ValueType);
157 static int ugGather_(
typename UG_NS<dim>::DDD_OBJ obj,
void* data)
162 UGMakeableEntity<0, dim, UGGrid<dim> > e(
reinterpret_cast<typename
UG_NS<dim>::Element*
>(obj),
nullptr);
164 if ((level == -1 && UG_NS<dim>::isLeaf(
reinterpret_cast<typename
UG_NS<dim>::Element*
>(obj))) || e.level() == level)
166 ThisType msgBuf(static_cast<DataType*>(data));
167 if (!duneDataHandle_->fixedsize(dim, codim))
168 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
169 duneDataHandle_->gather(msgBuf, e);
172 else if (codim == dim) {
175 UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > e(
reinterpret_cast<typename UG_NS<dim>::Node*
>(obj),
nullptr);
177 if ((level == -1 && UG_NS<dim>::isLeaf(
reinterpret_cast<typename UG_NS<dim>::Node*
>(obj))) || e.level() == level)
179 ThisType msgBuf(static_cast<DataType*>(data));
180 if (!duneDataHandle_->fixedsize(dim, codim))
181 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
182 duneDataHandle_->gather(msgBuf, e);
185 else if (codim == dim - 1) {
188 UGMakeableEntity<dim-1, dim,
Dune::UGGrid<dim> > e(
reinterpret_cast<typename UG_NS<dim>::Edge*
>(obj),
nullptr);
190 if ((level == -1 && UG_NS<dim>::isLeaf(
reinterpret_cast<typename UG_NS<dim>::Edge*
>(obj))) || e.level() == level)
192 ThisType msgBuf(static_cast<DataType*>(data));
193 if (!duneDataHandle_->fixedsize(dim, codim))
194 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
195 duneDataHandle_->gather(msgBuf, e);
199 DUNE_THROW(GridError,
200 "Only node and element wise "
201 "communication is currently "
202 "supported by UGGrid");
210 static int ugScatter_(
typename UG_NS<dim>::DDD_OBJ obj,
void* data)
213 typedef UGMakeableEntity<0, dim, UGGrid<dim> >
Entity;
218 if ((level == -1 && UG_NS<dim>::isLeaf(
reinterpret_cast<typename
UG_NS<dim>::Element*
>(obj))) || e.level() == level)
220 ThisType msgBuf(static_cast<DataType*>(data));
222 if (!duneDataHandle_->fixedsize(dim, codim))
225 n = duneDataHandle_->template size<Entity>(e);
227 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
230 else if (codim == dim) {
231 typedef UGMakeableEntity<dim, dim, Dune::UGGrid<dim> >
Entity;
234 Entity e(
reinterpret_cast<typename UG_NS<dim>::Node*
>(obj),
nullptr);
236 if ((level == -1 && UG_NS<dim>::isLeaf(
reinterpret_cast<typename UG_NS<dim>::Node*
>(obj))) || e.level() == level)
238 ThisType msgBuf(static_cast<DataType*>(data));
240 if (!duneDataHandle_->fixedsize(dim, codim))
243 n = duneDataHandle_->template size<Entity>(e);
245 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
248 else if (codim == dim - 1) {
252 Entity e(
reinterpret_cast<typename UG_NS<dim>::Edge*
>(obj),
nullptr);
254 if ((level == -1 && UG_NS<dim>::isLeaf(
reinterpret_cast<typename UG_NS<dim>::Edge*
>(obj))) || e.level() == level)
256 ThisType msgBuf(static_cast<DataType*>(data));
258 if (!duneDataHandle_->fixedsize(dim, codim))
261 n = duneDataHandle_->template size<Entity>(e);
263 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
267 DUNE_THROW(GridError,
268 "Only node and element wise "
269 "communication is currently "
270 "supported by UGGrid");
275 static DataHandle *duneDataHandle_;
282 template <
class DataHandle,
int Gr
idDim,
int codim>
283 class UGMessageBuffer
284 :
public UGMessageBufferBase<DataHandle, GridDim, codim>
286 typedef typename DataHandle::DataType DataType;
287 typedef UGMessageBufferBase<DataHandle, GridDim, codim> Base;
288 enum { dim = GridDim };
293 UGMessageBuffer(
void *ugData)
298 template <
class Gr
idView>
299 static unsigned ugBufferSize_(
const GridView &gv)
301 if (Base::duneDataHandle_->fixedsize(dim, codim)) {
302 return sizeof(DataType)
303 * Base::duneDataHandle_->
size(*gv.template begin<codim,InteriorBorder_Partition>());
306 typedef typename GridView::template Codim<codim>::Entity
Entity;
313 ::template Codim<codim>
314 ::template Partition<Dune::All_Partition>
316 Iterator it = gv.template begin<codim, Dune::All_Partition>();
317 const Iterator endIt = gv.template end<codim, Dune::All_Partition>();
318 for (; it != endIt; ++it) {
320 (
int) Base::duneDataHandle_->
size(*it));
324 maxSize = MPIHelper::getCollectiveCommunication().max(maxSize);
331 return sizeof(unsigned) +
sizeof(DataType)*maxSize;
335 template <
class DataHandle,
int Gr
idDim>
336 class UGEdgeMessageBuffer
337 :
public UGMessageBufferBase<DataHandle, GridDim, GridDim-1>
339 enum {codim = GridDim-1,
341 typedef typename DataHandle::DataType DataType;
342 typedef UGMessageBufferBase<DataHandle, GridDim, codim> Base;
346 UGEdgeMessageBuffer(
void *ugData)
351 template <
class Gr
idView>
352 static unsigned ugBufferSize_(
const GridView &gv)
354 if (Base::duneDataHandle_->fixedsize(dim, codim)) {
355 typedef typename GridView::template Codim<0>::Entity
Element;
356 const Element& element = *gv.template begin<0, InteriorBorder_Partition>();
357 return sizeof(DataType)
358 * Base::duneDataHandle_->
size(*element.template subEntity<codim>(0));
361 typedef typename GridView::template Codim<codim>::Entity
Entity;
369 ::template Partition<Dune::All_Partition>
371 Iterator it = gv.template begin<0, Dune::All_Partition>();
372 const Iterator endIt = gv.template end<0, Dune::All_Partition>();
373 for (; it != endIt; ++it) {
374 int numberOfEdges = it->template count<codim>();
375 for (
int k = 0; k < numberOfEdges; k++)
377 typedef typename GridView::template Codim<0>::Entity
Element;
378 typedef typename Element::template Codim<codim>::EntityPointer EdgePointer;
379 const EdgePointer edgePointer(it->template subEntity<codim>(k));
382 (
int) Base::duneDataHandle_->
size(*edgePointer));
387 maxSize = MPIHelper::getCollectiveCommunication().max(maxSize);
394 return sizeof(unsigned) +
sizeof(DataType)*maxSize;
398 template <
class DataHandle>
399 class UGMessageBuffer<DataHandle, 2, 1>
400 :
public UGEdgeMessageBuffer<DataHandle, 2>
403 template <
class DataHandle>
404 class UGMessageBuffer<DataHandle, 3, 2>
405 :
public UGEdgeMessageBuffer<DataHandle, 3>
410 template <
class DataHandle,
int Gr
idDim,
int codim>
411 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
413 template <
class DataHandle,
int Gr
idDim,
int codim>
414 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
421 class CollectiveCommunication<Dune::UGGrid<dim> > :
422 public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
424 typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
426 CollectiveCommunication()
427 : ParentType(MPIHelper::getCommunicator())
432 template<
int dim,
int dimworld>
440 UGGridLeafIntersection,
441 UGGridLevelIntersection,
442 UGGridLeafIntersectionIterator,
443 UGGridLevelIntersectionIterator,
444 UGGridHierarchicIterator,
446 UGGridLevelIndexSet< const UGGrid<dim> >,
447 UGGridLeafIndexSet< const UGGrid<dim> >,
448 UGGridIdSet< const UGGrid<dim>,
false >,
450 UGGridIdSet< const UGGrid<dim>,
true >,
452 CollectiveCommunication<Dune::UGGrid<dim> >,
506 friend class UGGridGeometry<0,dim,const
UGGrid<dim> >;
507 friend class UGGridGeometry<dim,dim,const
UGGrid<dim> >;
508 friend class UGGridGeometry<1,2,const
UGGrid<dim> >;
509 friend class UGGridGeometry<2,3,const
UGGrid<dim> >;
511 friend class UGGridEntity <0,dim,const
UGGrid<dim> >;
512 friend class UGGridEntity <dim,dim,const
UGGrid<dim> >;
513 friend class UGGridHierarchicIterator<const
UGGrid<dim> >;
514 friend class UGGridLeafIntersection<const
UGGrid<dim> >;
515 friend class UGGridLevelIntersection<const
UGGrid<dim> >;
516 friend class UGGridLeafIntersectionIterator<const
UGGrid<dim> >;
517 friend class UGGridLevelIntersectionIterator<const
UGGrid<dim> >;
519 friend class UGGridLevelIndexSet<const
UGGrid<dim> >;
520 friend class UGGridLeafIndexSet<const
UGGrid<dim> >;
521 friend class UGGridIdSet<const
UGGrid<dim>, false >;
522 friend class UGGridIdSet<const
UGGrid<dim>, true >;
526 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
528 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
530 template <
int codim_,
class Gr
idImp_>
534 dune_static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
577 template<
int codim, PartitionIteratorType PiType>
581 template<
int codim, PartitionIteratorType PiType>
593 return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
597 template<
int codim, PartitionIteratorType PiType>
603 template<
int codim, PartitionIteratorType PiType>
605 return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
610 typename Traits::template Codim<codim>::EntityPointer
618 int size (
int level,
int codim)
const;
642 return numBoundarySegments_;
661 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
662 return *levelIndexSets_[level];
668 return leafIndexSet_;
686 bool mark(
int refCount,
const typename Traits::template Codim<0>::Entity & e );
695 bool mark(
const typename Traits::template Codim<0>::Entity & e,
696 typename UG_NS<dim>::RefinementRule rule,
700 int getMark(
const typename Traits::template Codim<0>::Entity& e)
const;
720 return (codim==0) ? 1 : 0;
730 return (codim==0) ? 1 : 0;
746 template<
class DataHandle>
749 DUNE_THROW(NotImplemented,
"load balancing with data attached");
768 bool loadBalance(
int strategy,
int minlevel,
int depth,
int maxlevel,
int minelement);
780 template<
class DataHandle>
786 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
787 if (!dataHandle.contains(dim, curCodim))
791 communicateUG_<LevelGridView, DataHandle, 0>(this->
levelView(level), level, dataHandle, iftype, dir);
792 else if (curCodim == dim)
793 communicateUG_<LevelGridView, DataHandle, dim>(this->
levelView(level), level, dataHandle, iftype, dir);
794 else if (curCodim == dim - 1)
797 DUNE_THROW(NotImplemented,
798 className(*
this) <<
"::communicate(): Only "
799 "supported for codim=0 and "
800 "codim=dim(=" << dim <<
"), but "
801 "codim=" << curCodim <<
" was requested");
815 template<
class DataHandle>
823 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
824 if (!dataHandle.contains(dim, curCodim))
828 communicateUG_<LeafGridView, DataHandle, 0>(this->
leafView(), level, dataHandle, iftype, dir);
829 else if (curCodim == dim)
830 communicateUG_<LeafGridView, DataHandle, dim>(this->
leafView(), level, dataHandle, iftype, dir);
831 else if (curCodim == dim - 1)
833 communicateUG_<
LeafGridView, DataHandle, dim-1>(this->
leafView(), level, dataHandle, iftype, dir);
836 DUNE_THROW(NotImplemented,
837 className(*
this) <<
"::communicate(): Only "
838 "supported for codim=0 and "
839 "codim=dim(=" << dim <<
"), but "
840 "codim=" << curCodim <<
" was requested");
853 template <
class Gr
idView,
class DataHandle,
int codim>
854 void communicateUG_(
const GridView& gv,
int level,
855 DataHandle &dataHandle,
859 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
862 ugIfDir = UG_NS<dim>::IF_FORWARD();
864 ugIfDir = UG_NS<dim>::IF_BACKWARD();
866 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
867 UGMsgBuf::duneDataHandle_ = &dataHandle;
869 UGMsgBuf::level = level;
871 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
872 findDDDInterfaces_(ugIfs, iftype, codim);
874 unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
877 for (
unsigned i=0; i < ugIfs.size(); ++i)
878 UG_NS<dim>::DDD_IFOneway(ugIfs[i],
881 &UGMsgBuf::ugGather_,
882 &UGMsgBuf::ugScatter_);
885 void findDDDInterfaces_(std::vector<
typename UG_NS<dim>::DDD_IF > &dddIfaces,
902 dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
909 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
912 DUNE_THROW(GridError,
913 "Element communication not supported for "
914 "interfaces of type "
922 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
925 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
926 dddIfaces.push_back(UG_NS<dim>::NodeIF());
929 dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
932 DUNE_THROW(GridError,
933 "Node communication not supported for "
934 "interfaces of type "
942 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
945 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
950 DUNE_THROW(GridError,
951 "Edge communication not supported for "
952 "interfaces of type "
957 DUNE_THROW(GridError,
958 "Communication for codim "
960 <<
" entities is not yet supported "
961 <<
" by the DUNE UGGrid interface!");
991 std::vector<
typename Traits::template Codim<0>::EntityPointer>& childElements,
992 std::vector<unsigned char>& childElementSides)
const;
1010 refinementType_ = type;
1015 closureType_ = type;
1031 void setPosition(
const typename Traits::template Codim<dim>::EntityPointer& e,
1032 const FieldVector<double, dim>& pos);
1044 void saveState(
const std::string& filename)
const;
1050 void loadState(
const std::string& filename);
1054 typename UG_NS<dim>::MultiGrid* multigrid_;
1057 CollectiveCommunication<UGGrid> ccobj_;
1064 void setIndices(
bool setLevelZero,
1065 std::vector<unsigned int>* nodePermutation);
1072 std::vector<UGGridLevelIndexSet<const UGGrid<dim> >*> levelIndexSets_;
1074 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
1076 UGGridIdSet<const UGGrid<dim>,
false > globalIdSet_;
1078 UGGridIdSet<const UGGrid<dim>,
true > localIdSet_;
1093 static int numOfUGGrids;
1100 bool someElementHasBeenMarkedForRefinement_;
1107 bool someElementHasBeenMarkedForCoarsening_;
1113 static unsigned int heapSize_;
1116 std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
1123 unsigned int numBoundarySegments_;
1127 namespace Capabilities
1147 static const bool v =
true;
1156 static const bool v =
true;
1166 static const bool v =
true;
1168 static const bool v =
false;
1178 static const bool v =
true;
1187 static const bool v =
false;
1195 #endif // DUNE_UGGRID_HH