dune-localfunctions  2.2.1
p1localbasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_P1_LOCALBASIS_HH
3 #define DUNE_P1_LOCALBASIS_HH
4 
5 #include <dune/common/fmatrix.hh>
6 
8 
9 namespace Dune
10 {
22  template<class D, class R, int dim>
23  class P1LocalBasis
24  {
25  public:
27  typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
28  Dune::FieldMatrix<R,1,dim>, 2> Traits;
29 
31  unsigned int size () const
32  {
33  return dim+1;
34  }
35 
37  inline void evaluateFunction (const typename Traits::DomainType& in,
38  std::vector<typename Traits::RangeType>& out) const
39  {
40  out.resize(size());
41  out[0] = 1.0;
42  for (size_t i=0; i<dim; i++) {
43  out[0] -= in[i];
44  out[i+1] = in[i];
45  }
46  }
47 
49  inline void
50  evaluateJacobian (const typename Traits::DomainType& in, // position
51  std::vector<typename Traits::JacobianType>& out) const // return value
52  {
53  out.resize(size());
54 
55  for (int i=0; i<dim; i++)
56  out[0][0][i] = -1;
57 
58  for (int i=0; i<dim; i++)
59  for (int j=0; j<dim; j++)
60  out[i+1][0][j] = (i==j);
61 
62  }
63 
65  template<unsigned int k>
66  inline void evaluate (const typename Dune::array<int,k>& directions,
67  const typename Traits::DomainType& in,
68  std::vector<typename Traits::RangeType>& out) const
69  {
70  if (k==0)
71  evaluateFunction(in, out);
72  else if (k==1)
73  {
74  out.resize(size());
75 
76  out[0] = -1;
77  for (int i=0; i<dim; i++)
78  out[i+1] = (i==directions[0]);
79  }
80  else if (k==2)
81  {
82  out.resize(size());
83 
84  for (int i=0; i<dim+1; i++)
85  out[i] = 0;
86  }
87  }
88 
90  unsigned int order () const
91  {
92  return 1;
93  }
94  };
95 }
96 #endif