dune-localfunctions 3.0-git
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/fmatrix.hh>
10
12
13namespace Dune
14{
26 template<class D, class R, int dim>
28 {
29 public:
31 typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
32 Dune::FieldMatrix<R,1,dim>, 2> Traits;
33
35 unsigned int size () const
36 {
37 return dim+1;
38 }
39
41 inline void evaluateFunction (const typename Traits::DomainType& in,
42 std::vector<typename Traits::RangeType>& out) const
43 {
44 out.resize(size());
45 out[0] = 1.0;
46 for (size_t i=0; i<dim; i++) {
47 out[0] -= in[i];
48 out[i+1] = in[i];
49 }
50 }
51
53 inline void
54 evaluateJacobian (const typename Traits::DomainType& in, // position
55 std::vector<typename Traits::JacobianType>& out) const // return value
56 {
57 out.resize(size());
58
59 for (int i=0; i<dim; i++)
60 out[0][0][i] = -1;
61
62 for (int i=0; i<dim; i++)
63 for (int j=0; j<dim; j++)
64 out[i+1][0][j] = (i==j);
65
66 }
67
73 inline void partial(const std::array<unsigned int,dim>& order,
74 const typename Traits::DomainType& in,
75 std::vector<typename Traits::RangeType>& out) const
76 {
77 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
78
79 if (totalOrder==0)
80 evaluateFunction(in, out);
81 else if (totalOrder==1)
82 {
83 auto direction = std::find(order.begin(), order.end(), 1);
84 out.resize(size());
85
86 out[0] = -1;
87 for (int i=0; i<dim; i++)
88 out[i+1] = (i==(direction-order.begin()));
89 }
90 else // all higher order derivatives are zero
91 {
92 out.resize(size());
93
94 for (int i=0; i<dim+1; i++)
95 out[i] = 0;
96 }
97 }
98
100 template<unsigned int k>
101 inline void evaluate (const typename std::array<int,k>& directions,
102 const typename Traits::DomainType& in,
103 std::vector<typename Traits::RangeType>& out) const
104 {
105 if (k==0)
106 evaluateFunction(in, out);
107 else if (k==1)
108 {
109 out.resize(size());
110
111 out[0] = -1;
112 for (int i=0; i<dim; i++)
113 out[i+1] = (i==directions[0]);
114 }
115 else if (k==2)
116 {
117 out.resize(size());
118
119 for (int i=0; i<dim+1; i++)
120 out[i] = 0;
121 }
122 }
123
125 unsigned int order () const
126 {
127 return 1;
128 }
129 };
130}
131#endif
Definition brezzidouglasmarini1cube2d.hh:14
Type traits for LocalBasisVirtualInterface.
Definition localbasis.hh:38
D DomainType
domain type
Definition localbasis.hh:49
Linear Lagrange shape functions on the simplex.
Definition p1localbasis.hh:28
unsigned int order() const
Polynomial order of the shape functions.
Definition p1localbasis.hh:125
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition p1localbasis.hh:41
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition p1localbasis.hh:54
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:73
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:101
unsigned int size() const
number of shape functions
Definition p1localbasis.hh:35
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:32