dune-fem  2.8-git
blockmapper.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_BLOCKMAPPER_HH
2 #define DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_BLOCKMAPPER_HH
3 
4 #include <cassert>
5 #include <cstddef>
6 
7 #include <algorithm>
8 #include <functional>
9 #include <limits>
10 #include <type_traits>
11 #include <utility>
12 #include <vector>
13 
14 #include <dune/common/exceptions.hh>
15 
16 #include <dune/geometry/dimension.hh>
17 
18 #include <dune/grid/common/gridenums.hh>
19 #include <dune/grid/utility/persistentcontainer.hh>
20 
25 
27 
28 #include "datahandle.hh"
29 #include "localdofstorage.hh"
30 
31 namespace Dune
32 {
33 
34  namespace Fem
35  {
36 
37  namespace hpDG
38  {
39 
40  // Internal forward declaration
41  // ----------------------------
42 
43  template< class GridPart, class LocalKeys >
44  struct DiscontinuousGalerkinBlockMapper;
45 
46 
47 
48 #ifndef DOXYGEN
49 
50  // DiscontinuousGalerkinBlockMapperTraits
51  // --------------------------------------
52 
53  template< class GridPart, class LocalKeys >
54  struct DiscontinuousGalerkinBlockMapperTraits
55  {
56  using DofMapperType = DiscontinuousGalerkinBlockMapper< GridPart, LocalKeys >;
57  using ElementType = typename GridPart::template Codim< 0 >::EntityType;
58  using SizeType = std::size_t;
59  };
60 
61 #endif // #ifndef DOXYGEN
62 
63 
64 
65  // DiscontinuousGalerkinBlockMapper
66  // --------------------------------
67 
99  template< class GridPart, class LocalKeys >
101  : public AdaptiveDofMapper< DiscontinuousGalerkinBlockMapperTraits< GridPart, LocalKeys > >
102  {
105 
106  public:
108  using SizeType = typename BaseType::SizeType;
111 
113  using GridPartType = GridPart;
116 
118  using LocalKeysType = LocalKeys;
120  using KeyType = typename LocalKeysType::KeyType;
121 
122  private:
124  friend class DataHandle< DiscontinuousGalerkinBlockMapper< GridPart, LocalKeys > >;
125 
126  struct Reserve;
127  struct Resize;
128 
130 
131  using GridType = typename GridPartType::GridType;
132  using GridElementType = typename GridType::template Codim< 0 >::Entity;
133 
134  public:
137 
142  template< class Function >
144  const LocalKeysType &localKeys,
145  const KeyType &value, Function function )
146  : gridPart_( gridPart ),
147  localKeys_( localKeys ),
148  key_( value ),
149  keys_( gridPart_.get().grid(), 0, std::make_pair( key_, key_ ) ),
150  size_( 0u ),
151  dofs_( gridPart_.get().grid(), 0 ),
152  dofManager_( DofManagerType::instance( gridPart.grid() ) )
153  {
154  auto first = gridPart.template begin< 0, All_Partition >();
155  auto last = gridPart.template end< 0, All_Partition >();
156  for( ; first != last; ++first )
157  {
158  const ElementType &element = *first;
159  const KeyType key = function( element );
160 
161  const GridElementType &gridElement = gridEntity( element );
162  keys_[ gridElement ] = std::make_pair( key, key );
163 
164  auto &dofs = dofs_[ gridElement ];
165  resize( dofs, localKeys.blocks( gridElement.type(), key ) );
166  }
167 
168  dofManager_.get().addIndexSet( *this );
169  }
170 
172  const LocalKeysType &localKeys,
173  const KeyType &value )
174  : DiscontinuousGalerkinBlockMapper( gridPart, localKeys, value, [&value]( const ElementType & ){ return value; } )
175  {}
176 
185 
188 
190  ThisType &operator= ( const ThisType & ) = delete;
191 
193  ThisType &operator= ( ThisType && ) = default;
194 
197 #ifndef DOXYGEN
198 
200  {
201  dofManager_.get().removeIndexSet( *this );
202  }
203 
204 #endif // #ifndef DOXYGEN
205 
211  SizeType size () const { return size_; }
212 
214  int maxNumDofs () const { return localKeys().maxBlocks(); }
215 
217  SizeType numDofs ( const ElementType &element ) const
218  {
219  return numDofs( element, Codim< ElementType::codimension >() );
220  }
221 
223  template< class Entity >
224  SizeType numEntityDofs ( const Entity &entity ) const
225  {
226  return numDofs( entity, Codim< Entity::codimension >() );
227  }
228 
229  void onSubEntity ( const ElementType &element, int i, int c, std::vector< bool > &indices ) const
230  {
231  indices.resize( numDofs(element) );
232  if (c == 0)
233  std::fill(indices.begin(),indices.end(),true);
234  else
235  std::fill(indices.begin(),indices.end(),false);
236  }
237 
239  static constexpr bool contains ( const int codim ) { return (codim == 0); }
240 
242  static constexpr bool fixedDataSize ( int codim ) { return (codim != 0); }
243 
245  template< class Function >
246  void mapEach ( const ElementType &element, Function function ) const
247  {
248  mapEach( element, function, Codim< ElementType::codimension >() );
249  }
250 
252  template< class Entity, class Function >
253  void mapEachEntityDof ( const Entity &entity, Function function ) const
254  {
255  mapEach( entity, function, Codim< Entity::codimension >() );
256  }
257 
259  SizeType numberOfHoles ( const int block ) const
260  {
261  assert( block == 0 );
262  return indices_.size();
263  }
264 
266  GlobalKeyType oldIndex ( const int hole, const int block ) const
267  {
268  assert( block == 0 );
269  GlobalKeyType dof = indices_[ hole ].first;
271  return std::move( dof );
272  }
273 
275  GlobalKeyType newIndex ( const int hole, const int block ) const
276  {
277  assert( block == 0 );
278  return indices_[ hole ].second;
279  }
280 
282  static constexpr bool consecutive () { return true; }
283 
285  static constexpr SizeType oldOffSet ( const int block ) { return 0u; }
286 
288  static constexpr SizeType offSet ( const int block ) { return 0u; }
289 
291  static constexpr SizeType numBlocks () { return 1u; }
292 
300  template< class Element >
301  typename std::enable_if< (std::is_same< Element, ElementType >::value || std::is_same< Element, GridElementType >::value), const KeyType & >::type
302  key ( const Element &element ) const
303  {
304  return keys_[ gridEntity( element ) ].first;
305  }
306 
308  void mark ( const KeyType &key, const ElementType &element )
309  {
310  assert( gridPart().indexSet().contains( element ) );
311  keys_[ gridEntity( element ) ].second = key;
312  }
313 
315  KeyType getMark ( const ElementType &element ) const
316  {
317  assert( gridPart().indexSet().contains( element ) );
318  return keys_[ gridEntity( element ) ].second;
319  }
320 
322  template< class Function >
323  bool adapt ( Function function );
324 
326  bool adapt ()
327  {
328  auto function = []( const ElementType &,
329  const KeyType &, const KeyType &,
330  const std::vector< std::size_t > &,
331  const std::vector< std::size_t > & )
332  {};
333  return adapt( function );
334  }
335 
343  const DofManagerType &dofManager () const { return dofManager_.get(); }
344 
346  void insertEntity ( const GridElementType &gridElement )
347  {
348  resize( dofs_ );
349  resize( keys_, std::make_pair( key_, key_ ) );
350 
351  resize( dofs_[ gridElement ], localKeys().blocks( gridElement.type(), key( gridElement ) ) );
352  }
353 
355  void removeEntity ( const GridElementType &gridElement )
356  {
357  auto &dofs = dofs_[ gridElement ];
358 
359  size_ = dofs.reserve( 0u, Reserve( size_ ) );
360  Resize function( holes_ );
361  for( auto dof : dofs )
362  function( dof );
363  }
364 
366  void insertNewEntity ( const GridElementType &gridElement )
367  {
368  if( gridElement.isNew() )
369  {
370  resize( dofs_[ gridElement ], localKeys().blocks( gridElement.type(), key( gridElement ) ) );
371  assert( gridElement.hasFather() );
372  const GridElementType &father = gridElement.father();
373  insertNewEntity( father );
374  removeEntity( father );
375  }
376  else
377  {
378  auto &dofs = dofs_[ gridElement ];
379  resize( dofs, localKeys().blocks( gridElement.type(), key( gridElement ) ) );
380  for( auto dof : dofs )
381  {
382  auto iterator = std::lower_bound( holes_.begin(), holes_.end(), dof );
383  if( iterator != holes_.end() && *iterator == dof )
384  holes_.erase( iterator );
385  }
386  }
387  }
388 
389  void resize ()
390  {
391  resize( dofs_ );
392  resize( keys_, std::make_pair( key_, key_ ) );
393 
394  auto first = gridPart().template begin< 0, All_Partition >();
395  auto last = gridPart().template end< 0, All_Partition >();
396  for( ; first != last; ++first )
397  {
398  const ElementType &element = *first;
399  insertNewEntity( gridEntity( element ) );
400  }
401 #if 0
402  auto &dofs = dofs_[ gridEntity( element ) ];
403  resize( dofs, localKeys().blocks( element.type(), key( element ) ) );
404 #endif
405  }
406 
408  bool compress ();
409 
411  template< class Traits >
413  {
414  DUNE_THROW( NotImplemented, "Method write() not implemented yet" );
415  }
416 
418  template< class Traits >
420  {
421  DUNE_THROW( NotImplemented, "Method read() not implemented yet" );
422  }
423 
425  void backup () const
426  {
427  DUNE_THROW( NotImplemented, "Method backup() not implemented" );
428  }
429 
431  void restore ()
432  {
433  DUNE_THROW( NotImplemented, "Method restore() not implemented" );
434  }
435 
436  private:
437  bool compressed () const { return holes_.empty(); }
438 
439  template< class Element, class Function >
440  void mapEach ( const Element &element, Function function, Codim< 0 > ) const
441  {
442  const auto &dofs = dofs_[ gridEntity( element ) ];
443  std::size_t size = dofs.size();
444  for( std::size_t i = 0; i < size; ++i )
445  function( i, dofs[ i ] );
446  }
447 
448  template< class Entity, class Function, int codim >
449  void mapEach ( const Entity &entity, Function function, Codim< codim > ) const
450  {}
451 
452  template< class Element >
453  SizeType numDofs ( const Element &element, Codim< 0 > ) const
454  {
455  return dofs_[ gridEntity( element ) ].size();
456  }
457 
458  template< class Entity, int codim >
459  SizeType numDofs ( const Entity &entity, Codim< codim > ) const
460  {
461  return 0u;
462  }
463 
464  void resize ( LocalDofStorageType &dofs, SizeType n )
465  {
466  size_ = dofs.reserve( n, Reserve( size_ ) );
467  dofs.resize( Resize( holes_ ) );
468  }
469 
470  template< class T >
471  static void resize ( PersistentContainer< GridType, T > &container, const T &value = T() )
472  {
473  container.resize( value );
474  container.shrinkToFit();
475  }
476 
477  const GridPartType &gridPart () const { return gridPart_.get(); }
478 
479  const LocalKeysType &localKeys () const { return localKeys_.get(); }
480 
481  std::reference_wrapper< const GridPartType > gridPart_;
482  std::reference_wrapper< const LocalKeysType > localKeys_;
483 
484  const KeyType key_;
485  PersistentContainer< GridType, std::pair< KeyType, KeyType > > keys_;
486 
487  SizeType size_;
488  PersistentContainer< GridType, LocalDofStorageType > dofs_;
489 
490  std::vector< GlobalKeyType > holes_;
491  std::vector< std::pair< GlobalKeyType, GlobalKeyType > > indices_;
492 
493  std::reference_wrapper< DofManagerType > dofManager_;
494  };
495 
496 
497 
498  // DiscontinuousGalerkinBlockMapper::Reserve
499  // -----------------------------------------
500 
501  template< class GridPart, class LocalKeys >
502  struct DiscontinuousGalerkinBlockMapper< GridPart, LocalKeys >::Reserve
503  {
504  explicit Reserve ( SizeType &size ) : size_( size ) {}
505 
506  GlobalKeyType operator() () { return size_++; }
507 
508  operator SizeType () const { return size_; }
509 
510  private:
511  SizeType &size_;
512  };
513 
514 
515 
516  // DiscontinuousGalerkinBlockMapper::Resize
517  // ----------------------------------------
518 
519  template< class GridPart, class LocalKeys >
520  struct DiscontinuousGalerkinBlockMapper< GridPart, LocalKeys >::Resize
521  {
522  explicit Resize ( std::vector< GlobalKeyType > &holes )
523  : holes_( holes )
524  {
525  assert( std::is_sorted( holes_.begin(), holes_.end() ) );
526  }
527 
529  {
530  assert( std::is_sorted( holes_.begin(), holes_.end() ) );
531  }
532 
533  void operator() ( const GlobalKeyType &dof )
534  {
535  auto iterator = std::lower_bound( holes_.begin(), holes_.end(), dof );
536  if( iterator == holes_.end() || *iterator != dof )
537  holes_.insert( iterator, dof );
538  }
539 
540  private:
541  std::vector< GlobalKeyType > &holes_;
542  };
543 
544 
545 
546  // DiscontinuousGalerkinBlockMapper::adapt
547  // ---------------------------------------
548 
549  template< class GridPart, class LocalKeys >
550  template< class Function >
552  {
553  DataHandleType dataHandle( *this );
554  gridPart().communicate( dataHandle, InteriorBorder_All_Interface, ForwardCommunication );
555 
556  std::vector< std::size_t > origin, destination;
557  origin.reserve( maxNumDofs() );
558  destination.reserve( maxNumDofs() );
559 
560  std::size_t count = 0;
561 
562  auto first = gridPart().template begin< 0, All_Partition >();
563  auto last = gridPart().template end< 0, All_Partition >();
564  for( ; first != last; ++first )
565  {
566  const ElementType &element = *first;
567  const GridElementType &gridElement = gridEntity( element );
568  auto &keys = keys_[ gridElement ];
569  if( keys.first == keys.second )
570  continue;
571 
572  // remember old state
573  const KeyType prior = keys.first;
574  auto &dofs = dofs_[ gridElement ];
575  origin.resize( dofs.size() );
576  std::copy( dofs.begin(), dofs.end(), origin.begin() );
577 
578  // reset to new state
579  keys.first = keys.second;
580  resize( dofs, localKeys().blocks( gridElement.type(), keys.first ) );
581 
582  // call function
583  destination.resize( dofs.size() );
584  std::copy( dofs.begin(), dofs.end(), destination.begin() );
585  function( element, prior, keys.first, origin, destination );
586 
587  ++count;
588  }
589  return (count > 0u);
590  }
591 
592 
593 
594  // DiscontinuousGalerkinBlockMapper::compress
595  // ------------------------------------------
596 
597  template< class GridPart, class LocalKeys >
599  {
600  // resize persistent containers
601  resize( dofs_ );
602  resize( keys_, std::make_pair( key_, key_ ) );
603 
604  for( auto &dofs : dofs_ )
605  dofs.resize( Resize( holes_ ) );
606 
607  // check if compress is necessary
608  if( holes_.empty() )
609  {
610  indices_.clear();
611  return false;
612  }
613 
614  // get number of dofs after compress
615  assert( size_ >= holes_.size() );
616  size_ -= holes_.size();
617 
618  // remove trailing indices
619  auto iterator = std::lower_bound( holes_.begin(), holes_.end(), size_ );
620  holes_.erase( iterator, holes_.end() );
621 
622  // initialize vector of old and new indices
623  auto assign = []( GlobalKeyType newIndex ) -> std::pair< GlobalKeyType, GlobalKeyType >
624  {
626  return std::make_pair( oldIndex, newIndex );
627  };
628  indices_.resize( holes_.size() );
629  std::transform( holes_.begin(), holes_.end(), indices_.begin(), assign );
630 
631  // update local dof storages
632  std::size_t hole = 0u;
633  for( auto &dofs : dofs_ )
634  {
635  for( auto &dof : dofs )
636  {
637  if( dof >= size_ )
638  {
639  auto &indices = indices_.at( hole++ );
640  indices.first = dof;
641  dof = indices.second;
642  }
643  }
644  }
645 
646  assert( hole == holes_.size() );
647  holes_.clear();
648 
649  // replace all but markers currently used by default value
650  PersistentContainer< GridType, std::pair< KeyType, KeyType > > tmp( gridPart().grid(), 0, std::make_pair( key_, key_ ) );
651  auto first = gridPart().template begin< 0, All_Partition >();
652  auto last = gridPart().template end< 0, All_Partition >();
653  for( ; first != last; ++first )
654  {
655  const ElementType &element = *first;
656  const GridElementType &gridElement = gridEntity( element );
657  tmp[ gridElement ] = keys_[ gridElement ];
658  }
659  keys_.swap( tmp );
660 
661  return true;
662  }
663 
664  } // namespace hpDG
665 
666  namespace Capabilities
667  {
668  // isConsecutiveIndexSet
669  // ---------------------
670 
671  template< class GridPart, class LocalKeys >
672  struct isConsecutiveIndexSet< hpDG::DiscontinuousGalerkinBlockMapper< GridPart, LocalKeys > >
673  {
674  static const bool v = true;
675  };
676 
677  // isAdaptiveDofMapper
678  // -------------------
679 
680  template< class GridPart, class LocalKeys >
681  struct isAdaptiveDofMapper< hpDG::DiscontinuousGalerkinBlockMapper< GridPart, LocalKeys > >
682  {
683  static const bool v = true;
684  };
685 
686  } // namespace Capabilities
687 
688  } // namespace Fem
689 
690 } // namespace Dune
691 
692 #endif // #ifndef DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_BLOCKMAPPER_HH
Definition: bindguard.hh:11
std::tuple_element< i, Tuple >::type & get(Dune::TypeIndexedTuple< Tuple, Types > &tuple)
Definition: typeindexedtuple.hh:122
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
const GridEntityAccess< Entity >::GridEntityType & gridEntity(const Entity &entity)
Definition: gridpart.hh:412
static constexpr T max(T a)
Definition: utility.hh:77
Abstract class representing a function.
Definition: common/function.hh:50
specialize with true if index set implements the interface for consecutive index sets
Definition: common/indexset.hh:42
static const bool v
Definition: common/indexset.hh:49
abstract interface for an output stream
Definition: streams.hh:46
abstract interface for an input stream
Definition: streams.hh:179
Definition: dofmanager.hh:762
An -adaptive Dune::Fem::DofMapper.
Definition: blockmapper.hh:102
void backup() const
this mapper has no I/O capabilities
Definition: blockmapper.hh:425
static constexpr SizeType oldOffSet(const int block)
return 0 (this mapper has no offset)
Definition: blockmapper.hh:285
DiscontinuousGalerkinBlockMapper(const GridPartType &gridPart, const LocalKeysType &localKeys, const KeyType &value)
Definition: blockmapper.hh:171
SizeType size() const
return number of dofs
Definition: blockmapper.hh:211
GlobalKeyType newIndex(const int hole, const int block) const
return new index of given hole during compression
Definition: blockmapper.hh:275
typename LocalKeysType::KeyType KeyType
key type
Definition: blockmapper.hh:120
static constexpr bool consecutive()
return true (this mapper yields a consecutive DOF numbering)
Definition: blockmapper.hh:282
GridPart GridPartType
grid part type
Definition: blockmapper.hh:113
void write(OutStreamInterface< Traits > &)
this mapper has no I/O capabilities
Definition: blockmapper.hh:412
void mapEachEntityDof(const Entity &entity, Function function) const
map local dof to global key
Definition: blockmapper.hh:253
DiscontinuousGalerkinBlockMapper(const ThisType &)=delete
copy constructor
void mapEach(const ElementType &element, Function function) const
map local dof to global key
Definition: blockmapper.hh:246
LocalKeys LocalKeysType
basis function sets type
Definition: blockmapper.hh:118
ThisType & operator=(const ThisType &)=delete
assignment operator
SizeType numEntityDofs(const Entity &entity) const
return number of dofs for given element
Definition: blockmapper.hh:224
typename BaseType::ElementType ElementType
element type
Definition: blockmapper.hh:115
void insertNewEntity(const GridElementType &gridElement)
add DOFs for new element
Definition: blockmapper.hh:366
bool adapt()
please doc me
Definition: blockmapper.hh:326
void insertEntity(const GridElementType &gridElement)
add DOFs for element
Definition: blockmapper.hh:346
static constexpr bool fixedDataSize(int codim)
return true if number of dofs is fixed for given codimension
Definition: blockmapper.hh:242
const DofManagerType & dofManager() const
return DOF manager
Definition: blockmapper.hh:343
DiscontinuousGalerkinBlockMapper(const GridPartType &gridPart, const LocalKeysType &localKeys, const KeyType &value, Function function)
Definition: blockmapper.hh:143
KeyType getMark(const ElementType &element) const
get key to be assigned to an entity after next call to adapt()
Definition: blockmapper.hh:315
bool compress()
compress DOF mapping
Definition: blockmapper.hh:598
DiscontinuousGalerkinBlockMapper(ThisType &&)=default
move constructor
typename BaseType::GlobalKeyType GlobalKeyType
global key type
Definition: blockmapper.hh:110
void read(InStreamInterface< Traits > &)
this mapper has no I/O capabilities
Definition: blockmapper.hh:419
void onSubEntity(const ElementType &element, int i, int c, std::vector< bool > &indices) const
Definition: blockmapper.hh:229
SizeType numDofs(const ElementType &element) const
return number of dofs for given element
Definition: blockmapper.hh:217
int maxNumDofs() const
return upper bound for number of dofs
Definition: blockmapper.hh:214
SizeType numberOfHoles(const int block) const
return number of holes during compression
Definition: blockmapper.hh:259
static constexpr bool contains(const int codim)
return true if dofs are associated to codimension
Definition: blockmapper.hh:239
void removeEntity(const GridElementType &gridElement)
mark DOFs for removal
Definition: blockmapper.hh:355
std::enable_if<(std::is_same< Element, ElementType >::value||std::is_same< Element, GridElementType >::value), const KeyType & >::type key(const Element &element) const
get key currently assigned to an entity
Definition: blockmapper.hh:302
GlobalKeyType oldIndex(const int hole, const int block) const
return old index of given hole during compression
Definition: blockmapper.hh:266
static constexpr SizeType offSet(const int block)
return 0 (this mapper has no offset)
Definition: blockmapper.hh:288
void resize()
Definition: blockmapper.hh:389
void restore()
this mapper has no I/O capabilities
Definition: blockmapper.hh:431
void mark(const KeyType &key, const ElementType &element)
set key to be assigned to an entity after next call to adapt()
Definition: blockmapper.hh:308
static constexpr SizeType numBlocks()
return 1 (this mapper has one block)
Definition: blockmapper.hh:291
typename BaseType::SizeType SizeType
size type
Definition: blockmapper.hh:108
Reserve(SizeType &size)
Definition: blockmapper.hh:504
Resize(std::vector< GlobalKeyType > &holes)
Definition: blockmapper.hh:522
Definition: space/hpdg/datahandle.hh:32
Definition: localdofstorage.hh:25
Definition: space/mapper/capabilities.hh:22
static const bool v
Definition: space/mapper/capabilities.hh:23
Traits::ElementType ElementType
type of codimension 0 entities
Definition: mapper/dofmapper.hh:54
Extended interface for adaptive DoF mappers.
Definition: mapper/dofmapper.hh:219
SizeType GlobalKeyType
at the moment this should be similar to SizeType
Definition: mapper/dofmapper.hh:230
BaseType::SizeType SizeType
type of size integer
Definition: mapper/dofmapper.hh:227