11 #include <dune/common/classname.hh>
12 #include <dune/common/parallel/collectivecommunication.hh>
13 #include <dune/common/exceptions.hh>
14 #include <dune/common/parallel/mpihelper.hh>
15 #include <dune/common/static_assert.hh>
21 #if HAVE_UG || DOXYGEN
24 #include <dune/common/parallel/mpicollectivecommunication.hh>
46 #include "uggrid/ugincludes.hh"
51 #include "uggrid/ugwrapper.hh"
56 #include "uggrid/ug_undefs.hh"
73 #include "uggrid/ugincludes.hh"
78 #include "uggrid/ugwrapper.hh"
82 #include "uggrid/ug_undefs.hh"
88 #include "uggrid/uggridgeometry.hh"
89 #include "uggrid/uggridlocalgeometry.hh"
90 #include "uggrid/uggridentity.hh"
91 #include "uggrid/uggridentitypointer.hh"
92 #include "uggrid/uggridentityseed.hh"
93 #include "uggrid/uggridintersections.hh"
94 #include "uggrid/uggridintersectioniterators.hh"
95 #include "uggrid/uggridleveliterator.hh"
96 #include "uggrid/uggridleafiterator.hh"
97 #include "uggrid/uggridhieriterator.hh"
98 #include "uggrid/uggridindexsets.hh"
100 #include "uggrid/ugmessagebuffer.hh"
101 #include "uggrid/uglbgatherscatter.hh"
108 template <
class DataHandle,
int Gr
idDim,
int codim>
109 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
111 template <
class DataHandle,
int Gr
idDim,
int codim>
112 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
119 class CollectiveCommunication<Dune::UGGrid<dim> > :
120 public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
122 typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
124 CollectiveCommunication()
125 : ParentType(MPIHelper::getCommunicator())
138 UGGridLeafIntersection,
139 UGGridLevelIntersection,
140 UGGridLeafIntersectionIterator,
141 UGGridLevelIntersectionIterator,
142 UGGridHierarchicIterator,
144 UGGridLevelIndexSet< const UGGrid<dim> >,
145 UGGridLeafIndexSet< const UGGrid<dim> >,
146 UGGridIdSet< const UGGrid<dim> >,
147 typename UG_NS<dim>::UG_ID_TYPE,
148 UGGridIdSet< const UGGrid<dim> >,
149 typename UG_NS<dim>::UG_ID_TYPE,
150 CollectiveCommunication<Dune::UGGrid<dim> >,
207 friend class UGGridGeometry<0,dim,const
UGGrid<dim> >;
208 friend class UGGridGeometry<dim,dim,const
UGGrid<dim> >;
209 friend class UGGridGeometry<1,2,const
UGGrid<dim> >;
210 friend class UGGridGeometry<2,3,const
UGGrid<dim> >;
212 friend class UGGridEntity <0,dim,const
UGGrid<dim> >;
213 friend class UGGridEntity <dim,dim,const
UGGrid<dim> >;
214 friend class UGGridHierarchicIterator<const
UGGrid<dim> >;
215 friend class UGGridLeafIntersection<const
UGGrid<dim> >;
216 friend class UGGridLevelIntersection<const
UGGrid<dim> >;
217 friend class UGGridLeafIntersectionIterator<const
UGGrid<dim> >;
218 friend class UGGridLevelIntersectionIterator<const
UGGrid<dim> >;
220 friend class UGGridLevelIndexSet<const
UGGrid<dim> >;
221 friend class UGGridLeafIndexSet<const
UGGrid<dim> >;
222 friend class UGGridIdSet<const
UGGrid<dim> >;
227 friend class UGLBGatherScatter;
230 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
232 template <
int codim_, PartitionIteratorType PiType_,
class Gr
idImp_>
234 template <
int codim_,
class Gr
idImp_>
238 dune_static_assert(dim==2 || dim==3,
"Use UGGrid only for 2d and 3d!");
281 template<
int codim, PartitionIteratorType PiType>
285 template<
int codim, PartitionIteratorType PiType>
297 return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
301 template<
int codim, PartitionIteratorType PiType>
307 template<
int codim, PartitionIteratorType PiType>
309 return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
313 template <
typename Seed>
314 typename Traits::template Codim<Seed::codimension>::EntityPointer
317 enum {codim = Seed::codimension};
323 int size (
int level,
int codim)
const;
347 return numBoundarySegments_;
366 DUNE_THROW(
GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
367 return *levelIndexSets_[level];
373 return leafIndexSet_;
391 bool mark(
int refCount,
const typename Traits::template Codim<0>::Entity & e );
449 bool mark(
const typename Traits::template Codim<0>::Entity & e,
450 typename UG_NS<dim>::RefinementRule rule,
454 int getMark(
const typename Traits::template Codim<0>::Entity& e)
const;
474 return (codim==0) ? 1 : 0;
484 return (codim==0) ? 1 : 0;
502 template<
class DataHandle>
506 DUNE_THROW(NotImplemented,
"load balancing with data attached");
513 UGLBGatherScatter::template gather<dim>(this->
leafView(), dataHandle);
525 UGLBGatherScatter::template scatter<dim>(this->
leafView(), dataHandle);
529 #endif // HAVE_UG_PATCH10
548 bool loadBalance(
int strategy,
int minlevel,
int depth,
int maxlevel,
int minelement);
560 template<
class DataHandle>
566 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
567 if (!dataHandle.contains(dim, curCodim))
571 communicateUG_<LevelGridView, DataHandle, 0>(this->
levelGridView(level), level, dataHandle, iftype, dir);
572 else if (curCodim == dim)
573 communicateUG_<LevelGridView, DataHandle, dim>(this->
levelGridView(level), level, dataHandle, iftype, dir);
574 else if (curCodim == dim - 1)
576 else if (curCodim == 1)
577 communicateUG_<LevelGridView, DataHandle, 1>(this->
levelGridView(level), level, dataHandle, iftype, dir);
579 DUNE_THROW(NotImplemented,
580 className(*
this) <<
"::communicate(): Not "
581 "supported for dim " << dim <<
" and codim " << curCodim);
595 template<
class DataHandle>
601 for (
int curCodim = 0; curCodim <= dim; ++curCodim) {
602 if (!dataHandle.contains(dim, curCodim))
606 communicateUG_<LeafGridView, DataHandle, 0>(this->
leafView(), level, dataHandle, iftype, dir);
607 else if (curCodim == dim)
608 communicateUG_<LeafGridView, DataHandle, dim>(this->
leafView(), level, dataHandle, iftype, dir);
609 else if (curCodim == dim - 1)
610 communicateUG_<
LeafGridView, DataHandle, dim-1>(this->
leafView(), level, dataHandle, iftype, dir);
611 else if (curCodim == 1)
612 communicateUG_<LeafGridView, DataHandle, 1>(this->
leafView(), level, dataHandle, iftype, dir);
614 DUNE_THROW(NotImplemented,
615 className(*
this) <<
"::communicate(): Not "
616 "supported for dim " << dim <<
" and codim " << curCodim);
629 template <
class Gr
idView,
class DataHandle,
int codim>
630 void communicateUG_(
const GridView& gv,
int level,
631 DataHandle &dataHandle,
635 typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
638 ugIfDir = UG_NS<dim>::IF_FORWARD();
640 ugIfDir = UG_NS<dim>::IF_BACKWARD();
642 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
643 UGMsgBuf::duneDataHandle_ = &dataHandle;
645 UGMsgBuf::level = level;
647 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
648 findDDDInterfaces_(ugIfs, iftype, codim);
650 unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
653 for (
unsigned i=0; i < ugIfs.size(); ++i)
654 UG_NS<dim>::DDD_IFOneway(ugIfs[i],
657 &UGMsgBuf::ugGather_,
658 &UGMsgBuf::ugScatter_);
661 void findDDDInterfaces_(std::vector<
typename UG_NS<dim>::DDD_IF > &dddIfaces,
677 dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
684 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
687 DUNE_THROW(GridError,
688 "Element communication not supported for "
689 "interfaces of type "
693 else if (codim == dim)
698 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
701 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
702 dddIfaces.push_back(UG_NS<dim>::NodeIF());
705 dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
708 DUNE_THROW(GridError,
709 "Node communication not supported for "
710 "interfaces of type "
714 else if (codim == dim-1)
719 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
722 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
727 DUNE_THROW(GridError,
728 "Edge communication not supported for "
729 "interfaces of type "
739 dddIfaces.push_back(UG_NS<dim>::BorderVectorSymmIF());
742 DUNE_THROW(GridError,
743 "Face communication not supported for "
744 "interfaces of type "
750 DUNE_THROW(GridError,
751 "Communication for codim "
753 <<
" entities is not yet supported "
754 <<
" by the DUNE UGGrid interface!");
774 std::vector<
typename Traits::template Codim<0>::EntityPointer>& childElements,
775 std::vector<unsigned char>& childElementSides)
const;
795 refinementType_ = type;
816 void setPosition(
const typename Traits::template Codim<dim>::EntityPointer& e,
817 const FieldVector<double, dim>& pos);
829 void saveState(
const std::string& filename)
const;
835 void loadState(
const std::string& filename);
839 typename UG_NS<dim>::MultiGrid* multigrid_;
842 CollectiveCommunication<UGGrid> ccobj_;
849 void setIndices(
bool setLevelZero,
850 std::vector<unsigned int>* nodePermutation);
857 std::vector<shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_;
859 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
863 UGGridIdSet<const UGGrid<dim> > idSet_;
878 static int numOfUGGrids;
885 bool someElementHasBeenMarkedForRefinement_;
892 bool someElementHasBeenMarkedForCoarsening_;
898 static unsigned int heapSize_;
901 std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
908 unsigned int numBoundarySegments_;
912 namespace Capabilities
932 static const bool v =
true;
941 static const bool v =
true;
951 static const bool v =
true;
953 static const bool v =
false;
963 static const bool v =
true;
972 static const bool v =
false;
979 #endif // HAVE_UG || DOXYGEN
980 #endif // DUNE_UGGRID_HH