dune-localfunctions  2.5.1
p1localbasis.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_P1_LOCALBASIS_HH
4 #define DUNE_P1_LOCALBASIS_HH
5 
6 #include <array>
7 #include <numeric>
8 
9 #include <dune/common/deprecated.hh>
10 #include <dune/common/fmatrix.hh>
11 
13 
14 namespace Dune
15 {
27  template<class D, class R, int dim>
29  {
30  public:
32  typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
33  Dune::FieldMatrix<R,1,dim>, 2> Traits;
34 
36  unsigned int size () const
37  {
38  return dim+1;
39  }
40 
42  inline void evaluateFunction (const typename Traits::DomainType& in,
43  std::vector<typename Traits::RangeType>& out) const
44  {
45  out.resize(size());
46  out[0] = 1.0;
47  for (size_t i=0; i<dim; i++) {
48  out[0] -= in[i];
49  out[i+1] = in[i];
50  }
51  }
52 
54  inline void
55  evaluateJacobian (const typename Traits::DomainType& in, // position
56  std::vector<typename Traits::JacobianType>& out) const // return value
57  {
58  out.resize(size());
59 
60  for (int i=0; i<dim; i++)
61  out[0][0][i] = -1;
62 
63  for (int i=0; i<dim; i++)
64  for (int j=0; j<dim; j++)
65  out[i+1][0][j] = (i==j);
66 
67  }
68 
74  inline void partial(const std::array<unsigned int,dim>& order,
75  const typename Traits::DomainType& in,
76  std::vector<typename Traits::RangeType>& out) const
77  {
78  auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
79 
80  if (totalOrder==0)
81  evaluateFunction(in, out);
82  else if (totalOrder==1)
83  {
84  auto direction = std::find(order.begin(), order.end(), 1);
85  out.resize(size());
86 
87  out[0] = -1;
88  for (int i=0; i<dim; i++)
89  out[i+1] = (i==(direction-order.begin()));
90  }
91  else // all higher order derivatives are zero
92  {
93  out.resize(size());
94 
95  for (int i=0; i<dim+1; i++)
96  out[i] = 0;
97  }
98  }
99 
101  template<unsigned int k>
102  inline void DUNE_DEPRECATED_MSG("Use method 'partial' instead!")
103  evaluate (const typename std::array<int,k>& directions,
104  const typename Traits::DomainType& in,
105  std::vector<typename Traits::RangeType>& out) const
106  {
107  if (k==0)
108  evaluateFunction(in, out);
109  else if (k==1)
110  {
111  out.resize(size());
112 
113  out[0] = -1;
114  for (int i=0; i<dim; i++)
115  out[i+1] = (i==directions[0]);
116  }
117  else if (k==2)
118  {
119  out.resize(size());
120 
121  for (int i=0; i<dim+1; i++)
122  out[i] = 0;
123  }
124  }
125 
127  unsigned int order () const
128  {
129  return 1;
130  }
131  };
132 }
133 #endif
unsigned int order() const
Polynomial order of the shape functions.
Definition: p1localbasis.hh:127
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
Linear Lagrange shape functions on the simplex.
Definition: p1localbasis.hh:28
unsigned int size() const
number of shape functions
Definition: p1localbasis.hh:36
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p1localbasis.hh:42
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: p1localbasis.hh:55
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim >, 2 > Traits
export type traits for function signature
Definition: p1localbasis.hh:33
STL namespace.
void partial(const std::array< unsigned int, dim > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of any order of all shape functions.
Definition: p1localbasis.hh:74
void evaluate(const typename std::array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p1localbasis.hh:103
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37