1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_GEOGRID_CACHEDCOORDFUNCTION_HH
4 #define DUNE_GEOGRID_CACHEDCOORDFUNCTION_HH
5 
6 #include <cassert>
7 #include <memory>
8 
9 #include <dune/common/typetraits.hh>
10 
11 #include <dune/grid/common/gridenums.hh>
12 
13 #include <dune/grid/geometrygrid/capabilities.hh>
14 #include <dune/grid/geometrygrid/coordfunctioncaller.hh>
15 #include <dune/grid/utility/persistentcontainer.hh>
16 
17 namespace Dune
18 {
19 
20   // Internal Forward Declarations
21   // -----------------------------
22 
23   template< class HostGrid, class CoordFunction >
24   class CachedCoordFunction;
25 
26 
27 
28   // GeoGrid::CoordCache
29   // -------------------
30 
31   namespace GeoGrid
32   {
33 
34     template< class HostGrid, class Coordinate >
35     class CoordCache
36     {
37       typedef CoordCache< HostGrid, Coordinate > This;
38 
39       static const unsigned int dimension = HostGrid::dimension;
40 
41       typedef typename HostGrid::template Codim< dimension >::Entity Vertex;
42 
43       typedef PersistentContainer< HostGrid, Coordinate > DataCache;
44 
45     public:
CoordCache(const HostGrid & hostGrid)46       explicit CoordCache ( const HostGrid &hostGrid )
47         : data_( hostGrid, dimension )
48       {}
49 
50       template< class Entity >
operator ()(const Entity & entity,unsigned int corner) const51       const Coordinate &operator() ( const Entity &entity, unsigned int corner ) const
52       {
53         return data_( entity, corner );
54       }
55 
operator ()(const Vertex & vertex,unsigned int corner) const56       const Coordinate &operator() ( const Vertex &vertex, unsigned int corner ) const
57       {
58         assert( corner == 0 );
59         return data_[ vertex ];
60       }
61 
62       template< class Entity >
operator ()(const Entity & entity,unsigned int corner)63       Coordinate &operator() ( const Entity &entity, unsigned int corner )
64       {
65         return data_( entity,corner) ;
66       }
67 
operator ()(const Vertex & vertex,unsigned int corner)68       Coordinate &operator() ( const Vertex &vertex, unsigned int corner )
69       {
70         assert( corner == 0 );
71         return data_[ vertex ];
72       }
73 
adapt()74       void adapt ()
75       {
76         data_.resize();
77         data_.shrinkToFit();
78       }
79 
80     private:
81       CoordCache ( const This & );
82       This &operator= ( const This & );
83 
84       DataCache data_;
85     };
86 
87   } // namespace GeoGrid
88 
89 
90 
91   // CachedCoordFunction
92   // -------------------
93 
94   template< class HostGrid, class CoordFunction >
95   class CachedCoordFunction
96     : public DiscreteCoordFunction< typename CoordFunction::ctype, CoordFunction::dimRange, CachedCoordFunction< HostGrid, CoordFunction > >
97   {
98     typedef CachedCoordFunction< HostGrid, CoordFunction > This;
99     typedef DiscreteCoordFunction< typename CoordFunction::ctype, CoordFunction::dimRange, This > Base;
100 
101   public:
102     typedef typename Base::ctype ctype;
103 
104     typedef typename Base::RangeVector RangeVector;
105 
106   private:
107     typedef GeoGrid::CoordCache< HostGrid, RangeVector > Cache;
108 
109   public:
110     explicit
CachedCoordFunction(const HostGrid & hostGrid,const CoordFunction & coordFunction=CoordFunction ())111     CachedCoordFunction ( const HostGrid &hostGrid,
112                           const CoordFunction &coordFunction = CoordFunction() )
113       : hostGrid_( hostGrid ),
114         coordFunction_( coordFunction ),
115         cache_( hostGrid )
116     {
117       buildCache();
118     }
119 
adapt()120     void adapt ()
121     {
122       cache_.adapt();
123       buildCache();
124     }
125 
126     void buildCache ();
127 
128     template< class HostEntity >
129     void insertEntity ( const HostEntity &hostEntity );
130 
131     template< class HostEntity >
evaluate(const HostEntity & hostEntity,unsigned int corner,RangeVector & y) const132     void evaluate ( const HostEntity &hostEntity, unsigned int corner, RangeVector &y ) const
133     {
134       y = cache_( hostEntity, corner );
135 #ifndef NDEBUG
136       typedef GeoGrid::CoordFunctionCaller< HostEntity, typename CoordFunction::Interface >
137       CoordFunctionCaller;
138 
139       RangeVector z;
140       CoordFunctionCaller coordFunctionCaller( hostEntity, coordFunction_ );
141       coordFunctionCaller.evaluate( corner, z );
142       assert( ((y - z).two_norm() < 1e-6) );
143 #endif
144     }
145 
146   private:
147     const HostGrid &hostGrid_;
148     const CoordFunction &coordFunction_;
149     Cache cache_;
150   };
151 
152 
153 
154   // Implementation of CachedCoordFunction
155   // -------------------------------------
156 
157   template< class HostGrid, class CoordFunction >
buildCache()158   inline void CachedCoordFunction< HostGrid, CoordFunction >::buildCache ()
159   {
160     typedef typename HostGrid::template Codim< 0 >::Entity Element;
161     typedef typename HostGrid::LevelGridView MacroView;
162     typedef typename HostGrid::HierarchicIterator HierarchicIterator;
163 
164     typedef typename MacroView::template Codim< 0 >::template Partition< All_Partition >::Iterator MacroIterator;
165 
166     const MacroView macroView = hostGrid_.levelGridView( 0 );
167     const int maxLevel = hostGrid_.maxLevel();
168 
169     const MacroIterator mend = macroView.template end< 0, All_Partition >();
170     for( MacroIterator mit = macroView.template begin< 0, All_Partition >(); mit != mend; ++mit )
171     {
172       const Element &macroElement = *mit;
173       insertEntity( macroElement );
174 
175       const HierarchicIterator hend = macroElement.hend( maxLevel );
176       for( HierarchicIterator hit = macroElement.hbegin( maxLevel ); hit != hend; ++hit )
177         insertEntity( *hit );
178     }
179   }
180 
181 
182   template< class HostGrid, class CoordFunction >
183   template< class HostEntity >
184   inline void CachedCoordFunction< HostGrid, CoordFunction >
insertEntity(const HostEntity & hostEntity)185     ::insertEntity ( const HostEntity &hostEntity )
186   {
187     typedef GeoGrid::CoordFunctionCaller< HostEntity, typename CoordFunction::Interface >
188     CoordFunctionCaller;
189 
190     CoordFunctionCaller coordFunctionCaller( hostEntity, coordFunction_ );
191     auto refElement = referenceElement< ctype, HostEntity::dimension >( hostEntity.type() );
192 
193     const unsigned int numCorners = refElement.size( HostEntity::dimension );
194     for( unsigned int i = 0; i < numCorners; ++i )
195       coordFunctionCaller.evaluate( i, cache_( hostEntity, i ) );
196   }
197 
198 } // namespace Dune
199 
200 #endif // #ifndef DUNE_GEOGRID_CACHEDCOORDFUNCTION_HH
201