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