DUNE-FEM (unstable)

treeiterator.hh
1// SPDX-FileCopyrightText: Copyright © DUNE Project contributors, see file LICENSE.md in module root
2// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4// vi: set et ts=4 sw=2 sts=2:
5
6#ifndef DUNE_ALBERTA_TREEITERATOR_HH
7#define DUNE_ALBERTA_TREEITERATOR_HH
8
9#include <utility>
10
11#include <dune/common/hybridutilities.hh>
13
16
17#if HAVE_ALBERTA
18
19namespace Dune
20{
21
22 // AlbertaMarkerVector
23 // -------------------
24
33 template< int dim, int dimworld >
35 {
37
39
40 //friend class AlbertaGrid< dim, dimworld >;
41
42 static const int dimension = Grid::dimension;
43
44 typedef Alberta::HierarchyDofNumbering< dimension > DofNumbering;
45 typedef Alberta::ElementInfo< dimension > ElementInfo;
46
47 template< bool >
48 struct NoMarkSubEntities;
49 template< bool >
50 struct MarkSubEntities;
51
52 public:
54 explicit AlbertaMarkerVector ( const DofNumbering &dofNumbering )
55 : dofNumbering_( dofNumbering )
56 {
57 for( int codim = 0; codim <= dimension; ++codim )
58 marker_[ codim ] = 0;
59 }
60
61 AlbertaMarkerVector ( const This &other )
62 : dofNumbering_( other.dofNumbering_ )
63 {
64 for( int codim = 0; codim <= dimension; ++codim )
65 marker_[ codim ] = 0;
66 }
67
68 ~AlbertaMarkerVector ()
69 {
70 clear();
71 }
72
73 private:
74 This &operator= ( const This & );
75
76 public:
78 template< int codim >
79 bool subEntityOnElement ( const ElementInfo &elementInfo, int subEntity ) const;
80
81 template< int firstCodim, class Iterator >
82 void markSubEntities ( const Iterator &begin, const Iterator &end );
83
84 void clear ()
85 {
86 for( int codim = 0; codim <= dimension; ++codim )
87 {
88 if( marker_[ codim ] != 0 )
89 delete[] marker_[ codim ];
90 marker_[ codim ] = 0;
91 }
92 }
93
95 bool up2Date () const
96 {
97 return (marker_[ dimension ] != 0);
98 }
99
101 void print ( std::ostream &out = std::cout ) const;
102
103 private:
104 const DofNumbering &dofNumbering_;
105 int *marker_[ dimension+1 ];
106 };
107
108
109
110 // AlbertaMarkerVector::NoMarkSubEntities
111 // --------------------------------------
112
113 template< int dim, int dimworld >
114 template< bool >
115 struct AlbertaMarkerVector< dim, dimworld >::NoMarkSubEntities
116 {
117 template< int firstCodim, class Iterator >
118 static void mark ( [[maybe_unused]] const DofNumbering & dofNumbering,
119 [[maybe_unused]] int *(&marker)[ dimension + 1 ],
120 [[maybe_unused]] const Iterator &begin,
121 [[maybe_unused]] const Iterator &end )
122 {}
123 };
124
125
126
127 // AlbertaMarkerVector::MarkSubEntities
128 // ------------------------------------
129
130 template< int dim, int dimworld >
131 template< bool >
132 struct AlbertaMarkerVector< dim, dimworld >::MarkSubEntities
133 {
134 template< int codim >
135 struct Codim
136 {
137 static const int numSubEntities = Alberta::NumSubEntities< dimension, codim >::value;
138
139 typedef Alberta::ElementInfo< dimension > ElementInfo;
140
141 static void apply ( const DofNumbering &dofNumbering,
142 int *(&marker)[ dimension + 1 ],
143 const ElementInfo &elementInfo )
144 {
145 int *array = marker[ codim ];
146
147 const int index = dofNumbering( elementInfo, 0, 0 );
148 for( int i = 0; i < numSubEntities; ++i )
149 {
150 int &mark = array[ dofNumbering( elementInfo, codim, i ) ];
151 mark = std::max( index, mark );
152 }
153 }
154 };
155
156 template< int firstCodim, class Iterator >
157 static void mark ( const DofNumbering &dofNumbering, int *(&marker)[ dimension + 1 ],
158 const Iterator &begin, const Iterator &end )
159 {
160 for( int codim = firstCodim; codim <= dimension; ++codim )
161 {
162 const int size = dofNumbering.size( codim );
163 marker[ codim ] = new int[ size ];
164
165 int *array = marker[ codim ];
166 for( int i = 0; i < size; ++i )
167 array[ i ] = -1;
168 }
169
170 for( Iterator it = begin; it != end; ++it )
171 {
172 const ElementInfo &elementInfo = it->impl().elementInfo();
173 Hybrid::forEach( std::make_index_sequence< dimension+1-firstCodim >{},
174 [ & ]( auto i ){ Codim< i+firstCodim >::apply( dofNumbering, marker, elementInfo ); } );
175 }
176 }
177 };
178
179
180
181 // AlbertaGridTreeIterator
182 // -----------------------
183
187 template< int codim, class GridImp, bool leafIterator >
189 {
191
192 public:
193 static const int dimension = GridImp::dimension;
194 static const int codimension = codim;
195 static const int dimensionworld = GridImp::dimensionworld;
196
197 private:
198 friend class AlbertaGrid< dimension, dimensionworld >;
199
200 static const int numSubEntities
201 = Alberta::NumSubEntities< dimension, codimension >::value;
202
203 public:
204 typedef Alberta::MeshPointer< dimension > MeshPointer;
205 typedef typename MeshPointer::MacroIterator MacroIterator;
206
207 typedef typename GridImp::template Codim< codim >::Entity Entity;
209 typedef typename EntityObject::ImplementationType EntityImp;
210 typedef typename EntityImp::ElementInfo ElementInfo;
211
213
215
218
220 This &operator= ( const This &other );
221
223 AlbertaGridTreeIterator ( const GridImp &grid, int travLevel );
224
227 const MarkerVector *marker,
228 int travLevel );
229
231 bool equals ( const This &other ) const
232 {
233 return entity_.impl().equals( other.entity_.impl() );
234 }
235
237 Entity &dereference () const
238 {
239 return entity_;
240 }
241
243 int level () const
244 {
245 return entity_.impl().level();
246 }
247
249 void increment();
250
251 protected:
253 const GridImp &grid () const
254 {
255 return entity_.impl().grid();
256 }
257
258 private:
259 void nextElement ( ElementInfo &elementInfo );
260 void nextElementStop (ElementInfo &elementInfo );
261 bool stopAtElement ( const ElementInfo &elementInfo ) const;
262
263 void goNext ( ElementInfo &elementInfo );
264 void goNext ( const std::integral_constant< int, 0 > cdVariable,
265 ElementInfo &elementInfo );
266 void goNext ( const std::integral_constant< int, 1 > cdVariable,
267 ElementInfo &elementInfo );
268 template< int cd >
269 void goNext ( const std::integral_constant< int, cd > cdVariable,
270 ElementInfo &elementInfo );
271
272 mutable Entity entity_;
273
275 int level_;
276
278 int subEntity_;
279
280 MacroIterator macroIterator_;
281
282 // knows on which element a point,edge,face is viewed
283 const MarkerVector *marker_;
284 };
285
286
287
288 // Implementation of AlbertaMarkerVector
289 // -------------------------------------
290
291 template< int dim, int dimworld >
292 template< int codim >
294 ::subEntityOnElement ( const ElementInfo &elementInfo, int subEntity ) const
295 {
296 assert( marker_[ codim ] != 0 );
297
298 const int subIndex = dofNumbering_( elementInfo, codim, subEntity );
299 const int markIndex = marker_[ codim ][ subIndex ];
300 assert( (markIndex >= 0) );
301
302 const int index = dofNumbering_( elementInfo, 0, 0 );
303 return (markIndex == index);
304 }
305
306
307 template< int dim, int dimworld >
308 template< int firstCodim, class Iterator >
310 ::markSubEntities ( const Iterator &begin, const Iterator &end )
311 {
312 clear();
313 std::conditional< (firstCodim <= dimension), MarkSubEntities<true>, NoMarkSubEntities<false> >::type
314 ::template mark< firstCodim, Iterator >( dofNumbering_, marker_, begin, end );
315 }
316
317
318 template< int dim, int dimworld >
319 inline void AlbertaMarkerVector< dim, dimworld >::print ( std::ostream &out ) const
320 {
321 for( int codim = 1; codim <= dimension; ++codim )
322 {
323 int *marker = marker_[ codim ];
324 if( marker != 0 )
325 {
326 const int size = dofNumbering_.size( codim );
327 out << std::endl;
328 out << "Codimension " << codim << " (" << size << " entries)" << std::endl;
329 for( int i = 0; i < size; ++i )
330 out << "subentity " << i << " visited on Element " << marker[ i ] << std::endl;
331 }
332 }
333 }
334
335
336
337 // Implementation of AlbertaGridTreeIterator
338 // -----------------------------------------
339
340 template< int codim, class GridImp, bool leafIterator >
343 : entity_(),
344 level_( -1 ),
345 subEntity_( -1 ),
346 macroIterator_(),
347 marker_( NULL )
348 {}
349
350 template< int codim, class GridImp, bool leafIterator >
351 inline AlbertaGridTreeIterator< codim, GridImp, leafIterator >
352 ::AlbertaGridTreeIterator ( const GridImp &grid,
353 const MarkerVector *marker,
354 int travLevel )
355 : entity_( EntityImp( grid ) ),
356 level_( travLevel ),
357 subEntity_( (codim == 0 ? 0 : -1) ),
358 macroIterator_( grid.meshPointer().begin() ),
359 marker_( marker )
360 {
361 ElementInfo elementInfo = *macroIterator_;
362 nextElementStop( elementInfo );
363 if( codim > 0 )
364 goNext( elementInfo );
365 // it is ok to set the invalid ElementInfo
366 entity_.impl().setElement( elementInfo, subEntity_ );
367 }
368
369
370 // Make LevelIterator with point to element from previous iterations
371 template< int codim, class GridImp, bool leafIterator >
373 ::AlbertaGridTreeIterator ( const GridImp &grid,
374 int travLevel )
375 : entity_( EntityImp( grid ) ),
376 level_( travLevel ),
377 subEntity_( -1 ),
378 macroIterator_( grid.meshPointer().end() ),
379 marker_( 0 )
380 {}
381
382
383 // Make LevelIterator with point to element from previous iterations
384 template< int codim, class GridImp, bool leafIterator >
387 : entity_( other.entity_ ),
388 level_( other.level_ ),
389 subEntity_( other.subEntity_ ),
390 macroIterator_( other.macroIterator_ ),
391 marker_( other.marker_ )
392 {}
393
394
395 // Make LevelIterator with point to element from previous iterations
396 template< int codim, class GridImp, bool leafIterator >
399 {
400 entity_ = other.entity_;
401 level_ = other.level_;
402 subEntity_ = other.subEntity_;
403 macroIterator_ = other.macroIterator_;
404 marker_ = other.marker_;
405
406 return *this;
407 }
408
409
410 template< int codim, class GridImp, bool leafIterator >
412 {
413 ElementInfo elementInfo = entity_.impl().elementInfo_;
414 goNext ( elementInfo );
415 // it is ok to set the invalid ElementInfo
416 entity_.impl().setElement( elementInfo, subEntity_ );
417 }
418
419
420 template< int codim, class GridImp, bool leafIterator >
422 ::nextElement ( ElementInfo &elementInfo )
423 {
424 if( elementInfo.isLeaf() || (elementInfo.level() >= level_) )
425 {
426 while( (elementInfo.level() > 0) && (elementInfo.indexInFather() == 1) )
427 elementInfo = elementInfo.father();
428 if( elementInfo.level() == 0 )
429 {
430 ++macroIterator_;
431 elementInfo = *macroIterator_;
432 }
433 else
434 elementInfo = elementInfo.father().child( 1 );
435 }
436 else
437 elementInfo = elementInfo.child( 0 );
438 }
439
440
441 template< int codim, class GridImp, bool leafIterator >
442 inline void AlbertaGridTreeIterator< codim, GridImp, leafIterator >
443 ::nextElementStop ( ElementInfo &elementInfo )
444 {
445 while( !(!elementInfo || stopAtElement( elementInfo )) )
446 nextElement( elementInfo );
447 }
448
449
450 template< int codim, class GridImp, bool leafIterator >
451 inline bool AlbertaGridTreeIterator< codim, GridImp, leafIterator >
452 ::stopAtElement ( const ElementInfo &elementInfo ) const
453 {
454 if( !elementInfo )
455 return true;
456 return (leafIterator ? elementInfo.isLeaf() : (level_ == elementInfo.level()));
457 }
458
459
460 template< int codim, class GridImp, bool leafIterator >
461 inline void AlbertaGridTreeIterator< codim, GridImp, leafIterator >
462 ::goNext ( ElementInfo &elementInfo )
463 {
464 std::integral_constant< int, codim > codimVariable;
465 goNext( codimVariable, elementInfo );
466 }
467
468 template< int codim, class GridImp, bool leafIterator >
469 inline void AlbertaGridTreeIterator< codim, GridImp, leafIterator >
470 ::goNext ( const std::integral_constant< int, 0 > /* cdVariable */,
471 ElementInfo &elementInfo )
472 {
473 assert( stopAtElement( elementInfo ) );
474
475 nextElement( elementInfo );
476 nextElementStop( elementInfo );
477 }
478
479 template< int codim, class GridImp, bool leafIterator >
480 inline void AlbertaGridTreeIterator< codim, GridImp, leafIterator >
481 ::goNext ( const std::integral_constant< int, 1 > cdVariable,
482 ElementInfo &elementInfo )
483 {
484 assert( stopAtElement( elementInfo ) );
485
486 ++subEntity_;
487 if( subEntity_ >= numSubEntities )
488 {
489 subEntity_ = 0;
490 nextElement( elementInfo );
491 nextElementStop( elementInfo );
492 if( !elementInfo )
493 return;
494 }
495
496 if( leafIterator )
497 {
498 const int face = (dimension == 1 ? (numSubEntities-1)-subEntity_ : subEntity_);
499
500 const ALBERTA EL *neighbor = elementInfo.elInfo().neigh[ face ];
501 if( (neighbor != NULL) && !elementInfo.isBoundary( face ) )
502 {
503 // face is reached from element with largest number
504 const int elIndex = grid().dofNumbering() ( elementInfo, 0, 0 );
505 const int nbIndex = grid().dofNumbering() ( neighbor, 0, 0 );
506 if( elIndex < nbIndex )
507 goNext( cdVariable, elementInfo );
508 }
509 // uncomment this assertion only if codimension 1 entities are marked
510 // assert( marker_->template subEntityOnElement< 1 >( elementInfo, subEntity_ ) );
511 }
512 else
513 {
514 assert( marker_ != 0 );
515 if( !marker_->template subEntityOnElement< 1 >( elementInfo, subEntity_ ) )
516 goNext( cdVariable, elementInfo );
517 }
518 }
519
520 template< int codim, class GridImp, bool leafIterator >
521 template< int cd >
522 inline void AlbertaGridTreeIterator< codim, GridImp, leafIterator >
523 ::goNext ( const std::integral_constant< int, cd > cdVariable,
524 ElementInfo &elementInfo )
525 {
526 assert( stopAtElement( elementInfo ) );
527
528 ++subEntity_;
529 if( subEntity_ >= numSubEntities )
530 {
531 subEntity_ = 0;
532 nextElement( elementInfo );
533 nextElementStop( elementInfo );
534 if( !elementInfo )
535 return;
536 }
537
538 assert( marker_ != 0 );
539 if( !marker_->template subEntityOnElement< cd >( elementInfo, subEntity_ ) )
540 goNext( cdVariable, elementInfo );
541 }
542
543}
544
545#endif // #if HAVE_ALBERTA
546
547#endif // #ifndef DUNE_ALBERTA_TREEITERATOR_HH
Definition: treeiterator.hh:189
bool equals(const This &other) const
equality
Definition: treeiterator.hh:231
AlbertaGridTreeIterator(const GridImp &grid, const MarkerVector *marker, int travLevel)
Constructor making begin iterator.
Definition: treeiterator.hh:352
AlbertaGridTreeIterator(const GridImp &grid, int travLevel)
Constructor making end iterator.
Definition: treeiterator.hh:373
This & operator=(const This &other)
Constructor making end iterator.
Definition: treeiterator.hh:398
Entity & dereference() const
dereferencing
Definition: treeiterator.hh:237
int level() const
ask for level of entities
Definition: treeiterator.hh:243
void increment()
increment
Definition: treeiterator.hh:411
const GridImp & grid() const
obtain a reference to the grid
Definition: treeiterator.hh:253
AlbertaGridTreeIterator(const This &other)
Constructor making end iterator.
Definition: treeiterator.hh:386
marker assigning subentities to one element containing them
Definition: treeiterator.hh:35
AlbertaMarkerVector(const DofNumbering &dofNumbering)
create AlbertaMarkerVector with empty vectors
Definition: treeiterator.hh:54
bool up2Date() const
return true if marking is up to date
Definition: treeiterator.hh:95
bool subEntityOnElement(const ElementInfo &elementInfo, int subEntity) const
visit subentity on this element?
Definition: treeiterator.hh:294
void print(std::ostream &out=std::cout) const
print for debugin' only
Definition: treeiterator.hh:319
Wrapper class for entities.
Definition: entity.hh:66
provides a wrapper for ALBERTA's el_info structure
constexpr void forEach(Range &&range, F &&f)
Range based for loop.
Definition: hybridutilities.hh:256
provides a wrapper for ALBERTA's mesh structure
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
Static tag representing a codimension.
Definition: dimension.hh:24
Traits for type conversions and type information.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)