dune-fem  2.6-git
compile.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_DOFMAPPER_COMPILE_HH
2 #define DUNE_FEM_DOFMAPPER_COMPILE_HH
3 
4 #include <algorithm>
5 #include <type_traits>
6 #include <utility>
7 
8 #include <dune/geometry/referenceelements.hh>
9 #include <dune/geometry/typeindex.hh>
10 
13 
14 namespace Dune
15 {
16 
17  namespace Fem
18  {
19 
20  // generateCodimensionCode
21  // -----------------------
22 
23  template< class RefElement,
24  std::enable_if_t< std::is_same< std::decay_t< decltype( std::declval< const RefElement & >().size( 0 ) ) >, int >::value, int > = 0,
25  std::enable_if_t< std::is_same< std::decay_t< decltype( std::declval< const RefElement & >().type( 0, 0 ) ) >, GeometryType >::value, int > = 0 >
26  inline DofMapperCode generateCodimensionCode ( const RefElement &refElement, int codim, unsigned int blockSize = 1 )
27  {
28  unsigned int count = refElement.size( codim );
29  DofMapperCodeWriter code( count, count*blockSize );
30  unsigned int pos = 0;
31  for( unsigned int i = 0; i < count; ++i )
32  {
33  code[ pos++ ] = GlobalGeometryTypeIndex::index( refElement.type( i, codim ) );
34  code[ pos++ ] = i;
35  code[ pos++ ] = blockSize;
36  for( unsigned int j = 0; j < blockSize; ++j )
37  code[ pos++ ] = i*blockSize + j;
38  }
39  return code;
40  }
41 
42 
43 
44  // compile (for LocalCoefficients)
45  // -------------------------------
46 
47  template< class RefElement, class LocalCoefficients,
48  std::enable_if_t< std::is_same< std::decay_t< decltype( std::declval< const RefElement & >().size( 0 ) ) >, int >::value, int > = 0,
49  std::enable_if_t< std::is_same< std::decay_t< decltype( std::declval< const RefElement & >().type( 0, 0 ) ) >, GeometryType >::value, int > = 0 >
50  inline DofMapperCode compile ( const RefElement &refElement, const LocalCoefficients &localCoefficients )
51  {
52  const int dim = RefElement::dimension;
53 
54  const std::size_t numDofs = localCoefficients.size(); // total number of DoFs
55 
56  // count number of keys per subentity
57 
58  // total number of all sub-entities
59  unsigned int numSubEntities = 0;
60  for( int codim = 0; codim <= dim; ++codim )
61  numSubEntities += refElement.size( codim );
62  assert( numSubEntities > 0 );
63 
64  // form a "matrix" with variable lenght rows. This is the usual
65  // approach: pre-allocate the needed storage once and then
66  // insert the proper offsets into the row-pointer. After
67  // completion count[codim] is an array with one entry for each
68  // sub-entity for the given codim. It is initialized with zeros.
69  unsigned int *count[ dim+1 ];
70  count[ 0 ] = new unsigned int[ numSubEntities ];
71  assert( count[ 0 ] );
72  std::fill( count[ 0 ], count[ 0 ] + numSubEntities, (unsigned int)0 );
73  for( int codim = 0; codim < dim; ++codim )
74  count[ codim+1 ] = count[ codim ] + refElement.size( codim );
75 
76  // Now do the actual counting. After completion
77  // cound[codim][subEntity] will contain the number of DoFs
78  // attached to the particular sub-entity.
79  //
80  // numBlocks is the actual number of __USED__
81  // sub-entities. E.g. for continuous Lagrange-1 on a triangle numBlocks
82  // would be 3, after counting (only the vertices carry DoFs).
83  unsigned int numBlocks = 0;
84  for( std::size_t i = 0; i < numDofs; ++i )
85  {
86  const LocalKey &key = localCoefficients.localKey( i );
87 
88  const int codim = key.codim();
89  const int subEntity = key.subEntity();
90 
91  assert( (codim >= 0) && (codim <= dim) );
92  assert( (subEntity >= 0) && (subEntity < refElement.size( codim )) );
93 
94  if( count[ codim ][ subEntity ] == 0 )
95  ++numBlocks;
96  ++count[ codim ][ subEntity ];
97  }
98 
99  // format the code into subentity blocks
100  // result: count will hold the first local index in the block (0 = unused)
101  //
102  // I.e.: count[cd][subEntIdx] = local index offset for start of
103  // DoFs attached to sub entity
104 
105  DofMapperCodeWriter code( numBlocks, numDofs );
106 
107  unsigned int next = 0;
108  for( int codim = 0; codim <= dim; ++codim )
109  {
110  for( int i = 0; i < refElement.size( codim ); ++i )
111  {
112  const unsigned int cnt = count[ codim ][ i ];
113  if( cnt == 0 )
114  continue;
115 
116  code[ next++ ] = GlobalGeometryTypeIndex::index( refElement.type( i, codim ) );
117  code[ next++ ] = i;
118  code[ next++ ] = cnt;
119 
120  count[ codim ][ i ] = next;
121  next += cnt;
122  }
123  }
124 
125  // fill in the local indices
126  //
127  // Format of the code-array is described in code.hh
128  for( std::size_t i = 0; i < numDofs; ++i )
129  {
130  const LocalKey &key = localCoefficients.localKey( i );
131  const unsigned int block = count[ key.codim() ][ key.subEntity() ];
132  assert( block > 0 );
133  assert( (key.index() >= 0) && (key.index() < code[ block-1 ]) );
134  code[ block + key.index() ] = i;
135  }
136 
137  // clean up and return
138 
139  delete[] count[ 0 ];
140 
141  return code;
142  }
143 
144  } // namespace Fem
145 
146 } // namespace Dune
147 
148 #endif // #ifndef DUNE_FEM_DOFMAPPER_COMPILE_HH
Definition: bindguard.hh:11
DofMapperCode generateCodimensionCode(const RefElement &refElement, int codim, unsigned int blockSize=1)
Definition: compile.hh:26
DofMapperCode compile(const RefElement &refElement, const LocalCoefficients &localCoefficients)
Definition: compile.hh:50
Definition: code.hh:18
std::size_t size() const
Definition: code.hh:108
Definition: code.hh:152
Definition: localkey.hh:21
unsigned int subEntity() const
Definition: localkey.hh:26
unsigned int index() const
Definition: localkey.hh:28
unsigned int codim() const
Definition: localkey.hh:27