dune-pdelab  2.5-dev
l2orthonormal.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 
4 #ifndef DUNE_PDELAB_FINITEELEMENT_L2ORTHONORMAL_HH
5 #define DUNE_PDELAB_FINITEELEMENT_L2ORTHONORMAL_HH
6 
7 #include<array>
8 #include<iostream>
9 #include<algorithm>
10 #include<memory>
11 
12 #include<dune/common/fvector.hh>
13 #include<dune/common/fmatrix.hh>
14 #include<dune/common/gmpfield.hh>
15 #include<dune/common/exceptions.hh>
16 #include<dune/common/fvector.hh>
17 
18 #include<dune/geometry/referenceelements.hh>
19 #include<dune/geometry/quadraturerules.hh>
20 #include<dune/geometry/type.hh>
21 
22 #include<dune/localfunctions/common/localbasis.hh>
23 #include<dune/localfunctions/common/localkey.hh>
24 #include<dune/localfunctions/common/localfiniteelementtraits.hh>
25 
36 namespace Dune {
37  namespace PB {
38 
39 #if HAVE_GMP
40  typedef GMPField<512> DefaultComputationalFieldType;
41 #else
43 #endif
44 
45  //=====================================================
46  // Type to represent a multiindex in d dimensions
47  //=====================================================
48 
49  template<int d>
50  class MultiIndex : public std::array<int,d>
51  {
52  public:
53 
54  MultiIndex () : std::array<int,d>()
55  {
56  }
57 
58  };
59 
60  // the number of polynomials of at most degree k in space dimension d
61  constexpr int pk_size (int k, int d)
62  {
63  if (k==0) return 1;
64  if (d==1) return k+1;
65  int n=0;
66  for (int j=0; j<=k; j++)
67  n += pk_size(k-j,d-1);
68  return n;
69  }
70 
71  // the number of polynomials of exactly degree k in space dimension d (as run-time function)
72  inline int pk_size_exact (int k, int d)
73  {
74  if (k==0)
75  return 1;
76  else
77  return pk_size(k,d)-pk_size(k-1,d);
78  }
79 
80  // k is the polynomial degree and d is the space dimension
81  // then PkSize<k,d> is the number of polynomials of at most total degree k.
82  template<int k, int d>
83  struct PkSize
84  {
85  enum{
87  };
88  };
89 
90  // enumerate all multiindices of degree k and find the i'th
91  template<int d>
92  void pk_enumerate_multiindex (MultiIndex<d>& alpha, int& count, int norm, int dim, int k, int i)
93  {
94  // check if dim is valid
95  if (dim>=d) return;
96 
97  // put all k to current dimension and check
98  alpha[dim]=k; count++; // alpha has correct norm
99  // std::cout << "dadi alpha=" << alpha << " count=" << count << " norm=" << norm+k << " dim=" << dim << " k=" << k << " i=" << i << std::endl;
100  if (count==i) return; // found the index
101 
102  // search recursively
103  for (int m=k-1; m>=0; m--)
104  {
105  alpha[dim]=m;
106  //std::cout << "dada alpha=" << alpha << " count=" << count << " norm=" << norm+m << " dim=" << dim << " k=" << k << " i=" << i << std::endl;
107  pk_enumerate_multiindex(alpha,count,norm+m,dim+1,k-m,i);
108  if (count==i) return;
109  }
110 
111  // reset if index is not yet found
112  alpha[dim]=0;
113  }
114 
115  // map integer 0<=i<pk_size(k,d) to multiindex
116  template<int d>
117  void pk_multiindex (int i, MultiIndex<d>& alpha)
118  {
119  for (int j=0; j<d; j++) alpha[j] = 0; // set alpha to 0
120  if (i==0) return; // handle case i==0 and k==0
121  for (int k=1; k<25; k++)
122  if (i>=pk_size(k-1,d) && i<pk_size(k,d))
123  {
124  int count = -1;
125  pk_enumerate_multiindex(alpha,count,0,0,k,i-pk_size(k-1,d));
126  return;
127  }
128  DUNE_THROW(Dune::NotImplemented,"Polynomial degree greater than 24 in pk_multiindex");
129  }
130 
131  // the number of polynomials of at most degree k in space dimension d (as run-time function)
132  constexpr int qk_size (int k, int d)
133  {
134  int n = 1;
135  for (int i=0; i<d; ++i)
136  n *= (k+1);
137  return n;
138  }
139 
140  // map integer 0<=i<qk_size(k,d) to multiindex
141  template<int d>
142  void qk_multiindex (int i, int k, MultiIndex<d>& alpha)
143  {
144  for (int j = 0; j<d; ++j) {
145  alpha[j] = i % (k+1);
146  i /= (k+1);
147  }
148  }
149 
150  //=====================================================
151  // Traits classes to group Pk and Qk specifics
152  //=====================================================
153  enum BasisType {
155  };
156 
157  template <BasisType basisType>
158  struct BasisTraits;
159 
160  template <>
162  {
163  template <int k, int d>
164  struct Size
165  {
166  enum{
167  value = pk_size(k,d)
168  };
169  };
170 
171  template <int k, int d>
172  struct Order
173  {
174  enum{
175  value = k
176  };
177  };
178 
179  static int size(int k, int d)
180  {
181  return pk_size(k,d);
182  }
183 
184  template <int d>
185  static void multiindex(int i, int k, MultiIndex<d>& alpha)
186  {
187  pk_multiindex(i,alpha);
188  }
189  };
190 
191  template <>
193  {
194  template <int k, int d>
195  struct Size
196  {
197  enum{
198  value = qk_size(k,d)
199  };
200  };
201 
202  template <int k, int d>
203  struct Order
204  {
205  enum{
206  // value = d*k
207  value = k
208  };
209  };
210 
211  static int size(int k, int d)
212  {
213  return qk_size(k,d);
214  }
215 
216  template <int d>
217  static void multiindex(int i, int k, MultiIndex<d>& alpha)
218  {
219  return qk_multiindex(i,k,alpha);
220  }
221  };
222 
223  //=====================================================
224  // Integration kernels for monomials
225  //=====================================================
226 
228  inline long binomial (long n, long k)
229  {
230  // pick the shorter version of
231  // n*(n-1)*...*(k+1)/((n-k)*(n-k-1)*...*1)
232  // and
233  // n*(n-1)*...*(n-k+1)/(k*(k-1)*...*1)
234  if (2*k>=n)
235  {
236  long nominator=1;
237  for (long i=k+1; i<=n; i++) nominator *= i;
238  long denominator=1;
239  for (long i=2; i<=n-k; i++) denominator *= i;
240  return nominator/denominator;
241  }
242  else
243  {
244  long nominator=1;
245  for (long i=n-k+1; i<=n; i++) nominator *= i;
246  long denominator=1;
247  for (long i=2; i<=k; i++) denominator *= i;
248  return nominator/denominator;
249  }
250  }
251 
258  template<typename ComputationFieldType, Dune::GeometryType::BasicType bt, int d>
260  {
261  public:
263  ComputationFieldType integrate (const MultiIndex<d>& a) const
264  {
265  DUNE_THROW(Dune::NotImplemented,"non-specialized version of MonomalIntegrator called. Please implement.");
266  }
267  };
268 
271  template<typename ComputationFieldType, int d>
272  class MonomialIntegrator<ComputationFieldType,Dune::GeometryType::cube,d>
273  {
274  public:
276  ComputationFieldType integrate (const MultiIndex<d>& a) const
277  {
278  ComputationFieldType result(1.0);
279  for (int i=0; i<d; i++)
280  {
281  ComputationFieldType exponent(a[i]+1);
282  result = result/exponent;
283  }
284  return result;
285  }
286  };
287 
290  template<typename ComputationFieldType>
291  class MonomialIntegrator<ComputationFieldType,Dune::GeometryType::simplex,1>
292  {
293  public:
295  ComputationFieldType integrate (const MultiIndex<1>& a) const
296  {
297  ComputationFieldType one(1.0);
298  ComputationFieldType exponent0(a[0]+1);
299  return one/exponent0;
300  }
301  };
302 
305  template<typename ComputationFieldType>
306  class MonomialIntegrator<ComputationFieldType,Dune::GeometryType::simplex,2>
307  {
308  public:
310  ComputationFieldType integrate (const MultiIndex<2>& a) const
311  {
312  ComputationFieldType sum(0.0);
313  for (int k=0; k<=a[1]+1; k++)
314  {
315  int sign=1;
316  if (k%2==1) sign=-1;
317  ComputationFieldType nom(sign*binomial(a[1]+1,k));
318  ComputationFieldType denom(a[0]+k+1);
319  sum = sum + (nom/denom);
320  }
321  ComputationFieldType denom(a[1]+1);
322  return sum/denom;
323  }
324  };
325 
328  template<typename ComputationFieldType>
329  class MonomialIntegrator<ComputationFieldType,Dune::GeometryType::simplex,3>
330  {
331  public:
333  ComputationFieldType integrate (const MultiIndex<3>& a) const
334  {
335  ComputationFieldType sumk(0.0);
336  for (int k=0; k<=a[2]+1; k++)
337  {
338  int sign=1;
339  if (k%2==1) sign=-1;
340  ComputationFieldType nom(sign*binomial(a[2]+1,k));
341  ComputationFieldType denom(a[1]+k+1);
342  sumk = sumk + (nom/denom);
343  }
344  ComputationFieldType sumj(0.0);
345  for (int j=0; j<=a[1]+a[2]+2; j++)
346  {
347  int sign=1;
348  if (j%2==1) sign=-1;
349  ComputationFieldType nom(sign*binomial(a[1]+a[2]+2,j));
350  ComputationFieldType denom(a[0]+j+1);
351  sumj = sumj + (nom/denom);
352  }
353  ComputationFieldType denom(a[2]+1);
354  return sumk*sumj/denom;
355  }
356  };
357 
358  //=====================================================
359  // construct orthonormal basis
360  //=====================================================
361 
363  template<typename F, int d>
365  {
366  template<typename X, typename A>
367  static F compute (const X& x, const A& a)
368  {
369  F prod(1.0);
370  for (int i=0; i<a[d]; i++)
371  prod = prod*x[d];
372  return prod*MonomialEvaluate<F,d-1>::compute(x,a);
373  }
374  };
375 
376  template<typename F>
377  struct MonomialEvaluate<F,0>
378  {
379  template<typename X, typename A>
380  static F compute (const X& x, const A& a)
381  {
382  F prod(1.0);
383  for (int i=0; i<a[0]; i++)
384  prod = prod*x[0];
385  return prod;
386  }
387  };
388 
418  template<typename FieldType, int k, int d, Dune::GeometryType::BasicType bt, typename ComputationFieldType=FieldType, BasisType basisType = BasisType::Pk>
420  {
422  public:
423  enum{ n = Traits::template Size<k,d>::value };
424  typedef Dune::FieldMatrix<FieldType,n,n> LowprecMat;
425  typedef Dune::FieldMatrix<ComputationFieldType,n,n> HighprecMat;
426 
427  // construct orthonormal basis
429  : coeffs(new LowprecMat)
430  {
431  for (int i=0; i<d; ++i)
432  gradcoeffs[i].reset(new LowprecMat());
433  // compute index to multiindex map
434  for (int i=0; i<n; i++)
435  {
436  alpha[i].reset(new MultiIndex<d>());
437  Traits::multiindex(i,k,*alpha[i]);
438  //std::cout << "i=" << i << " alpha_i=" << alpha[i] << std::endl;
439  }
440 
441  orthonormalize();
442  }
443 
444  // construct orthonormal basis from an other basis
445  template<class LFE>
446  OrthonormalPolynomialBasis (const LFE & lfe)
447  : coeffs(new LowprecMat)
448  {
449  for (int i=0; i<d; ++i)
450  gradcoeffs[i].reset(new LowprecMat());
451  // compute index to multiindex map
452  for (int i=0; i<n; i++)
453  {
454  alpha[i].reset(new MultiIndex<d>());
455  Traits::multiindex(i,k,*alpha[i]);
456  //std::cout << "i=" << i << " alpha_i=" << alpha[i] << std::endl;
457  }
458 
459  orthonormalize();
460  }
461 
462  // return dimension of P_l
463  int size (int l)
464  {
465  return Traits::size(l,d);
466  }
467 
468  // evaluate all basis polynomials at given point
469  template<typename Point, typename Result>
470  void evaluateFunction (const Point& x, Result& r) const
471  {
472  std::fill(r.begin(),r.end(),0.0);
473  for (int j=0; j<n; ++j)
474  {
475  const FieldType monomial_value = MonomialEvaluate<FieldType,d-1>::compute(x,*alpha[j]);
476  for (int i=j; i<n; ++i)
477  r[i] += (*coeffs)[i][j] * monomial_value;
478  }
479  }
480 
481  // evaluate all basis polynomials at given point
482  template<typename Point, typename Result>
483  void evaluateJacobian (const Point& x, Result& r) const
484  {
485  std::fill(r.begin(),r.end(),0.0);
486 
487  for (int j=0; j<n; ++j)
488  {
489  const FieldType monomial_value = MonomialEvaluate<FieldType,d-1>::compute(x,*alpha[j]);
490  for (int i=j; i<n; ++i)
491  for (int s=0; s<d; ++s)
492  r[i][0][s] += (*gradcoeffs[s])[i][j]*monomial_value;
493  }
494  }
495 
496  // evaluate all basis polynomials at given point up to order l <= k
497  template<typename Point, typename Result>
498  void evaluateFunction (int l, const Point& x, Result& r) const
499  {
500  if (l>k)
501  DUNE_THROW(Dune::RangeError,"l>k in OrthonormalPolynomialBasis::evaluateFunction");
502 
503  for (int i=0; i<Traits::size(l,d); i++)
504  {
505  FieldType sum(0.0);
506  for (int j=0; j<=i; j++)
507  sum = sum + (*coeffs)[i][j]*MonomialEvaluate<FieldType,d-1>::compute(x,*alpha[j]);
508  r[i] = sum;
509  }
510  }
511 
512  // evaluate all basis polynomials at given point
513  template<typename Point, typename Result>
514  void evaluateJacobian (int l, const Point& x, Result& r) const
515  {
516  if (l>k)
517  DUNE_THROW(Dune::RangeError,"l>k in OrthonormalPolynomialBasis::evaluateFunction");
518 
519  for (int i=0; i<Traits::size(l,d); i++)
520  {
521  FieldType sum[d];
522  for (int s=0; s<d; s++)
523  {
524  sum[s] = 0.0;
525  for (int j=0; j<=i; j++)
526  sum[s] += (*gradcoeffs[s])[i][j]*MonomialEvaluate<FieldType,d-1>::compute(x,*alpha[j]);
527  }
528  for (int s=0; s<d; s++) r[i][0][s] = sum[s];
529  }
530  }
531 
532  private:
533  // store multiindices and coefficients on heap
534  std::array<std::shared_ptr<MultiIndex<d> >,n> alpha; // store index to multiindex map
535  std::shared_ptr<LowprecMat> coeffs; // coefficients with respect to monomials
536  std::array<std::shared_ptr<LowprecMat>,d > gradcoeffs; // coefficients of gradient
537 
538  // compute orthonormalized shapefunctions from a given set of coefficients
539  void orthonormalize()
540  {
541  // run Gram-Schmidt orthogonalization procedure in high precission
542  gram_schmidt();
543 
544  // std::cout << "orthogonal basis monomial representation" << std::endl;
545  // for (int i=0; i<n; i++)
546  // {
547  // std::cout << "phi_" << i << ":" ;
548  // for (int j=0; j<=i; j++)
549  // std::cout << " (" << alpha[j] << "," << coeffs[i][j] << ")";
550  // std::cout << std::endl;
551  // }
552 
553  // compute coefficients of gradient
554  for (int s=0; s<d; s++)
555  for (int i=0; i<n; i++)
556  for (int j=0; j<=i; j++)
557  (*gradcoeffs[s])[i][j] = 0;
558  for (int i=0; i<n; i++)
559  for (int j=0; j<=i; j++)
560  for (int s=0; s<d; s++)
561  if ((*alpha[j])[s]>0)
562  {
563  MultiIndex<d> beta = *alpha[j]; // get exponents
564  FieldType factor = beta[s];
565  beta[s] -= 1;
566  int l = invert_index(beta);
567  (*gradcoeffs[s])[i][l] += (*coeffs)[i][j]*factor;
568  }
569 
570  // for (int s=0; s<d; s++)
571  // {
572  // std::cout << "derivative in direction " << s << std::endl;
573  // for (int i=0; i<n; i++)
574  // {
575  // std::cout << "phi_" << i << ":" ;
576  // for (int j=0; j<=i; j++)
577  // std::cout << " (" << alpha[j] << "," << gradcoeffs[s][i][j] << ")";
578  // std::cout << std::endl;
579  // }
580  // }
581  }
582 
583  // get index from a given multiindex
584  int invert_index (MultiIndex<d>& a)
585  {
586  for (int i=0; i<n; i++)
587  {
588  bool found(true);
589  for (int j=0; j<d; j++)
590  if (a[j]!=(*alpha[i])[j]) found=false;
591  if (found) return i;
592  }
593  DUNE_THROW(Dune::RangeError,"index not found in invertindex");
594  }
595 
596  void gram_schmidt ()
597  {
598  // allocate a high precission matrix on the heap
599  HighprecMat *p = new HighprecMat();
600  HighprecMat& c = *p;
601 
602  // fill identity matrix
603  for (int i=0; i<n; i++)
604  for (int j=0; j<n; j++)
605  if (i==j)
606  c[i][j] = ComputationFieldType(1.0);
607  else
608  c[i][j] = ComputationFieldType(0.0);
609 
610  // the Gram-Schmidt loop
611  MonomialIntegrator<ComputationFieldType,bt,d> integrator;
612  for (int i=0; i<n; i++)
613  {
614  // store orthogonalization coefficients for scaling
615  ComputationFieldType bi[n];
616 
617  // make p_i orthogonal to previous polynomials p_j
618  for (int j=0; j<i; j++)
619  {
620  // p_j is represented with monomials
621  bi[j] = ComputationFieldType(0.0);
622  for (int l=0; l<=j; l++)
623  {
624  MultiIndex<d> a;
625  for (int m=0; m<d; m++) a[m] = (*alpha[i])[m] + (*alpha[l])[m];
626  bi[j] = bi[j] + c[j][l]*integrator.integrate(a);
627  }
628  for (int l=0; l<=j; l++)
629  c[i][l] = c[i][l] - bi[j]*c[j][l];
630  }
631 
632  // scale ith polynomial
633  ComputationFieldType s2(0.0);
634  MultiIndex<d> a;
635  for (int m=0; m<d; m++) a[m] = (*alpha[i])[m] + (*alpha[i])[m];
636  s2 = s2 + integrator.integrate(a);
637  for (int j=0; j<i; j++)
638  s2 = s2 - bi[j]*bi[j];
639  ComputationFieldType s(1.0);
640  using std::sqrt;
641  s = s/sqrt(s2);
642  for (int l=0; l<=i; l++)
643  c[i][l] = s*c[i][l];
644  }
645 
646  // store coefficients in low precission type
647  for (int i=0; i<n; i++)
648  for (int j=0; j<n; j++)
649  (*coeffs)[i][j] = c[i][j];
650 
651  delete p;
652 
653  //std::cout << coeffs << std::endl;
654  }
655  };
656 
657  } // PB namespace
658 
659  // define the local finite element here
660 
661  template<class D, class R, int k, int d, Dune::GeometryType::BasicType bt, typename ComputationFieldType=Dune::PB::DefaultComputationalFieldType, PB::BasisType basisType = PB::BasisType::Pk>
663  {
666  PolynomialBasis opb;
667  Dune::GeometryType gt;
668 
669  public:
670  typedef Dune::LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> > Traits;
671  enum{ n = BasisTraits::template Size<k,d>::value };
672 
673 DUNE_NO_DEPRECATED_BEGIN
674 
675  OPBLocalBasis (int order_) : opb(), gt(bt,d) {}
676 
677  template<class LFE>
678  OPBLocalBasis (int order_, const LFE & lfe) : opb(lfe), gt(bt,d) {}
679 
680 DUNE_NO_DEPRECATED_END
681 
682  unsigned int size () const { return n; }
683 
685  inline void evaluateFunction (const typename Traits::DomainType& in,
686  std::vector<typename Traits::RangeType>& out) const {
687  out.resize(n);
688  opb.evaluateFunction(in,out);
689  }
690 
692  inline void
693  evaluateJacobian (const typename Traits::DomainType& in,
694  std::vector<typename Traits::JacobianType>& out) const {
695  out.resize(n);
696  opb.evaluateJacobian(in,out);
697  }
698 
700  unsigned int order () const {
701  return BasisTraits::template Order<k,d>::value;
702  }
703 
704  Dune::GeometryType type () const { return gt; }
705  };
706 
707  template<int k, int d, PB::BasisType basisType = PB::BasisType::Pk>
709  {
711  public:
712  OPBLocalCoefficients (int order_) : li(n) {
713  for (int i=0; i<n; i++) li[i] = Dune::LocalKey(0,0,i);
714  }
715 
717  std::size_t size () const { return n; }
718 
720  const Dune::LocalKey& localKey (int i) const {
721  return li[i];
722  }
723 
724  private:
725  std::vector<Dune::LocalKey> li;
726  };
727 
728  template<class LB>
730  {
731  const LB& lb;
732 
733  public:
734  OPBLocalInterpolation (const LB& lb_, int order_) : lb(lb_) {}
735 
737  template<typename F, typename C>
738  void interpolate (const F& f, std::vector<C>& out) const
739  {
740  // select quadrature rule
741  typedef typename FieldTraits<typename LB::Traits::RangeFieldType>::real_type RealFieldType;
742 
743  typedef typename LB::Traits::RangeType RangeType;
744  const int d = LB::Traits::dimDomain;
745  const Dune::QuadratureRule<RealFieldType,d>&
746  rule = Dune::QuadratureRules<RealFieldType,d>::rule(lb.type(),2*lb.order());
747 
748  // prepare result
749  out.resize(LB::n);
750  for (int i=0; i<LB::n; i++) out[i] = 0.0;
751 
752  // loop over quadrature points
753  for (typename Dune::QuadratureRule<RealFieldType,d>::const_iterator
754  it=rule.begin(); it!=rule.end(); ++it)
755  {
756  // evaluate function at quadrature point
757  typename LB::Traits::DomainType x;
758  RangeType y;
759  for (int i=0; i<d; i++) x[i] = it->position()[i];
760  f.evaluate(x,y);
761 
762  // evaluate the basis
763  std::vector<RangeType> phi(LB::n);
764  lb.evaluateFunction(it->position(),phi);
765 
766  // do integration
767  for (int i=0; i<LB::n; i++)
768  out[i] += y*phi[i]*it->weight();
769  }
770  }
771  };
772 
773  template<class D, class R, int k, int d, Dune::GeometryType::BasicType bt, typename ComputationFieldType=Dune::PB::DefaultComputationalFieldType, PB::BasisType basisType = PB::BasisType::Pk>
775  {
776  Dune::GeometryType gt;
780  public:
781  typedef Dune::LocalFiniteElementTraits<OPBLocalBasis<D,R,k,d,bt,ComputationFieldType,basisType>,
784 
785 DUNE_NO_DEPRECATED_BEGIN
786 
788  : gt(bt,d), basis(k), coefficients(k), interpolation(basis,k)
789  {}
790 
791  template<class LFE>
792  explicit OPBLocalFiniteElement (const LFE & lfe)
793  : gt(bt,d), basis(k, lfe), coefficients(k), interpolation(basis,k)
794  {}
795 
796 DUNE_NO_DEPRECATED_END
797 
799  : gt(other.gt), basis(other.basis), coefficients(other.coefficients), interpolation(basis,k)
800  {}
801 
802  const typename Traits::LocalBasisType& localBasis () const
803  {
804  return basis;
805  }
806 
807  const typename Traits::LocalCoefficientsType& localCoefficients () const
808  {
809  return coefficients;
810  }
811 
812  const typename Traits::LocalInterpolationType& localInterpolation () const
813  {
814  return interpolation;
815  }
816 
817  Dune::GeometryType type () const { return gt; }
818 
820  return new OPBLocalFiniteElement(*this);
821  }
822  };
823 
824 }
825 
828 #endif // DUNE_PDELAB_FINITEELEMENT_L2ORTHONORMAL_HH
Dune::OPBLocalCoefficients::localKey
const Dune::LocalKey & localKey(int i) const
map index i to local key
Definition: l2orthonormal.hh:720
Dune::PB::BasisType
BasisType
Definition: l2orthonormal.hh:153
Dune::PB::MonomialEvaluate::compute
static F compute(const X &x, const A &a)
Definition: l2orthonormal.hh:367
Dune::OPBLocalBasis
Definition: l2orthonormal.hh:662
Dune::PB::MonomialEvaluate
compute
Definition: l2orthonormal.hh:364
Dune::OPBLocalFiniteElement::OPBLocalFiniteElement
OPBLocalFiniteElement(const LFE &lfe)
Definition: l2orthonormal.hh:792
Dune::PB::Pk
@ Pk
Definition: l2orthonormal.hh:154
Dune::OPBLocalFiniteElement::OPBLocalFiniteElement
DUNE_NO_DEPRECATED_END OPBLocalFiniteElement(const OPBLocalFiniteElement &other)
Definition: l2orthonormal.hh:798
Dune::OPBLocalBasis::evaluateFunction
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: l2orthonormal.hh:685
Dune::PB::OrthonormalPolynomialBasis::OrthonormalPolynomialBasis
OrthonormalPolynomialBasis()
Definition: l2orthonormal.hh:428
dim
static const int dim
Definition: adaptivity.hh:84
Dune::OPBLocalFiniteElement::type
Dune::GeometryType type() const
Definition: l2orthonormal.hh:817
Dune::PB::OrthonormalPolynomialBasis::HighprecMat
Dune::FieldMatrix< ComputationFieldType, n, n > HighprecMat
Definition: l2orthonormal.hh:425
Dune::PB::DefaultComputationalFieldType
double DefaultComputationalFieldType
Definition: l2orthonormal.hh:42
Dune::OPBLocalBasis::type
Dune::GeometryType type() const
Definition: l2orthonormal.hh:704
Dune::OPBLocalFiniteElement::clone
OPBLocalFiniteElement * clone() const
Definition: l2orthonormal.hh:819
Dune::PB::MonomialIntegrator< ComputationFieldType, Dune::GeometryType::simplex, 1 >::integrate
ComputationFieldType integrate(const MultiIndex< 1 > &a) const
integrate one monomial
Definition: l2orthonormal.hh:295
Dune::PB::OrthonormalPolynomialBasis
Integrate monomials over the reference element.
Definition: l2orthonormal.hh:419
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PB::MonomialIntegrator::integrate
ComputationFieldType integrate(const MultiIndex< d > &a) const
integrate one monomial
Definition: l2orthonormal.hh:263
Dune::PB::OrthonormalPolynomialBasis::OrthonormalPolynomialBasis
OrthonormalPolynomialBasis(const LFE &lfe)
Definition: l2orthonormal.hh:446
Dune::OPBLocalFiniteElement::Traits
Dune::LocalFiniteElementTraits< OPBLocalBasis< D, R, k, d, bt, ComputationFieldType, basisType >, OPBLocalCoefficients< k, d, basisType >, OPBLocalInterpolation< OPBLocalBasis< D, R, k, d, bt, ComputationFieldType, basisType > > > Traits
Definition: l2orthonormal.hh:783
Dune::OPBLocalBasis::evaluateJacobian
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: l2orthonormal.hh:693
Dune::PB::OrthonormalPolynomialBasis::size
int size(int l)
Definition: l2orthonormal.hh:463
Dune::PB::MonomialIntegrator< ComputationFieldType, Dune::GeometryType::cube, d >::integrate
ComputationFieldType integrate(const MultiIndex< d > &a) const
integrate one monomial
Definition: l2orthonormal.hh:276
Dune::OPBLocalBasis::size
DUNE_NO_DEPRECATED_END unsigned int size() const
Definition: l2orthonormal.hh:682
Dune::OPBLocalInterpolation::interpolate
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: l2orthonormal.hh:738
Dune::OPBLocalFiniteElement::localInterpolation
const Traits::LocalInterpolationType & localInterpolation() const
Definition: l2orthonormal.hh:812
Dune::PB::OrthonormalPolynomialBasis::evaluateJacobian
void evaluateJacobian(int l, const Point &x, Result &r) const
Definition: l2orthonormal.hh:514
Dune::OPBLocalFiniteElement::OPBLocalFiniteElement
DUNE_NO_DEPRECATED_BEGIN OPBLocalFiniteElement()
Definition: l2orthonormal.hh:787
Dune::PB::MultiIndex
Definition: l2orthonormal.hh:50
Dune::PB::OrthonormalPolynomialBasis::LowprecMat
Dune::FieldMatrix< FieldType, n, n > LowprecMat
Definition: l2orthonormal.hh:424
Dune::PB::pk_size
constexpr int pk_size(int k, int d)
Definition: l2orthonormal.hh:61
Dune::PB::MonomialIntegrator< ComputationFieldType, Dune::GeometryType::simplex, 3 >::integrate
ComputationFieldType integrate(const MultiIndex< 3 > &a) const
integrate one monomial
Definition: l2orthonormal.hh:333
Dune::PB::binomial
long binomial(long n, long k)
compute binomial coefficient "n over k"
Definition: l2orthonormal.hh:228
Dune::QkStuff::multiindex
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdglagrange.hh:117
Dune::PB::MultiIndex::MultiIndex
MultiIndex()
Definition: l2orthonormal.hh:54
Dune::PB::MonomialIntegrator
Integrate monomials over the reference element.
Definition: l2orthonormal.hh:259
Dune::OPBLocalCoefficients
Definition: l2orthonormal.hh:708
Dune::PB::qk_size
constexpr int qk_size(int k, int d)
Definition: l2orthonormal.hh:132
Dune::PB::pk_size_exact
int pk_size_exact(int k, int d)
Definition: l2orthonormal.hh:72
Dune::PB::BasisTraits< BasisType::Qk >::size
static int size(int k, int d)
Definition: l2orthonormal.hh:211
Dune::PB::pk_multiindex
void pk_multiindex(int i, MultiIndex< d > &alpha)
Definition: l2orthonormal.hh:117
Dune::OPBLocalBasis::OPBLocalBasis
OPBLocalBasis(int order_, const LFE &lfe)
Definition: l2orthonormal.hh:678
Dune::PB::BasisTraits< BasisType::Pk >::multiindex
static void multiindex(int i, int k, MultiIndex< d > &alpha)
Definition: l2orthonormal.hh:185
Dune::OPBLocalFiniteElement::localCoefficients
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: l2orthonormal.hh:807
Dune::OPBLocalFiniteElement
Definition: l2orthonormal.hh:774
Dune::OPBLocalCoefficients::size
std::size_t size() const
number of coefficients
Definition: l2orthonormal.hh:717
Dune::PB::OrthonormalPolynomialBasis::evaluateFunction
void evaluateFunction(const Point &x, Result &r) const
Definition: l2orthonormal.hh:470
Dune::OPBLocalInterpolation::OPBLocalInterpolation
OPBLocalInterpolation(const LB &lb_, int order_)
Definition: l2orthonormal.hh:734
Dune::PB::Qk
@ Qk
Definition: l2orthonormal.hh:154
Dune::PB::BasisTraits
Definition: l2orthonormal.hh:158
Dune::OPBLocalBasis::n
@ n
Definition: l2orthonormal.hh:671
Dune::PB::PkSize::value
@ value
Definition: l2orthonormal.hh:86
Dune::OPBLocalBasis::OPBLocalBasis
DUNE_NO_DEPRECATED_BEGIN OPBLocalBasis(int order_)
Definition: l2orthonormal.hh:675
Dune::PB::OrthonormalPolynomialBasis::evaluateFunction
void evaluateFunction(int l, const Point &x, Result &r) const
Definition: l2orthonormal.hh:498
s
const std::string s
Definition: function.hh:830
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PB::OrthonormalPolynomialBasis::n
@ n
Definition: l2orthonormal.hh:423
Dune::PB::PkSize
Definition: l2orthonormal.hh:83
Dune::PB::OrthonormalPolynomialBasis::evaluateJacobian
void evaluateJacobian(const Point &x, Result &r) const
Definition: l2orthonormal.hh:483
Dune::PB::pk_enumerate_multiindex
void pk_enumerate_multiindex(MultiIndex< d > &alpha, int &count, int norm, int dim, int k, int i)
Definition: l2orthonormal.hh:92
Dune::PB::BasisTraits< BasisType::Pk >::size
static int size(int k, int d)
Definition: l2orthonormal.hh:179
Dune::PB::BasisTraits< BasisType::Qk >::multiindex
static void multiindex(int i, int k, MultiIndex< d > &alpha)
Definition: l2orthonormal.hh:217
Dune::OPBLocalBasis::order
unsigned int order() const
Polynomial order of the shape functions.
Definition: l2orthonormal.hh:700
Dune::OPBLocalCoefficients::OPBLocalCoefficients
OPBLocalCoefficients(int order_)
Definition: l2orthonormal.hh:712
Dune::PB::MonomialEvaluate< F, 0 >::compute
static F compute(const X &x, const A &a)
Definition: l2orthonormal.hh:380
Dune::OPBLocalFiniteElement::localBasis
const Traits::LocalBasisType & localBasis() const
Definition: l2orthonormal.hh:802
Dune::PB::qk_multiindex
void qk_multiindex(int i, int k, MultiIndex< d > &alpha)
Definition: l2orthonormal.hh:142
Dune::PB::MonomialIntegrator< ComputationFieldType, Dune::GeometryType::simplex, 2 >::integrate
ComputationFieldType integrate(const MultiIndex< 2 > &a) const
integrate one monomial
Definition: l2orthonormal.hh:310
Dune::OPBLocalInterpolation
Definition: l2orthonormal.hh:729
p
const P & p
Definition: constraints.hh:147
Dune::OPBLocalBasis::Traits
Dune::LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: l2orthonormal.hh:670