dune-grid  2.3.1
meshpointer.hh
Go to the documentation of this file.
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_ALBERTA_MESHPOINTER_HH
4 #define DUNE_ALBERTA_MESHPOINTER_HH
5 
11 #include <limits>
12 #include <string>
13 
18 
19 #if HAVE_ALBERTA
20 
21 namespace Dune
22 {
23 
24  namespace Alberta
25  {
26 
27  // External Forward Declarations
28  // -----------------------------
29 
30  template< int dim >
31  class HierarchyDofNumbering;
32 
33 
34 
35  // MeshPointer
36  // -----------
37 
38  template< int dim >
40  {
42  typedef typename ElementInfo::MacroElement MacroElement;
43  typedef typename ElementInfo::FillFlags FillFlags;
44 
45  class BoundaryProvider;
46 
47  template< int dimWorld >
48  struct Library;
49 
50  public:
51  class MacroIterator;
52 
54  : mesh_( 0 )
55  {}
56 
57  explicit MeshPointer ( Mesh *mesh )
58  : mesh_( mesh )
59  {}
60 
61  operator Mesh * () const
62  {
63  return mesh_;
64  }
65 
66  operator bool () const
67  {
68  return (bool)mesh_;
69  }
70 
71  MacroIterator begin () const
72  {
73  return MacroIterator( *this, false );
74  }
75 
76  MacroIterator end () const
77  {
78  return MacroIterator( *this, true );
79  }
80 
81  int numMacroElements () const;
82  int size ( int codim ) const;
83 
84  // create a mesh from a macrodata structure
85  // params: macroData - macro data structure
86  // returns: number of boundary segments
87  unsigned int create ( const MacroData< dim > &macroData );
88 
89  // create a mesh from a macrodata structure, adding projections
90  // params: macroData - macro data structure
91  // projectionFactory - factory for the projections
92  // returns: number of boundary segments
93  template< class Proj, class Impl >
94  unsigned int create ( const MacroData< dim > &macroData,
95  const ProjectionFactoryInterface< Proj, Impl > &projectionFactory );
96 
97  // create a mesh from a file
98  // params: filename - file name of an Alberta macro triangulation
99  // binary - read binary?
100  // returns: number of boundary segments
101  unsigned int create ( const std::string &filename, bool binary = false );
102 
103  // read back a mesh from a file
104  // params: filename - file name of an Alberta save file
105  // time - variable to receive the time stored in the file
106  // returns: number of boundary segments
107  //
108  // notes: - projections are not preserved
109  // - we assume that projections are added in the same order they
110  // inserted in when the grid was created (otherwise the boundary
111  // indices change)
112  unsigned int read ( const std::string &filename, Real &time );
113 
114  bool write ( const std::string &filename, Real time ) const;
115 
116  void release ();
117 
118  template< class Functor >
119  void hierarchicTraverse ( Functor &functor,
120  typename FillFlags::Flags fillFlags = FillFlags::standard ) const;
121 
122  template< class Functor >
123  void leafTraverse ( Functor &functor,
124  typename FillFlags::Flags fillFlags = FillFlags::standard ) const;
125 
126  bool coarsen ( typename FillFlags::Flags fillFlags = FillFlags::nothing );
127 
128  bool refine ( typename FillFlags::Flags fillFlags = FillFlags::nothing );
129 
130  private:
131  static ALBERTA NODE_PROJECTION *
132  initNodeProjection ( Mesh *mesh, ALBERTA MACRO_EL *macroElement, int n );
133  template< class ProjectionProvider >
134  static ALBERTA NODE_PROJECTION *
135  initNodeProjection ( Mesh *mesh, ALBERTA MACRO_EL *macroElement, int n );
136 
137  Mesh *mesh_;
138  };
139 
140 
141 
142  // MeshPointer::Library
143  // --------------------
144 
145  template< int dim >
146  template< int dimWorld >
147  struct MeshPointer< dim >::Library
148  {
150 
151  static unsigned int boundaryCount;
152  static const void *projectionFactory;
153 
154  static void
155  create ( MeshPointer &ptr, const MacroData< dim > &macroData,
156  ALBERTA NODE_PROJECTION *(*initNodeProjection)( Mesh *, ALBERTA MACRO_EL *, int ) );
157  static void release ( MeshPointer &ptr );
158  };
159 
160 
161 
162  // MeshPointer::MacroIterator
163  // --------------------------
164 
165  template< int dim >
167  {
168  typedef MacroIterator This;
169 
170  friend class MeshPointer< dim >;
171 
172  public:
175 
176  private:
177  explicit MacroIterator ( const MeshPointer &mesh, bool end = false )
178  : mesh_( mesh ),
179  index_( end ? mesh.numMacroElements() : 0 )
180  {}
181 
182  public:
183  bool done () const
184  {
185  return (index_ >= mesh().numMacroElements());
186  }
187 
188  bool equals ( const MacroIterator &other ) const
189  {
190  return (index_ == other.index_);
191  }
192 
193  void increment ()
194  {
195  assert( !done() );
196  ++index_;
197  }
198 
199  const MacroElement &macroElement () const
200  {
201  assert( !done() );
202  return static_cast< const MacroElement & >( mesh().mesh_->macro_els[ index_ ] );
203  }
204 
205  const MeshPointer &mesh () const
206  {
207  return mesh_;
208  }
209 
210  This &operator++ ()
211  {
212  increment();
213  return *this;
214  }
215 
216  ElementInfo operator* () const
217  {
218  return elementInfo();
219  }
220 
221  bool operator== ( const MacroIterator &other ) const
222  {
223  return equals( other );
224  }
225 
226  bool operator!= ( const MacroIterator &other ) const
227  {
228  return !equals( other );
229  }
230 
232  elementInfo ( typename FillFlags::Flags fillFlags = FillFlags::standard ) const
233  {
234  if( done() )
235  return ElementInfo();
236  else
237  return ElementInfo( mesh(), macroElement(), fillFlags );
238  }
239 
240  private:
241  MeshPointer mesh_;
242  int index_;
243  };
244 
245 
246 
247  // Implementation of MeshPointer
248  // -----------------------------
249 
250  template< int dim >
252  {
253  return (mesh_ ? mesh_->n_macro_el : 0);
254  }
255 
256 
257  template<>
258  inline int MeshPointer< 1 >::size( int codim ) const
259  {
260  assert( (codim >= 0) && (codim <= 1) );
261  return (codim == 0 ? mesh_->n_elements : mesh_->n_vertices);
262  }
263 
264  template<>
265  inline int MeshPointer< 2 >::size( int codim ) const
266  {
267  assert( (codim >= 0) && (codim <= 2) );
268  if( codim == 0 )
269  return mesh_->n_elements;
270  else if( codim == 2 )
271  return mesh_->n_vertices;
272  else
273  return mesh_->n_edges;
274  }
275 
276  template<>
277  inline int MeshPointer< 3 >::size( int codim ) const
278  {
279  assert( (codim >= 0) && (codim <= 3) );
280  if( codim == 0 )
281  return mesh_->n_elements;
282  else if( codim == 3 )
283  return mesh_->n_vertices;
284  else if( codim == 1 )
285  return mesh_->n_faces;
286  else
287  return mesh_->n_edges;
288  }
289 
290 
291  template< int dim >
292  inline unsigned int MeshPointer< dim >
293  ::create ( const MacroData< dim > &macroData )
294  {
295  release();
296 
298  Library< dimWorld >::create( *this, macroData, &initNodeProjection );
300  }
301 
302 
303  template< int dim >
304  template< class Proj, class Impl >
305  inline unsigned int MeshPointer< dim >
306  ::create ( const MacroData< dim > &macroData,
307  const ProjectionFactoryInterface< Proj, Impl > &projectionFactory )
308  {
310 
311  release();
312 
314  Library< dimWorld >::projectionFactory = &projectionFactory;
315  Library< dimWorld >::create( *this, macroData, &initNodeProjection< ProjectionFactory > );
318  }
319 
320 
321 
322 
323  template< int dim >
324  inline unsigned int MeshPointer< dim >
325  ::create ( const std::string &filename, bool binary )
326  {
327  MacroData< dim > macroData;
328  macroData.read( filename, binary );
329  const unsigned int boundaryCount = create( macroData );
330  macroData.release();
331  return boundaryCount;
332  }
333 
334 
335  template< int dim >
336  inline unsigned int MeshPointer< dim >::read ( const std::string &filename, Real &time )
337  {
338  release();
339 
341 #if DUNE_ALBERTA_VERSION >= 0x300
342  mesh_ = ALBERTA read_mesh_xdr( filename.c_str(), &time, NULL, NULL );
343 #else
344  mesh_ = ALBERTA read_mesh_xdr( filename.c_str(), &time, NULL );
345 #endif
347  }
348 
349 
350  template< int dim >
351  inline bool MeshPointer< dim >::write ( const std::string &filename, Real time ) const
352  {
353  int success = ALBERTA write_mesh_xdr( mesh_, filename.c_str(), time );
354  return (success == 0);
355  }
356 
357 
358  template< int dim >
360  {
362  }
363 
364 
365  template< int dim >
366  template< class Functor >
367  inline void MeshPointer< dim >
368  ::hierarchicTraverse ( Functor &functor,
369  typename FillFlags::Flags fillFlags ) const
370  {
371  const MacroIterator eit = end();
372  for( MacroIterator it = begin(); it != eit; ++it )
373  {
374  const ElementInfo info = it.elementInfo( fillFlags );
375  info.hierarchicTraverse( functor );
376  }
377  }
378 
379 
380  template< int dim >
381  template< class Functor >
382  inline void MeshPointer< dim >
383  ::leafTraverse ( Functor &functor,
384  typename FillFlags::Flags fillFlags ) const
385  {
386  const MacroIterator eit = end();
387  for( MacroIterator it = begin(); it != eit; ++it )
388  {
389  const ElementInfo info = it.elementInfo( fillFlags );
390  info.leafTraverse( functor );
391  }
392  }
393 
394 
395 #if DUNE_ALBERTA_VERSION >= 0x300
396  template< int dim >
397  inline bool MeshPointer< dim >::coarsen ( typename FillFlags::Flags fillFlags )
398  {
399  const bool coarsened = (ALBERTA coarsen( mesh_, fillFlags ) == meshCoarsened);
400  if( coarsened )
401  ALBERTA dof_compress( mesh_ );
402  return coarsened;
403  }
404 #endif // #if DUNE_ALBERTA_VERSION >= 0x300
405 
406 #if DUNE_ALBERTA_VERSION <= 0x200
407  template< int dim >
408  inline bool MeshPointer< dim >::coarsen ( typename FillFlags::Flags fillFlags )
409  {
410  assert( fillFlags == FillFlags::nothing );
411  const bool coarsened = (ALBERTA coarsen( mesh_ ) == meshCoarsened);
412  if( coarsened )
413  ALBERTA dof_compress( mesh_ );
414  return coarsened;
415  }
416 #endif // #if DUNE_ALBERTA_VERSION <= 0x200
417 
418 
419 #if DUNE_ALBERTA_VERSION >= 0x300
420  template< int dim >
421  inline bool MeshPointer< dim >::refine ( typename FillFlags::Flags fillFlags )
422  {
423  return (ALBERTA refine( mesh_, fillFlags ) == meshRefined);
424  }
425 #endif // #if DUNE_ALBERTA_VERSION >= 0x300
426 
427 #if DUNE_ALBERTA_VERSION <= 0x200
428  template< int dim >
429  inline bool MeshPointer< dim >::refine ( typename FillFlags::Flags fillFlags )
430  {
431  assert( fillFlags == FillFlags::nothing );
432  return (ALBERTA refine( mesh_ ) == meshRefined);
433  }
434 #endif // #if DUNE_ALBERTA_VERSION <= 0x200
435 
436 
437  template< int dim >
438  inline ALBERTA NODE_PROJECTION *
439  MeshPointer< dim >::initNodeProjection ( Mesh *mesh, ALBERTA MACRO_EL *macroEl, int n )
440  {
441  const MacroElement &macroElement = static_cast< const MacroElement & >( *macroEl );
442  if( (n > 0) && macroElement.isBoundary( n-1 ) )
443  return new BasicNodeProjection( Library< dimWorld >::boundaryCount++ );
444  else
445  return 0;
446  }
447 
448 
449  template< int dim >
450  template< class ProjectionFactory >
451  inline ALBERTA NODE_PROJECTION *
452  MeshPointer< dim >::initNodeProjection ( Mesh *mesh, ALBERTA MACRO_EL *macroEl, int n )
453  {
454  typedef typename ProjectionFactory::Projection Projection;
455 
456  const MacroElement &macroElement = static_cast< const MacroElement & >( *macroEl );
457 
458  MeshPointer< dim > meshPointer( mesh );
459  ElementInfo elementInfo( meshPointer, macroElement, FillFlags::standard );
460  const ProjectionFactory &projectionFactory = *static_cast< const ProjectionFactory * >( Library< dimWorld >::projectionFactory );
461  if( (n > 0) && macroElement.isBoundary( n-1 ) )
462  {
463  const unsigned int boundaryIndex = Library< dimWorld >::boundaryCount++;
464  if( projectionFactory.hasProjection( elementInfo, n-1 ) )
465  {
466  Projection projection = projectionFactory.projection( elementInfo, n-1 );
467  return new NodeProjection< dim, Projection >( boundaryIndex, projection );
468  }
469  else
470  return new BasicNodeProjection( boundaryIndex );
471  }
472  else if( (dim < dimWorld) && (n == 0) )
473  {
474  const unsigned int boundaryIndex = std::numeric_limits< unsigned int >::max();
475  if( projectionFactory.hasProjection( elementInfo ) )
476  {
477  Projection projection = projectionFactory.projection( elementInfo );
478  return new NodeProjection< dim, Projection >( boundaryIndex, projection );
479  }
480  else
481  return 0;
482  }
483  else
484  return 0;
485  }
486 
487  } // namespace Alberta
488 
489 } // namespace Dune
490 
491 #endif // #if HAVE_ALBERTA
492 
493 #endif // #ifndef DUNE_ALBERTA_MESHPOINTER_HH