dune-localfunctions  2.5-git
brezzidouglasmarini1simplex2dlocalbasis.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_BREZZIDOUGLASMARINI1_SIMPLEX2D_LOCALBASIS_HH
4 #define DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1_SIMPLEX2D_LOCALBASIS_HH
5 
6 #include <array>
7 #include <bitset>
8 #include <numeric>
9 #include <vector>
10 
11 #include <dune/common/fmatrix.hh>
12 
13 #include "../../common/localbasis.hh"
14 
15 namespace Dune
16 {
26  template<class D, class R>
28  {
29 
30  public:
31  typedef LocalBasisTraits<D,2,Dune::FieldVector<D,2>,R,2,Dune::FieldVector<R,2>,
32  Dune::FieldMatrix<R,2,2> > Traits;
33 
36  {
37  for (size_t i=0; i<3; i++)
38  sign_[i] = 1.0;
39  }
40 
46  BDM1Simplex2DLocalBasis (std::bitset<3> s)
47  {
48  for (size_t i=0; i<3; i++)
49  sign_[i] = s[i] ? -1.0 : 1.0;
50  }
51 
53  unsigned int size () const
54  {
55  return 6;
56  }
57 
64  inline void evaluateFunction (const typename Traits::DomainType& in,
65  std::vector<typename Traits::RangeType>& out) const
66  {
67  out.resize(6);
68 
69  out[0][0] = sign_[0]*in[0];
70  out[0][1] = sign_[0]*(in[1] - 1.0);
71  out[1][0] = sign_[1]*(in[0] - 1.0);
72  out[1][1] = sign_[1]*in[1];
73  out[2][0] = sign_[2]*in[0];
74  out[2][1] = sign_[2]*in[1];
75  out[3][0] = 3.0*in[0];
76  out[3][1] = 3.0 - 6.0*in[0] - 3.0*in[1];
77  out[4][0] = -3.0 + 3.0*in[0] + 6.0*in[1];
78  out[4][1] = -3.0*in[1];
79  out[5][0] = -3.0*in[0];
80  out[5][1] = 3.0*in[1];
81  }
82 
89  inline void evaluateJacobian (const typename Traits::DomainType& in,
90  std::vector<typename Traits::JacobianType>& out) const
91  {
92  out.resize(6);
93 
94  out[0][0][0] = sign_[0];
95  out[0][0][1] = 0.0;
96  out[0][1][0] = 0.0;
97  out[0][1][1] = sign_[0];
98 
99  out[1][0][0] = sign_[1];
100  out[1][0][1] = 0.0;
101  out[1][1][0] = 0.0;
102  out[1][1][1] = sign_[1];
103 
104  out[2][0][0] = sign_[2];
105  out[2][0][1] = 0.0;
106  out[2][1][0] = 0.0;
107  out[2][1][1] = sign_[2];
108 
109  out[3][0][0] = 3.0;
110  out[3][0][1] = 0.0;
111  out[3][1][0] = -6.0;
112  out[3][1][1] = -3.0;
113 
114  out[4][0][0] = 3.0;
115  out[4][0][1] = 6.0;
116  out[4][1][0] = 0.0;
117  out[4][1][1] = -3.0;
118 
119  out[5][0][0] = -3.0;
120  out[5][0][1] = 0.0;
121  out[5][1][0] = 0.0;
122  out[5][1][1] = 3.0;
123  }
124 
126  void partial (const std::array<unsigned int, 2>& order,
127  const typename Traits::DomainType& in, // position
128  std::vector<typename Traits::RangeType>& out) const // return value
129  {
130  auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
131  if (totalOrder == 0) {
132  evaluateFunction(in, out);
133  } else if (totalOrder == 1) {
134  out.resize(size());
135  auto const direction = std::distance(order.begin(), std::find(order.begin(), order.end(), 1));
136 
137  switch (direction) {
138  case 0:
139  out[0][0] = sign_[0];
140  out[0][1] = 0.0;
141 
142  out[1][0] = sign_[1];
143  out[1][1] = 0.0;
144 
145  out[2][0] = sign_[2];
146  out[2][1] = 0.0;
147 
148  out[3][0] = 3.0;
149  out[3][1] = -6.0;
150 
151  out[4][0] = 3.0;
152  out[4][1] = 0.0;
153 
154  out[5][0] = -3.0;
155  out[5][1] = 0.0;
156  break;
157  case 1:
158  out[0][0] = 0.0;
159  out[0][1] = sign_[0];
160 
161  out[1][0] = 0.0;
162  out[1][1] = sign_[1];
163 
164  out[2][0] = 0.0;
165  out[2][1] = sign_[2];
166 
167  out[3][0] = 0.0;
168  out[3][1] = -3.0;
169 
170  out[4][0] = 6.0;
171  out[4][1] = -3.0;
172 
173  out[5][0] = 0.0;
174  out[5][1] = 3.0;
175  break;
176  default:
177  DUNE_THROW(RangeError, "Component out of range.");
178  }
179  } else {
180  DUNE_THROW(NotImplemented, "Desired derivative order is not implemented");
181  }
182  }
183 
185  unsigned int order () const
186  {
187  return 1;
188  }
189 
190  private:
191  std::array<R,3> sign_;
192  };
193 }
194 #endif // DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1_SIMPLEX2D_LOCALBASIS_HH
BDM1Simplex2DLocalBasis(std::bitset< 3 > s)
Make set number s, where 0 <= s < 8.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:46
LocalBasisTraits< D, 2, Dune::FieldVector< D, 2 >, R, 2, Dune::FieldVector< R, 2 >, Dune::FieldMatrix< R, 2, 2 > > Traits
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:32
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:64
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
BDM1Simplex2DLocalBasis()
Standard constructor.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:35
D DomainType
domain type
Definition: localbasis.hh:49
unsigned int size() const
number of shape functions
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:53
First order Brezzi-Douglas-Marini shape functions on the reference triangle.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:27
unsigned int order() const
Polynomial order of the shape functions.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:185
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
void partial(const std::array< unsigned int, 2 > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of all shape functions.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:126
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: brezzidouglasmarini1simplex2dlocalbasis.hh:89