dune-localfunctions  2.4.1
dualq1localbasis.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_DUAL_Q1_LOCALBASIS_HH
4 #define DUNE_DUAL_Q1_LOCALBASIS_HH
5 
6 #include <array>
7 
8 #include <dune/common/fvector.hh>
9 #include <dune/common/fmatrix.hh>
10 
12 
13 namespace Dune
14 {
24  template<class D, class R, int dim>
26  {
27  public:
28  typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
29  Dune::FieldMatrix<R,1,dim> > Traits;
30 
31  void setCoefficients(const std::array<Dune::FieldVector<R, (1<<dim)> ,(1<<dim)>& coefficients)
32  {
33  coefficients_ = coefficients;
34  }
35 
37  unsigned int size () const
38  {
39  return 1<<dim;
40  }
41 
43  inline void evaluateFunction (const typename Traits::DomainType& in,
44  std::vector<typename Traits::RangeType>& out) const
45  {
46  // compute q1 values
47  std::vector<typename Traits::RangeType> q1Values(size());
48 
49  for (size_t i=0; i<size(); i++) {
50 
51  q1Values[i] = 1;
52 
53  for (int j=0; j<dim; j++)
54  // if j-th bit of i is set multiply with in[j], else with 1-in[j]
55  q1Values[i] *= (i & (1<<j)) ? in[j] : 1-in[j];
56 
57  }
58 
59  // compute the dual values by using that they are linear combinations of q1 functions
60  out.resize(size());
61  for (size_t i=0; i<size(); i++)
62  out[i] = 0;
63 
64  for (size_t i=0; i<size(); i++)
65  for (size_t j=0; j<size(); j++)
66  out[i] += coefficients_[i][j]*q1Values[j];
67 
68 
69  }
70 
72  inline void
73  evaluateJacobian (const typename Traits::DomainType& in, // position
74  std::vector<typename Traits::JacobianType>& out) const // return value
75  {
76  // compute q1 jacobians
77  std::vector<typename Traits::JacobianType> q1Jacs(size());
78 
79  // Loop over all shape functions
80  for (size_t i=0; i<size(); i++) {
81 
82  // Loop over all coordinate directions
83  for (int j=0; j<dim; j++) {
84 
85  // Initialize: the overall expression is a product
86  // if j-th bit of i is set to -1, else 1
87  q1Jacs[i][0][j] = (i & (1<<j)) ? 1 : -1;
88 
89  for (int k=0; k<dim; k++) {
90 
91  if (j!=k)
92  // if k-th bit of i is set multiply with in[j], else with 1-in[j]
93  q1Jacs[i][0][j] *= (i & (1<<k)) ? in[k] : 1-in[k];
94 
95  }
96 
97  }
98 
99  }
100 
101  // compute the dual jacobians by using that they are linear combinations of q1 functions
102  out.resize(size());
103  for (size_t i=0; i<size(); i++)
104  out[i] = 0;
105 
106  for (size_t i=0; i<size(); i++)
107  for (size_t j=0; j<size(); j++)
108  out[i].axpy(coefficients_[i][j],q1Jacs[j]);
109 
110  }
111 
113  unsigned int order () const
114  {
115  return 1;
116  }
117 
118  private:
119  std::array<Dune::FieldVector<R, (1<<dim)> ,(1<<dim)> coefficients_;
120  };
121 }
122 #endif
unsigned int order() const
Polynomial order of the shape functions.
Definition: dualq1localbasis.hh:113
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: dualq1localbasis.hh:43
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: dualq1localbasis.hh:73
unsigned int size() const
number of shape functions
Definition: dualq1localbasis.hh:37
void setCoefficients(const std::array< Dune::FieldVector< R,(1<< dim)>,(1<< dim)> &coefficients)
Definition: dualq1localbasis.hh:31
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
D DomainType
domain type
Definition: localbasis.hh:49
Dual Lagrange shape functions of order 1 on the reference cube.
Definition: dualq1localbasis.hh:25
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim > > Traits
Definition: dualq1localbasis.hh:29