dune-localfunctions  2.8.0
nedelecsimplexinterpolation.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_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXINTERPOLATION_HH
4 #define DUNE_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXINTERPOLATION_HH
5 
6 #include <fstream>
7 #include <utility>
8 #include <numeric>
9 
10 #include <dune/common/exceptions.hh>
11 
12 #include <dune/geometry/quadraturerules.hh>
13 #include <dune/geometry/referenceelements.hh>
14 #include <dune/geometry/type.hh>
15 
20 
21 namespace Dune
22 {
23 
24  // Internal Forward Declarations
25  // -----------------------------
26 
27  template < unsigned int dim, class Field >
28  struct NedelecL2InterpolationFactory;
29 
30 
31 
32  // LocalCoefficientsContainer
33  // --------------------------
34 
36  {
38 
39  public:
40  template <class Setter>
41  LocalCoefficientsContainer ( const Setter &setter )
42  {
43  setter.setLocalKeys(localKey_);
44  }
45 
46  const LocalKey &localKey ( const unsigned int i ) const
47  {
48  assert( i < size() );
49  return localKey_[ i ];
50  }
51 
52  std::size_t size () const
53  {
54  return localKey_.size();
55  }
56 
57  private:
58  std::vector< LocalKey > localKey_;
59  };
60 
61 
62 
63  // NedelecCoefficientsFactory
64  // --------------------------------
65 
66  template < unsigned int dim >
68  {
69  typedef std::size_t Key;
71 
72  template< GeometryType::Id geometryId >
73  static Object *create( const Key &key )
74  {
75  typedef NedelecL2InterpolationFactory< dim, double > InterpolationFactory;
76  if( !supports< geometryId >( key ) )
77  return nullptr;
78  typename InterpolationFactory::Object *interpolation = InterpolationFactory::template create< geometryId >( key );
79  Object *localKeys = new Object( *interpolation );
80  InterpolationFactory::release( interpolation );
81  return localKeys;
82  }
83 
84  template< GeometryType::Id geometryId >
85  static bool supports ( const Key &key )
86  {
87  GeometryType gt = geometryId;
88  return gt.isTriangle() || gt.isTetrahedron() ;
89  }
90  static void release( Object *object ) { delete object; }
91  };
92 
93 
94 
95  // NedelecL2InterpolationBuilder
96  // ------------------------
97 
98  // L2 Interpolation requires:
99  // - for element
100  // - test basis
101  // - for each face (dynamic)
102  // - test basis
103  // - tangents
104  // - for each edge (dynamic)
105  // - test basis
106  // - tangent
107  template< unsigned int dim, class Field >
109  {
110  static const unsigned int dimension = dim;
111 
112  // for the dofs associated to the element
115 
116  // for the dofs associated to the faces
119 
120  // for the dofs associated to the edges
123 
124  // the tangent of the edges
125  typedef FieldVector< Field, dimension > Tangent;
126 
127  // the normal and the tangents of the faces
128  typedef FieldVector< Field, dimension > Normal;
129  typedef std::array<FieldVector< Field, dimension >,dim-1> FaceTangents;
130 
132 
135 
137  {
138  TestBasisFactory::release( testBasis_ );
139  for( FaceStructure &f : faceStructure_ )
140  TestFaceBasisFactory::release( f.basis_ );
141  for( EdgeStructure& e : edgeStructure_ )
142  TestEdgeBasisFactory::release( e.basis_ );
143  }
144 
145  unsigned int topologyId () const
146  {
147  return geometry_.id();
148  }
149 
150  GeometryType type () const
151  {
152  return geometry_;
153  }
154 
155  std::size_t order () const
156  {
157  return order_;
158  }
159 
160  // number of faces
161  unsigned int faceSize () const
162  {
163  return numberOfFaces_;
164  }
165 
166  // number of edges
167  unsigned int edgeSize () const
168  {
169  return numberOfEdges_;
170  }
171 
172  // basis associated to the element
174  {
175  return testBasis_;
176  }
177 
178  // basis associated to face f
179  TestFaceBasis *testFaceBasis ( unsigned int f ) const
180  {
181  assert( f < faceSize() );
182  return faceStructure_[ f ].basis_;
183  }
184 
185  // basis associated to edge e
186  TestEdgeBasis *testEdgeBasis ( unsigned int e ) const
187  {
188  assert( e < edgeSize() );
189  return edgeStructure_[ e ].basis_;
190  }
191 
192  const Tangent& edgeTangent ( unsigned int e ) const
193  {
194  assert( e < edgeSize() );
195  return edgeStructure_[ e ].tangent_;
196  }
197 
198  const FaceTangents& faceTangents ( unsigned int f ) const
199  {
200  assert( f < faceSize() );
201  return faceStructure_[ f ].faceTangents_;
202  }
203 
204  const Normal &normal ( unsigned int f ) const
205  {
206  assert( f < faceSize() );
207  return faceStructure_[ f ].normal_;
208  }
209 
210  template< GeometryType::Id geometryId >
211  void build ( std::size_t order )
212  {
213  constexpr GeometryType geometry = geometryId;
214  order_ = order;
215  geometry_ = geometry;
216 
217  /*
218  * The Nedelec parameter begins at 1.
219  * This is the numbering used by J.C. Nedelec himself.
220  * See "Mixed Finite Elements in \R^3" published in 1980.
221  *
222  * This construction is based on the construction of Raviart-Thomas elements.
223  * There the numbering starts at 0.
224  * Because of this we reduce the order internally by 1.
225  */
226  order--;
227 
228  // if dimension == 2: order-1 on element
229  // if dimension == 3: order-2 on element
230  int requiredOrder = static_cast<int>(dimension==3);
231  testBasis_ = (order > requiredOrder ? TestBasisFactory::template create< geometry >( order-1-requiredOrder ) : nullptr);
232 
233  const auto &refElement = ReferenceElements< Field, dimension >::general( type() );
234 
235  numberOfFaces_ = refElement.size( 1 );
236  faceStructure_.reserve( numberOfFaces_ );
237 
238  // compute the basis, tangents and normals of each face
239  for (std::size_t i=0; i<numberOfFaces_; i++)
240  {
241  FieldVector<Field,dimension> zero(0);
243  faceTangents.fill(zero);
244 
245  // use the first dim-1 vertices of a face to compute the tangents
246  auto vertices = refElement.subEntities(i,1,dim).begin();
247  auto vertex1 = *vertices;
248  for(int j=1; j<dim;j++)
249  {
250  auto vertex2 = vertices[j];
251 
252  faceTangents[j-1] = refElement.position(vertex2,dim) - refElement.position(vertex1,dim);
253 
254  // By default, edges point from the vertex with the smaller index
255  // to the vertex with the larger index.
256  if (vertex1>vertex2)
257  faceTangents[j-1] *=-1;
258 
259  vertex1 = vertex2;
260  }
261 
262  /* For simplices or cubes of arbitrary dimension you could just use
263  *
264  * ```
265  * GeometryType faceGeometry = Impl::getBase(geometry_);
266  * TestFaceBasis *faceBasis = ( dim == 3 && order > 0 ? TestFaceBasisFactory::template create< faceGeometry >( order-1 ) : nullptr);
267  * ```
268  *
269  * For i.e. Prisms and Pyramids in 3d this does not work because they contain squares and triangles as faces.
270  * And depending on the dynamic face index a different face geometry is needed.
271  *
272  */
273  TestFaceBasis *faceBasis = ( dim == 3 && order > 0 ? Impl::IfGeometryType< CreateFaceBasis, dimension-1 >::apply( refElement.type( i, 1 ), order-1 ) : nullptr);
274  faceStructure_.emplace_back( faceBasis, refElement.integrationOuterNormal(i), faceTangents );
275  }
276  assert( faceStructure_.size() == numberOfFaces_ );
277 
278  numberOfEdges_ = refElement.size( dim-1 );
279  edgeStructure_.reserve( numberOfEdges_ );
280 
281  // compute the basis and tangent of each edge
282  for (std::size_t i=0; i<numberOfEdges_; i++)
283  {
284  auto vertexIterator = refElement.subEntities(i,dim-1,dim).begin();
285  auto v0 = *vertexIterator;
286  auto v1 = *(++vertexIterator);
287 
288  // By default, edges point from the vertex with the smaller index
289  // to the vertex with the larger index.
290  if (v0>v1)
291  std::swap(v0,v1);
292  auto tangent = std::move(refElement.position(v1,dim) - refElement.position(v0,dim));
293 
294  TestEdgeBasis *edgeBasis = Impl::IfGeometryType< CreateEdgeBasis, 1 >::apply( refElement.type( i, dim-1 ), order );
295  edgeStructure_.emplace_back( edgeBasis, tangent );
296  }
297  assert( edgeStructure_.size() == numberOfEdges_ );
298  }
299 
300  private:
301 
302  // helper struct for edges
303  struct EdgeStructure
304  {
305  EdgeStructure( TestEdgeBasis *teb, const Tangent &t )
306  : basis_( teb ), tangent_( t )
307  {}
308 
309  TestEdgeBasis *basis_;
310  const Dune::FieldVector< Field, dimension > tangent_;
311  };
312 
313  template< GeometryType::Id edgeGeometryId >
314  struct CreateEdgeBasis
315  {
316  static TestEdgeBasis *apply ( std::size_t order ) { return TestEdgeBasisFactory::template create< edgeGeometryId >( order ); }
317  };
318 
319  // helper struct for faces
320  struct FaceStructure
321  {
322  FaceStructure( TestFaceBasis *tfb, const Normal& normal, const FaceTangents& faceTangents )
323  : basis_( tfb ), normal_(normal), faceTangents_( faceTangents )
324  {}
325 
326  TestFaceBasis *basis_;
327  const Dune::FieldVector< Field, dimension > normal_;
328  const FaceTangents faceTangents_;
329  };
330 
331  template< GeometryType::Id faceGeometryId >
332  struct CreateFaceBasis
333  {
334  static TestFaceBasis *apply ( std::size_t order ) { return TestFaceBasisFactory::template create< faceGeometryId >( order ); }
335  };
336 
337  TestBasis *testBasis_ = nullptr;
338  std::vector< FaceStructure > faceStructure_;
339  unsigned int numberOfFaces_;
340  std::vector< EdgeStructure > edgeStructure_;
341  unsigned int numberOfEdges_;
342  GeometryType geometry_;
343  std::size_t order_;
344  };
345 
346 
347 
348  // NedelecL2Interpolation
349  // ----------------------------
350 
356  template< unsigned int dimension, class F>
358  : public InterpolationHelper< F ,dimension >
359  {
362 
363  public:
364  typedef F Field;
367 
369  : order_(0),
370  size_(0)
371  {}
372 
373  template< class Function, class Vector >
374  auto interpolate ( const Function &function, Vector &coefficients ) const
375  -> std::enable_if_t< std::is_same< decltype(std::declval<Vector>().resize(1) ),void >::value,void>
376  {
377  coefficients.resize(size());
378  typename Base::template Helper<Function,Vector,true> func( function,coefficients );
379  interpolate(func);
380  }
381 
382  template< class Basis, class Matrix >
383  auto interpolate ( const Basis &basis, Matrix &matrix ) const
384  -> std::enable_if_t< std::is_same<
385  decltype(std::declval<Matrix>().rowPtr(0)), typename Matrix::Field* >::value,void>
386  {
387  matrix.resize( size(), basis.size() );
388  typename Base::template Helper<Basis,Matrix,false> func( basis,matrix );
389  interpolate(func);
390  }
391 
392  std::size_t order() const
393  {
394  return order_;
395  }
396  std::size_t size() const
397  {
398  return size_;
399  }
400 
401  template <GeometryType::Id geometryId>
402  void build( std::size_t order )
403  {
404  size_ = 0;
405  order_ = order;
406  builder_.template build<geometryId>(order_);
407  if (builder_.testBasis())
408  size_ += dimension*builder_.testBasis()->size();
409 
410  for ( unsigned int f=0; f<builder_.faceSize(); ++f )
411  if (builder_.testFaceBasis(f))
412  size_ += (dimension-1)*builder_.testFaceBasis(f)->size();
413 
414  for ( unsigned int e=0; e<builder_.edgeSize(); ++e )
415  if (builder_.testEdgeBasis(e))
416  size_ += builder_.testEdgeBasis(e)->size();
417  }
418 
419  void setLocalKeys(std::vector< LocalKey > &keys) const
420  {
421  keys.resize(size());
422  unsigned int row = 0;
423  for (unsigned int e=0; e<builder_.edgeSize(); ++e)
424  {
425  if (builder_.edgeSize())
426  for (unsigned int i=0; i<builder_.testEdgeBasis(e)->size(); ++i,++row)
427  keys[row] = LocalKey(e,dimension-1,i);
428  }
429  for (unsigned int f=0; f<builder_.faceSize(); ++f)
430  {
431  if (builder_.faceSize())
432  for (unsigned int i=0; i<builder_.testFaceBasis(f)->size()*(dimension-1); ++i,++row)
433  keys[row] = LocalKey(f,1,i);
434  }
435 
436  if (builder_.testBasis())
437  for (unsigned int i=0; i<builder_.testBasis()->size()*dimension; ++i,++row)
438  keys[row] = LocalKey(0,0,i);
439  assert( row == size() );
440  }
441 
442  protected:
443  template< class Func, class Container, bool type >
444  void interpolate ( typename Base::template Helper<Func,Container,type> &func ) const
445  {
446  const Dune::GeometryType geoType( builder_.topologyId(), dimension );
447 
448  std::vector<Field> testBasisVal;
449 
450  for (unsigned int i=0; i<size(); ++i)
451  for (unsigned int j=0; j<func.size(); ++j)
452  func.set(i,j,0);
453 
454  unsigned int row = 0;
455 
456  // edge dofs:
457  typedef Dune::QuadratureRule<Field, 1> EdgeQuadrature;
458  typedef Dune::QuadratureRules<Field, 1> EdgeQuadratureRules;
459 
460  const auto &refElement = Dune::ReferenceElements< Field, dimension >::general( geoType );
461 
462  for (unsigned int e=0; e<builder_.edgeSize(); ++e)
463  {
464  if (!builder_.testEdgeBasis(e))
465  continue;
466  testBasisVal.resize(builder_.testEdgeBasis(e)->size());
467 
468  const auto &geometry = refElement.template geometry< dimension-1 >( e );
469  const Dune::GeometryType subGeoType( geometry.type().id(), 1 );
470  const EdgeQuadrature &edgeQuad = EdgeQuadratureRules::rule( subGeoType, 2*order_+2 );
471 
472  const unsigned int quadratureSize = edgeQuad.size();
473  for( unsigned int qi = 0; qi < quadratureSize; ++qi )
474  {
475  if (dimension>1)
476  builder_.testEdgeBasis(e)->template evaluate<0>(edgeQuad[qi].position(),testBasisVal);
477  else
478  testBasisVal[0] = 1.;
479  computeEdgeDofs(row,
480  testBasisVal,
481  func.evaluate( geometry.global( edgeQuad[qi].position() ) ),
483  edgeQuad[qi].weight(),
484  func);
485  }
486 
487  row += builder_.testEdgeBasis(e)->size();
488  }
489 
490  // face dofs:
491  typedef Dune::QuadratureRule<Field, dimension-1> FaceQuadrature;
492  typedef Dune::QuadratureRules<Field, dimension-1> FaceQuadratureRules;
493 
494  for (unsigned int f=0; f<builder_.faceSize(); ++f)
495  {
496  if (builder_.testFaceBasis(f))
497  {
498  testBasisVal.resize(builder_.testFaceBasis(f)->size());
499 
500  const auto &geometry = refElement.template geometry< 1 >( f );
501  const Dune::GeometryType subGeoType( geometry.type().id(), dimension-1 );
502  const FaceQuadrature &faceQuad = FaceQuadratureRules::rule( subGeoType, 2*order_+2 );
503 
504  const unsigned int quadratureSize = faceQuad.size();
505  for( unsigned int qi = 0; qi < quadratureSize; ++qi )
506  {
507  if (dimension>1)
508  builder_.testFaceBasis(f)->template evaluate<0>(faceQuad[qi].position(),testBasisVal);
509  else
510  testBasisVal[0] = 1.;
511 
512  computeFaceDofs( row,
513  testBasisVal,
514  func.evaluate( geometry.global( faceQuad[qi].position() ) ),
516  builder_.normal(f),
517  faceQuad[qi].weight(),
518  func);
519  }
520 
521  row += builder_.testFaceBasis(f)->size()*(dimension-1);
522  }
523  }
524 
525  // element dofs
526  if (builder_.testBasis())
527  {
528  testBasisVal.resize(builder_.testBasis()->size());
529 
530  typedef Dune::QuadratureRule<Field, dimension> Quadrature;
531  typedef Dune::QuadratureRules<Field, dimension> QuadratureRules;
532  const Quadrature &elemQuad = QuadratureRules::rule( geoType, 2*order_+1 );
533 
534  const unsigned int quadratureSize = elemQuad.size();
535  for( unsigned int qi = 0; qi < quadratureSize; ++qi )
536  {
537  builder_.testBasis()->template evaluate<0>(elemQuad[qi].position(),testBasisVal);
538  computeInteriorDofs(row,
539  testBasisVal,
540  func.evaluate(elemQuad[qi].position()),
541  elemQuad[qi].weight(),
542  func );
543  }
544 
545  row += builder_.testBasis()->size()*dimension;
546  }
547  assert(row==size());
548  }
549 
550  private:
560  template <class MVal, class NedVal,class Matrix>
561  void computeEdgeDofs (unsigned int startRow,
562  const MVal &mVal,
563  const NedVal &nedVal,
564  const FieldVector<Field,dimension> &tangent,
565  const Field &weight,
566  Matrix &matrix) const
567  {
568  const unsigned int endRow = startRow+mVal.size();
569  typename NedVal::const_iterator nedIter = nedVal.begin();
570  for ( unsigned int col = 0; col < nedVal.size() ; ++nedIter,++col)
571  {
572  Field cFactor = (*nedIter)*tangent;
573  typename MVal::const_iterator mIter = mVal.begin();
574  for (unsigned int row = startRow; row!=endRow; ++mIter, ++row )
575  matrix.add(row,col, (weight*cFactor)*(*mIter) );
576 
577  assert( mIter == mVal.end() );
578  }
579  }
580 
591  template <class MVal, class NedVal,class Matrix>
592  void computeFaceDofs (unsigned int startRow,
593  const MVal &mVal,
594  const NedVal &nedVal,
595  const FaceTangents& faceTangents,
596  const FieldVector<Field,dimension> &normal,
597  const Field &weight,
598  Matrix &matrix) const
599  {
600  const unsigned int endRow = startRow+mVal.size()*(dimension-1);
601  typename NedVal::const_iterator nedIter = nedVal.begin();
602  for ( unsigned int col = 0; col < nedVal.size() ; ++nedIter,++col)
603  {
604  auto const& u=*nedIter;
605  auto const& n=normal;
606  FieldVector<Field,dimension> nedTimesNormal = { u[1]*n[2]-u[2]*n[1],
607  u[2]*n[0]-u[0]*n[2],
608  u[0]*n[1]-u[1]*n[0]};
609  typename MVal::const_iterator mIter = mVal.begin();
610  for (unsigned int row = startRow; row!=endRow; ++mIter)
611  {
612  for(int i=0; i<dimension-1;i++)
613  {
614  auto test = *mIter*faceTangents[i];
615  matrix.add(row+i,col, weight*(nedTimesNormal*test) );
616  }
617  row += dimension-1;
618  }
619 
620  assert( mIter == mVal.end() );
621  }
622  }
623 
632  template <class MVal, class NedVal,class Matrix>
633  void computeInteriorDofs (unsigned int startRow,
634  const MVal &mVal,
635  const NedVal &nedVal,
636  Field weight,
637  Matrix &matrix) const
638  {
639  const unsigned int endRow = startRow+mVal.size()*dimension;
640  typename NedVal::const_iterator nedIter = nedVal.begin();
641  for ( unsigned int col = 0; col < nedVal.size() ; ++nedIter,++col)
642  {
643  typename MVal::const_iterator mIter = mVal.begin();
644  for (unsigned int row = startRow; row!=endRow; ++mIter,row+=dimension )
645  for (unsigned int i=0; i<dimension; ++i)
646  matrix.add(row+i,col, (weight*(*mIter))*(*nedIter)[i] );
647 
648  assert( mIter == mVal.end() );
649  }
650  }
651 
652  public:
654  std::size_t order_;
655  std::size_t size_;
656  };
657 
658  template < unsigned int dim, class Field >
660  {
663  typedef std::size_t Key;
664  typedef typename std::remove_const<Object>::type NonConstObject;
665 
666  template <GeometryType::Id geometryId>
667  static Object *create( const Key &key )
668  {
669  if ( !supports<geometryId>(key) )
670  return 0;
671  NonConstObject *interpol = new NonConstObject();
672  interpol->template build<geometryId>(key);
673  return interpol;
674  }
675 
676  template <GeometryType::Id geometryId>
677  static bool supports( const Key &key )
678  {
679  GeometryType gt = geometryId;
680  return gt.isTriangle() || gt.isTetrahedron() ;
681  }
682  static void release( Object *object ) { delete object; }
683  };
684 
685 } // namespace Dune
686 
687 #endif // #ifndef DUNE_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXINTERPOLATION_HH
Definition: bdfmcube.hh:16
@ value
Definition: tensor.hh:166
Describe position of one degree of freedom.
Definition: localkey.hh:21
Definition: nedelecsimplexinterpolation.hh:660
static Object * create(const Key &key)
Definition: nedelecsimplexinterpolation.hh:667
const NedelecL2Interpolation< dim, Field > Object
Definition: nedelecsimplexinterpolation.hh:662
NedelecL2InterpolationBuilder< dim, Field > Builder
Definition: nedelecsimplexinterpolation.hh:661
std::size_t Key
Definition: nedelecsimplexinterpolation.hh:663
static bool supports(const Key &key)
Definition: nedelecsimplexinterpolation.hh:677
std::remove_const< Object >::type NonConstObject
Definition: nedelecsimplexinterpolation.hh:664
static void release(Object *object)
Definition: nedelecsimplexinterpolation.hh:682
Definition: nedelecsimplexinterpolation.hh:36
const LocalKey & localKey(const unsigned int i) const
Definition: nedelecsimplexinterpolation.hh:46
LocalCoefficientsContainer(const Setter &setter)
Definition: nedelecsimplexinterpolation.hh:41
std::size_t size() const
Definition: nedelecsimplexinterpolation.hh:52
Definition: nedelecsimplexinterpolation.hh:68
static bool supports(const Key &key)
Definition: nedelecsimplexinterpolation.hh:85
const LocalCoefficientsContainer Object
Definition: nedelecsimplexinterpolation.hh:70
std::size_t Key
Definition: nedelecsimplexinterpolation.hh:69
static Object * create(const Key &key)
Definition: nedelecsimplexinterpolation.hh:73
static void release(Object *object)
Definition: nedelecsimplexinterpolation.hh:90
Definition: nedelecsimplexinterpolation.hh:109
~NedelecL2InterpolationBuilder()
Definition: nedelecsimplexinterpolation.hh:136
TestBasis * testBasis() const
Definition: nedelecsimplexinterpolation.hh:173
GeometryType type() const
Definition: nedelecsimplexinterpolation.hh:150
TestBasisFactory::Object TestBasis
Definition: nedelecsimplexinterpolation.hh:114
FieldVector< Field, dimension > Tangent
Definition: nedelecsimplexinterpolation.hh:125
TestFaceBasisFactory::Object TestFaceBasis
Definition: nedelecsimplexinterpolation.hh:118
TestEdgeBasisFactory::Object TestEdgeBasis
Definition: nedelecsimplexinterpolation.hh:122
FieldVector< Field, dimension > Normal
Definition: nedelecsimplexinterpolation.hh:128
const Tangent & edgeTangent(unsigned int e) const
Definition: nedelecsimplexinterpolation.hh:192
void build(std::size_t order)
Definition: nedelecsimplexinterpolation.hh:211
OrthonormalBasisFactory< dimension, Field > TestBasisFactory
Definition: nedelecsimplexinterpolation.hh:113
TestEdgeBasis * testEdgeBasis(unsigned int e) const
Definition: nedelecsimplexinterpolation.hh:186
OrthonormalBasisFactory< dimension-1, Field > TestFaceBasisFactory
Definition: nedelecsimplexinterpolation.hh:117
unsigned int faceSize() const
Definition: nedelecsimplexinterpolation.hh:161
const FaceTangents & faceTangents(unsigned int f) const
Definition: nedelecsimplexinterpolation.hh:198
std::array< FieldVector< Field, dimension >, dim-1 > FaceTangents
Definition: nedelecsimplexinterpolation.hh:129
OrthonormalBasisFactory< 1, Field > TestEdgeBasisFactory
Definition: nedelecsimplexinterpolation.hh:121
NedelecL2InterpolationBuilder(NedelecL2InterpolationBuilder &&)=delete
TestFaceBasis * testFaceBasis(unsigned int f) const
Definition: nedelecsimplexinterpolation.hh:179
std::size_t order() const
Definition: nedelecsimplexinterpolation.hh:155
unsigned int edgeSize() const
Definition: nedelecsimplexinterpolation.hh:167
const Normal & normal(unsigned int f) const
Definition: nedelecsimplexinterpolation.hh:204
unsigned int topologyId() const
Definition: nedelecsimplexinterpolation.hh:145
NedelecL2InterpolationBuilder(const NedelecL2InterpolationBuilder &)=delete
static const unsigned int dimension
Definition: nedelecsimplexinterpolation.hh:110
An L2-based interpolation for Nedelec.
Definition: nedelecsimplexinterpolation.hh:359
Builder::FaceTangents FaceTangents
Definition: nedelecsimplexinterpolation.hh:366
F Field
Definition: nedelecsimplexinterpolation.hh:364
auto interpolate(const Function &function, Vector &coefficients) const -> std::enable_if_t< std::is_same< decltype(std::declval< Vector >().resize(1)), void >::value, void >
Definition: nedelecsimplexinterpolation.hh:374
std::size_t size() const
Definition: nedelecsimplexinterpolation.hh:396
void interpolate(typename Base::template Helper< Func, Container, type > &func) const
Definition: nedelecsimplexinterpolation.hh:444
std::size_t order_
Definition: nedelecsimplexinterpolation.hh:654
NedelecL2InterpolationBuilder< dimension, Field > Builder
Definition: nedelecsimplexinterpolation.hh:365
auto interpolate(const Basis &basis, Matrix &matrix) const -> std::enable_if_t< std::is_same< decltype(std::declval< Matrix >().rowPtr(0)), typename Matrix::Field * >::value, void >
Definition: nedelecsimplexinterpolation.hh:383
std::size_t size_
Definition: nedelecsimplexinterpolation.hh:655
NedelecL2Interpolation()
Definition: nedelecsimplexinterpolation.hh:368
void build(std::size_t order)
Definition: nedelecsimplexinterpolation.hh:402
std::size_t order() const
Definition: nedelecsimplexinterpolation.hh:392
void setLocalKeys(std::vector< LocalKey > &keys) const
Definition: nedelecsimplexinterpolation.hh:419
Builder builder_
Definition: nedelecsimplexinterpolation.hh:653
Definition: orthonormalbasis.hh:18
static void release(Object *object)
Definition: orthonormalbasis.hh:55
Definition: interpolationhelper.hh:20
Definition: interpolationhelper.hh:22
Definition: polynomialbasis.hh:63
unsigned int size() const
Definition: polynomialbasis.hh:111