dune-grid  2.2.0
albertagrid/refinement.hh
Go to the documentation of this file.
00001 #ifndef DUNE_ALBERTA_REFINEMENT_HH
00002 #define DUNE_ALBERTA_REFINEMENT_HH
00003 
00010 #include <cassert>
00011 
00012 #include <dune/grid/albertagrid/misc.hh>
00013 #include <dune/grid/albertagrid/elementinfo.hh>
00014 
00015 #if HAVE_ALBERTA
00016 
00017 namespace Dune
00018 {
00019 
00020   namespace Alberta
00021   {
00022 
00023     // Internal Forward Declarations
00024     // -----------------------------
00025 
00026     template< int dim, int codim >
00027     struct ForEachInteriorSubChild;
00028 
00029 
00030 
00031     // Patch
00032     // -----
00033 
00034     template< int dim >
00035     class Patch
00036     {
00037       typedef Patch< dim > This;
00038 
00039       dune_static_assert( ((dim >= 1) && (dim <= 3)),
00040                           "Alberta supports only dimensions 1, 2, 3" );
00041 
00042     public:
00043       static const int dimension = dim;
00044 
00045       typedef Alberta::ElementInfo< dimension > ElementInfo;
00046 
00047       typedef ALBERTA RC_LIST_EL ElementList;
00048 
00049     private:
00050       ElementList *list_;
00051       int count_;
00052 
00053     public:
00054       Patch ( ElementList *list, int count )
00055       : list_( list ),
00056         count_( count )
00057       {
00058         assert( count > 0 );
00059       }
00060 
00061       Element *operator[] ( int i ) const;
00062 
00063       int count () const
00064       {
00065         return count_;
00066       }
00067 
00068       template< class LevelProvider >
00069       ElementInfo elementInfo ( int i, const LevelProvider &levelProvider ) const;
00070 
00071       int elementType ( int i ) const;
00072       bool hasNeighbor ( int i, int neighbor ) const;
00073       int neighborIndex ( int i, int neighbor ) const;
00074 
00075       template< class Functor >
00076       void forEach ( Functor &functor ) const
00077       {
00078         for( int i = 0; i < count(); ++i )
00079           functor( (*this)[ i ] );
00080       }
00081 
00082       template< int codim, class Functor >
00083       void forEachInteriorSubChild ( Functor &functor ) const
00084       {
00085         ForEachInteriorSubChild< dim, codim >::apply( functor, *this );
00086       }
00087     };
00088 
00089 
00090     template< int dim >
00091     inline Element *Patch< dim >::operator[] ( int i ) const
00092     {
00093       assert( (i >= 0) && (i < count()) );
00094       return list_[ i ].el_info.el;
00095     }
00096 
00097 
00098     template< int dim >
00099     template< class LevelProvider >
00100     inline typename Patch< dim >::ElementInfo
00101     Patch< dim >::elementInfo ( int i, const LevelProvider &levelProvider ) const
00102     {
00103       assert( (i >= 0) && (i < count()) );
00104       return ElementInfo::createFake( list_[ i ].el_info );
00105     }
00106 
00107     template<>
00108     template< class LevelProvider >
00109     inline typename Patch< 2 >::ElementInfo
00110     Patch< 2 >::elementInfo ( int i, const LevelProvider &levelProvider ) const
00111     {
00112       assert( (i >= 0) && (i < count()) );
00113       const MeshPointer< 2 > &mesh = levelProvider.mesh();
00114       const Element *element = (*this)[ i ];
00115       const int level = levelProvider( element );
00116       return ElementInfo::createFake( mesh, element, level );
00117     }
00118 
00119 
00120     template< int dim >
00121     inline int Patch< dim >::elementType ( int i ) const
00122     {
00123       assert( (i >= 0) && (i < count()) );
00124       return list_[ i ].el_info.el_type;
00125     }
00126 
00127 
00128     template< int dim >
00129     inline bool Patch< dim >::hasNeighbor ( int i, int neighbor ) const
00130     {
00131       return (list_[ i ].neigh[ neighbor ] != NULL);
00132     }
00133 
00134     template< int dim >
00135     inline int Patch< dim >::neighborIndex ( int i, int neighbor ) const
00136     {
00137       assert( hasNeighbor( i, neighbor ) );
00138       return (list_[ i ].neigh[ neighbor ]->no);
00139     }
00140 
00141 
00142 
00143     // ForEachInteriorSubChild
00144     // -----------------------
00145 
00146     template< int dim >
00147     struct ForEachInteriorSubChild< dim, 0 >
00148     {
00149       template< class Functor >
00150       static void apply ( Functor &functor, const Patch< dim > &patch )
00151       {
00152         for( int i = 0; i < patch.count(); ++i )
00153         {
00154           Element *const father = patch[ i ];
00155           functor( father->child[ 0 ], 0 );
00156           functor( father->child[ 1 ], 0 );
00157         }
00158       }
00159     };
00160 
00161     template< int dim >
00162     struct ForEachInteriorSubChild< dim, dim >
00163     {
00164       template< class Functor >
00165       static void apply ( Functor &functor, const Patch< dim > &patch )
00166       {
00167         functor( patch[ 0 ]->child[ 0 ], dim );
00168       }
00169     };
00170 
00171     template<>
00172     struct ForEachInteriorSubChild< 2, 1 >
00173     {
00174       template< class Functor >
00175       static void apply ( Functor &functor, const Patch< 2 > &patch )
00176       {
00177         // see alberta/src/2d/lagrange_2_2d.c for details
00178         Element *const firstFather = patch[ 0 ];
00179 
00180         Element *const firstChild = firstFather->child[ 0 ];
00181         functor( firstChild, 0 );
00182         functor( firstChild, 1 );
00183 
00184         functor( firstFather->child[ 1 ], 1 );
00185 
00186         if( patch.count() > 1 )
00187         {
00188           Element *const father = patch[ 1 ];
00189           functor( father->child[ 0 ], 1 );
00190         }
00191       }
00192     };
00193 
00194     template<>
00195     struct ForEachInteriorSubChild< 3, 1 >
00196     {
00197       template< class Functor >
00198       static void apply ( Functor &functor, const Patch< 3 > &patch )
00199       {
00200         // see alberta/src/3d/lagrange_3_3d.c for details
00201         Element *const firstFather = patch[ 0 ];
00202 
00203         Element *const firstChild = firstFather->child[ 0 ];
00204         functor( firstChild, 0 );
00205         functor( firstChild, 1 );
00206         functor( firstChild, 2 );
00207 
00208         Element *const secondChild = firstFather->child[ 1 ];
00209         functor( secondChild, 1 );
00210         functor( secondChild, 2 );
00211 
00212         for( int i = 1; i < patch.count(); ++i )
00213         {
00214           Element *const father = patch[ i ];
00215           const int type = patch.elementType( i );
00216 
00217           int lr_set = 0;
00218           if( patch.hasNeighbor( i, 0 ) && (patch.neighborIndex( i, 0 ) < i) )
00219             lr_set |= 1;
00220           if( patch.hasNeighbor( i, 1 ) && (patch.neighborIndex( i, 1 ) < i) )
00221             lr_set |= 2;
00222           assert( lr_set != 0 );
00223 
00224           functor( father->child[ 0 ], 0 );
00225           switch( lr_set )
00226           {
00227           case 1:
00228             functor( father->child[ 0 ], 2 );
00229             functor( father->child[ 1 ], (type == 0 ? 1 : 2) );
00230             break;
00231 
00232           case 2:
00233             functor( father->child[ 0 ], 1 );
00234             functor( father->child[ 1 ], (type == 0 ? 2 : 1) );
00235             break;
00236           }
00237         }
00238       }
00239     };
00240 
00241     template<>
00242     struct ForEachInteriorSubChild< 3, 2 >
00243     {
00244       template< class Functor >
00245       static void apply ( Functor &functor, const Patch< 3 > &patch )
00246       {
00247         // see alberta/src/3d/lagrange_2_3d.c for details
00248         Element *const firstFather = patch[ 0 ];
00249 
00250         Element *const firstChild = firstFather->child[ 0 ];
00251         functor( firstChild, 2 );
00252         functor( firstChild, 4 );
00253         functor( firstChild, 5 );
00254 
00255         functor( firstFather->child[ 1 ], 2 );
00256 
00257         for( int i = 1; i < patch.count(); ++i )
00258         {
00259           Element *const father = patch[ i ];
00260 
00261           int lr_set = 0;
00262           if( patch.hasNeighbor( i, 0 ) && (patch.neighborIndex( i, 0 ) < i) )
00263             lr_set = 1;
00264           if( patch.hasNeighbor( i, 1 ) && (patch.neighborIndex( i, 1 ) < i) )
00265             lr_set += 2;
00266           assert( lr_set != 0 );
00267 
00268           switch( lr_set )
00269           {
00270           case 1:
00271             functor( father->child[ 0 ], 4 );
00272             break;
00273 
00274           case 2:
00275             functor( father->child[ 0 ], 5 );
00276             break;
00277           }
00278         }
00279       }
00280     };
00281 
00282 
00283 
00284     // GeometryInFather
00285     // ----------------
00286 
00287     template< int dim >
00288     struct GeometryInFather;
00289 
00290     template<>
00291     struct GeometryInFather< 1 >
00292     {
00293       static const int dim = 1;
00294 
00295       typedef Real LocalVector[ dim ];
00296 
00297       static const LocalVector &
00298       coordinate ( int child, int orientation, int i )
00299       {
00300         static const Real coords[ 2 ][ dim+1 ][ dim ]
00301           = { { {0.0}, {0.5} }, { {0.5}, {1.0} } };
00302         assert( (i >= 0) && (i <= dim) );
00303         return coords[ child ][ i ];
00304       }
00305     };
00306     
00307     template<>
00308     struct GeometryInFather< 2 >
00309     {
00310       static const int dim = 2;
00311 
00312       typedef Real LocalVector[ dim ];
00313 
00314       static const LocalVector &
00315       coordinate ( int child, int orientation, int i )
00316       {
00317         static const Real coords[ 2 ][ dim+1 ][ dim ]
00318           = { { {0.0, 1.0}, {0.0, 0.0}, {0.5, 0.0} },
00319               { {1.0, 0.0}, {0.0, 1.0}, {0.5, 0.0} } };
00320         assert( (i >= 0) && (i <= dim) );
00321         return coords[ child ][ i ];
00322       }
00323     };
00324 
00325     template<>
00326     struct GeometryInFather< 3 >
00327     {
00328       static const int dim = 3;
00329 
00330       typedef Real LocalVector[ dim ];
00331 
00332       static const LocalVector &
00333       coordinate ( int child, int orientation, int i )
00334       {
00335         static const Real coords[ 2 ][ dim+1 ][ dim ]
00336           = { { {0.0, 0.0, 0.0}, {0.0, 1.0, 0.0}, {0.0, 0.0, 1.0}, {0.5, 0.0, 0.0} },
00337               { {1.0, 0.0, 0.0}, {0.0, 1.0, 0.0}, {0.0, 0.0, 1.0}, {0.5, 0.0, 0.0} } };
00338         static const int flip[ 2 ][ 2 ][ dim+1 ]
00339           = { { {0, 1, 2, 3}, {0, 1, 2, 3} }, { {0, 2, 1, 3}, {0, 1, 2, 3} } };
00340         assert( (i >= 0) && (i <= dim) );
00341         i = flip[ child ][ orientation ][ i ];
00342         return coords[ child ][ i ];
00343       }
00344     };
00345 
00346   }
00347 
00348 }
00349 
00350 #endif // #if HAVE_ALBERTA
00351 
00352 #endif