dune-localfunctions  2.4.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 
8 #include <dune/common/fmatrix.hh>
9 
11 
12 namespace Dune
13 {
25  template<class D, class R, int dim>
27  {
28  public:
30  typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
31  Dune::FieldMatrix<R,1,dim>, 2> Traits;
32 
34  unsigned int size () const
35  {
36  return dim+1;
37  }
38 
40  inline void evaluateFunction (const typename Traits::DomainType& in,
41  std::vector<typename Traits::RangeType>& out) const
42  {
43  out.resize(size());
44  out[0] = 1.0;
45  for (size_t i=0; i<dim; i++) {
46  out[0] -= in[i];
47  out[i+1] = in[i];
48  }
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(size());
57 
58  for (int i=0; i<dim; i++)
59  out[0][0][i] = -1;
60 
61  for (int i=0; i<dim; i++)
62  for (int j=0; j<dim; j++)
63  out[i+1][0][j] = (i==j);
64 
65  }
66 
68  template<unsigned int k>
69  inline void evaluate (const typename std::array<int,k>& directions,
70  const typename Traits::DomainType& in,
71  std::vector<typename Traits::RangeType>& out) const
72  {
73  if (k==0)
74  evaluateFunction(in, out);
75  else if (k==1)
76  {
77  out.resize(size());
78 
79  out[0] = -1;
80  for (int i=0; i<dim; i++)
81  out[i+1] = (i==directions[0]);
82  }
83  else if (k==2)
84  {
85  out.resize(size());
86 
87  for (int i=0; i<dim+1; i++)
88  out[i] = 0;
89  }
90  }
91 
93  unsigned int order () const
94  {
95  return 1;
96  }
97  };
98 }
99 #endif
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
Linear Lagrange shape functions on the simplex.
Definition: p1localbasis.hh:26
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p1localbasis.hh:40
unsigned int order() const
Polynomial order of the shape functions.
Definition: p1localbasis.hh:93
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: p1localbasis.hh:53
unsigned int size() const
number of shape functions
Definition: p1localbasis.hh:34
D DomainType
domain type
Definition: localbasis.hh:49
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:31
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:69