dune-localfunctions  2.7.1
monomialbasis.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_MONOMIALBASIS_HH
4 #define DUNE_MONOMIALBASIS_HH
5 
6 #include <vector>
7 
8 #include <dune/common/fvector.hh>
9 #include <dune/common/fmatrix.hh>
10 
11 #include <dune/geometry/topologyfactory.hh>
12 #include <dune/geometry/type.hh>
13 
17 
18 namespace Dune
19 {
20  /************************************************
21  * Classes for evaluating ''Monomials'' on any order
22  * for all reference element type.
23  * For a simplex topology these are the normal
24  * monomials for cube topologies the bimonomials.
25  * The construction follows the construction of the
26  * generic geometries using tensor products for
27  * prism generation and duffy transform for pyramid
28  * construction.
29  * A derivative argument can be applied, in which case
30  * all derivatives up to the desired order are
31  * evaluated. Note that for higher order derivatives
32  * only the ''lower'' part of the symmetric tensor
33  * is evaluated, e.g., passing derivative equal to 2
34  * to the class will provide the vector
35  * (d/dxdx p, d/dxydx p, d/dydy p,
36  * d/dx p, d/dy p, p)
37  * Important:
38  * So far the computation of the derivatives has not
39  * been fully implemented for general pyramid
40  * construction, i.e., in the case where a pyramid is
41  * build over a non simplex base geometry.
42  *
43  * Central classes:
44  * 1) template< class Topology, class F >
45  * class MonomialBasisImpl;
46  * Implementation of the monomial evaluation for
47  * a given topology and field type.
48  * The method evaluate fills a F* vector
49  * 2) template< class Topology, class F >
50  * class MonomialBasis
51  * The base class for the static monomial evaluation
52  * providing addiional evaluate methods including
53  * one taking std::vector<F>.
54  * 3) template< int dim, class F >
55  * class VirtualMonomialBasis
56  * Virtualization of the MonomialBasis.
57  * 4) template< int dim, class F >
58  * struct MonomialBasisFactory;
59  * A factory class for the VirtualMonomialBasis
60  * 5) template< int dim, class F >
61  * struct MonomialBasisProvider
62  * A singleton container for the virtual monomial
63  * basis
64  ************************************************/
65 
66  // Internal Forward Declarations
67  // -----------------------------
68 
69  template< class Topology >
70  class MonomialBasisSize;
71 
72  template< class Topology, class F >
73  class MonomialBasis;
74 
75 
76 
77  // MonomialBasisSize
78  // -----------------
79 
80  template< class TopologyType >
82  {
84 
85  public:
86  static This &instance ()
87  {
88  static This _instance;
89  return _instance;
90  }
91 
92  unsigned int maxOrder_;
93 
94  // sizes_[ k ]: number of basis functions of exactly order k
95  mutable unsigned int *sizes_;
96 
97  // numBaseFunctions_[ k ] = sizes_[ 0 ] + ... + sizes_[ k ]
98  mutable unsigned int *numBaseFunctions_;
99 
101  : maxOrder_( 0 ),
102  sizes_( 0 ),
103  numBaseFunctions_( 0 )
104  {
105  computeSizes( 2 );
106  }
107 
109  {
110  delete[] sizes_;
111  delete[] numBaseFunctions_;
112  }
113 
114  unsigned int operator() ( const unsigned int order ) const
115  {
116  return numBaseFunctions_[ order ];
117  }
118 
119  unsigned int maxOrder() const
120  {
121  return maxOrder_;
122  }
123 
124  void computeSizes ( unsigned int order )
125  {
126  if (order <= maxOrder_)
127  return;
128 
129  maxOrder_ = order;
130 
131  delete[] sizes_;
132  delete[] numBaseFunctions_;
133  sizes_ = new unsigned int[ order+1 ];
134  numBaseFunctions_ = new unsigned int[ order+1 ];
135 
136  constexpr auto dim = TopologyType::dimension;
137 
138  sizes_[ 0 ] = 1;
139  for( unsigned int k = 1; k <= order; ++k )
140  sizes_[ k ] = 0;
141 
142  std::fill(numBaseFunctions_, numBaseFunctions_+order+1, 1);
143 
144  for( int codim=dim-1; codim>=0; codim--)
145  {
146  if (Impl::isPrism(TopologyType::id,dim,codim))
147  {
148  for( unsigned int k = 1; k <= order; ++k )
149  {
150  sizes_[ k ] = numBaseFunctions_[ k ] + k*sizes_[ k ];
151  numBaseFunctions_[ k ] = numBaseFunctions_[ k-1 ] + sizes_[ k ];
152  }
153  }
154  else
155  {
156  for( unsigned int k = 1; k <= order; ++k )
157  {
158  sizes_[ k ] = numBaseFunctions_[ k ];
159  numBaseFunctions_[ k ] = numBaseFunctions_[ k-1 ] + sizes_[ k ];
160  }
161  }
162  }
163  }
164  };
165 
166 
167 
168  // MonomialBasisHelper
169  // -------------------
170 
171 
172  template< int mydim, int dim, class F >
174  {
177 
178  static void copy ( const unsigned int deriv, F *&wit, F *&rit,
179  const unsigned int numBaseFunctions, const F &z )
180  {
181  // n(d,k) = size<k>[d];
182  MySize &mySize = MySize::instance();
183  Size &size = Size::instance();
184 
185  const F *const rend = rit + size( deriv )*numBaseFunctions;
186  for( ; rit != rend; )
187  {
188  F *prit = rit;
189 
190  *wit = z * *rit;
191  ++rit, ++wit;
192 
193  for( unsigned d = 1; d <= deriv; ++d )
194  {
195  #ifndef NDEBUG
196  const F *const derivEnd = rit + mySize.sizes_[ d ];
197  #endif
198 
199  {
200  const F *const drend = rit + mySize.sizes_[ d ] - mySize.sizes_[ d-1 ];
201  for( ; rit != drend ; ++rit, ++wit )
202  *wit = z * *rit;
203  }
204 
205  for (unsigned int j=1; j<d; ++j)
206  {
207  const F *const drend = rit + mySize.sizes_[ d-j ] - mySize.sizes_[ d-j-1 ];
208  for( ; rit != drend ; ++prit, ++rit, ++wit )
209  *wit = F(j) * *prit + z * *rit;
210  }
211  *wit = F(d) * *prit + z * *rit;
212  ++prit, ++rit, ++wit;
213  assert(derivEnd == rit);
214  rit += size.sizes_[d] - mySize.sizes_[d];
215  prit += size.sizes_[d-1] - mySize.sizes_[d-1];
216  const F *const emptyWitEnd = wit + size.sizes_[d] - mySize.sizes_[d];
217  for ( ; wit != emptyWitEnd; ++wit )
218  *wit = Zero<F>();
219  }
220  }
221  }
222  };
223 
224 
225 
226  // MonomialBasisImpl
227  // -----------------
228 
229  template< class Topology, class F >
231 
232  template< class F >
233  class MonomialBasisImpl< Impl::Point, F >
234  {
236 
237  public:
238  typedef Impl::Point Topology;
239  typedef F Field;
240 
241  static const unsigned int dimDomain = Topology::dimension;
242 
243  typedef FieldVector< Field, dimDomain > DomainVector;
244 
245  private:
246  friend class MonomialBasis< Topology, Field >;
247  friend class MonomialBasisImpl< Impl::Prism< Topology >, Field >;
248  friend class MonomialBasisImpl< Impl::Pyramid< Topology >, Field >;
249 
250  template< int dimD >
251  void evaluate ( const unsigned int deriv, const unsigned int order,
252  const FieldVector< Field, dimD > &x,
253  const unsigned int block, const unsigned int *const offsets,
254  Field *const values ) const
255  {
256  *values = Unity< F >();
257  F *const end = values + block;
258  for( Field *it = values+1 ; it != end; ++it )
259  *it = Zero< F >();
260  }
261 
262  void integrate ( const unsigned int order,
263  const unsigned int *const offsets,
264  Field *const values ) const
265  {
266  values[ 0 ] = Unity< Field >();
267  }
268  };
269 
270  template< class BaseTopology, class F >
271  class MonomialBasisImpl< Impl::Prism< BaseTopology >, F >
272  {
274 
275  public:
276  typedef Impl::Prism< BaseTopology > Topology;
277  typedef F Field;
278 
279  static const unsigned int dimDomain = Topology::dimension;
280 
281  typedef FieldVector< Field, dimDomain > DomainVector;
282 
283  private:
284  friend class MonomialBasis< Topology, Field >;
285  friend class MonomialBasisImpl< Impl::Prism< Topology >, Field >;
286  friend class MonomialBasisImpl< Impl::Pyramid< Topology >, Field >;
287 
288  typedef MonomialBasisSize< BaseTopology > BaseSize;
289  typedef MonomialBasisSize< Topology > Size;
290 
291  MonomialBasisImpl< BaseTopology, Field > baseBasis_;
292 
294  {}
295 
296  template< int dimD >
297  void evaluate ( const unsigned int deriv, const unsigned int order,
298  const FieldVector< Field, dimD > &x,
299  const unsigned int block, const unsigned int *const offsets,
300  Field *const values ) const
301  {
303  const BaseSize &size = BaseSize::instance();
304  const_cast<BaseSize&>(size).computeSizes(order);
305 
306  const Field &z = x[ dimDomain-1 ];
307 
308  // fill first column
309  baseBasis_.evaluate( deriv, order, x, block, offsets, values );
310 
311  Field *row0 = values;
312  for( unsigned int k = 1; k <= order; ++k )
313  {
314  Field *row1 = values + block*offsets[ k-1 ];
315  Field *wit = row1 + block*size.sizes_[ k ];
316  Helper::copy( deriv, wit, row1, k*size.sizes_[ k ], z );
317  Helper::copy( deriv, wit, row0, size( k-1 ), z );
318  row0 = row1;
319  }
320  }
321 
322  void integrate ( const unsigned int order,
323  const unsigned int *const offsets,
324  Field *const values ) const
325  {
326  const BaseSize &size = BaseSize::instance();
327  const Size &mySize = Size::instance();
328  // fill first column
329  baseBasis_.integrate( order, offsets, values );
330  const unsigned int *const baseSizes = size.sizes_;
331 
332  Field *row0 = values;
333  for( unsigned int k = 1; k <= order; ++k )
334  {
335  Field *const row1begin = values + offsets[ k-1 ];
336  Field *const row1End = row1begin + mySize.sizes_[ k ];
337  assert( (unsigned int)(row1End - values) <= offsets[ k ] );
338 
339  Field *row1 = row1begin;
340  Field *it = row1begin + baseSizes[ k ];
341  for( unsigned int j = 1; j <= k; ++j )
342  {
343  Field *const end = it + baseSizes[ k ];
344  assert( (unsigned int)(end - values) <= offsets[ k ] );
345  for( ; it != end; ++row1, ++it )
346  *it = (Field( j ) / Field( j+1 )) * (*row1);
347  }
348  for( ; it != row1End; ++row0, ++it )
349  *it = (Field( k ) / Field( k+1 )) * (*row0);
350  row0 = row1;
351  }
352  }
353  };
354 
355  template< class BaseTopology, class F >
356  class MonomialBasisImpl< Impl::Pyramid< BaseTopology >, F >
357  {
359 
360  public:
361  typedef Impl::Pyramid< BaseTopology > Topology;
362  typedef F Field;
363 
364  static const unsigned int dimDomain = Topology::dimension;
365 
366  typedef FieldVector< Field, dimDomain > DomainVector;
367 
368  private:
369  friend class MonomialBasis< Topology, Field >;
370  friend class MonomialBasisImpl< Impl::Prism< Topology >, Field >;
371  friend class MonomialBasisImpl< Impl::Pyramid< Topology >, Field >;
372 
373  typedef MonomialBasisSize< BaseTopology > BaseSize;
374  typedef MonomialBasisSize< Topology > Size;
375 
376  MonomialBasisImpl< BaseTopology, Field > baseBasis_;
377 
379  {}
380 
381  template< int dimD >
382  void evaluateSimplexBase ( const unsigned int deriv, const unsigned int order,
383  const FieldVector< Field, dimD > &x,
384  const unsigned int block, const unsigned int *const offsets,
385  Field *const values,
386  const BaseSize &size ) const
387  {
388  baseBasis_.evaluate( deriv, order, x, block, offsets, values );
389  }
390 
391  template< int dimD >
392  void evaluatePyramidBase ( const unsigned int deriv, const unsigned int order,
393  const FieldVector< Field, dimD > &x,
394  const unsigned int block, const unsigned int *const offsets,
395  Field *const values,
396  const BaseSize &size ) const
397  {
398  Field omz = Unity< Field >() - x[ dimDomain-1 ];
399 
400  if( Zero< Field >() < omz )
401  {
402  const Field invomz = Unity< Field >() / omz;
403  FieldVector< Field, dimD > y;
404  for( unsigned int i = 0; i < dimDomain-1; ++i )
405  y[ i ] = x[ i ] * invomz;
406 
407  // fill first column
408  baseBasis_.evaluate( deriv, order, y, block, offsets, values );
409 
410  Field omzk = omz;
411  for( unsigned int k = 1; k <= order; ++k )
412  {
413  Field *it = values + block*offsets[ k-1 ];
414  Field *const end = it + block*size.sizes_[ k ];
415  for( ; it != end; ++it )
416  *it *= omzk;
417  omzk *= omz;
418  }
419  }
420  else
421  {
422  assert( deriv==0 );
423  *values = Unity< Field >();
424  for( unsigned int k = 1; k <= order; ++k )
425  {
426  Field *it = values + block*offsets[ k-1 ];
427  Field *const end = it + block*size.sizes_[ k ];
428  for( ; it != end; ++it )
429  *it = Zero< Field >();
430  }
431  }
432  }
433 
434  template< int dimD >
435  void evaluate ( const unsigned int deriv, const unsigned int order,
436  const FieldVector< Field, dimD > &x,
437  const unsigned int block, const unsigned int *const offsets,
438  Field *const values ) const
439  {
440  typedef MonomialBasisHelper< dimDomain, dimD, Field > Helper;
441  const BaseSize &size = BaseSize::instance();
442  const_cast<BaseSize&>(size).computeSizes(order);
443 
445  evaluateSimplexBase( deriv, order, x, block, offsets, values, size );
446  else
447  evaluatePyramidBase( deriv, order, x, block, offsets, values, size );
448 
449  Field *row0 = values;
450  for( unsigned int k = 1; k <= order; ++k )
451  {
452  Field *row1 = values + block*offsets[ k-1 ];
453  Field *wit = row1 + block*size.sizes_[ k ];
454  Helper::copy( deriv, wit, row0, size( k-1 ), x[ dimDomain-1 ] );
455  row0 = row1;
456  }
457  }
458 
459  void integrate ( const unsigned int order,
460  const unsigned int *const offsets,
461  Field *const values ) const
462  {
463  const BaseSize &size = BaseSize::instance();
464 
465  // fill first column
466  baseBasis_.integrate( order, offsets, values );
467 
468  const unsigned int *const baseSizes = size.sizes_;
469 
470  {
471  Field *const col0End = values + baseSizes[ 0 ];
472  for( Field *it = values; it != col0End; ++it )
473  *it *= Field( 1 ) / Field( int(dimDomain) );
474  }
475 
476  Field *row0 = values;
477  for( unsigned int k = 1; k <= order; ++k )
478  {
479  const Field factor = (Field( 1 ) / Field( k + dimDomain ));
480 
481  Field *const row1 = values+offsets[ k-1 ];
482  Field *const col0End = row1 + baseSizes[ k ];
483  Field *it = row1;
484  for( ; it != col0End; ++it )
485  *it *= factor;
486  for( unsigned int i = 1; i <= k; ++i )
487  {
488  Field *const end = it + baseSizes[ k-i ];
489  assert( (unsigned int)(end - values) <= offsets[ k ] );
490  for( ; it != end; ++row0, ++it )
491  *it = (*row0) * (Field( i ) * factor);
492  }
493  row0 = row1;
494  }
495  }
496  };
497 
498 
499 
500  // MonomialBasis
501  // -------------
502 
503  template< class Topology, class F >
505  : public MonomialBasisImpl< Topology, F >
506  {
509 
510  public:
511  static const unsigned int dimension = Base::dimDomain;
512  static const unsigned int dimRange = 1;
513 
514  typedef typename Base::Field Field;
515 
516  typedef typename Base::DomainVector DomainVector;
517 
518  typedef Dune::FieldVector<Field,dimRange> RangeVector;
519 
521 
522  MonomialBasis (unsigned int order)
523  : Base(),
524  order_(order),
525  size_(Size::instance())
526  {
527  assert(order<=1024); // avoid wrapping of unsigned int (0-1) order=1024 is quite hight...)
528  }
529 
530  const unsigned int *sizes ( unsigned int order ) const
531  {
532  size_.computeSizes( order );
533  return size_.numBaseFunctions_;
534  }
535 
536  const unsigned int *sizes () const
537  {
538  return sizes( order_ );
539  }
540 
541  unsigned int size () const
542  {
543  size_.computeSizes( order_ );
544  return size_( order_ );
545  }
546 
547  unsigned int derivSize ( const unsigned int deriv ) const
548  {
549  typedef typename Impl::SimplexTopology< dimension >::type SimplexTopology;
550  MonomialBasisSize< SimplexTopology >::instance().computeSizes( deriv );
552  }
553 
554  unsigned int order () const
555  {
556  return order_ ;
557  }
558 
559  unsigned int topologyId ( ) const
560  {
561  return Topology::id;
562  }
563 
564  void evaluate ( const unsigned int deriv, const DomainVector &x,
565  Field *const values ) const
566  {
567  Base::evaluate( deriv, order_, x, derivSize( deriv ), sizes( order_ ), values );
568  }
569 
570  template <unsigned int deriv>
571  void evaluate ( const DomainVector &x,
572  Field *const values ) const
573  {
574  evaluate( deriv, x, values );
575  }
576 
577  template<unsigned int deriv, class Vector >
578  void evaluate ( const DomainVector &x,
579  Vector &values ) const
580  {
581  evaluate<deriv>(x,&(values[0]));
582  }
583  template<unsigned int deriv, DerivativeLayoutNS::DerivativeLayout layout >
584  void evaluate ( const DomainVector &x,
586  {
587  evaluate<deriv>(x,&(values->block()));
588  }
589  template< unsigned int deriv >
590  void evaluate ( const DomainVector &x,
592  {
593  evaluate(0,x,&(values[0][0]));
594  }
595 
596  template<class Vector >
597  void evaluate ( const DomainVector &x,
598  Vector &values ) const
599  {
600  evaluate<0>(x,&(values[0]));
601  }
602 
603  template< class DVector, class RVector >
604  void evaluate ( const DVector &x, RVector &values ) const
605  {
606  assert( DVector::dimension == dimension);
607  DomainVector bx;
608  for( int d = 0; d < dimension; ++d )
609  field_cast( x[ d ], bx[ d ] );
610  evaluate<0>( bx, values );
611  }
612 
613  void integrate ( Field *const values ) const
614  {
615  Base::integrate( order_, sizes( order_ ), values );
616  }
617  template <class Vector>
618  void integrate ( Vector &values ) const
619  {
620  integrate( &(values[ 0 ]) );
621  }
622  private:
623  MonomialBasis(const This&);
624  This& operator=(const This&);
625  unsigned int order_;
626  Size &size_;
627  };
628 
629 
630 
631  // StdMonomialBasis
632  // ----------------
633 
634  template< int dim,class F >
636  : public MonomialBasis< typename Impl::SimplexTopology< dim >::type, F >
637  {
640 
641  public:
642  typedef typename Impl::SimplexTopology< dim >::type Topology;
643  static const int dimension = dim;
644 
645  StandardMonomialBasis ( unsigned int order )
646  : Base( order )
647  {}
648  };
649 
650 
651 
652  // StandardBiMonomialBasis
653  // -----------------------
654 
655  template< int dim, class F >
657  : public MonomialBasis< typename Impl::CubeTopology< dim >::type, F >
658  {
661 
662  public:
663  typedef typename Impl::CubeTopology< dim >::type Topology;
664  static const int dimension = dim;
665 
667  : Base( order )
668  {}
669  };
670 
671  // -----------------------------------------------------------
672  // -----------------------------------------------------------
673  // VirtualMonomialBasis
674  // -------------------
675 
676  template< int dim, class F >
678  {
680 
681  public:
682  typedef F Field;
683  typedef F StorageField;
684  static const int dimension = dim;
685  static const unsigned int dimRange = 1;
686 
687  typedef FieldVector<Field,dimension> DomainVector;
688  typedef FieldVector<Field,dimRange> RangeVector;
689 
690  explicit VirtualMonomialBasis(unsigned int topologyId,
691  unsigned int order)
693 
695 
696  virtual const unsigned int *sizes ( ) const = 0;
697 
698  unsigned int size ( ) const
699  {
700  return sizes( )[ order_ ];
701  }
702 
703  unsigned int order () const
704  {
705  return order_;
706  }
707 
708  unsigned int topologyId ( ) const
709  {
710  return topologyId_;
711  }
712 
713  virtual void evaluate ( const unsigned int deriv, const DomainVector &x,
714  Field *const values ) const = 0;
715  template < unsigned int deriv >
716  void evaluate ( const DomainVector &x,
717  Field *const values ) const
718  {
719  evaluate( deriv, x, values );
720  }
721  template < unsigned int deriv, int size >
722  void evaluate ( const DomainVector &x,
723  Dune::FieldVector<Field,size> *const values ) const
724  {
725  evaluate( deriv, x, &(values[0][0]) );
726  }
727  template<unsigned int deriv, DerivativeLayoutNS::DerivativeLayout layout >
728  void evaluate ( const DomainVector &x,
730  {
731  evaluate<deriv>(x,&(values->block()));
732  }
733  template <unsigned int deriv, class Vector>
734  void evaluate ( const DomainVector &x,
735  Vector &values ) const
736  {
737  evaluate<deriv>( x, &(values[ 0 ]) );
738  }
739  template< class Vector >
740  void evaluate ( const DomainVector &x,
741  Vector &values ) const
742  {
743  evaluate<0>(x,values);
744  }
745  template< class DVector, class RVector >
746  void evaluate ( const DVector &x, RVector &values ) const
747  {
748  assert( DVector::dimension == dimension);
749  DomainVector bx;
750  for( int d = 0; d < dimension; ++d )
751  field_cast( x[ d ], bx[ d ] );
752  evaluate<0>( bx, values );
753  }
754  template< unsigned int deriv, class DVector, class RVector >
755  void evaluate ( const DVector &x, RVector &values ) const
756  {
757  assert( DVector::dimension == dimension);
758  DomainVector bx;
759  for( int d = 0; d < dimension; ++d )
760  field_cast( x[ d ], bx[ d ] );
761  evaluate<deriv>( bx, values );
762  }
763 
764  virtual void integrate ( Field *const values ) const = 0;
765  template <class Vector>
766  void integrate ( Vector &values ) const
767  {
768  integrate( &(values[ 0 ]) );
769  }
770  protected:
771  unsigned int order_;
772  unsigned int topologyId_;
773  };
774 
775  template< class Topology, class F >
777  : public VirtualMonomialBasis< Topology::dimension, F >
778  {
781 
782  public:
783  typedef typename Base::Field Field;
785 
787  : Base(Topology::id,order), basis_(order)
788  {}
789 
790  const unsigned int *sizes ( ) const
791  {
792  return basis_.sizes(order_);
793  }
794 
795  void evaluate ( const unsigned int deriv, const DomainVector &x,
796  Field *const values ) const
797  {
798  basis_.evaluate(deriv,x,values);
799  }
800 
801  void integrate ( Field *const values ) const
802  {
803  basis_.integrate(values);
804  }
805 
806  private:
808  using Base::order_;
809  };
810 
811  // MonomialBasisFactory
812  // --------------------
813 
814  template< int dim, class F >
816  {
817  static const unsigned int dimension = dim;
818  typedef F StorageField;
819 
820  typedef unsigned int Key;
822 
823  template < int dd, class FF >
825  {
827  };
828 
829  template< class Topology >
830  static Object* create ( const Key &order )
831  {
833  }
834  static void release( Object *object ) { delete object; }
835  };
836 
837 
838 
839  // MonomialBasisProvider
840  // ---------------------
841 
842  template< int dim, class SF >
844  : public TopologySingletonFactory< MonomialBasisFactory< dim, SF > >
845  {
846  static const unsigned int dimension = dim;
847  typedef SF StorageField;
848  template < int dd, class FF >
850  {
852  };
853  };
854 
855 }
856 
857 #endif
Definition: bdfmcube.hh:16
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition: field.hh:157
@ value
Definition: tensor.hh:166
A class representing the unit of a given Field.
Definition: field.hh:28
A class representing the zero of a given Field.
Definition: field.hh:77
Definition: monomialbasis.hh:82
unsigned int * sizes_
Definition: monomialbasis.hh:95
~MonomialBasisSize()
Definition: monomialbasis.hh:108
unsigned int * numBaseFunctions_
Definition: monomialbasis.hh:98
MonomialBasisSize()
Definition: monomialbasis.hh:100
static This & instance()
Definition: monomialbasis.hh:86
void computeSizes(unsigned int order)
Definition: monomialbasis.hh:124
unsigned int operator()(const unsigned int order) const
Definition: monomialbasis.hh:114
unsigned int maxOrder() const
Definition: monomialbasis.hh:119
unsigned int maxOrder_
Definition: monomialbasis.hh:92
Definition: monomialbasis.hh:506
Base::Field Field
Definition: monomialbasis.hh:514
static const unsigned int dimension
Definition: monomialbasis.hh:511
unsigned int size() const
Definition: monomialbasis.hh:541
unsigned int derivSize(const unsigned int deriv) const
Definition: monomialbasis.hh:547
void evaluate(const DomainVector &x, Derivatives< Field, dimension, 1, deriv, layout > *values) const
Definition: monomialbasis.hh:584
void evaluate(const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:571
void integrate(Vector &values) const
Definition: monomialbasis.hh:618
unsigned int topologyId() const
Definition: monomialbasis.hh:559
Dune::FieldVector< Field, dimRange > RangeVector
Definition: monomialbasis.hh:518
void evaluate(const DomainVector &x, FieldVector< Field, Derivatives< Field, dimension, 1, deriv, DerivativeLayoutNS::value >::size > *values) const
Definition: monomialbasis.hh:590
unsigned int order() const
Definition: monomialbasis.hh:554
MonomialBasis(unsigned int order)
Definition: monomialbasis.hh:522
void evaluate(const DomainVector &x, Vector &values) const
Definition: monomialbasis.hh:597
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:564
Base::DomainVector DomainVector
Definition: monomialbasis.hh:516
void integrate(Field *const values) const
Definition: monomialbasis.hh:613
static const unsigned int dimRange
Definition: monomialbasis.hh:512
const unsigned int * sizes() const
Definition: monomialbasis.hh:536
MonomialBasisSize< Topology > Size
Definition: monomialbasis.hh:520
void evaluate(const DVector &x, RVector &values) const
Definition: monomialbasis.hh:604
const unsigned int * sizes(unsigned int order) const
Definition: monomialbasis.hh:530
void evaluate(const DomainVector &x, Vector &values) const
Definition: monomialbasis.hh:578
Definition: monomialbasis.hh:174
MonomialBasisSize< typename Impl::SimplexTopology< mydim >::type > MySize
Definition: monomialbasis.hh:175
MonomialBasisSize< typename Impl::SimplexTopology< dim >::type > Size
Definition: monomialbasis.hh:176
static void copy(const unsigned int deriv, F *&wit, F *&rit, const unsigned int numBaseFunctions, const F &z)
Definition: monomialbasis.hh:178
Definition: monomialbasis.hh:230
Definition: monomialbasis.hh:234
F Field
Definition: monomialbasis.hh:239
Impl::Point Topology
Definition: monomialbasis.hh:238
FieldVector< Field, dimDomain > DomainVector
Definition: monomialbasis.hh:243
Impl::Prism< BaseTopology > Topology
Definition: monomialbasis.hh:276
FieldVector< Field, dimDomain > DomainVector
Definition: monomialbasis.hh:281
Impl::Pyramid< BaseTopology > Topology
Definition: monomialbasis.hh:361
FieldVector< Field, dimDomain > DomainVector
Definition: monomialbasis.hh:366
Definition: monomialbasis.hh:637
StandardMonomialBasis(unsigned int order)
Definition: monomialbasis.hh:645
Impl::SimplexTopology< dim >::type Topology
Definition: monomialbasis.hh:642
static const int dimension
Definition: monomialbasis.hh:643
Definition: monomialbasis.hh:658
static const int dimension
Definition: monomialbasis.hh:664
Impl::CubeTopology< dim >::type Topology
Definition: monomialbasis.hh:663
StandardBiMonomialBasis(unsigned int order)
Definition: monomialbasis.hh:666
Definition: monomialbasis.hh:678
unsigned int topologyId() const
Definition: monomialbasis.hh:708
FieldVector< Field, dimension > DomainVector
Definition: monomialbasis.hh:687
unsigned int order_
Definition: monomialbasis.hh:771
void evaluate(const DomainVector &x, Vector &values) const
Definition: monomialbasis.hh:740
void evaluate(const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:716
F Field
Definition: monomialbasis.hh:682
void evaluate(const DomainVector &x, Vector &values) const
Definition: monomialbasis.hh:734
void evaluate(const DVector &x, RVector &values) const
Definition: monomialbasis.hh:746
unsigned int order() const
Definition: monomialbasis.hh:703
virtual const unsigned int * sizes() const =0
unsigned int topologyId_
Definition: monomialbasis.hh:772
static const unsigned int dimRange
Definition: monomialbasis.hh:685
F StorageField
Definition: monomialbasis.hh:683
static const int dimension
Definition: monomialbasis.hh:684
VirtualMonomialBasis(unsigned int topologyId, unsigned int order)
Definition: monomialbasis.hh:690
unsigned int size() const
Definition: monomialbasis.hh:698
FieldVector< Field, dimRange > RangeVector
Definition: monomialbasis.hh:688
virtual ~VirtualMonomialBasis()
Definition: monomialbasis.hh:694
virtual void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const =0
virtual void integrate(Field *const values) const =0
void evaluate(const DomainVector &x, Dune::FieldVector< Field, size > *const values) const
Definition: monomialbasis.hh:722
void evaluate(const DVector &x, RVector &values) const
Definition: monomialbasis.hh:755
void evaluate(const DomainVector &x, Derivatives< Field, dimension, 1, deriv, layout > *values) const
Definition: monomialbasis.hh:728
void integrate(Vector &values) const
Definition: monomialbasis.hh:766
Definition: monomialbasis.hh:778
VirtualMonomialBasisImpl(unsigned int order)
Definition: monomialbasis.hh:786
Base::DomainVector DomainVector
Definition: monomialbasis.hh:784
void integrate(Field *const values) const
Definition: monomialbasis.hh:801
const unsigned int * sizes() const
Definition: monomialbasis.hh:790
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:795
Base::Field Field
Definition: monomialbasis.hh:783
Definition: monomialbasis.hh:816
static void release(Object *object)
Definition: monomialbasis.hh:834
const VirtualMonomialBasis< dimension, F > Object
Definition: monomialbasis.hh:821
F StorageField
Definition: monomialbasis.hh:818
static const unsigned int dimension
Definition: monomialbasis.hh:817
static Object * create(const Key &order)
Definition: monomialbasis.hh:830
unsigned int Key
Definition: monomialbasis.hh:820
MonomialBasisFactory< dd, FF > Type
Definition: monomialbasis.hh:826
Definition: monomialbasis.hh:845
static const unsigned int dimension
Definition: monomialbasis.hh:846
SF StorageField
Definition: monomialbasis.hh:847
MonomialBasisProvider< dd, FF > Type
Definition: monomialbasis.hh:851
Definition: tensor.hh:170