uggrid.hh

Go to the documentation of this file.
00001 #ifndef DUNE_UGGRID_HH
00002 #define DUNE_UGGRID_HH
00003 
00008 #include <dune/grid/common/capabilities.hh>
00009 #include <dune/grid/common/grid.hh>
00010 #include <dune/grid/common/boundarysegment.hh>
00011 #include <dune/common/collectivecommunication.hh>
00012 #include <dune/common/deprecated.hh>
00013 #include <dune/common/static_assert.hh>
00014 
00015 /* The following lines including the necessary UG headers are somewhat
00016    tricky.  Here's what's happening:
00017    UG can support two- and three-dimensional grids.  You choose be setting
00018    either _2 oder _3 while compiling.  This changes all sorts of stuff, in
00019    particular data structures in the headers.
00020    UG was never supposed to provide 2d and 3d grids at the same time. 
00021    However, when compiling it as c++, the dimension-dependent parts are
00022    wrapped up cleanly in the namespaces UG::D2 and UG::D3, respectively.  That
00023    way it is possible to link together the UG lib for 2d and the one for 3d.
00024    But we also need the headers twice!  Once with _2 set and once with _3!
00025    So here we go:*/
00026 
00027 /* The following define tells the UG headers that we want access to a few
00028    special fields, for example the extra index fields in the element data structures. */
00029 #define FOR_DUNE
00030 
00031 // Set UG's space-dimension flag to 2d
00032 #define _2
00033 // And include all necessary UG headers
00034 #include "uggrid/ugincludes.hh"
00035 
00036 // Wrap a few large UG macros by functions before they get undef'ed away.
00037 // Here: The 2d-version of the macros
00038 #define UG_DIM 2
00039 #include "uggrid/ugwrapper.hh"
00040 #undef UG_DIM
00041 
00042 // UG defines a whole load of preprocessor macros.  ug_undefs.hh undefines
00043 // them all, so we don't get name clashes.
00044 #ifdef UG_LGMDOMAIN
00045 #include "uggrid/ug_undefs_lgm_seq.hh"
00046 #else
00047 #include "uggrid/ug_undefs.hh"
00048 #endif
00049 #undef _2
00050 
00051 /* Now we're done with 2d, and we can do the whole thing over again for 3d */
00052 
00053 /* All macros set by UG have been unset.  This includes the macros that ensure
00054    single inclusion of headers.  We can thus include them again.  However, we
00055    only want to include those headers again that contain dimension-dependent stuff.
00056    Therefore, we set a few single-inclusion defines manually before including
00057    ugincludes.hh again.
00058 */
00059 #define UGTYPES_H
00060 #define __HEAPS__
00061 #define __UGENV__
00062 #define __PARGM_H__
00063 #define __DEVICESH__
00064 #define __SM__
00065 
00066 #define _3
00067 #include "uggrid/ugincludes.hh"
00068 
00069 // Wrap a few large UG macros by functions before they get undef'ed away.
00070 // This time it's the 3d-versions.
00071 #define UG_DIM 3
00072 #include "uggrid/ugwrapper.hh"
00073 #undef UG_DIM
00074 
00075 // undef all macros defined by UG
00076 #ifdef UG_LGMDOMAIN
00077 #include "uggrid/ug_undefs_lgm_seq.hh"
00078 #else
00079 #include "uggrid/ug_undefs.hh"
00080 #endif
00081 
00082 #undef _3
00083 #undef FOR_DUNE
00084 
00085 // The components of the UGGrid interface
00086 #include "uggrid/uggridgeometry.hh"
00087 #include "uggrid/uggridentity.hh"
00088 #include "uggrid/uggridentitypointer.hh"
00089 #include "uggrid/ugintersectionit.hh"
00090 #include "uggrid/uggridleveliterator.hh"
00091 #include "uggrid/uggridleafiterator.hh"
00092 #include "uggrid/uggridhieriterator.hh"
00093 #include "uggrid/uggridindexsets.hh"
00094 
00095 // Not needed here, but included for user convenience
00096 #include "uggrid/uggridfactory.hh"
00097 
00098 
00099 #ifdef ModelP
00100 namespace Dune {
00101 
00102 // converts the UG speak message buffers to DUNE speak and vince-versa
00103     template <class DataHandle, int GridDim, int codim>
00104     class UGMessageBuffer {
00105     protected:
00106         typedef UGMessageBuffer<DataHandle, GridDim, codim>  ThisType;
00107         typedef typename DataHandle::DataType                DataType;
00108 
00109         enum {
00110             dim = GridDim
00111         };
00112 
00113         UGMessageBuffer(DataType *ugData)
00114         {
00115             ugData_ = ugData;
00116         };
00117         
00118     public:
00119         
00120         void write(const DataType &t) 
00121         {
00122             *ugData_ = t;
00123         }
00124         
00125         void read(DataType &t)
00126         {
00127             t = *ugData_;
00128             ++ ugData_;
00129         }
00130 
00131     protected:
00132         friend class Dune::UGGrid<dim>;
00133 
00134         // called by DDD_IFOneway to serialize the data structure to
00135         // be send
00136         static int ugGather(DDD_OBJ obj, void* data)
00137         {
00138             if (codim == 0) {
00139                 std::cout << "ugGather element index: " << UG_NS<dim>::levelIndex((typename UG_NS<dim>::Element*)obj) << "\n";
00140 
00141                 UGMakeableEntity<0, dim, UGGrid<dim> > e((typename UG_NS<dim>::Element*) obj);
00142                 ThisType msgBuf(static_cast<DataType*>(data));
00143                 duneDataHandle_->gather(msgBuf, e);
00144             }
00145             else if (codim == dim) {
00146                 std::cout << "ugGather node index: " << UG_NS<dim>::levelIndex((typename UG_NS<dim>::Node*)obj) << "\n";
00147 
00148                 UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > e((typename UG_NS<dim>::Node*)obj);
00149                 ThisType msgBuf(static_cast<DataType*>(data));
00150                 duneDataHandle_->gather(msgBuf, e);
00151             }
00152             else {
00153                 DUNE_THROW(GridError, 
00154                            "Only node and element wise "
00155                            "communication is currently "
00156                            "supported by UGGrid");
00157             }
00158 
00159             return 0;
00160         }
00161 
00162         // called by DDD_IFOneway to deserialize the data structure
00163         // which has been received
00164         static int ugScatter(DDD_OBJ obj, void* data)
00165         {
00166             
00167             if (codim == 0) {
00168                 std::cout << "ugScatter element index: " << UG_NS<dim>::levelIndex((typename UG_NS<dim>::Element*)obj) << "\n";
00169 
00170                 UGMakeableEntity<0, dim, UGGrid<dim> > e((typename UG_NS<dim>::Element*) obj);
00171                 ThisType msgBuf(static_cast<DataType*>(data));
00172                 duneDataHandle_->scatter(msgBuf, e, 1);
00173             }
00174             else if (codim == dim) {
00175                 std::cout << "ugScatter node index: " << UG_NS<dim>::levelIndex((typename UG_NS<dim>::Node*)obj) << "\n";
00176 
00177                 UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > e((typename UG_NS<dim>::Node*)obj);
00178                 ThisType msgBuf(static_cast<DataType*>(data));
00179                 duneDataHandle_->scatter(msgBuf, e, 1);
00180             }
00181             else {
00182                 DUNE_THROW(GridError, 
00183                            "Only node and element wise "
00184                            "communication is currently "
00185                            "supported by UGGrid");
00186             }
00187             
00188             return 0;
00189         }
00190         static DataHandle *duneDataHandle_;
00191 
00192         DataType          *ugData_;
00193     };
00194 
00195 }   // end namespace Dune
00196 
00197 template <class DataHandle, int GridDim, int codim>
00198 DataHandle *Dune::UGMessageBuffer<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
00199 #endif
00200 
00201 namespace Dune {
00202 
00203 template<int dim, int dimworld>
00204 struct UGGridFamily
00205 {
00206   typedef GridTraits<dim,dimworld,Dune::UGGrid<dim>,
00207                      UGGridGeometry,
00208                      UGGridEntity,
00209                      UGGridEntityPointer,
00210                      UGGridLevelIterator,
00211                      UGGridLeafIntersectionIterator, // leaf  intersection
00212                      UGGridLevelIntersectionIterator, // level intersection
00213                      UGGridLeafIntersectionIterator, // leaf  intersection iterator
00214                      UGGridLevelIntersectionIterator, // level intersection iterator
00215                      UGGridHierarchicIterator,
00216                      UGGridLeafIterator,
00217                      UGGridLevelIndexSet< const UGGrid<dim> >,
00218                      UGGridLevelIndexSetTypes< const UGGrid<dim> >,
00219                      UGGridLeafIndexSet< const UGGrid<dim> >,
00220                      UGGridLeafIndexSetTypes< const UGGrid<dim> >,
00221                      UGGridIdSet< const UGGrid<dim>, false >,
00222                      unsigned int,
00223                      UGGridIdSet< const UGGrid<dim>, true >,
00224                      unsigned int,
00225                      CollectiveCommunication<Dune::UGGrid<dim> > > 
00226   Traits;
00227 };
00228 
00229 
00230 //**********************************************************************
00231 //
00232 // --UGGrid
00233 //
00234 //**********************************************************************
00235 
00272 template <int dim>
00273 class UGGrid : public GridDefaultImplementation  <dim, dim, double, UGGridFamily<dim,dim> >
00274 {
00275     friend class UGGridGeometry<0,dim,const UGGrid<dim> >;
00276     friend class UGGridGeometry<dim,dim,const UGGrid<dim> >;
00277     friend class UGGridGeometry<1,2,const UGGrid<dim> >;
00278     friend class UGGridGeometry<2,3,const UGGrid<dim> >;
00279 
00280     friend class UGGridEntity <0,dim,const UGGrid<dim> >;
00281     friend class UGGridEntity <dim,dim,const UGGrid<dim> >;
00282     friend class UGGridHierarchicIterator<const UGGrid<dim> >;
00283     friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >;
00284     friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >;
00285 
00286     friend class UGGridLevelIndexSet<const UGGrid<dim> >;
00287     friend class UGGridLeafIndexSet<const UGGrid<dim> >;
00288     friend class UGGridIdSet<const UGGrid<dim>, false >;
00289     friend class UGGridIdSet<const UGGrid<dim>, true >;
00290 
00291     template <int codim_, PartitionIteratorType PiType_, class GridImp_>
00292     friend class UGGridLeafIterator;
00293     template <int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
00294     friend class Entity;
00295 
00297     dune_static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!");
00298 
00299     // The different instantiations are mutual friends so they can access
00300     // each others numOfUGGrids field
00301     friend class UGGrid<2>;
00302     friend class UGGrid<3>;
00303     //**********************************************************
00304       // The Interface Methods
00305       //**********************************************************
00306 public:  
00308       typedef UGGridFamily<dim,dim>  GridFamily;
00309 
00310     // the Traits 
00311     typedef typename UGGridFamily<dim,dim>::Traits Traits;
00312 
00314     typedef UG::DOUBLE ctype;
00315 
00320     UGGrid(unsigned int heapSize=500);
00321 
00323     ~UGGrid();
00324    
00327      int maxLevel() const;
00328      
00330     template<int codim>
00331     typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
00332 
00334     template<int codim>
00335     typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
00336 
00338     template<int codim, PartitionIteratorType PiType>
00339     typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
00340 
00342     template<int codim, PartitionIteratorType PiType>
00343     typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
00344 
00346     template<int codim>
00347     typename Traits::template Codim<codim>::LeafIterator leafbegin() const {
00348         return typename Traits::template Codim<codim>::template Partition<All_Partition>::LeafIterator(*this);
00349     }
00350 
00352     template<int codim>
00353     typename Traits::template Codim<codim>::LeafIterator leafend() const {
00354         return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
00355     }
00356 
00358     template<int codim, PartitionIteratorType PiType>
00359     typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const {
00360         return typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator(*this);
00361     }
00362 
00364     template<int codim, PartitionIteratorType PiType>
00365     typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const {
00366         return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
00367     }
00368 
00371     int size (int level, int codim) const;
00372 
00374   int size (int codim) const
00375   {
00376       return leafIndexSet().size(codim);
00377   }
00378 
00380   int size (int level, GeometryType type) const
00381   {
00382         return this->levelIndexSet(level).size(type);
00383   }
00384 
00386   int size (GeometryType type) const
00387   {
00388         return this->leafIndexSet().size(type);
00389   }
00390 
00392     const typename Traits::GlobalIdSet& globalIdSet() const
00393     {
00394         return globalIdSet_;
00395     }
00396     
00398     const typename Traits::LocalIdSet& localIdSet() const
00399     {
00400         return localIdSet_;
00401     }
00402     
00404     const typename Traits::LevelIndexSet& levelIndexSet(int level) const
00405     {
00406         if (level<0 || level>maxLevel())
00407             DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
00408         return *levelIndexSets_[level];
00409     }
00410     
00412     const typename Traits::LeafIndexSet& leafIndexSet() const
00413     {
00414         return leafIndexSet_;
00415     }
00416 
00419 
00432     bool mark(int refCount, const typename Traits::template Codim<0>::EntityPointer & e ) DUNE_DEPRECATED;
00433 
00446     bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e );
00447 
00455     bool mark(const typename Traits::template Codim<0>::EntityPointer & e, 
00456               typename UG_NS<dim>::RefinementRule rule,
00457               int side=0) DUNE_DEPRECATED ;
00458 
00466     bool mark(const typename Traits::template Codim<0>::Entity & e, 
00467               typename UG_NS<dim>::RefinementRule rule,
00468               int side=0);
00469 
00471     int getMark(const typename Traits::template Codim<0>::EntityPointer& e) const DUNE_DEPRECATED ;
00472 
00474     int getMark(const typename Traits::template Codim<0>::Entity& e) const;
00475 
00478     bool preAdapt();
00479     
00481     bool adapt();
00482 
00484     void postAdapt();
00488     std::string name () const { return "UGGrid"; };
00489 
00491     unsigned int overlapSize(int codim) const {
00492         return 0;
00493     }
00494 
00496     unsigned int ghostSize(int codim) const {
00497         return (codim==0) ? 1 : 0;
00498     }
00499 
00501     unsigned int overlapSize(int level, int codim) const {
00502         return 0;
00503     }
00504 
00506     unsigned int ghostSize(int level, int codim) const {
00507         return (codim==0) ? 1 : 0;
00508     }
00509 
00514     bool loadBalance() {
00515         return loadBalance(0,0,2,32,1);
00516     }
00517 
00523     template<class DataHandle>
00524     bool loadBalance (DataHandle& data)
00525     {
00526         DUNE_THROW(NotImplemented, "load balancing with data attached");
00527     }
00528 
00545     bool loadBalance(int strategy, int minlevel, int depth, int maxlevel, int minelement);
00546 
00547 #if 0
00548 
00564     template<class T, template<class> class P, int codim>
00565     void communicate (T& t, InterfaceType iftype, CommunicationDirection dir, int level);
00566 #endif
00567 
00578     template<class DataHandle>
00579     void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int level) const
00580     {
00581         DUNE_THROW(NotImplemented, "Level communication has not been implemented yet!");
00582     }
00583 
00593     template<class DataHandle>
00594     void communicate (DataHandle& dataHandle, 
00595                       InterfaceType iftype,
00596                       CommunicationDirection dir) const
00597     {
00598 #ifdef ModelP  
00599         // Translate the communication direction from Dune-Speak to UG-Speak
00600         DDD_IF_DIR UGIfDir = (dir==ForwardCommunication) ? IF_FORWARD : IF_BACKWARD;
00601         
00602         for (int curCodim = 0; curCodim <= dim; ++curCodim) {
00603             if (!dataHandle.contains(dim, curCodim))
00604                 continue;
00605             else if (!dataHandle.fixedsize(dim, curCodim))
00606                 DUNE_THROW(GridError, "Currently UG supports supports communication of fixed-size data types!");
00607 
00608             if (curCodim == 0 && 
00609                 iftype == InteriorBorder_All_Interface)
00610             {
00611                 typedef UGMessageBuffer<DataHandle,dim,0> UGMsgBuf;
00612                 UGMsgBuf::duneDataHandle_ = &dataHandle;
00613         
00614                 DDD_IFOneway(UG_NS<dim>::ElementVHIF(),
00615                              UGIfDir,
00616                              sizeof(typename DataHandle::DataType),
00617                              &UGMsgBuf::ugGather,
00618                              &UGMsgBuf::ugScatter);
00619             }
00620             else if (curCodim == dim && 
00621                      iftype == InteriorBorder_All_Interface)
00622             {
00623                 typedef UGMessageBuffer<DataHandle,dim,dim> UGMsgBuf;
00624                 UGMsgBuf::duneDataHandle_ = &dataHandle;
00625 
00626                 DDD_IFOneway(UG_NS<dim>::NodeIF(),
00627                              UGIfDir,
00628                              sizeof(typename DataHandle::DataType),
00629                              &UGMsgBuf::ugGather,
00630                              &UGMsgBuf::ugScatter);
00631             }
00632             else
00633             {
00634                 DUNE_THROW(GridError, 
00635                            "Communication for codim "
00636                            << curCodim
00637                            << " entities is not yet supported "
00638                            << "for communication interfaces of type " 
00639                            << iftype 
00640                            << " by the DUNE UGGrid interface!");
00641             }
00642         }
00643 #endif
00644     }
00645 
00646 
00648         const CollectiveCommunication<UGGrid>& comm () const
00649         {
00650           return ccobj;
00651         }
00652 
00653     // **********************************************************
00654     // End of Interface Methods
00655     // **********************************************************
00656     
00659     
00663     void createBegin() DUNE_DEPRECATED;
00664 
00666     void createEnd() DUNE_DEPRECATED;
00667 
00671     void createLGMGrid(const std::string& name);
00672 
00678     void insertBoundarySegment(const std::vector<unsigned int> vertices,
00679                                const BoundarySegment<dim>* boundarySegment) DUNE_DEPRECATED;
00680 
00682     void insertVertex(const FieldVector<double,dim>& pos) DUNE_DEPRECATED;
00683 
00688     void insertElement(GeometryType type,
00689                        const std::vector<unsigned int>& vertices) DUNE_DEPRECATED;
00690     
00698     void getChildrenOfSubface(typename Traits::template Codim<0>::EntityPointer & e,
00699                               int elementSide,
00700                               int maxl, 
00701                               std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements,
00702                               std::vector<unsigned char>& childElementSides) const;
00703     
00705     enum RefinementType {
00707         LOCAL, 
00709         COPY};
00710 
00712     enum ClosureType {
00714         GREEN,
00716         NONE};
00717 
00719     void setRefinementType(RefinementType type) {
00720         refinementType_ = type;
00721     }
00722 
00724     void setClosureType(ClosureType type) {
00725         closureType_ = type;
00726     }
00727 
00729     void collapse() {
00730         if (Collapse(multigrid_))
00731             DUNE_THROW(GridError, "UG" << dim << "d::Collapse() returned error code!");
00732 
00733         setIndices();
00734     }
00735 
00739     void setPosition(typename Traits::template Codim<dim>::EntityPointer& e,
00740                      const FieldVector<double, dim>& pos);
00741 
00743     FieldVector<ctype,dim> getBoundaryPosition(const typename Traits::LevelIntersectionIterator& iIt,
00744                                                const FieldVector<ctype,dim-1>& localPos) const;
00745 
00750     void globalRefine(int n);
00751 
00752     void saveState(const std::string& filename) const;
00753 
00754     void loadState(const std::string& filename);
00755 
00756 private:
00758     typename UG_NS<dim>::MultiGrid* multigrid_;
00759 
00761     std::vector<const BoundarySegment<dim>*> boundarySegments_;
00762 
00764     std::vector<array<unsigned int, dim*2-2> > boundarySegmentVertices_;
00765 
00766   CollectiveCommunication<UGGrid> ccobj;
00767 
00768     // Recomputes entity indices after the grid was changed
00769     void setIndices(std::vector<unsigned int>* nodePermutation=0);
00770 
00771     // Each UGGrid object has a unique name to identify it in the
00772     // UG environment structure
00773     std::string name_;
00774 
00775     // Our set of level indices
00776     std::vector<UGGridLevelIndexSet<const UGGrid<dim> >*> levelIndexSets_;
00777 
00778     UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
00779 
00780     UGGridIdSet<const UGGrid<dim>, false > globalIdSet_;
00781 
00782     UGGridIdSet<const UGGrid<dim>, true > localIdSet_;
00783 
00785     RefinementType refinementType_;
00786 
00788     ClosureType closureType_;
00789 
00792     std::vector<unsigned char> elementTypes_;
00793 
00796     std::vector<unsigned int> elementVertices_;
00797 
00799     std::vector<FieldVector<double, dim> > vertexPositions_;
00800 
00808     static int numOfUGGrids;
00809 
00815     bool someElementHasBeenMarkedForRefinement_;
00816 
00822     bool someElementHasBeenMarkedForCoarsening_;
00823 
00828     unsigned int heapsize;
00829 
00830 
00831 }; // end Class UGGrid
00832 
00833 namespace Capabilities
00834 {
00850   template<int dim>
00851   struct hasEntity< UGGrid<dim>, 0>
00852   {
00853     static const bool v = true;
00854   };
00855 
00859   template<int dim>
00860   struct hasEntity< UGGrid<dim>, dim>
00861   {
00862     static const bool v = true;
00863   };
00864   
00868   template<int dim>
00869   struct isParallel< UGGrid<dim> >
00870   {
00871 #ifdef ModelP
00872       static const bool v = true;
00873 #else
00874       static const bool v = false;
00875 #endif
00876   };
00877 
00881   template<int dim>
00882   struct isLevelwiseConforming< UGGrid<dim> >
00883   {
00884     static const bool v = true;
00885   };
00886 
00890   template<int dim>
00891   struct isLeafwiseConforming< UGGrid<dim> >
00892   {
00893     static const bool v = false;
00894   };
00895 
00899   template<int dim>
00900   struct hasHangingNodes< UGGrid<dim> >
00901   {
00902     static const bool v = true;
00903   };
00904   
00905 }
00906 
00907 } // namespace Dune
00908 
00909 #endif

Generated on Tue Jul 28 22:28:21 2009 for dune-grid by  doxygen 1.5.6