dune-localfunctions  2.5-git
prismp1localbasis.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_PRISM_P1_LOCALBASIS_HH
4 #define DUNE_PRISM_P1_LOCALBASIS_HH
5 
6 #include <numeric>
7 
8 #include <dune/common/fmatrix.hh>
9 
11 
12 namespace Dune
13 {
24  template<class D, class R>
26  {
27  public:
29  typedef LocalBasisTraits<D,3,Dune::FieldVector<D,3>,R,1,Dune::FieldVector<R,1>,
30  Dune::FieldMatrix<R,1,3> > Traits;
31 
33  unsigned int size () const
34  {
35  return 6;
36  }
37 
39  inline void evaluateFunction (const typename Traits::DomainType& in,
40  std::vector<typename Traits::RangeType>& out) const
41  {
42  out.resize(6);
43  out[0] = (1.0-in[0]-in[1])*(1.0-in[2]);
44  out[1] = in[0]*(1-in[2]);
45  out[2] = in[1]*(1-in[2]);
46  out[3] = in[2]*(1.0-in[0]-in[1]);
47  out[4] = in[0]*in[2];
48  out[5] = in[1]*in[2];
49  }
50 
52  inline void
53  evaluateJacobian (const typename Traits::DomainType& in, // position
54  std::vector<typename Traits::JacobianType>& out) const // return value
55  {
56  out.resize(6);
57  out[0][0][0] = in[2]-1; out[0][0][1] = in[2]-1; out[0][0][2] = in[0]+in[1]-1; // basis function 0
58  out[1][0][0] = 1-in[2]; out[1][0][1] = 0; out[1][0][2] = -in[0]; // basis function 1
59  out[2][0][0] = 0; out[2][0][1] = 1-in[2]; out[2][0][2] = -in[1]; // basis function 2
60  out[3][0][0] = -in[2]; out[3][0][1] = -in[2]; out[3][0][2] = 1-in[0]-in[1]; // basis function 3
61  out[4][0][0] = in[2]; out[4][0][1] = 0; out[4][0][2] = in[0]; // basis function 4
62  out[5][0][0] = 0; out[5][0][1] = in[2]; out[5][0][2] = in[1]; // basis function 5
63  }
64 
66  void partial (const std::array<unsigned int, 3>& order,
67  const typename Traits::DomainType& in, // position
68  std::vector<typename Traits::RangeType>& out) const // return value
69  {
70  auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
71  if (totalOrder == 0) {
72  evaluateFunction(in, out);
73  } else if (totalOrder == 1) {
74  auto direction = std::distance(order.begin(), std::find(order.begin(), order.end(), 1));
75  out.resize(size());
76 
77  switch (direction) {
78  case 0:
79  out[0] = in[2]-1;
80  out[1] = 1-in[2];
81  out[2] = 0;
82  out[3] = -in[2];
83  out[4] = in[2];
84  out[5] = 0;
85  break;
86  case 1:
87  out[0] = in[2]-1;
88  out[1] = 0;
89  out[2] = 1-in[2];
90  out[3] = -in[2];
91  out[4] = 0;
92  out[5] = in[2];
93  break;
94  case 2:
95  out[0] = in[0]+in[1]-1;
96  out[1] = -in[0];
97  out[2] = -in[1];
98  out[3] = 1-in[0]-in[1];
99  out[4] = in[0];
100  out[5] = in[1];
101  break;
102  default:
103  DUNE_THROW(RangeError, "Component out of range.");
104  }
105  } else if (totalOrder == 2) {
106  out.resize(size());
107  if (order[0] == 1 && order[2] == 1) {
108  out[0] = 1;
109  out[1] =-1;
110  out[2] = 0;
111  out[3] =-1;
112  out[4] = 1;
113  out[5] = 0;
114  } else if (order[1] == 1 && order[2] == 1) {
115  out[0] = 1;
116  out[1] = 0;
117  out[2] =-1;
118  out[3] =-1;
119  out[4] = 0;
120  out[5] = 1;
121  } else {
122  for (std::size_t i = 0; i < size(); ++i)
123  out[i] = 0;
124  }
125  } else {
126  out.resize(size());
127  for (std::size_t i = 0; i < size(); ++i)
128  out[i] = 0;
129  }
130  }
131 
133  unsigned int order () const
134  {
135  return 1;
136  }
137  };
138 }
139 #endif
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 3 > > Traits
export type traits for function signature
Definition: prismp1localbasis.hh:30
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: prismp1localbasis.hh:39
void partial(const std::array< unsigned int, 3 > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of all shape functions.
Definition: prismp1localbasis.hh:66
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: prismp1localbasis.hh:53
D DomainType
domain type
Definition: localbasis.hh:49
unsigned int order() const
Polynomial order of the shape functions.
Definition: prismp1localbasis.hh:133
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
unsigned int size() const
number of shape functions
Definition: prismp1localbasis.hh:33
Linear Lagrange shape functions on the prism.
Definition: prismp1localbasis.hh:25