- Home
- About DUNE
- Download
- Documentation
- Community
- Development
00001 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- 00002 // vi: set et ts=4 sw=4 sts=4: 00003 00004 #ifndef DUNE_UGGRID_HH 00005 #define DUNE_UGGRID_HH 00006 00011 #include <dune/common/classname.hh> 00012 #include <dune/common/collectivecommunication.hh> 00013 #include <dune/common/exceptions.hh> 00014 #include <dune/common/mpihelper.hh> 00015 #include <dune/common/static_assert.hh> 00016 00017 #include <dune/grid/common/boundarysegment.hh> 00018 #include <dune/grid/common/capabilities.hh> 00019 #include <dune/grid/common/grid.hh> 00020 00021 #if HAVE_UG 00022 00023 #ifdef ModelP 00024 #include <dune/common/mpicollectivecommunication.hh> 00025 #endif 00026 00027 /* The following lines including the necessary UG headers are somewhat 00028 tricky. Here's what's happening: 00029 UG can support two- and three-dimensional grids. You choose be setting 00030 either _2 oder _3 while compiling. This changes all sorts of stuff, in 00031 particular data structures in the headers. 00032 UG was never supposed to provide 2d and 3d grids at the same time. 00033 However, when compiling it as c++, the dimension-dependent parts are 00034 wrapped up cleanly in the namespaces UG::D2 and UG::D3, respectively. That 00035 way it is possible to link together the UG lib for 2d and the one for 3d. 00036 But we also need the headers twice! Once with _2 set and once with _3! 00037 So here we go:*/ 00038 00039 /* The following define tells the UG headers that we want access to a few 00040 special fields, for example the extra index fields in the element data structures. */ 00041 #define FOR_DUNE 00042 00043 // Set UG's space-dimension flag to 2d 00044 #define _2 00045 // And include all necessary UG headers 00046 #include "uggrid/ugincludes.hh" 00047 00048 // Wrap a few large UG macros by functions before they get undef'ed away. 00049 // Here: The 2d-version of the macros 00050 #define UG_DIM 2 00051 #include "uggrid/ugwrapper.hh" 00052 #undef UG_DIM 00053 00054 // UG defines a whole load of preprocessor macros. ug_undefs.hh undefines 00055 // them all, so we don't get name clashes. 00056 #ifdef UG_LGMDOMAIN 00057 #include "uggrid/ug_undefs_lgm_seq.hh" 00058 #else 00059 #include "uggrid/ug_undefs.hh" 00060 #endif 00061 #undef _2 00062 00063 /* Now we're done with 2d, and we can do the whole thing over again for 3d */ 00064 00065 /* All macros set by UG have been unset. This includes the macros that ensure 00066 single inclusion of headers. We can thus include them again. However, we 00067 only want to include those headers again that contain dimension-dependent stuff. 00068 Therefore, we set a few single-inclusion defines manually before including 00069 ugincludes.hh again. 00070 */ 00071 #define UGTYPES_H 00072 #define __HEAPS__ 00073 #define __UGENV__ 00074 #define __DEVICESH__ 00075 00076 #define _3 00077 #include "uggrid/ugincludes.hh" 00078 00079 // Wrap a few large UG macros by functions before they get undef'ed away. 00080 // This time it's the 3d-versions. 00081 #define UG_DIM 3 00082 #include "uggrid/ugwrapper.hh" 00083 #undef UG_DIM 00084 00085 // undef all macros defined by UG 00086 #ifdef UG_LGMDOMAIN 00087 #include "uggrid/ug_undefs_lgm_seq.hh" 00088 #else 00089 #include "uggrid/ug_undefs.hh" 00090 #endif 00091 00092 #undef _3 00093 #undef FOR_DUNE 00094 00095 // The components of the UGGrid interface 00096 #include "uggrid/uggridgeometry.hh" 00097 #include "uggrid/uggridentity.hh" 00098 #include "uggrid/uggridentitypointer.hh" 00099 #include "uggrid/uggridintersections.hh" 00100 #include "uggrid/uggridintersectioniterators.hh" 00101 #include "uggrid/uggridleveliterator.hh" 00102 #include "uggrid/uggridleafiterator.hh" 00103 #include "uggrid/uggridhieriterator.hh" 00104 #include "uggrid/uggridindexsets.hh" 00105 00106 // Not needed here, but included for user convenience 00107 #include "uggrid/uggridfactory.hh" 00108 00109 #ifdef ModelP 00110 namespace Dune { 00111 00112 // converts the UG speak message buffers to DUNE speak and vice-versa 00113 template <class DataHandle, int GridDim, int codim> 00114 class UGMessageBuffer { 00115 protected: 00116 typedef UGMessageBuffer<DataHandle, GridDim, codim> ThisType; 00117 typedef UGGrid<GridDim> GridType; 00118 typedef typename DataHandle::DataType DataType; 00119 00120 enum { 00121 dim = GridDim 00122 }; 00123 00124 UGMessageBuffer(void *ugData) 00125 { 00126 ugData_ = static_cast<char*>(ugData); 00127 }; 00128 00129 public: 00130 void write(const DataType &t) 00131 { this->writeRaw_<DataType>(t); } 00132 00133 void read(DataType &t) 00134 { this->readRaw_<DataType>(t); } 00135 00136 protected: 00137 friend class Dune::UGGrid<dim>; 00138 00139 template <class ValueType> 00140 void writeRaw_(const ValueType &v) 00141 { 00142 *reinterpret_cast<ValueType*>(ugData_) = v; 00143 ugData_ += sizeof(ValueType); 00144 } 00145 00146 template <class ValueType> 00147 void readRaw_(ValueType &v) 00148 { 00149 v = *reinterpret_cast<ValueType*>(ugData_); 00150 ugData_ += sizeof(ValueType); 00151 } 00152 00153 // returns number of bytes required for the UG message buffer 00154 template <class GridView> 00155 static unsigned ugBufferSize_(const GridView &gv) 00156 { 00157 if (duneDataHandle_->fixedsize(dim, codim)) { 00158 return sizeof(DataType) 00159 * duneDataHandle_->size(*gv.template begin<codim,InteriorBorder_Partition>()); 00160 } 00161 00162 typedef typename GridView::template Codim<codim>::Entity Entity; 00163 00164 // iterate over all entities, find the maximum size for 00165 // the current rank 00166 int maxSize = 0; 00167 typedef typename 00168 GridView 00169 ::template Codim<codim> 00170 ::template Partition<Dune::All_Partition> 00171 ::Iterator Iterator; 00172 Iterator it = gv.template begin<codim, Dune::All_Partition>(); 00173 const Iterator endIt = gv.template end<codim, Dune::All_Partition>(); 00174 for (; it != endIt; ++it) { 00175 maxSize = std::max((int) maxSize, 00176 (int) duneDataHandle_->size(*it)); 00177 } 00178 00179 // find maximum size for all ranks 00180 maxSize = MPIHelper::getCollectiveCommunication().max(maxSize); 00181 if (!maxSize) 00182 return 0; 00183 00184 // add the size of an unsigned integer to the actual 00185 // buffer size. (we somewhere have to store the actual 00186 // number of objects for each entity.) 00187 return sizeof(unsigned) + sizeof(DataType)*maxSize; 00188 } 00189 00190 // called by DDD_IFOneway to serialize the data structure to 00191 // be send 00192 static int ugGather_(typename UG_NS<dim>::DDD_OBJ obj, void* data) 00193 { 00194 if (codim == 0) { 00195 UGMakeableEntity<0, dim, UGGrid<dim> > e((typename UG_NS<dim>::Element*) obj); 00196 // safety check to only communicate what is needed 00197 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Element*) obj)) || e.level() == level) 00198 { 00199 ThisType msgBuf(static_cast<DataType*>(data)); 00200 if (!duneDataHandle_->fixedsize(dim, codim)) 00201 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e)); 00202 duneDataHandle_->gather(msgBuf, e); 00203 } 00204 } 00205 else if (codim == dim) { 00206 UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > e((typename UG_NS<dim>::Node*)obj); 00207 // safety check to only communicate what is needed 00208 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Node*)obj)) || e.level() == level) 00209 { 00210 ThisType msgBuf(static_cast<DataType*>(data)); 00211 if (!duneDataHandle_->fixedsize(dim, codim)) 00212 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e)); 00213 duneDataHandle_->gather(msgBuf, e); 00214 } 00215 } 00216 else { 00217 DUNE_THROW(GridError, 00218 "Only node and element wise " 00219 "communication is currently " 00220 "supported by UGGrid"); 00221 } 00222 00223 return 0; 00224 } 00225 00226 // called by DDD_IFOneway to deserialize the data structure 00227 // that has been received 00228 static int ugScatter_(typename UG_NS<dim>::DDD_OBJ obj, void* data) 00229 { 00230 00231 if (codim == 0) { 00232 typedef UGMakeableEntity<0, dim, UGGrid<dim> > Entity; 00233 Entity e((typename UG_NS<dim>::Element*) obj); 00234 // safety check to only communicate what is needed 00235 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Element*) obj)) || e.level() == level) 00236 { 00237 ThisType msgBuf(static_cast<DataType*>(data)); 00238 int n; 00239 if (!duneDataHandle_->fixedsize(dim, codim)) 00240 msgBuf.readRaw_(n); 00241 else 00242 n = duneDataHandle_->template size<Entity>(e); 00243 if (n > 0) 00244 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n); 00245 } 00246 } 00247 else if (codim == dim) { 00248 typedef UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > Entity; 00249 Entity e((typename UG_NS<dim>::Node*)obj); 00250 // safety check to only communicate what is needed 00251 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Node*)obj)) || e.level() == level) 00252 { 00253 ThisType msgBuf(static_cast<DataType*>(data)); 00254 int n; 00255 if (!duneDataHandle_->fixedsize(dim, codim)) 00256 msgBuf.readRaw_(n); 00257 else 00258 n = duneDataHandle_->template size<Entity>(e); 00259 if (n > 0) 00260 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n); 00261 } 00262 } 00263 else { 00264 DUNE_THROW(GridError, 00265 "Only node and element wise " 00266 "communication is currently " 00267 "supported by UGGrid"); 00268 } 00269 00270 return 0; 00271 } 00272 static DataHandle *duneDataHandle_; 00273 00274 static int level; 00275 00276 char *ugData_; 00277 }; 00278 00279 } // end namespace Dune 00280 00281 template <class DataHandle, int GridDim, int codim> 00282 DataHandle *Dune::UGMessageBuffer<DataHandle,GridDim,codim>::duneDataHandle_ = 0; 00283 00284 template <class DataHandle, int GridDim, int codim> 00285 int Dune::UGMessageBuffer<DataHandle,GridDim,codim>::level = -1; 00286 #endif // ModelP 00287 00288 namespace Dune { 00289 00290 #ifdef ModelP 00291 template <int dim> 00292 class CollectiveCommunication<Dune::UGGrid<dim> > : 00293 public CollectiveCommunication< Dune::MPIHelper::MPICommunicator > 00294 { 00295 typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType; 00296 public: 00297 CollectiveCommunication() 00298 : ParentType(MPIHelper::getCommunicator()) 00299 {} 00300 }; 00301 #endif // ModelP 00302 00303 template<int dim, int dimworld> 00304 struct UGGridFamily 00305 { 00306 typedef GridTraits<dim,dimworld,Dune::UGGrid<dim>, 00307 UGGridGeometry, 00308 UGGridEntity, 00309 UGGridEntityPointer, 00310 UGGridLevelIterator, 00311 UGGridLeafIntersection, 00312 UGGridLevelIntersection, 00313 UGGridLeafIntersectionIterator, 00314 UGGridLevelIntersectionIterator, 00315 UGGridHierarchicIterator, 00316 UGGridLeafIterator, 00317 UGGridLevelIndexSet< const UGGrid<dim> >, 00318 UGGridLeafIndexSet< const UGGrid<dim> >, 00319 UGGridIdSet< const UGGrid<dim>, false >, 00320 unsigned int, 00321 UGGridIdSet< const UGGrid<dim>, true >, 00322 unsigned int, 00323 CollectiveCommunication<Dune::UGGrid<dim> > > 00324 Traits; 00325 }; 00326 00327 00328 //********************************************************************** 00329 // 00330 // --UGGrid 00331 // 00332 //********************************************************************** 00333 00370 template <int dim> 00371 class UGGrid : public GridDefaultImplementation <dim, dim, double, UGGridFamily<dim,dim> > 00372 { 00373 friend class UGGridGeometry<0,dim,const UGGrid<dim> >; 00374 friend class UGGridGeometry<dim,dim,const UGGrid<dim> >; 00375 friend class UGGridGeometry<1,2,const UGGrid<dim> >; 00376 friend class UGGridGeometry<2,3,const UGGrid<dim> >; 00377 00378 friend class UGGridEntity <0,dim,const UGGrid<dim> >; 00379 friend class UGGridEntity <dim,dim,const UGGrid<dim> >; 00380 friend class UGGridHierarchicIterator<const UGGrid<dim> >; 00381 friend class UGGridLeafIntersection<const UGGrid<dim> >; 00382 friend class UGGridLevelIntersection<const UGGrid<dim> >; 00383 friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >; 00384 friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >; 00385 00386 friend class UGGridLevelIndexSet<const UGGrid<dim> >; 00387 friend class UGGridLeafIndexSet<const UGGrid<dim> >; 00388 friend class UGGridIdSet<const UGGrid<dim>, false >; 00389 friend class UGGridIdSet<const UGGrid<dim>, true >; 00390 00391 friend class GridFactory<UGGrid<dim> >; 00392 00393 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 00394 friend class UGGridLeafIterator; 00395 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 00396 friend class UGGridLevelIterator; 00397 template <int codim_, class GridImp_> 00398 friend class UGGridEntityPointer; 00399 00401 dune_static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!"); 00402 00403 // The different instantiations are mutual friends so they can access 00404 // each others numOfUGGrids field 00405 friend class UGGrid<2>; 00406 friend class UGGrid<3>; 00407 00408 //********************************************************** 00409 // The Interface Methods 00410 //********************************************************** 00411 public: 00413 typedef UGGridFamily<dim,dim> GridFamily; 00414 00415 // the Traits 00416 typedef typename UGGridFamily<dim,dim>::Traits Traits; 00417 00419 typedef UG::DOUBLE ctype; 00420 00427 UGGrid(unsigned int heapSize) DUNE_DEPRECATED; 00428 00435 UGGrid(); 00436 00438 ~UGGrid(); 00439 00440 private: 00443 void init(); 00444 00445 public: 00448 int maxLevel() const; 00449 00451 template<int codim> 00452 typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const; 00453 00455 template<int codim> 00456 typename Traits::template Codim<codim>::LevelIterator lend (int level) const; 00457 00459 template<int codim, PartitionIteratorType PiType> 00460 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const; 00461 00463 template<int codim, PartitionIteratorType PiType> 00464 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const; 00465 00467 template<int codim> 00468 typename Traits::template Codim<codim>::LeafIterator leafbegin() const { 00469 return typename Traits::template Codim<codim>::template Partition<All_Partition>::LeafIterator(*this); 00470 } 00471 00473 template<int codim> 00474 typename Traits::template Codim<codim>::LeafIterator leafend() const { 00475 return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >(); 00476 } 00477 00479 template<int codim, PartitionIteratorType PiType> 00480 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const { 00481 return typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator(*this); 00482 } 00483 00485 template<int codim, PartitionIteratorType PiType> 00486 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const { 00487 return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >(); 00488 } 00489 00492 int size (int level, int codim) const; 00493 00495 int size (int codim) const 00496 { 00497 return leafIndexSet().size(codim); 00498 } 00499 00501 int size (int level, GeometryType type) const 00502 { 00503 return this->levelIndexSet(level).size(type); 00504 } 00505 00507 int size (GeometryType type) const 00508 { 00509 return this->leafIndexSet().size(type); 00510 } 00511 00513 size_t numBoundarySegments() const { 00514 // The number is stored as a member of UGGrid upon grid creation. 00515 // The corresponding data structure is not exported by UG. (It is in ug/dom/std/std_internal.h) 00516 return numBoundarySegments_; 00517 } 00518 00520 const typename Traits::GlobalIdSet& globalIdSet() const 00521 { 00522 return globalIdSet_; 00523 } 00524 00526 const typename Traits::LocalIdSet& localIdSet() const 00527 { 00528 return localIdSet_; 00529 } 00530 00532 const typename Traits::LevelIndexSet& levelIndexSet(int level) const 00533 { 00534 if (level<0 || level>maxLevel()) 00535 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!"); 00536 return *levelIndexSets_[level]; 00537 } 00538 00540 const typename Traits::LeafIndexSet& leafIndexSet() const 00541 { 00542 return leafIndexSet_; 00543 } 00544 00547 00560 bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e ); 00561 00569 bool mark(const typename Traits::template Codim<0>::Entity & e, 00570 typename UG_NS<dim>::RefinementRule rule, 00571 int side=0); 00572 00574 int getMark(const typename Traits::template Codim<0>::Entity& e) const; 00575 00578 bool preAdapt(); 00579 00581 bool adapt(); 00582 00584 void postAdapt(); 00588 unsigned int overlapSize(int codim) const { 00589 return 0; 00590 } 00591 00593 unsigned int ghostSize(int codim) const { 00594 return (codim==0) ? 1 : 0; 00595 } 00596 00598 unsigned int overlapSize(int level, int codim) const { 00599 return 0; 00600 } 00601 00603 unsigned int ghostSize(int level, int codim) const { 00604 return (codim==0) ? 1 : 0; 00605 } 00606 00611 bool loadBalance() { 00612 return loadBalance(0,0,2,32,1); 00613 } 00614 00620 template<class DataHandle> 00621 bool loadBalance (DataHandle& data) 00622 { 00623 DUNE_THROW(NotImplemented, "load balancing with data attached"); 00624 } 00625 00642 bool loadBalance(int strategy, int minlevel, int depth, int maxlevel, int minelement); 00643 00654 template<class DataHandle> 00655 void communicate (DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const 00656 { 00657 #ifdef ModelP 00658 typedef typename UGGrid::LevelGridView LevelGridView; 00659 00660 for (int curCodim = 0; curCodim <= dim; ++curCodim) { 00661 if (!dataHandle.contains(dim, curCodim)) 00662 continue; 00663 00664 if (curCodim == 0) 00665 communicateUG_<LevelGridView, DataHandle, 0>(this->levelView(level), level, dataHandle, iftype, dir); 00666 else if (curCodim == dim) 00667 communicateUG_<LevelGridView, DataHandle, dim>(this->levelView(level), level, dataHandle, iftype, dir); 00668 else 00669 DUNE_THROW(NotImplemented, 00670 className(*this) << "::communicate(): Only " 00671 "supported for codim=0 and " 00672 "codim=dim(=" << dim << "), but " 00673 "codim=" << curCodim << " was requested"); 00674 } 00675 #endif // ModelP 00676 } 00677 00687 template<class DataHandle> 00688 void communicate(DataHandle& dataHandle, 00689 InterfaceType iftype, 00690 CommunicationDirection dir) const 00691 { 00692 #ifdef ModelP 00693 typedef typename UGGrid::LeafGridView LeafGridView; 00694 00695 for (int curCodim = 0; curCodim <= dim; ++curCodim) { 00696 if (!dataHandle.contains(dim, curCodim)) 00697 continue; 00698 00699 int level = -1; 00700 if (curCodim == 0) 00701 communicateUG_<LeafGridView, DataHandle, 0>(this->leafView(), level, dataHandle, iftype, dir); 00702 else if (curCodim == dim) 00703 communicateUG_<LeafGridView, DataHandle, dim>(this->leafView(), level, dataHandle, iftype, dir); 00704 else 00705 DUNE_THROW(NotImplemented, 00706 className(*this) << "::communicate(): Only " 00707 "supported for codim=0 and " 00708 "codim=dim(=" << dim << "), but " 00709 "codim=" << curCodim << " was requested"); 00710 } 00711 #endif // ModelP 00712 } 00713 00715 const CollectiveCommunication<UGGrid>& comm () const 00716 { 00717 return ccobj_; 00718 } 00719 00720 protected: 00721 #ifdef ModelP 00722 template <class GridView, class DataHandle, int codim> 00723 void communicateUG_(const GridView& gv, int level, 00724 DataHandle &dataHandle, 00725 InterfaceType iftype, 00726 CommunicationDirection dir) const 00727 { 00728 typename UG_NS<dim>::DDD_IF_DIR ugIfDir; 00729 // Translate the communication direction from Dune-Speak to UG-Speak 00730 if (dir==ForwardCommunication) 00731 ugIfDir = UG_NS<dim>::IF_FORWARD(); 00732 else 00733 ugIfDir = UG_NS<dim>::IF_BACKWARD(); 00734 00735 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf; 00736 UGMsgBuf::duneDataHandle_ = &dataHandle; 00737 00738 UGMsgBuf::level = level; 00739 00740 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs; 00741 findDDDInterfaces_(ugIfs, iftype, codim); 00742 00743 unsigned bufSize = UGMsgBuf::ugBufferSize_(gv); 00744 if (!bufSize) 00745 return; // we don't need to communicate if we don't have any data! 00746 for (unsigned i=0; i < ugIfs.size(); ++i) 00747 UG_NS<dim>::DDD_IFOneway(ugIfs[i], 00748 ugIfDir, 00749 bufSize, 00750 &UGMsgBuf::ugGather_, 00751 &UGMsgBuf::ugScatter_); 00752 } 00753 00754 void findDDDInterfaces_(std::vector<typename UG_NS<dim>::DDD_IF > &dddIfaces, 00755 InterfaceType iftype, 00756 int codim) const 00757 { 00758 dddIfaces.clear(); 00759 switch (codim) 00760 { 00761 case 0: 00762 switch (iftype) { 00763 case InteriorBorder_InteriorBorder_Interface: 00764 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF()); 00765 return; 00766 case InteriorBorder_All_Interface: 00767 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF()); 00768 return; 00769 case All_All_Interface: 00770 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF()); 00771 return; 00772 default: 00773 DUNE_THROW(GridError, 00774 "Element communication not supported for " 00775 "interfaces of type " 00776 << iftype); 00777 } 00778 00779 case dim: 00780 switch (iftype) 00781 { 00782 case InteriorBorder_InteriorBorder_Interface: 00783 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF()); 00784 return; 00785 case InteriorBorder_All_Interface: 00786 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF()); 00787 dddIfaces.push_back(UG_NS<dim>::NodeIF()); 00788 return; 00789 case All_All_Interface: 00790 dddIfaces.push_back(UG_NS<dim>::NodeAllIF()); 00791 return; 00792 default: 00793 DUNE_THROW(GridError, 00794 "Node communication not supported for " 00795 "interfaces of type " 00796 << iftype); 00797 } 00798 00799 default: 00800 DUNE_THROW(GridError, 00801 "Communication for codim " 00802 << codim 00803 << " entities is not yet supported " 00804 << " by the DUNE UGGrid interface!"); 00805 } 00806 }; 00807 #endif // ModelP 00808 00809 public: 00810 // ********************************************************** 00811 // End of Interface Methods 00812 // ********************************************************** 00813 00816 00820 void createLGMGrid(const std::string& name); 00821 00829 void getChildrenOfSubface(const typename Traits::template Codim<0>::EntityPointer & e, 00830 int elementSide, 00831 int maxl, 00832 std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements, 00833 std::vector<unsigned char>& childElementSides) const; 00834 00836 enum RefinementType { 00838 LOCAL, 00840 COPY}; 00841 00843 enum ClosureType { 00845 GREEN, 00847 NONE}; 00848 00850 void setRefinementType(RefinementType type) { 00851 refinementType_ = type; 00852 } 00853 00855 void setClosureType(ClosureType type) { 00856 closureType_ = type; 00857 } 00858 00865 static void setDefaultHeapSize(unsigned size) { 00866 heapSize_ = size; 00867 } 00868 00872 void setPosition(const typename Traits::template Codim<dim>::EntityPointer& e, 00873 const FieldVector<double, dim>& pos); 00874 00879 void globalRefine(int n); 00880 00885 void saveState(const std::string& filename) const; 00886 00891 void loadState(const std::string& filename); 00892 00893 private: 00895 typename UG_NS<dim>::MultiGrid* multigrid_; 00896 00898 CollectiveCommunication<UGGrid> ccobj_; 00899 00905 void setIndices(bool setLevelZero, 00906 std::vector<unsigned int>* nodePermutation); 00907 00908 // Each UGGrid object has a unique name to identify it in the 00909 // UG environment structure 00910 std::string name_; 00911 00912 // Our set of level indices 00913 std::vector<UGGridLevelIndexSet<const UGGrid<dim> >*> levelIndexSets_; 00914 00915 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_; 00916 00917 UGGridIdSet<const UGGrid<dim>, false > globalIdSet_; 00918 00919 UGGridIdSet<const UGGrid<dim>, true > localIdSet_; 00920 00922 RefinementType refinementType_; 00923 00925 ClosureType closureType_; 00926 00934 static int numOfUGGrids; 00935 00941 bool someElementHasBeenMarkedForRefinement_; 00942 00948 bool someElementHasBeenMarkedForCoarsening_; 00949 00954 static unsigned int heapSize_; 00955 00957 std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_; 00958 00964 unsigned int numBoundarySegments_; 00965 00966 }; // end Class UGGrid 00967 00968 namespace Capabilities 00969 { 00985 template<int dim> 00986 struct hasEntity< UGGrid<dim>, 0> 00987 { 00988 static const bool v = true; 00989 }; 00990 00994 template<int dim> 00995 struct hasEntity< UGGrid<dim>, dim> 00996 { 00997 static const bool v = true; 00998 }; 00999 01003 template<int dim> 01004 struct isParallel< UGGrid<dim> > 01005 { 01006 #ifdef ModelP 01007 static const bool v = true; 01008 #else // !ModelP 01009 static const bool v = false; 01010 #endif // !ModelP 01011 }; 01012 01016 template<int dim> 01017 struct isLevelwiseConforming< UGGrid<dim> > 01018 { 01019 static const bool v = true; 01020 }; 01021 01025 template<int dim> 01026 struct isLeafwiseConforming< UGGrid<dim> > 01027 { 01028 static const bool v = false; 01029 }; 01030 01031 } 01032 01033 } // namespace Dune 01034 01035 #endif // HAVE_UG 01036 #endif // DUNE_UGGRID_HH
Generated on Fri Apr 29 2011 with Doxygen (ver 1.7.1) [doxygen-log,error-log].