dune-grid  2.2.0
elementinfo.hh
Go to the documentation of this file.
00001 #ifndef DUNE_ALBERTA_ELEMENTINFO_HH
00002 #define DUNE_ALBERTA_ELEMENTINFO_HH
00003 
00009 #include <cassert>
00010 #include <vector>
00011 
00012 #include <dune/grid/albertagrid/geometrycache.hh>
00013 #include <dune/grid/albertagrid/macroelement.hh>
00014 
00015 #if HAVE_ALBERTA
00016 
00017 namespace Dune
00018 {
00019 
00020   namespace Alberta
00021   {
00022 
00023     // External Forward Declarations
00024     // -----------------------------
00025 
00026     template< int dim >
00027     class MeshPointer;
00028 
00029     class BasicNodeProjection;
00030 
00031 
00032 
00033     // ElementInfo
00034     // -----------
00035 
00036     template< int dim >
00037     class ElementInfo
00038     {
00039       class Instance;
00040       class Stack;
00041 
00042       template< int >
00043       struct Library;
00044 
00045       typedef Instance *InstancePtr;
00046 
00047     public:
00048       static const int dimension = dim;
00049 
00050       static const int numVertices = NumSubEntities< dimension, dimension >::value;
00051       static const int numFaces = NumSubEntities< dimension, 1 >::value;
00052 
00053       typedef Alberta::MacroElement< dimension > MacroElement;
00054       typedef Alberta::MeshPointer< dimension > MeshPointer;
00055       typedef Alberta::FillFlags< dimension > FillFlags;
00056 
00057       static const int maxNeighbors = N_NEIGH_MAX;
00058 
00059       static const int maxLevelNeighbors = Library< dimWorld >::maxLevelNeighbors;
00060 
00061 #if !DUNE_ALBERTA_CACHE_COORDINATES
00062       typedef GeometryCacheProxy< dim > GeometryCache;
00063 #endif
00064 
00065       struct Seed;
00066 
00067     private:
00068       explicit ElementInfo ( const InstancePtr &instance );
00069 
00070     public:
00071       ElementInfo ();
00072       ElementInfo ( const MeshPointer &mesh, const MacroElement &macroElement,
00073                     typename FillFlags::Flags fillFlags = FillFlags::standard );
00074       ElementInfo ( const MeshPointer &mesh, const Seed &seed,
00075                     typename FillFlags::Flags fillFlags = FillFlags::standard );
00076       ElementInfo ( const ElementInfo &other );
00077 
00078       ~ElementInfo ();
00079 
00080       ElementInfo &operator= ( const ElementInfo &other );
00081 
00082       operator bool () const { return (instance_ != null()); }
00083 
00084       bool operator== ( const ElementInfo &other ) const;
00085       bool operator!= ( const ElementInfo &other ) const;
00086 
00087       const MacroElement &macroElement () const;
00088       ElementInfo father () const;
00089       int indexInFather () const;
00090       ElementInfo child ( int i ) const;
00091       bool isLeaf () const;
00092 
00093       Seed seed () const;
00094 
00095       MeshPointer mesh () const;
00096 
00097       bool mightVanish () const;
00098 
00099       int level () const;
00100       // see ALBERTA documentation for definition of element type
00101       // values are 0, 1, 2
00102       int type () const;
00103 
00104       int getMark () const;
00105       void setMark ( int refCount ) const;
00106 
00107       bool hasLeafNeighbor ( const int face ) const;
00108       ElementInfo leafNeighbor ( const int face ) const;
00109 
00110       /* obtain all level neighbors of a face
00111        *
00112        * param[in]  face            face for which the neighbors are desired
00113        * param[out] neighbor        array storing the neighbors
00114        * param[out] faceInNeighbor  array storing the faces in neighbor
00115        *                            (-1, if this neighbor does not exist)
00116        *
00117        * returns (potential) number of neighbors (i.e., the number of valid
00118        *         entries in the output arrays
00119        */
00120       int levelNeighbors ( const int face, ElementInfo (&neighbor)[ maxLevelNeighbors ], int (&faceInNeighbor)[ maxLevelNeighbors ] ) const;
00121 
00122       template< int codim >
00123       int twist ( int subEntity ) const;
00124       int twistInNeighbor ( int face ) const;
00125       bool isBoundary ( int face ) const;
00126       int boundaryId ( int face ) const;
00127       AffineTransformation *transformation ( int face ) const;
00128       BasicNodeProjection *boundaryProjection ( int face ) const;
00129 
00130       bool hasCoordinates () const;
00131       const GlobalVector &coordinate ( int vertex ) const;
00132 #if !DUNE_ALBERTA_CACHE_COORDINATES
00133       GeometryCache geometryCache () const
00134       {
00135         return GeometryCache( instance_->geometryCache, instance_->elInfo );
00136       }
00137 #endif
00138 
00139       template< class Functor >
00140       void hierarchicTraverse ( Functor &functor ) const;
00141 
00142       template< class Functor >
00143       void leafTraverse ( Functor &functor ) const;
00144 
00145       const Element *element () const;
00146       const Element *neighbor ( int face ) const;
00147       Element *el () const;
00148       ALBERTA EL_INFO &elInfo () const;
00149 
00150       static ElementInfo
00151       createFake ( const MeshPointer &mesh,
00152                    const Element *element, int level, int type = 0 );
00153       static ElementInfo createFake ( const ALBERTA EL_INFO &elInfo );
00154 
00155     private:
00156       static bool isLeaf ( Element *element );
00157       static bool mightVanish ( Element *element, int depth );
00158 
00159       static void fill ( Mesh *mesh, const ALBERTA MACRO_EL *mel, ALBERTA EL_INFO &elInfo );
00160       static void fill ( int ichild, const ALBERTA EL_INFO &parentInfo, ALBERTA EL_INFO &elInfo );
00161 
00162       void addReference () const;
00163       void removeReference () const;
00164 
00165       static InstancePtr null ();
00166       static Stack &stack ();
00167 
00168       InstancePtr instance_;
00169     };
00170 
00171 
00172 
00173     // ElementInfo::Instance
00174     // ---------------------
00175 
00176     template< int dim >
00177     struct ElementInfo< dim >::Instance
00178     {
00179       ALBERTA EL_INFO elInfo;
00180       unsigned int refCount;
00181 
00182       InstancePtr &parent ()
00183       {
00184         return parent_;
00185       }
00186 
00187     private:
00188       InstancePtr parent_;
00189 
00190 #if !DUNE_ALBERTA_CACHE_COORDINATES
00191     public:
00192       Alberta::GeometryCache< dim > geometryCache;
00193 #endif
00194     };
00195 
00196 
00197 
00198     // ElementInfo::Stack
00199     // ------------------
00200 
00201     template< int dim >
00202     class ElementInfo< dim >::Stack
00203     {
00204       InstancePtr top_;
00205       Instance null_;
00206 
00207     public:
00208       Stack ();
00209       ~Stack ();
00210 
00211       InstancePtr allocate ();
00212       void release ( InstancePtr &p );
00213       InstancePtr null ();
00214     };
00215 
00216 
00217 
00218     // ElementInfo::Library
00219     // --------------------
00220 
00221     template< int dim >
00222     template< int >
00223     struct ElementInfo< dim >::Library
00224     {
00225       typedef Alberta::ElementInfo< dim > ElementInfo;
00226 
00227       static const int maxLevelNeighbors = (1 << (dim-1));
00228 
00229       static int
00230       leafNeighbor ( const ElementInfo &element, const int face, ElementInfo &neighbor );
00231 
00232       static int
00233       levelNeighbors ( const ElementInfo &element, const int face,
00234                        ElementInfo (&neighbor)[ maxLevelNeighbors ], int (&faceInNeighbor)[ maxLevelNeighbors ] );
00235 
00236     private:
00237       static int
00238       macroNeighbor ( const ElementInfo &element, const int face, ElementInfo &neighbor );
00239     };
00240 
00241 
00242 
00243     // ElementInfo::Seed
00244     // -----------------
00245 
00246     template< int dim >
00247     struct ElementInfo< dim >::Seed
00248     {
00249       Seed ( const int macroIndex, const int level, const unsigned long path )
00250       : macroIndex_( macroIndex ), level_( level ), path_( path )
00251       {}
00252 
00253       bool operator== ( const Seed &other ) const
00254       {
00255         return (macroIndex() == other.macroIndex()) && (level() == other.level()) && (path() == other.path());
00256       }
00257 
00258       bool operator< ( const Seed &other ) const
00259       {
00260         const bool ml = (macroIndex() < other.macroIndex());
00261         const bool me = (macroIndex() == other.macroIndex());
00262         const bool ll = (level() < other.level());
00263         const bool le = (level() == other.level());
00264         const bool pl = (path < other.path());
00265         return ml | (me & (ll | (le & pl)));
00266       }
00267 
00268       bool operator!= ( const Seed &other ) const { return !(*this == other); }
00269       bool operator<= ( const Seed &other ) const { return !(other < *this); }
00270       bool operator> ( const Seed &other ) const { return (other < *this); }
00271       bool operator>= ( const Seed &other ) const { return !(*this < other); }
00272 
00273       int macroIndex () const { return macroIndex_; }
00274       int level () const { return level_; }
00275       unsigned long path () const { return path_; }
00276 
00277     private:
00278       int macroIndex_;
00279       int level_;
00280       unsigned long path_;
00281     };
00282 
00283 
00284 
00285     // Implementation of ElementInfo
00286     // -----------------------------
00287 
00288     template< int dim >
00289     inline ElementInfo< dim >::ElementInfo ( const InstancePtr &instance )
00290     : instance_( instance )
00291     {
00292       addReference();
00293     }
00294 
00295 
00296     template< int dim >
00297     inline ElementInfo< dim >::ElementInfo ()
00298     : instance_( null() )
00299     {
00300       addReference();
00301     }
00302 
00303 
00304     template< int dim >
00305     inline ElementInfo< dim >
00306       ::ElementInfo ( const MeshPointer &mesh, const MacroElement &macroElement,
00307                       typename FillFlags::Flags fillFlags )
00308     {
00309       instance_ = stack().allocate();
00310       instance_->parent() = null();
00311       ++(instance_->parent()->refCount);
00312 
00313       addReference();
00314 
00315       elInfo().fill_flag = fillFlags;
00316 
00317       // Alberta fills opp_vertex only if there is a neighbor
00318       for( int k = 0; k < maxNeighbors; ++k )
00319         elInfo().opp_vertex[ k ] = -1;
00320 
00321       fill( mesh, &macroElement, elInfo() );
00322     }
00323 
00324 
00325     template< int dim >
00326     inline ElementInfo< dim >
00327       ::ElementInfo ( const MeshPointer &mesh, const Seed &seed,
00328                       typename FillFlags::Flags fillFlags )
00329     {
00330       instance_ = stack().allocate();
00331       instance_->parent() = null();
00332       ++(instance_->parent()->refCount);
00333 
00334       addReference();
00335 
00336       // fill in macro element info
00337       elInfo().fill_flag = fillFlags;
00338 
00339       // Alberta fills opp_vertex only if there is a neighbor
00340       for( int k = 0; k < maxNeighbors; ++k )
00341         elInfo().opp_vertex[ k ] = -1;
00342 
00343       fill( mesh, ((Mesh *)mesh)->macro_els + seed.macroIndex(), elInfo() );
00344 
00345       // traverse the seed's path
00346       unsigned long path = seed.path();
00347       for( int i = 0; i < seed.level(); ++i )
00348       {
00349         InstancePtr child = stack().allocate();
00350         child->parent() = instance_;
00351 
00352         // Alberta fills opp_vertex only if there is a neighbor
00353         for( int k = 0; k < maxNeighbors; ++k )
00354           child->elInfo.opp_vertex[ k ] = -2;
00355 
00356         fill( path & 1, elInfo(), child->elInfo );
00357 
00358         instance_ = child;
00359         addReference();
00360 
00361         path = path >> 1;
00362       }
00363 
00364       assert( this->seed() == seed );
00365     }
00366 
00367 
00368     template< int dim >
00369     inline ElementInfo< dim >::ElementInfo ( const ElementInfo &other )
00370     : instance_( other.instance_ )
00371     {
00372       addReference();
00373     }
00374 
00375 
00376     template< int dim >
00377     inline ElementInfo< dim >::~ElementInfo ()
00378     {
00379       removeReference();
00380     }
00381 
00382 
00383     template< int dim >
00384     inline ElementInfo< dim > &
00385     ElementInfo< dim >::operator= ( const ElementInfo< dim > &other )
00386     {
00387       other.addReference();
00388       removeReference();
00389       instance_ = other.instance_;
00390       return *this;
00391     }
00392 
00393 
00394     template< int dim >
00395     inline bool
00396     ElementInfo< dim >::operator== ( const ElementInfo< dim > &other ) const
00397     {
00398       return (instance_->elInfo.el == other.instance_->elInfo.el);
00399     }
00400 
00401 
00402     template< int dim >
00403     inline bool
00404     ElementInfo< dim >::operator!= ( const ElementInfo< dim > &other ) const
00405     {
00406       return (instance_->elInfo.el != other.instance_->elInfo.el);
00407     }
00408 
00409 
00410     template< int dim >
00411     inline const typename ElementInfo< dim >::MacroElement &
00412     ElementInfo< dim >::macroElement () const
00413     {
00414       assert( !!(*this) );
00415       assert( elInfo().macro_el != NULL );
00416       return static_cast< const MacroElement & >( *(elInfo().macro_el) );
00417     }
00418 
00419 
00420     template< int dim >
00421     inline ElementInfo< dim > ElementInfo< dim >::father () const
00422     {
00423       assert( !!(*this) );
00424       return ElementInfo< dim >( instance_->parent() );
00425     }
00426 
00427 
00428     template< int dim >
00429     inline int ElementInfo< dim >::indexInFather () const
00430     {
00431       const Element *element = elInfo().el;
00432 #if DUNE_ALBERTA_VERSION >= 0x300
00433       const Element *father = elInfo().parent->el;
00434 #else
00435       const Element *father = elInfo().parent;
00436 #endif
00437       assert( father != NULL );
00438 
00439       const int index = (father->child[ 0 ] == element ? 0 : 1);
00440       assert( father->child[ index ] == element );
00441       return index;
00442     }
00443 
00444 
00445     template< int dim >
00446     inline ElementInfo< dim > ElementInfo< dim >::child ( int i ) const
00447     {
00448       assert( !isLeaf() );
00449 
00450       InstancePtr child = stack().allocate();
00451       child->parent() = instance_;
00452       addReference();
00453 
00454       // Alberta fills opp_vertex only if there is a neighbor
00455       for( int k = 0; k < maxNeighbors; ++k )
00456         child->elInfo.opp_vertex[ k ] = -2;
00457 
00458       fill( i, elInfo(), child->elInfo );
00459       return ElementInfo< dim >( child );
00460     }
00461 
00462 
00463     template< int dim >
00464     inline bool ElementInfo< dim >::isLeaf () const
00465     {
00466       assert( !(*this) == false );
00467       return isLeaf( el() );
00468     }
00469 
00470 
00471     template< int dim >
00472     inline typename ElementInfo< dim >::Seed ElementInfo< dim >::seed () const
00473     {
00474       assert( !!(*this) );
00475 
00476       int level = 0;
00477       unsigned long path = 0;
00478       for( InstancePtr p = instance_; p->parent() != null(); p = p->parent() )
00479       {
00480         const Element *element = p->elInfo.el;
00481         const Element *father = p->parent()->elInfo.el;
00482         const unsigned long child = static_cast< unsigned long >( father->child[ 1 ] == element );
00483         path = (path << 1) | child;
00484         ++level;
00485       }
00486 
00487       if( level != elInfo().level )
00488         DUNE_THROW( NotImplemented, "Seed for fake elements not implemented." );
00489       
00490       return Seed( macroElement().index, level, path );
00491     };
00492 
00493 
00494     template< int dim >
00495     inline typename ElementInfo< dim >::MeshPointer ElementInfo< dim >::mesh () const
00496     {
00497       return MeshPointer( elInfo().mesh );
00498     }
00499 
00500 
00501     template< int dim >
00502     inline bool ElementInfo< dim >::mightVanish () const
00503     {
00504       return mightVanish( el(), 0 );
00505     }
00506 
00507 
00508     template< int dim >
00509     inline int ElementInfo< dim >::level () const
00510     {
00511       return elInfo().level;
00512     }
00513 
00514 
00515     template< int dim >
00516     inline int ElementInfo< dim >::type () const
00517     {
00518       return 0;
00519     }
00520 
00521 
00522     template<>
00523     inline int ElementInfo< 3 >::type () const
00524     {
00525       return instance_->elInfo.el_type;
00526     }
00527 
00528 
00529     template< int dim >
00530     inline int ElementInfo< dim >::getMark () const
00531     {
00532       return el()->mark;
00533     }
00534 
00535 
00536     template< int dim >
00537     inline void ElementInfo< dim >::setMark ( int refCount ) const
00538     {
00539       assert( isLeaf() );
00540       assert( (refCount >= -128) && (refCount < 127) );
00541       el()->mark = refCount;
00542     }
00543 
00544 
00545 #if DUNE_ALBERTA_VERSION >= 0x300
00546     template< int dim >
00547     inline bool ElementInfo< dim >::hasLeafNeighbor ( const int face ) const
00548     {
00549       assert( !!(*this) );
00550       assert( (face >= 0) && (face < maxNeighbors) );
00551 
00552       assert( (elInfo().fill_flag & FillFlags::boundaryId) != 0 );
00553       const int macroFace = elInfo().macro_wall[ face ];
00554       if( macroFace >= 0 )
00555         return (macroElement().neighbor( macroFace ) != NULL);
00556       else
00557         return true;
00558     }
00559 #endif // DUNE_ALBERTA_VERSION >= 0x300
00560 
00561 #if DUNE_ALBERTA_VERSION < 0x300
00562     template< int dim >
00563     inline bool ElementInfo< dim >::hasLeafNeighbor ( const int face ) const
00564     {
00565       return (neighbor( face ) != NULL);
00566     }
00567 #endif // DUNE_ALBERTA_VERSION < 0x300
00568 
00569     
00570     template< int dim >
00571     inline ElementInfo< dim > ElementInfo< dim >::leafNeighbor ( const int face ) const
00572     {
00573       assert( (face >= 0) && (face < numFaces) );
00574       ElementInfo neighbor;
00575       Library< dimWorld >::leafNeighbor( *this, face, neighbor );
00576       return neighbor;
00577     }
00578 
00579 
00580     template< int dim >
00581     inline int ElementInfo< dim >
00582       ::levelNeighbors ( const int face, ElementInfo (&neighbor)[ maxLevelNeighbors ], int (&faceInNeighbor)[ maxLevelNeighbors ] ) const
00583     {
00584       assert( (face >= 0) && (face < numFaces) );
00585       return Library< dimWorld >::levelNeighbors( *this, face, neighbor, faceInNeighbor );
00586     }
00587 
00588 
00589     template< int dim >
00590     template< int codim >
00591     inline int ElementInfo< dim >::twist ( int subEntity ) const
00592     {
00593       return Twist< dim, dim-codim >::twist( element(), subEntity );
00594     }
00595 
00596 
00597     template< int dim >
00598     inline int ElementInfo< dim >::twistInNeighbor ( const int face ) const
00599     {
00600       assert( neighbor( face ) != NULL );
00601       return Twist< dim, dim-1 >::twist( neighbor( face ), elInfo().opp_vertex[ face ] );
00602     }
00603 
00604 
00605 #if DUNE_ALBERTA_VERSION >= 0x300
00606     template< int dim >
00607     inline bool ElementInfo< dim >::isBoundary ( int face ) const
00608     {
00609       assert( !!(*this) );
00610       assert( (face >= 0) && (face < maxNeighbors) );
00611 
00612       assert( (elInfo().fill_flag & FillFlags::boundaryId) != 0 );
00613       const int macroFace = elInfo().macro_wall[ face ];
00614       if( macroFace >= 0 )
00615         return macroElement().isBoundary( macroFace );
00616       else
00617         return false;
00618     }
00619 #endif // DUNE_ALBERTA_VERSION >= 0x300
00620 
00621 #if DUNE_ALBERTA_VERSION <= 0x200
00622     template< int dim >
00623     inline bool ElementInfo< dim >::isBoundary ( int face ) const
00624     {
00625       assert( !!(*this) );
00626       assert( (face >= 0) && (face < maxNeighbors) );
00627       return (elInfo().neigh[ face ] == 0);
00628     }
00629 #endif // DUNE_ALBERTA_VERSION <= 0x200
00630 
00631 
00632 #if DUNE_ALBERTA_VERSION >= 0x300
00633     template< int dim >
00634     inline int ElementInfo< dim >::boundaryId ( int face ) const
00635     {
00636       assert( !!(*this) );
00637       assert( (face >= 0) && (face < N_WALLS_MAX) );
00638       
00639       assert( (elInfo().fill_flag & FillFlags::boundaryId) != 0 );
00640       const int macroFace = elInfo().macro_wall[ face ];
00641       const int id = macroElement().boundaryId( macroFace );
00642       // this assertion is only allowed, if FILL_BOUND is set
00643       // assert( id == elInfo().wall_bound[ face ] );
00644       return id;
00645     }
00646 #endif // #if DUNE_ALBERTA_VERSION >= 0x300
00647 
00648 
00649 #if DUNE_ALBERTA_VERSION == 0x200
00650     template<>
00651     inline int ElementInfo< 1 >::boundaryId ( int face ) const
00652     {
00653       assert( !!(*this) );
00654       assert( (face >= 0) && (face < N_VERTICES_MAX) );
00655       return elInfo().vertex_bound[ 1-face ];
00656     }
00657 
00658     template<>
00659     inline int ElementInfo< 2 >::boundaryId ( int face ) const
00660     {
00661       assert( !!(*this) );
00662       assert( (face >= 0) && (face < N_EDGES_MAX) );
00663       return elInfo().edge_bound[ face ];
00664     }
00665 
00666     template<>
00667     inline int ElementInfo< 3 >::boundaryId ( int face ) const
00668     {
00669       assert( !!(*this) );
00670       assert( (face >= 0) && (face < N_FACES_MAX) );
00671       return elInfo().face_bound[ face ];
00672     }
00673 #endif // #if DUNE_ALBERTA_VERSION == 0x200
00674 
00675 
00676 #if DUNE_ALBERTA_VERSION >= 0x300
00677     template< int dim >
00678     inline AffineTransformation *
00679     ElementInfo< dim >::transformation ( int face ) const
00680     {
00681       assert( !!(*this) );
00682       assert( (face >= 0) && (face < N_WALLS_MAX) );
00683       
00684       assert( (elInfo().fill_flag & FillFlags::boundaryId) != 0 );
00685       const int macroFace = elInfo().macro_wall[ face ];
00686       return (macroFace < 0 ? NULL : macroElement().wall_trafo[ macroFace ]);
00687     }
00688 #endif // #if DUNE_ALBERTA_VERSION >= 0x300
00689 
00690 #if DUNE_ALBERTA_VERSION <= 0x200
00691     template< int dim >
00692     inline AffineTransformation *
00693     ElementInfo< dim >::transformation ( int face ) const
00694     {
00695       return NULL;
00696     }
00697 #endif // #if DUNE_ALBERTA_VERSION <= 0x200
00698 
00699 
00700 #if DUNE_ALBERTA_VERSION >= 0x300
00701     template< int dim >
00702     inline BasicNodeProjection *
00703     ElementInfo< dim >::boundaryProjection ( int face ) const
00704     {
00705       assert( !!(*this) );
00706       assert( (face >= 0) && (face < N_WALLS_MAX) );
00707       
00708       assert( (elInfo().fill_flag & FillFlags::boundaryId) != 0 );
00709       const int macroFace = elInfo().macro_wall[ face ];
00710       if( macroFace >= 0 )
00711         return static_cast< BasicNodeProjection * >( macroElement().projection[ macroFace+1 ] );
00712       else
00713         return 0;
00714     }
00715 #endif // #if DUNE_ALBERTA_VERSION >= 0x300
00716 
00717 #if DUNE_ALBERTA_VERSION <= 0x200
00718     template< int dim >
00719     inline BasicNodeProjection *
00720     ElementInfo< dim >::boundaryProjection ( int face ) const
00721     {
00722       assert( !!(*this) );
00723       assert( (face >= 0) && (face < maxNeighbors) );
00724       const int idx = (dim == 1 ? 2-face : 1+face);
00725       return static_cast< BasicNodeProjection * >( elInfo().projections[ idx ] );
00726     }
00727 #endif // #if DUNE_ALBERTA_VERSION <= 0x200
00728 
00729 
00730     template< int dim >
00731     inline bool ElementInfo< dim >::hasCoordinates () const
00732     {
00733       return ((elInfo().fill_flag & FillFlags::coords) != 0);
00734     }
00735 
00736     template< int dim >
00737     inline const GlobalVector &ElementInfo< dim >::coordinate ( int vertex ) const
00738     {
00739       assert( hasCoordinates() );
00740       assert( (vertex >= 0) && (vertex < numVertices) );
00741       return elInfo().coord[ vertex ];
00742     }
00743 
00744 
00745     template< int dim >
00746     template< class Functor >
00747     inline void ElementInfo< dim >::hierarchicTraverse ( Functor &functor ) const
00748     {
00749       functor( *this );
00750       if( !isLeaf() )
00751       {
00752         child( 0 ).hierarchicTraverse( functor );
00753         child( 1 ).hierarchicTraverse( functor );
00754       }
00755     }
00756 
00757 
00758     template< int dim >
00759     template< class Functor >
00760     inline void ElementInfo< dim >::leafTraverse ( Functor &functor ) const
00761     {
00762       if( !isLeaf() )
00763       {
00764         child( 0 ).leafTraverse( functor );
00765         child( 1 ).leafTraverse( functor );
00766       }
00767       else
00768         functor( *this );
00769     }
00770 
00771 
00772     template< int dim >
00773     inline const Element *ElementInfo< dim >::element () const
00774     {
00775       return elInfo().el;
00776     }
00777 
00778 
00779     template< int dim >
00780     inline const Element *ElementInfo< dim >::neighbor ( int face ) const
00781     {
00782       assert( (face >= 0) && (face < numFaces) );
00783       assert( (elInfo().fill_flag & FillFlags::neighbor) != 0 );
00784       return elInfo().neigh[ face ];
00785     }
00786 
00787 
00788     template< int dim >
00789     inline Element *ElementInfo< dim >::el () const
00790     {
00791       return elInfo().el;
00792     }
00793 
00794 
00795     template< int dim >
00796     inline ALBERTA EL_INFO &ElementInfo< dim >::elInfo () const
00797     {
00798       return (instance_->elInfo);
00799     }
00800 
00801 
00802     template< int dim >
00803     inline ElementInfo< dim >
00804     ElementInfo< dim >::createFake ( const MeshPointer &mesh,
00805                                      const Element *element, int level, int type )
00806     {
00807       InstancePtr instance = stack().allocate();
00808       instance->parent() = null();
00809       ++(instance->parent()->refCount);
00810 
00811       instance->elInfo.mesh = mesh;
00812       instance->elInfo.macro_el = NULL;
00813       instance->elInfo.el = const_cast< Element * >( element );
00814       instance->elInfo.parent = NULL;
00815       instance->elInfo.fill_flag = FillFlags::nothing;
00816       instance->elInfo.level = level;
00817       instance->elInfo.el_type = type;
00818 
00819       return ElementInfo< dim >( instance );
00820     }
00821 
00822 
00823     template< int dim >
00824     inline ElementInfo< dim >
00825     ElementInfo< dim >::createFake ( const ALBERTA EL_INFO &elInfo )
00826     {
00827       InstancePtr instance = stack().allocate();
00828       instance->parent() = null();
00829       ++(instance->parent()->refCount);
00830 
00831       instance->elInfo = elInfo;
00832       return ElementInfo< dim >( instance );
00833     }
00834 
00835 
00836     template< int dim >
00837     inline bool ElementInfo< dim >::isLeaf ( Element *element )
00838     {
00839       return IS_LEAF_EL( element );
00840     }
00841 
00842 
00843     template< int dim >
00844     inline bool ElementInfo< dim >::mightVanish ( Alberta::Element *element, int depth )
00845     {
00846       if( isLeaf( element ) )
00847         return (element->mark < depth);
00848       else
00849         return (mightVanish( element->child[ 0 ], depth-1 ) && mightVanish( element->child[ 1 ], depth-1 ));
00850     }
00851 
00852 
00853     template< int dim >
00854     inline void ElementInfo< dim >
00855       ::fill ( Mesh *mesh, const ALBERTA MACRO_EL *mel, ALBERTA EL_INFO &elInfo )
00856     {
00857       ALBERTA fill_macro_info( mesh, mel, &elInfo );
00858 
00859 #if DUNE_ALBERTA_VERSION < 0x300
00860       // The 1d grid does not fill in projections, so we do it here
00861       if( (dim == 1) && (elInfo.fill_flag & FILL_PROJECTION) )
00862       {
00863         for( int i = 0; i <= N_VERTICES_1D; ++i )
00864           elInfo.projections[ i ] = mel->projection[ i ];
00865       }
00866 #endif
00867     }
00868 
00869     template< int dim >
00870     inline void ElementInfo< dim >
00871       ::fill ( int ichild, const ALBERTA EL_INFO &parentInfo, ALBERTA EL_INFO &elInfo )
00872     {
00873 #if DUNE_ALBERTA_VERSION >= 0x300
00874       ALBERTA fill_elinfo( ichild, FILL_ANY, &parentInfo, &elInfo );
00875 #else
00876       ALBERTA fill_elinfo( ichild, &parentInfo, &elInfo );
00877 
00878       // The 1d grid does not fill in projections, so we do it here
00879       if( (dim == 1) && (elInfo.fill_flag & FILL_PROJECTION) )
00880       {
00881         elInfo.projections[ 0 ] = parentInfo.projections[ 0 ];
00882         if( ichild == 0 )
00883         {
00884           elInfo.projections[ 1 ] = parentInfo.projections[ 0 ];
00885           elInfo.projections[ 2 ] = parentInfo.projections[ 2 ];
00886         }
00887         else
00888         {
00889           elInfo.projections[ 1 ] = parentInfo.projections[ 1 ];
00890           elInfo.projections[ 2 ] = parentInfo.projections[ 0 ];
00891         }
00892       }
00893 #endif
00894     }
00895 
00896 
00897     template< int dim >
00898     inline void ElementInfo< dim >::addReference () const
00899     {
00900       ++(instance_->refCount);
00901     }
00902 
00903 
00904     template< int dim >
00905     inline void ElementInfo< dim >::removeReference () const
00906     {
00907       // this loop breaks when instance becomes null()
00908       for( InstancePtr instance = instance_; --(instance->refCount) == 0; )
00909       {
00910         const InstancePtr parent = instance->parent();
00911         stack().release( instance );
00912         instance = parent;
00913       }
00914     }
00915 
00916 
00917     template< int dim >
00918     inline typename ElementInfo< dim >::InstancePtr
00919     ElementInfo< dim >::null ()
00920     {
00921       return stack().null();
00922     }
00923 
00924 
00925     template< int dim >
00926     inline typename ElementInfo< dim >::Stack &
00927     ElementInfo< dim >::stack ()
00928     {
00929       static Stack s;
00930       return s;
00931     }
00932 
00933 
00934 
00935     // Implementation of ElementInfo::Stack
00936     // ------------------------------------
00937 
00938     template< int dim >
00939     inline ElementInfo< dim >::Stack::Stack ()
00940     : top_( 0 )
00941     {
00942       null_.elInfo.el = NULL;
00943       null_.refCount = 1;
00944       null_.parent() = 0;
00945     }
00946 
00947 
00948     template< int dim >
00949     inline ElementInfo< dim >::Stack::~Stack ()
00950     {
00951       while( top_ != 0 )
00952       {
00953         InstancePtr p = top_;
00954         top_ = p->parent();
00955         delete p;
00956       }
00957     }
00958 
00959 
00960     template< int dim >
00961     inline typename ElementInfo< dim >::InstancePtr
00962     ElementInfo< dim >::Stack::allocate ()
00963     {
00964       InstancePtr p = top_;
00965       if( p != 0 )
00966         top_ = p->parent();
00967       else
00968         p = new Instance;
00969       p->refCount = 0;
00970       return p;
00971     }
00972 
00973 
00974     template< int dim >
00975     inline void ElementInfo< dim >::Stack::release ( InstancePtr &p )
00976     {
00977       assert( (p != null()) && (p->refCount == 0) );
00978       p->parent() = top_;
00979       top_ = p;
00980     }
00981 
00982 
00983     template< int dim >
00984     inline typename ElementInfo< dim >::InstancePtr
00985     ElementInfo< dim >::Stack::null ()
00986     {
00987       return &null_;
00988     }
00989 
00990   } // namespace Alberta
00991 
00992 } // namespace Dune
00993 
00994 #endif // #if HAVE_ALBERTA
00995 
00996 #endif // #ifndef DUNE_ALBERTA_ELEMENTINFO_HH