dune-localfunctions  2.4.1-rc2
q1localbasis.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_Q1_LOCALBASIS_HH
4 #define DUNE_Q1_LOCALBASIS_HH
5 
6 #include <dune/common/fmatrix.hh>
7 
9 
10 namespace Dune
11 {
23  template<class D, class R, int dim>
25  {
26  public:
27  typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
28  Dune::FieldMatrix<R,1,dim> > Traits;
29 
31  unsigned int size () const
32  {
33  return 1<<dim;
34  }
35 
37  inline void evaluateFunction (const typename Traits::DomainType& in,
38  std::vector<typename Traits::RangeType>& out) const
39  {
40  out.resize(size());
41 
42  for (size_t i=0; i<size(); i++) {
43 
44  out[i] = 1;
45 
46  for (int j=0; j<dim; j++)
47  // if j-th bit of i is set multiply with in[j], else with 1-in[j]
48  out[i] *= (i & (1<<j)) ? in[j] : 1-in[j];
49 
50  }
51 
52  }
53 
55  inline void
56  evaluateJacobian (const typename Traits::DomainType& in, // position
57  std::vector<typename Traits::JacobianType>& out) const // return value
58  {
59  out.resize(size());
60 
61  // Loop over all shape functions
62  for (size_t i=0; i<size(); i++) {
63 
64  // Loop over all coordinate directions
65  for (int j=0; j<dim; j++) {
66 
67  // Initialize: the overall expression is a product
68  // if j-th bit of i is set to -1, else 1
69  out[i][0][j] = (i & (1<<j)) ? 1 : -1;
70 
71  for (int k=0; k<dim; k++) {
72 
73  if (j!=k)
74  // if k-th bit of i is set multiply with in[j], else with 1-in[j]
75  out[i][0][j] *= (i & (1<<k)) ? in[k] : 1-in[k];
76 
77  }
78 
79  }
80 
81  }
82 
83  }
84 
86  unsigned int order () const
87  {
88  return 1;
89  }
90  };
91 }
92 #endif
unsigned int order() const
Polynomial order of the shape functions.
Definition: q1localbasis.hh:86
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
unsigned int size() const
number of shape functions
Definition: q1localbasis.hh:31
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim > > Traits
Definition: q1localbasis.hh:28
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
Lagrange shape functions of order 1 on the reference cube.
Definition: q1localbasis.hh:24
D DomainType
domain type
Definition: localbasis.hh:49
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: q1localbasis.hh:56
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: q1localbasis.hh:37