1#ifndef DUNE_FEM_DOFMAPPER_INDEXSETDOFMAPPER_HH
2#define DUNE_FEM_DOFMAPPER_INDEXSETDOFMAPPER_HH
9#include <dune/geometry/referenceelements.hh>
13#include <dune/fem/gridpart/common/gridpart.hh>
14#include <dune/fem/gridpart/common/indexset.hh>
15#include <dune/fem/misc/functor.hh>
16#include <dune/fem/space/common/dofmanager.hh>
17#include <dune/fem/space/mapper/code.hh>
18#include <dune/fem/space/mapper/exceptions.hh>
19#include <dune/fem/space/mapper/dofmapper.hh>
22#pragma GCC diagnostic push
24#pragma GCC diagnostic ignored "-Wattributes"
35 template<
class Gr
idPart >
36 class DefaultLocalDofMapping
40 template<
class Iterator,
class Functor >
41 void operator() ( std::size_t index,
unsigned int numDofs, Iterator begin, Iterator end, Functor functor )
const
44 functor( *(begin++), index++ );
49 DefaultLocalDofMapping () {}
50 DefaultLocalDofMapping (
const GridPart & ) {}
52 Mapping operator() (
const typename GridPart::template Codim< 0 >::EntityType &element,
unsigned int subEntity,
unsigned int codim )
const {
return {}; }
57 namespace __IndexSetDofMapper
63 template<
class Gr
idPart,
class LocalDofMapping >
66 typedef DofMapperBase< GridPart, LocalDofMapping > ThisType;
69 typedef std::size_t SizeType;
80 SizeType offset, oldOffset;
83 enum CodimType { CodimEmpty, CodimFixedSize, CodimVariableSize };
84 static const int dimension = GridPart::dimension;
90 struct SubEntityFilter
92 SubEntityFilter(
const RefElementType &refElement,
int subEntity,
int codim)
93 : active_(dimension+1), size_(0)
95 for (
int c=0;c<=dimension;++c)
97 std::vector<bool> &a = active_[c];
98 a.resize( refElement.size( c ),
false );
99 if (c<codim)
continue;
100 if (c==codim) { a[subEntity]=
true; ++size_;
continue; }
101 for (
int i=0;i<refElement.size(subEntity, codim, c);++i)
103 a[refElement.subEntity(subEntity, codim, i, c)] =
true;
108 bool operator()(
int i,
int c)
const {
return active_[c][i]; }
110 std::vector< std::vector< bool > > active_;
114 template<
class Functor >
118 typedef SizeType GlobalKeyType;
120 typedef GridPart GridPartType;
121 typedef LocalDofMapping LocalDofMappingType;
122 typedef typename GridPart::GridType GridType;
124 typedef typename GridPartType::template Codim< 0 >::EntityType ElementType;
126 template<
class CodeFactory >
127 DofMapperBase (
const GridPartType &gridPart, LocalDofMappingType localDofMapping,
const CodeFactory &codeFactory );
150 template<
class Functor >
151 void mapEach (
const ElementType &element, Functor f )
const;
153 void map (
const ElementType &element, std::vector< GlobalKeyType > &indices )
const;
155 [[deprecated(
"Use onSubEntity method with char vector instead")]]
156 void onSubEntity (
const ElementType &element,
int i,
int c, std::vector< bool > &indices )
const
158 std::vector< char > _idx;
159 onSubEntity(element, i, c, _idx);
160 indices.resize( _idx.size() );
161 for (std::size_t i=0; i<_idx.size();++i)
162 _idx[i] = indices[i] > 0;
180 void onSubEntity(
const ElementType &element,
int i,
int c, std::vector< char > &indices )
const;
182 unsigned int maxNumDofs ()
const {
return maxNumDofs_; }
183 unsigned int numDofs (
const ElementType &element )
const {
return code( element ).numDofs(); }
186 template<
class Entity,
class Functor >
187 void mapEachEntityDof (
const Entity &entity, Functor f )
const;
189 template<
class Entity >
190 void mapEntityDofs (
const Entity &entity, std::vector< GlobalKeyType > &indices )
const;
192 template<
class Entity >
193 unsigned int numEntityDofs (
const Entity &entity )
const;
197 bool contains (
int codim )
const {
return (codimType_[ codim ] != CodimEmpty); }
199 bool fixedDataSize (
int codim )
const {
return (codimType_[ codim ] == CodimFixedSize); }
201 SizeType
size ()
const {
return size_; }
222 static constexpr bool consecutive () noexcept {
return false; }
224 SizeType numBlocks ()
const
226 DUNE_THROW( NotImplemented,
"Method numBlocks() called on non-adaptive block mapper" );
229 SizeType numberOfHoles (
int )
const
231 DUNE_THROW( NotImplemented,
"Method numberOfHoles() called on non-adaptive block mapper" );
234 GlobalKeyType oldIndex (
int hole,
int )
const
236 DUNE_THROW( NotImplemented,
"Method oldIndex() called on non-adaptive block mapper" );
239 GlobalKeyType newIndex (
int hole,
int )
const
241 DUNE_THROW( NotImplemented,
"Method newIndex() called on non-adaptive block mapper" );
244 SizeType oldOffSet (
int )
const
246 DUNE_THROW( NotImplemented,
"Method oldOffSet() called on non-adaptive block mapper" );
249 SizeType offSet (
int )
const
251 DUNE_THROW( NotImplemented,
"Method offSet() called on non-adaptive block mapper" );
254 template<
class Entity >
255 void insertEntity (
const Entity &entity )
257 DUNE_THROW( NotImplemented,
"Method insertEntity(entity) called on non-adaptive block mapper" );
259 template<
class Entity >
260 void removeEntity (
const Entity &entity )
262 DUNE_THROW( NotImplemented,
"Method removeEntity(entity) called on non-adaptive block mapper" );
265 void resize () { update(); }
266 bool compress () { update();
return true;}
267 void backup ()
const {}
269 template <
class StreamTraits>
270 void write( OutStreamInterface< StreamTraits >& out )
const {}
271 template <
class StreamTraits>
272 void read( InStreamInterface< StreamTraits >& in ) { update(); }
278 void requestCodimensions ();
281 typedef std::vector< GeometryType > BlockMapType;
283 const DofMapperCode &code (
const GeometryType &
gt )
const;
284 const DofMapperCode &code (
const ElementType &element )
const {
return code( element.type() ); }
286 template<
class Entity >
287 const SubEntityInfo &subEntityInfo (
const Entity &entity )
const;
289 const IndexSetType &indexSet ()
const {
return indexSet_; }
291 const IndexSetType &indexSet_;
292 LocalDofMapping localDofMapping_;
293 std::vector< DofMapperCode > code_;
294 unsigned int maxNumDofs_;
296 std::vector< SubEntityInfo > subEntityInfo_;
297 BlockMapType blockMap_;
298 CodimType codimType_[ dimension+1 ];
306 template<
class Gr
idPart,
class LocalDofMapping >
307 struct DofMapperBase< GridPart, LocalDofMapping >::BuildFunctor
309 explicit BuildFunctor ( std::vector< SubEntityInfo > &subEntityInfo )
310 : subEntityInfo_( subEntityInfo )
313 template<
class Iterator >
314 void operator() (
unsigned int gtIndex,
unsigned int subEntity, Iterator it, Iterator end )
316 SubEntityInfo &info = subEntityInfo_[ gtIndex ];
317 const unsigned int numDofs = end - it;
318 if( info.numDofs == 0 )
319 info.numDofs = numDofs;
320 else if( info.numDofs != numDofs )
321 DUNE_THROW( DofMapperError,
"Inconsistent number of DoFs on subEntity (codim = " << info.codim <<
")." );
325 std::vector< SubEntityInfo > &subEntityInfo_;
334 template<
class Gr
idPart,
class LocalDofMapping >
335 const int DofMapperBase< GridPart, LocalDofMapping >::dimension;
337 template<
class Gr
idPart,
class LocalDofMapping >
338 template<
class CodeFactory >
339 inline DofMapperBase< GridPart, LocalDofMapping >
340 ::DofMapperBase (
const GridPartType &gridPart, LocalDofMappingType localDofMapping,
const CodeFactory &codeFactory )
344 : indexSet_( gridPart.indexSet() ),
345 localDofMapping_(
std::move( localDofMapping ) ),
346 code_( LocalGeometryTypeIndex::
size( dimension ) ),
348 subEntityInfo_( GlobalGeometryTypeIndex::
size( dimension ) )
350 std::vector< GeometryType >
gt( GlobalGeometryTypeIndex::size( dimension ) );
352 const typename RefElementsType::Iterator end = RefElementsType::end();
353 for(
typename RefElementsType::Iterator it = RefElementsType::begin(); it != end; ++it )
355 const RefElementType refElement = *it;
357 for(
int codim = 0; codim <= dimension; ++codim )
359 for(
int i = 0; i < refElement.size( codim ); ++i )
361 const unsigned int gtIdx = GlobalGeometryTypeIndex::index( refElement.type( i, codim ) );
362 gt[ gtIdx ] = refElement.type( i, codim );
363 subEntityInfo_[ gtIdx ].codim = codim;
367 DofMapperCode &code = code_[ LocalGeometryTypeIndex::index( refElement.type() ) ];
368 code = codeFactory( refElement );
369 maxNumDofs_ = std::max( code.numDofs(), maxNumDofs_ );
370 code( BuildFunctor( subEntityInfo_ ) );
373 for(
int codim = 0; codim <= dimension; ++codim )
374 codimType_[ codim ] = CodimEmpty;
376 unsigned int codimDofs[ dimension+1 ];
377 for(
unsigned int i = 0; i < subEntityInfo_.size(); ++i )
379 const SubEntityInfo &info = subEntityInfo_[ i ];
380 if( info.numDofs == 0 )
388 const auto & geomTypes = indexSet().types(info.codim);
390 if (hasSingleGeometryType && geomTypes[0] !=
gt[i])
395 if( codimType_[ info.codim ] == CodimEmpty )
396 codimType_[ info.codim ] = CodimFixedSize;
397 else if( codimDofs[ info.codim ] != info.numDofs )
398 codimType_[ info.codim ] = CodimVariableSize;
400 codimDofs[ info.codim ] = info.numDofs;
401 blockMap_.push_back(
gt[ i ] );
405 requestCodimensions ();
412 template<
class Gr
idPart,
class LocalDofMapping >
413 template<
class Functor >
414 inline void DofMapperBase< GridPart, LocalDofMapping >
415 ::mapEach (
const ElementType &element, Functor f )
const
417 const auto &idxSet = indexSet();
419 code( element )( [
this, &idxSet, &element, f ] (
unsigned int gtIndex,
unsigned int subEntity,
auto begin,
auto end ) {
420 const SubEntityInfo &info = subEntityInfo_[ gtIndex ];
421 const SizeType subIndex = idxSet.subIndex( element, subEntity, info.codim );
422 SizeType index = info.offset + SizeType( info.numDofs ) * subIndex;
423 localDofMapping_( element, subEntity, info.codim )( index, info.numDofs, begin, end, f );
428 template<
class Gr
idPart,
class LocalDofMapping >
429 inline void DofMapperBase< GridPart, LocalDofMapping >
430 ::map (
const ElementType &element, std::vector< SizeType > &indices )
const
432 indices.resize( numDofs( element ) );
433 mapEach( element, AssignFunctor< std::vector< SizeType > >( indices ) );
436 template<
class Gr
idPart,
class LocalDofMapping >
437 template<
class Entity,
class Functor >
438 inline void DofMapperBase< GridPart, LocalDofMapping >
439 ::mapEachEntityDof (
const Entity &entity, Functor f )
const
441 const SubEntityInfo &info = subEntityInfo( entity );
442 const unsigned int numDofs = info.numDofs;
443 SizeType index = info.offset + numDofs * SizeType( indexSet().index( entity ) );
444 for(
unsigned int i = 0; i < info.numDofs; ++i )
449 template<
class Gr
idPart,
class LocalDofMapping >
450 template<
class Entity >
451 inline void DofMapperBase< GridPart, LocalDofMapping >
452 ::mapEntityDofs (
const Entity &entity, std::vector< SizeType > &indices )
const
454 indices.resize( numEntityDofs( entity ) );
455 mapEachEntityDof( entity, AssignFunctor< std::vector< SizeType > >( indices ) );
458 template<
class Gr
idPart,
class LocalDofMapping >
459 inline void DofMapperBase< GridPart, LocalDofMapping >
460 ::onSubEntity(
const ElementType &element,
int i,
int c, std::vector< char > &indices )
const
462 const SubEntityFilter filter( RefElementsType::general( element.type() ), i, c );
463 indices.resize( numDofs( element ) );
464 code( element )( [
this, &indices, &filter ] (
unsigned int gtIndex,
unsigned int subEntity,
auto begin,
auto end ) {
465 const bool active = filter( subEntity, subEntityInfo_[ gtIndex ].codim );
466 while( begin != end )
467 indices[ *(begin++) ] = active? 1:0;
471 template<
class Gr
idPart,
class LocalDofMapping >
472 template<
class Entity >
474 DofMapperBase< GridPart, LocalDofMapping >
475 ::numEntityDofs (
const Entity &entity )
const
477 return subEntityInfo( entity ).numDofs;
481 template<
class Gr
idPart,
class LocalDofMapping >
482 inline void DofMapperBase< GridPart, LocalDofMapping >::requestCodimensions ()
485 if constexpr ( Capabilities::isDuneFemIndexSet< IndexSetType >::v )
488 std::vector< int > codimensions;
489 codimensions.reserve( dimension+1 );
491 for(
typename BlockMapType::const_iterator it = blockMap_.begin(); it != blockMap_.end(); ++it )
493 SubEntityInfo &info = subEntityInfo_[ GlobalGeometryTypeIndex::index( *it ) ];
494 codimensions.push_back( info.codim );
498 indexSet().requestCodimensions( codimensions );
502 template<
class Gr
idPart,
class LocalDofMapping >
503 inline void DofMapperBase< GridPart, LocalDofMapping >::update ()
506 for(
const auto& geomType : blockMap_ )
508 SubEntityInfo &info = subEntityInfo_[ GlobalGeometryTypeIndex::index( geomType ) ];
509 info.oldOffset = info.offset;
511 size_ += SizeType( info.numDofs ) * SizeType( indexSet().size( geomType ) );
516 template<
class Gr
idPart,
class LocalDofMapping >
517 inline const DofMapperCode &DofMapperBase< GridPart, LocalDofMapping >
518 ::code (
const GeometryType &
gt )
const
520 return code_[ LocalGeometryTypeIndex::index(
gt ) ];
524 template<
class Gr
idPart,
class LocalDofMapping >
525 template<
class Entity >
526 inline const typename DofMapperBase< GridPart, LocalDofMapping >::SubEntityInfo &
527 DofMapperBase< GridPart, LocalDofMapping >::subEntityInfo (
const Entity &entity )
const
529 return subEntityInfo_[ GlobalGeometryTypeIndex::index( entity.type() ) ];
537 template<
class Gr
idPart,
class LocalDofMapping >
539 :
public DofMapperBase< GridPart, LocalDofMapping >
541 typedef DofMapper< GridPart, LocalDofMapping > ThisType;
542 typedef DofMapperBase< GridPart, LocalDofMapping > BaseType;
545 typedef typename BaseType::SubEntityInfo SubEntityInfo;
546 typedef typename BaseType::GridPartType::GridType GridType;
547 typedef DofManager< GridType > DofManagerType;
550 typedef typename BaseType::GridPartType GridPartType;
551 typedef typename BaseType::LocalDofMappingType LocalDofMappingType;
552 typedef typename BaseType::SizeType SizeType;
554 template<
class CodeFactory >
555 DofMapper (
const GridPartType &gridPart, LocalDofMappingType localDofMapping,
const CodeFactory &codeFactory )
556 : BaseType( gridPart,
std::move( localDofMapping ), codeFactory ),
557 dofManager_( DofManagerType::instance( gridPart.grid() ) )
559 dofManager_.addIndexSet( *
this );
562 DofMapper (
const ThisType & ) =
delete;
566 dofManager_.removeIndexSet( *
this );
569 ThisType &operator= (
const ThisType & ) =
delete;
572 DofManagerType& dofManager_;
579 template<
class Gr
idPart,
class LocalDofMapping >
580 class AdaptiveDofMapper
581 :
public DofMapperBase< GridPart, LocalDofMapping >
583 typedef AdaptiveDofMapper< GridPart, LocalDofMapping > ThisType;
584 typedef DofMapperBase< GridPart, LocalDofMapping > BaseType;
587 typedef typename BaseType::SubEntityInfo SubEntityInfo;
588 typedef typename BaseType::GridPartType::GridType GridType;
589 typedef DofManager< GridType > DofManagerType;
592 typedef typename BaseType::GridPartType GridPartType;
593 typedef typename BaseType::LocalDofMappingType LocalDofMappingType;
594 typedef typename BaseType::SizeType SizeType;
596 template<
class CodeFactory >
597 AdaptiveDofMapper (
const GridPartType &gridPart, LocalDofMappingType localDofMapping,
const CodeFactory &codeFactory )
598 : BaseType( gridPart,
std::move( localDofMapping ), codeFactory ),
599 dofManager_( DofManagerType::instance( gridPart.grid() ) )
601 dofManager_.addIndexSet( *
this );
604 AdaptiveDofMapper (
const ThisType & ) =
delete;
606 ~AdaptiveDofMapper ()
608 dofManager_.removeIndexSet( *
this );
611 ThisType &operator= (
const ThisType & ) =
delete;
619 int numBlocks ()
const {
return blockMap_.size(); }
620 SizeType offSet (
int blk )
const;
621 SizeType oldOffSet (
int blk )
const;
623 SizeType numberOfHoles (
int blk )
const;
625 SizeType oldIndex ( SizeType hole,
int blk )
const;
626 SizeType newIndex ( SizeType hole,
int blk )
const;
630 bool consecutive ()
const {
return true; }
632 template<
class Entity >
633 void insertEntity (
const Entity &entity ) { BaseType::update(); }
635 template<
class Entity >
636 void removeEntity (
const Entity &entity ) {}
638 void resize () { BaseType::update(); }
640 bool compress () { BaseType::update();
return true; }
642 template <
class StreamTraits>
643 void write( OutStreamInterface< StreamTraits >& out )
const {}
645 template <
class StreamTraits>
646 void read( InStreamInterface< StreamTraits >& in )
651 void backup ()
const {}
655 using BaseType::indexSet;
657 using BaseType::blockMap_;
658 using BaseType::subEntityInfo_;
660 DofManagerType& dofManager_;
668 template<
class Gr
idPart,
class LocalDofMapping >
669 inline typename AdaptiveDofMapper< GridPart, LocalDofMapping >::SizeType
670 AdaptiveDofMapper< GridPart, LocalDofMapping >::offSet (
int blk )
const
672 assert( (blk >= 0) && (blk < numBlocks()) );
673 const unsigned int gtIdx = GlobalGeometryTypeIndex::index( blockMap_[ blk ] );
674 return subEntityInfo_[ gtIdx ].offset;
678 template<
class Gr
idPart,
class LocalDofMapping >
679 inline typename AdaptiveDofMapper< GridPart, LocalDofMapping >::SizeType
680 AdaptiveDofMapper< GridPart, LocalDofMapping >::oldOffSet (
int blk )
const
682 assert( (blk >= 0) && (blk < numBlocks()) );
683 const unsigned int gtIdx = GlobalGeometryTypeIndex::index( blockMap_[ blk ] );
684 return subEntityInfo_[ gtIdx ].oldOffset;
688 template<
class Gr
idPart,
class LocalDofMapping >
689 inline typename AdaptiveDofMapper< GridPart, LocalDofMapping >::SizeType
690 AdaptiveDofMapper< GridPart, LocalDofMapping >::numberOfHoles (
int blk )
const
692 assert( (blk >= 0) && (blk < numBlocks()) );
693 const unsigned int gtIdx = GlobalGeometryTypeIndex::index( blockMap_[ blk ] );
694 const SubEntityInfo &info = subEntityInfo_[ gtIdx ];
695 return SizeType( info.numDofs ) * SizeType( indexSet().numberOfHoles( blockMap_[ blk ] ) );
699 template<
class Gr
idPart,
class LocalDofMapping >
700 inline typename AdaptiveDofMapper< GridPart, LocalDofMapping >::SizeType
701 AdaptiveDofMapper< GridPart, LocalDofMapping >::oldIndex ( SizeType hole,
int blk )
const
703 assert( (hole >= 0) && (hole < numberOfHoles( blk )) );
704 const unsigned int gtIdx = GlobalGeometryTypeIndex::index( blockMap_[ blk ] );
705 const SubEntityInfo &info = subEntityInfo_[ gtIdx ];
706 const unsigned int numDofs = info.numDofs;
707 const SizeType index = indexSet().oldIndex( hole / numDofs, blockMap_[ blk ] );
708 return info.offset + numDofs * index + (hole % numDofs);
712 template<
class Gr
idPart,
class LocalDofMapping >
713 inline typename AdaptiveDofMapper< GridPart, LocalDofMapping >::SizeType
714 AdaptiveDofMapper< GridPart, LocalDofMapping >::newIndex ( SizeType hole,
int blk )
const
716 assert( (hole >= 0) && (hole < numberOfHoles( blk )) );
717 const unsigned int gtIdx = GlobalGeometryTypeIndex::index( blockMap_[ blk ] );
718 const SubEntityInfo &info = subEntityInfo_[ gtIdx ];
719 const unsigned int numDofs = info.numDofs;
720 const SizeType index = indexSet().newIndex( hole / numDofs, blockMap_[ blk ] );
721 return info.offset + numDofs * index + (hole % numDofs);
729 template< class GridPart, class LocalDofMapping, bool adaptive = Capabilities::isAdaptiveIndexSet< typename GridPart::IndexSetType >::v >
730 struct Implementation
732 typedef typename std::conditional< adaptive, AdaptiveDofMapper< GridPart, LocalDofMapping >, DofMapper< GridPart, LocalDofMapping > >::type Type;
742 template<
class Gr
idPart,
class LocalDofMapping = DefaultLocalDofMapping< Gr
idPart > >
743 class IndexSetDofMapper
744 :
public __IndexSetDofMapper::template Implementation< GridPart, LocalDofMapping >::Type
746 typedef typename __IndexSetDofMapper::template Implementation< GridPart, LocalDofMapping >::Type BaseType;
749 typedef typename BaseType::GridPartType GridPartType;
750 typedef typename BaseType::LocalDofMappingType LocalDofMappingType;
752 template<
class CodeFactory >
753 IndexSetDofMapper (
const GridPartType &gridPart, LocalDofMappingType localDofMapping,
const CodeFactory &codeFactory )
754 : BaseType( gridPart,
std::move( localDofMapping ), codeFactory )
757 template<
class CodeFactory >
758 IndexSetDofMapper (
const GridPartType &gridPart,
const CodeFactory &codeFactory )
759 : BaseType( gridPart, LocalDofMappingType( gridPart ), codeFactory )
767 namespace Capabilities
772 template<
class Gr
idPart,
class LocalDofMapping >
773 struct isAdaptiveDofMapper< IndexSetDofMapper< GridPart, LocalDofMapping > >
775 static const bool v = Capabilities::isAdaptiveIndexSet< typename GridPart::IndexSetType >::v;
782 template<
class Gr
idPart,
class LocalDofMapping >
783 struct isConsecutiveIndexSet< __IndexSetDofMapper::AdaptiveDofMapper< GridPart, LocalDofMapping > >
785 static const bool v =
true;
795#pragma GCC diagnostic pop
consecutive, persistent index set for the leaf level based on the grid's hierarchy index set
Definition: adaptiveleafindexset.hh:1351
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
bool gt(const T &first, const T &second, typename EpsilonType< T >::Type epsilon)
test if first greater than second
Definition: float_cmp.cc:158
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
constexpr std::bool_constant<((II==value)||...)> contains(std::integer_sequence< T, II... >, std::integral_constant< T, value >)
Checks whether or not a given sequence contains a value.
Definition: integersequence.hh:137
Specialize with 'true' for if the codimension 0 entity of the grid has only one possible geometry typ...
Definition: capabilities.hh:27
Class providing access to the singletons of the reference elements.
Definition: referenceelements.hh:128
typename Container::ReferenceElement ReferenceElement
The reference element type.
Definition: referenceelements.hh:146
A unique label for each type of element that can occur in a grid.
Helper classes to provide indices for geometrytypes for use in a vector.