dune-grid  2.2.0
uggrid.hh
Go to the documentation of this file.
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/uggridlocalgeometry.hh"
00098 #include "uggrid/uggridentity.hh"
00099 #include "uggrid/uggridentitypointer.hh"
00100 #include "uggrid/uggridentityseed.hh"
00101 #include "uggrid/uggridintersections.hh"
00102 #include "uggrid/uggridintersectioniterators.hh"
00103 #include "uggrid/uggridleveliterator.hh"
00104 #include "uggrid/uggridleafiterator.hh"
00105 #include "uggrid/uggridhieriterator.hh"
00106 #include "uggrid/uggridindexsets.hh"
00107 
00108 // Not needed here, but included for user convenience
00109 #include "uggrid/uggridfactory.hh"
00110 
00111 #ifdef ModelP
00112 namespace Dune {
00113 
00114 // converts the UG speak message buffers to DUNE speak and vice-versa
00115     template <class DataHandle, int GridDim, int codim>
00116     class UGMessageBufferBase {
00117     protected:
00118         typedef UGMessageBufferBase<DataHandle, GridDim, codim>  ThisType;
00119         typedef UGGrid<GridDim>                              GridType;
00120         typedef typename DataHandle::DataType                DataType;
00121 
00122         enum {
00123             dim = GridDim
00124         };
00125 
00126         UGMessageBufferBase(void *ugData)
00127         {
00128             ugData_ = static_cast<char*>(ugData);
00129         };
00130         
00131     public:
00132         void write(const DataType &t) 
00133         { this->writeRaw_<DataType>(t);  }
00134         
00135         void read(DataType &t)
00136         { this->readRaw_<DataType>(t);  }
00137 
00138     protected:
00139         friend class Dune::UGGrid<dim>;
00140 
00141         template <class ValueType>
00142         void writeRaw_(const ValueType &v) 
00143         {
00144             *reinterpret_cast<ValueType*>(ugData_) = v;
00145             ugData_ += sizeof(ValueType);
00146         }
00147         
00148         template <class ValueType>
00149         void readRaw_(ValueType &v)
00150         {
00151             v = *reinterpret_cast<ValueType*>(ugData_);
00152             ugData_ += sizeof(ValueType);
00153         }
00154 
00155        // called by DDD_IFOneway to serialize the data structure to
00156         // be send
00157         static int ugGather_(typename UG_NS<dim>::DDD_OBJ obj, void* data)
00158         {
00159             if (codim == 0) {
00162                 UGMakeableEntity<0, dim, UGGrid<dim> > e((typename UG_NS<dim>::Element*) obj,nullptr);
00163                 // safety check to only communicate what is needed
00164                 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Element*) obj)) || e.level() == level)
00165                 {
00166                         ThisType msgBuf(static_cast<DataType*>(data));
00167                         if (!duneDataHandle_->fixedsize(dim, codim))
00168                                 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
00169                         duneDataHandle_->gather(msgBuf, e);
00170                 }
00171             }
00172             else if (codim == dim) {
00175                 UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > e((typename UG_NS<dim>::Node*)obj,nullptr);
00176                 // safety check to only communicate what is needed
00177                 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Node*)obj)) || e.level() == level)
00178                 {
00179                         ThisType msgBuf(static_cast<DataType*>(data));
00180                         if (!duneDataHandle_->fixedsize(dim, codim))
00181                                 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
00182                         duneDataHandle_->gather(msgBuf, e);
00183                 }
00184             }
00185             else if (codim == dim - 1) {
00188                 UGMakeableEntity<dim-1, dim, Dune::UGGrid<dim> > e((typename UG_NS<dim>::Edge*)obj,nullptr);
00189                 // safety check to only communicate what is needed
00190                 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Edge*)obj)) || e.level() == level)
00191                 {
00192                         ThisType msgBuf(static_cast<DataType*>(data));
00193                         if (!duneDataHandle_->fixedsize(dim, codim))
00194                                 msgBuf.template writeRaw_<unsigned>(duneDataHandle_->size(e));
00195                         duneDataHandle_->gather(msgBuf, e);
00196                 }
00197             }
00198             else {
00199                 DUNE_THROW(GridError, 
00200                            "Only node and element wise "
00201                            "communication is currently "
00202                            "supported by UGGrid");
00203             }
00204 
00205             return 0;
00206         }
00207 
00208         // called by DDD_IFOneway to deserialize the data structure
00209         // that has been received
00210         static int ugScatter_(typename UG_NS<dim>::DDD_OBJ obj, void* data)
00211         {
00212             if (codim == 0) {
00213                 typedef UGMakeableEntity<0, dim, UGGrid<dim> > Entity;
00216                 Entity e((typename UG_NS<dim>::Element*) obj,nullptr);
00217                 // safety check to only communicate what is needed
00218                 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Element*) obj)) || e.level() == level)
00219                 {
00220                         ThisType msgBuf(static_cast<DataType*>(data));
00221                         int n;
00222                         if (!duneDataHandle_->fixedsize(dim, codim))
00223                                 msgBuf.readRaw_(n);
00224                         else
00225                                 n = duneDataHandle_->template size<Entity>(e);
00226                         if (n > 0)
00227                                 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
00228                 }
00229             }
00230             else if (codim == dim) {
00231                 typedef UGMakeableEntity<dim, dim, Dune::UGGrid<dim> > Entity;
00234                 Entity e((typename UG_NS<dim>::Node*)obj,nullptr);
00235                 // safety check to only communicate what is needed
00236                 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Node*)obj)) || e.level() == level)
00237                 {
00238                         ThisType msgBuf(static_cast<DataType*>(data));
00239                         int n;
00240                         if (!duneDataHandle_->fixedsize(dim, codim))
00241                                 msgBuf.readRaw_(n);
00242                         else
00243                                 n = duneDataHandle_->template size<Entity>(e);
00244                         if (n > 0)
00245                                 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
00246                 }
00247             }
00248             else if (codim == dim - 1) {  // !!!ALEX!!! Is it possible to send codim 1 in UG<2>?
00249                 typedef UGMakeableEntity<dim-1, dim, Dune::UGGrid<dim> > Entity;
00252                 Entity e((typename UG_NS<dim>::Edge*)obj,nullptr);
00253                 // safety check to only communicate what is needed
00254                 if ((level == -1 && UG_NS<dim>::isLeaf((typename UG_NS<dim>::Edge*)obj)) || e.level() == level)
00255                 {
00256                         ThisType msgBuf(static_cast<DataType*>(data));
00257                         int n;
00258                         if (!duneDataHandle_->fixedsize(dim, codim))
00259                                 msgBuf.readRaw_(n);
00260                         else
00261                                 n = duneDataHandle_->template size<Entity>(e);
00262                         if (n > 0)
00263                                 duneDataHandle_->template scatter<ThisType, Entity>(msgBuf, e, n);
00264                 }
00265             }
00266             else {
00267                 DUNE_THROW(GridError, 
00268                            "Only node and element wise "
00269                            "communication is currently "
00270                            "supported by UGGrid");
00271             }
00272             
00273             return 0;
00274         }
00275         static DataHandle *duneDataHandle_;
00276 
00277         static int level;
00278 
00279         char              *ugData_;
00280     };
00281     
00282     template <class DataHandle, int GridDim, int codim>
00283     class UGMessageBuffer 
00284     : public UGMessageBufferBase<DataHandle, GridDim, codim>
00285     {
00286       typedef typename DataHandle::DataType                DataType;
00287       typedef UGMessageBufferBase<DataHandle, GridDim, codim> Base;
00288       enum { dim = GridDim };
00289  
00290       protected:
00291        friend class Dune::UGGrid<dim>;
00292 
00293        UGMessageBuffer(void *ugData)
00294         : Base(ugData)
00295         {}
00296         
00297         // returns number of bytes required for the UG message buffer
00298         template <class GridView>
00299         static unsigned ugBufferSize_(const GridView &gv)
00300         {
00301             if (Base::duneDataHandle_->fixedsize(dim, codim)) {
00302                 return sizeof(DataType)
00303                     * Base::duneDataHandle_->size(*gv.template begin<codim,InteriorBorder_Partition>());
00304             }
00305 
00306             typedef typename GridView::template Codim<codim>::Entity Entity;
00307 
00308             // iterate over all entities, find the maximum size for
00309             // the current rank
00310             int maxSize = 0;
00311             typedef typename 
00312                 GridView
00313                 ::template Codim<codim>
00314                 ::template Partition<Dune::All_Partition>
00315                 ::Iterator Iterator;
00316             Iterator it = gv.template begin<codim, Dune::All_Partition>();
00317             const Iterator endIt = gv.template end<codim, Dune::All_Partition>();
00318             for (; it != endIt; ++it) {
00319                 maxSize = std::max((int) maxSize, 
00320                                    (int) Base::duneDataHandle_->size(*it));
00321             }
00322             
00323             // find maximum size for all ranks
00324             maxSize = MPIHelper::getCollectiveCommunication().max(maxSize);
00325             if (!maxSize)
00326                 return 0;
00327             
00328             // add the size of an unsigned integer to the actual
00329             // buffer size. (we somewhere have to store the actual
00330             // number of objects for each entity.)
00331             return sizeof(unsigned) + sizeof(DataType)*maxSize;
00332         }
00333     };
00334     
00335     template <class DataHandle, int GridDim>
00336     class UGEdgeMessageBuffer
00337     : public UGMessageBufferBase<DataHandle, GridDim, GridDim-1>
00338     {
00339       enum {codim = GridDim-1, 
00340         dim = GridDim};
00341       typedef typename DataHandle::DataType                DataType;
00342       typedef UGMessageBufferBase<DataHandle, GridDim, codim> Base;
00343     protected:
00344        friend class Dune::UGGrid<dim>;
00345 
00346        UGEdgeMessageBuffer(void *ugData)
00347         : Base(ugData)
00348         {}
00349         
00350         // returns number of bytes required for the UG message buffer
00351         template <class GridView>
00352         static unsigned ugBufferSize_(const GridView &gv)
00353         {
00354              if (Base::duneDataHandle_->fixedsize(dim, codim)) {
00355                typedef typename GridView::template Codim<0>::Entity Element;
00356                const Element& element = *gv.template begin<0, InteriorBorder_Partition>();
00357                 return sizeof(DataType)
00358                     * Base::duneDataHandle_->size(*element.template subEntity<codim>(0));
00359             }
00360 
00361             typedef typename GridView::template Codim<codim>::Entity Entity;
00362 
00363             // iterate over all entities, find the maximum size for
00364             // the current rank
00365             int maxSize = 0;
00366             typedef typename 
00367                 GridView
00368                 ::template Codim<0>
00369                 ::template Partition<Dune::All_Partition>
00370                 ::Iterator Iterator;
00371             Iterator it = gv.template begin<0, Dune::All_Partition>();
00372             const Iterator endIt = gv.template end<0, Dune::All_Partition>();
00373             for (; it != endIt; ++it) {
00374                 int numberOfEdges = it->template count<codim>();
00375                 for (int k = 0; k < numberOfEdges; k++)
00376                     {
00377                         typedef typename GridView::template Codim<0>::Entity Element;
00378                         typedef typename Element::template Codim<codim>::EntityPointer EdgePointer;
00379                         const EdgePointer edgePointer(it->template subEntity<codim>(k));
00380                         
00381                         maxSize = std::max((int) maxSize, 
00382                                            (int) Base::duneDataHandle_->size(*edgePointer));
00383                     }
00384             }
00385             
00386             // find maximum size for all ranks
00387             maxSize = MPIHelper::getCollectiveCommunication().max(maxSize);
00388             if (!maxSize)
00389                 return 0;
00390             
00391             // add the size of an unsigned integer to the actual
00392             // buffer size. (we somewhere have to store the actual
00393             // number of objects for each entity.)
00394             return sizeof(unsigned) + sizeof(DataType)*maxSize;
00395        }
00396     };
00397 
00398     template <class DataHandle>
00399     class UGMessageBuffer<DataHandle, 2, 1>
00400     : public UGEdgeMessageBuffer<DataHandle, 2>
00401     {};
00402 
00403     template <class DataHandle>
00404     class UGMessageBuffer<DataHandle, 3, 2>
00405     : public UGEdgeMessageBuffer<DataHandle, 3>
00406     {};
00407 
00408 }   // end namespace Dune
00409 
00410 template <class DataHandle, int GridDim, int codim>
00411 DataHandle *Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::duneDataHandle_ = 0;
00412 
00413 template <class DataHandle, int GridDim, int codim>
00414 int Dune::UGMessageBufferBase<DataHandle,GridDim,codim>::level = -1;
00415 #endif // ModelP
00416 
00417 namespace Dune {
00418 
00419 #ifdef ModelP
00420 template <int dim>
00421 class CollectiveCommunication<Dune::UGGrid<dim> > :
00422     public CollectiveCommunication< Dune::MPIHelper::MPICommunicator >
00423 {
00424     typedef CollectiveCommunication< Dune::MPIHelper::MPICommunicator > ParentType;
00425 public:
00426     CollectiveCommunication()
00427         : ParentType(MPIHelper::getCommunicator())
00428     {}
00429 };
00430 #endif // ModelP
00431 
00432 template<int dim, int dimworld>
00433 struct UGGridFamily
00434 {
00435   typedef GridTraits<dim,dimworld,Dune::UGGrid<dim>,
00436                      UGGridGeometry,
00437                      UGGridEntity,
00438                      UGGridEntityPointer,
00439                      UGGridLevelIterator,
00440                      UGGridLeafIntersection,
00441                      UGGridLevelIntersection,
00442                      UGGridLeafIntersectionIterator,
00443                      UGGridLevelIntersectionIterator,
00444                      UGGridHierarchicIterator,
00445                      UGGridLeafIterator,
00446                      UGGridLevelIndexSet< const UGGrid<dim> >,
00447                      UGGridLeafIndexSet< const UGGrid<dim> >,
00448                      UGGridIdSet< const UGGrid<dim>, false >,
00449                      unsigned int,
00450                      UGGridIdSet< const UGGrid<dim>, true >,
00451                      unsigned int,
00452                      CollectiveCommunication<Dune::UGGrid<dim> >,
00453                      DefaultLevelGridViewTraits,
00454                      DefaultLeafGridViewTraits,
00455                      UGGridEntitySeed,
00456                      UGGridLocalGeometry>
00457   Traits;
00458 };
00459 
00460 
00461 //**********************************************************************
00462 //
00463 // --UGGrid
00464 //
00465 //**********************************************************************
00466 
00503 template <int dim>
00504 class UGGrid : public GridDefaultImplementation  <dim, dim, double, UGGridFamily<dim,dim> >
00505 {
00506     friend class UGGridGeometry<0,dim,const UGGrid<dim> >;
00507     friend class UGGridGeometry<dim,dim,const UGGrid<dim> >;
00508     friend class UGGridGeometry<1,2,const UGGrid<dim> >;
00509     friend class UGGridGeometry<2,3,const UGGrid<dim> >;
00510 
00511     friend class UGGridEntity <0,dim,const UGGrid<dim> >;
00512     friend class UGGridEntity <dim,dim,const UGGrid<dim> >;
00513     friend class UGGridHierarchicIterator<const UGGrid<dim> >;
00514     friend class UGGridLeafIntersection<const UGGrid<dim> >;
00515     friend class UGGridLevelIntersection<const UGGrid<dim> >;
00516     friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >;
00517     friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >;
00518 
00519     friend class UGGridLevelIndexSet<const UGGrid<dim> >;
00520     friend class UGGridLeafIndexSet<const UGGrid<dim> >;
00521     friend class UGGridIdSet<const UGGrid<dim>, false >;
00522     friend class UGGridIdSet<const UGGrid<dim>, true >;
00523 
00524     friend class GridFactory<UGGrid<dim> >;
00525 
00526     template <int codim_, PartitionIteratorType PiType_, class GridImp_>
00527     friend class UGGridLeafIterator;
00528     template <int codim_, PartitionIteratorType PiType_, class GridImp_>
00529     friend class UGGridLevelIterator;
00530     template <int codim_, class GridImp_>
00531     friend class UGGridEntityPointer;
00532 
00534     dune_static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!");
00535 
00536     // The different instantiations are mutual friends so they can access
00537     // each others numOfUGGrids field
00538     friend class UGGrid<2>;
00539     friend class UGGrid<3>;
00540 
00541     //**********************************************************
00542       // The Interface Methods
00543       //**********************************************************
00544 public:  
00546       typedef UGGridFamily<dim,dim>  GridFamily;
00547 
00548     // the Traits 
00549     typedef typename UGGridFamily<dim,dim>::Traits Traits;
00550 
00552     typedef UG::DOUBLE ctype;
00553 
00559     UGGrid();
00560 
00562     ~UGGrid();
00563    
00566      int maxLevel() const;
00567      
00569     template<int codim>
00570     typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
00571 
00573     template<int codim>
00574     typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
00575 
00577     template<int codim, PartitionIteratorType PiType>
00578     typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
00579 
00581     template<int codim, PartitionIteratorType PiType>
00582     typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
00583 
00585     template<int codim>
00586     typename Traits::template Codim<codim>::LeafIterator leafbegin() const {
00587         return typename Traits::template Codim<codim>::template Partition<All_Partition>::LeafIterator(*this);
00588     }
00589 
00591     template<int codim>
00592     typename Traits::template Codim<codim>::LeafIterator leafend() const {
00593         return UGGridLeafIterator<codim,All_Partition, const UGGrid<dim> >();
00594     }
00595 
00597     template<int codim, PartitionIteratorType PiType>
00598     typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const {
00599         return typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator(*this);
00600     }
00601 
00603     template<int codim, PartitionIteratorType PiType>
00604     typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const {
00605         return UGGridLeafIterator<codim,PiType, const UGGrid<dim> >();
00606     }
00607 
00609     template <int codim>
00610     typename Traits::template Codim<codim>::EntityPointer
00611     entityPointer(const UGGridEntitySeed<codim, const UGGrid<dim> >& seed) const
00612     {
00613         return typename Traits::template Codim<codim>::EntityPointer(UGGridEntityPointer<codim,const UGGrid<dim> >(seed.target(),this));
00614     }
00615 
00618     int size (int level, int codim) const;
00619 
00621   int size (int codim) const
00622   {
00623       return leafIndexSet().size(codim);
00624   }
00625 
00627   int size (int level, GeometryType type) const
00628   {
00629         return this->levelIndexSet(level).size(type);
00630   }
00631 
00633   int size (GeometryType type) const
00634   {
00635         return this->leafIndexSet().size(type);
00636   }
00637 
00639     size_t numBoundarySegments() const {
00640         // The number is stored as a member of UGGrid upon grid creation.
00641         // The corresponding data structure is not exported by UG.  (It is in ug/dom/std/std_internal.h)
00642         return numBoundarySegments_;
00643     }
00644 
00646     const typename Traits::GlobalIdSet& globalIdSet() const
00647     {
00648         return globalIdSet_;
00649     }
00650     
00652     const typename Traits::LocalIdSet& localIdSet() const
00653     {
00654         return localIdSet_;
00655     }
00656     
00658     const typename Traits::LevelIndexSet& levelIndexSet(int level) const
00659     {
00660         if (level<0 || level>maxLevel())
00661             DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
00662         return *levelIndexSets_[level];
00663     }
00664     
00666     const typename Traits::LeafIndexSet& leafIndexSet() const
00667     {
00668         return leafIndexSet_;
00669     }
00670 
00673 
00686     bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e );
00687 
00695     bool mark(const typename Traits::template Codim<0>::Entity & e, 
00696               typename UG_NS<dim>::RefinementRule rule,
00697               int side=0);
00698 
00700     int getMark(const typename Traits::template Codim<0>::Entity& e) const;
00701 
00704     bool preAdapt();
00705     
00707     bool adapt();
00708 
00710     void postAdapt();
00714     unsigned int overlapSize(int codim) const {
00715         return 0;
00716     }
00717 
00719     unsigned int ghostSize(int codim) const {
00720         return (codim==0) ? 1 : 0;
00721     }
00722 
00724     unsigned int overlapSize(int level, int codim) const {
00725         return 0;
00726     }
00727 
00729     unsigned int ghostSize(int level, int codim) const {
00730         return (codim==0) ? 1 : 0;
00731     }
00732 
00737     bool loadBalance() {
00738         return loadBalance(0,0,2,32,1);
00739     }
00740 
00746     template<class DataHandle>
00747     bool loadBalance (DataHandle& data)
00748     {
00749         DUNE_THROW(NotImplemented, "load balancing with data attached");
00750     }
00751 
00768     bool loadBalance(int strategy, int minlevel, int depth, int maxlevel, int minelement);
00769 
00780     template<class DataHandle>
00781     void communicate (DataHandle& dataHandle, InterfaceType iftype, CommunicationDirection dir, int level) const
00782     {
00783 #ifdef ModelP  
00784         typedef typename UGGrid::LevelGridView LevelGridView;
00785 
00786         for (int curCodim = 0; curCodim <= dim; ++curCodim) {
00787                 if (!dataHandle.contains(dim, curCodim))
00788                         continue;
00789 
00790                 if (curCodim == 0)
00791                         communicateUG_<LevelGridView, DataHandle, 0>(this->levelView(level), level, dataHandle, iftype, dir);
00792                 else if (curCodim == dim)
00793                         communicateUG_<LevelGridView, DataHandle, dim>(this->levelView(level), level, dataHandle, iftype, dir);
00794         else if (curCodim == dim - 1)   
00795                         communicateUG_<LevelGridView, DataHandle, dim-1>(this->levelView(level), level, dataHandle, iftype, dir);
00796             else
00797                 DUNE_THROW(NotImplemented,
00798                            className(*this) << "::communicate(): Only "
00799                            "supported for codim=0 and "
00800                            "codim=dim(=" << dim << "), but "
00801                            "codim=" << curCodim << " was requested");
00802         }
00803 #endif // ModelP
00804     }
00805 
00815     template<class DataHandle>
00816     void communicate(DataHandle& dataHandle, 
00817                 InterfaceType iftype,
00818                 CommunicationDirection dir) const
00819     {
00820 #ifdef ModelP  
00821         typedef typename UGGrid::LeafGridView LeafGridView;
00822 
00823         for (int curCodim = 0; curCodim <= dim; ++curCodim) {
00824                 if (!dataHandle.contains(dim, curCodim))
00825                         continue;
00826                 int level = -1;
00827                 if (curCodim == 0)
00828                         communicateUG_<LeafGridView, DataHandle, 0>(this->leafView(), level, dataHandle, iftype, dir);
00829                 else if (curCodim == dim)
00830                         communicateUG_<LeafGridView, DataHandle, dim>(this->leafView(), level, dataHandle, iftype, dir);
00831         else if (curCodim == dim - 1)   // !!!ALEX!!! Is it possible to send codim 1 in UG<2>?
00832         {
00833                         communicateUG_<LeafGridView, DataHandle, dim-1>(this->leafView(), level, dataHandle, iftype, dir);
00834         }
00835         else
00836           DUNE_THROW(NotImplemented,
00837                      className(*this) << "::communicate(): Only "
00838                      "supported for codim=0 and "
00839                      "codim=dim(=" << dim << "), but "
00840                      "codim=" << curCodim << " was requested");
00841         }
00842 #endif // ModelP
00843     }
00844 
00846         const CollectiveCommunication<UGGrid>& comm () const
00847         {
00848         return ccobj_;
00849         }
00850 
00851 protected:
00852 #ifdef ModelP
00853     template <class GridView, class DataHandle, int codim>
00854     void communicateUG_(const GridView& gv, int level,
00855                         DataHandle &dataHandle, 
00856                         InterfaceType iftype,
00857                         CommunicationDirection dir) const
00858     {
00859         typename UG_NS<dim>::DDD_IF_DIR ugIfDir;
00860         // Translate the communication direction from Dune-Speak to UG-Speak
00861         if (dir==ForwardCommunication) 
00862             ugIfDir = UG_NS<dim>::IF_FORWARD();
00863         else 
00864             ugIfDir = UG_NS<dim>::IF_BACKWARD();
00865 
00866         typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf;
00867         UGMsgBuf::duneDataHandle_ = &dataHandle;              
00868         
00869         UGMsgBuf::level = level;
00870 
00871         std::vector<typename UG_NS<dim>::DDD_IF> ugIfs;
00872         findDDDInterfaces_(ugIfs, iftype, codim);
00873 
00874         unsigned bufSize = UGMsgBuf::ugBufferSize_(gv);
00875         if (!bufSize)
00876             return; // we don't need to communicate if we don't have any data!
00877         for (unsigned i=0; i < ugIfs.size(); ++i)
00878             UG_NS<dim>::DDD_IFOneway(ugIfs[i],
00879                                      ugIfDir,
00880                                      bufSize,
00881                                      &UGMsgBuf::ugGather_,
00882                                      &UGMsgBuf::ugScatter_);
00883     }
00884 
00885     void findDDDInterfaces_(std::vector<typename UG_NS<dim>::DDD_IF > &dddIfaces,
00886                             InterfaceType iftype, 
00887                             int codim) const
00888     {
00889         dddIfaces.clear();
00890         switch (codim)
00891         {
00892         case 0:
00893             switch (iftype) {
00894             case InteriorBorder_InteriorBorder_Interface:
00895                 // do not communicate anything: Elements can not be in
00896                 // the interior of two processes at the same time
00897                 return;
00898             case InteriorBorder_All_Interface:
00899                 // The communicated elements are in the sender's
00900                 // interior and it does not matter what they are for
00901                 // the receiver
00902                 dddIfaces.push_back(UG_NS<dim>::ElementVHIF());
00903                 return;
00904             case All_All_Interface:
00905                 // It does neither matter what the communicated
00906                 // elements are for sender nor for the receiver. If
00907                 // they are seen by these two processes, data is send
00908                 // and received.
00909                 dddIfaces.push_back(UG_NS<dim>::ElementSymmVHIF());
00910                 return;
00911             default:
00912                 DUNE_THROW(GridError, 
00913                            "Element communication not supported for "
00914                            "interfaces of type  "
00915                            << iftype);
00916             }
00917 
00918         case dim:
00919             switch (iftype)
00920             {
00921             case InteriorBorder_InteriorBorder_Interface:
00922                 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
00923                 return;
00924             case InteriorBorder_All_Interface:
00925                 dddIfaces.push_back(UG_NS<dim>::BorderNodeSymmIF());
00926                 dddIfaces.push_back(UG_NS<dim>::NodeIF());
00927                 return;
00928             case All_All_Interface:
00929                 dddIfaces.push_back(UG_NS<dim>::NodeAllIF());
00930                 return;
00931             default:
00932                 DUNE_THROW(GridError, 
00933                            "Node communication not supported for "
00934                            "interfaces of type  "
00935                            << iftype);
00936             }
00937 
00938         case dim-1:  
00939             switch (iftype)
00940             {   
00941             case InteriorBorder_InteriorBorder_Interface:
00942                 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
00943                 return;
00944             case InteriorBorder_All_Interface:
00945                 dddIfaces.push_back(UG_NS<dim>::BorderEdgeSymmIF());
00946                 // Is the following line needed or not?
00947                 // dddIfaces.push_back(UG_NS<dim>::EdgeIF());
00948                 return;
00949             default:
00950               DUNE_THROW(GridError,
00951                            "Edge communication not supported for "
00952                            "interfaces of type  "
00953                            << iftype);
00954             }
00955 
00956         default:
00957             DUNE_THROW(GridError, 
00958                        "Communication for codim "
00959                        << codim
00960                        << " entities is not yet supported "
00961                        << " by the DUNE UGGrid interface!");
00962         }
00963     };
00964 #endif // ModelP
00965 
00966 public:
00967     // **********************************************************
00968     // End of Interface Methods
00969     // **********************************************************
00970     
00973     
00977     void createLGMGrid(const std::string& name);
00978 
00988     void getChildrenOfSubface(const typename Traits::template Codim<0>::EntityPointer & e,
00989                               int elementSide,
00990                               int maxl, 
00991                               std::vector<typename Traits::template Codim<0>::EntityPointer>& childElements,
00992                               std::vector<unsigned char>& childElementSides) const;
00993     
00995     enum RefinementType {
00997         LOCAL, 
00999         COPY};
01000 
01002     enum ClosureType {
01004         GREEN,
01006         NONE};
01007 
01009     void setRefinementType(RefinementType type) {
01010         refinementType_ = type;
01011     }
01012 
01014     void setClosureType(ClosureType type) {
01015         closureType_ = type;
01016     }
01017 
01024     static void setDefaultHeapSize(unsigned size) {
01025         heapSize_ = size;
01026     }
01027 
01031     void setPosition(const typename Traits::template Codim<dim>::EntityPointer& e,
01032                      const FieldVector<double, dim>& pos);
01033 
01038     void globalRefine(int n);
01039 
01044     void saveState(const std::string& filename) const;
01045 
01050     void loadState(const std::string& filename);
01051 
01052 private:
01054     typename UG_NS<dim>::MultiGrid* multigrid_;
01055 
01057     CollectiveCommunication<UGGrid> ccobj_;
01058 
01064     void setIndices(bool setLevelZero,
01065                     std::vector<unsigned int>* nodePermutation);
01066 
01067     // Each UGGrid object has a unique name to identify it in the
01068     // UG environment structure
01069     std::string name_;
01070 
01071     // Our set of level indices
01072     std::vector<UGGridLevelIndexSet<const UGGrid<dim> >*> levelIndexSets_;
01073 
01074     UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_;
01075 
01076     UGGridIdSet<const UGGrid<dim>, false > globalIdSet_;
01077 
01078     UGGridIdSet<const UGGrid<dim>, true > localIdSet_;
01079 
01081     RefinementType refinementType_;
01082 
01084     ClosureType closureType_;
01085 
01093     static int numOfUGGrids;
01094 
01100     bool someElementHasBeenMarkedForRefinement_;
01101 
01107     bool someElementHasBeenMarkedForCoarsening_;
01108 
01113     static unsigned int heapSize_;
01114 
01116     std::vector<shared_ptr<BoundarySegment<dim> > > boundarySegments_;
01117 
01123     unsigned int numBoundarySegments_;
01124 
01125 }; // end Class UGGrid
01126 
01127 namespace Capabilities
01128 {
01144   template<int dim>
01145   struct hasEntity< UGGrid<dim>, 0>
01146   {
01147     static const bool v = true;
01148   };
01149 
01153   template<int dim>
01154   struct hasEntity< UGGrid<dim>, dim>
01155   {
01156     static const bool v = true;
01157   };
01158   
01162   template<int dim>
01163   struct isParallel< UGGrid<dim> >
01164   {
01165 #ifdef ModelP
01166       static const bool v = true;
01167 #else // !ModelP
01168       static const bool v = false;
01169 #endif // !ModelP
01170   };
01171 
01175   template<int dim>
01176   struct isLevelwiseConforming< UGGrid<dim> >
01177   {
01178     static const bool v = true;
01179   };
01180 
01184   template<int dim>
01185   struct isLeafwiseConforming< UGGrid<dim> >
01186   {
01187     static const bool v = false;
01188   };
01189 
01190 }
01191 
01192 } // namespace Dune
01193 
01194 #endif   // HAVE_UG
01195 #endif   // DUNE_UGGRID_HH