3 #ifndef DUNE_LOCALFUNCTIONS_MONOMIAL_MONOMIALLOCALBASIS_HH 4 #define DUNE_LOCALFUNCTIONS_MONOMIAL_MONOMIALLOCALBASIS_HH 9 #include <dune/common/fmatrix.hh> 11 #include "../common/localbasis.hh" 19 template<
int d,
int k>
51 template <
typename Traits>
53 std::vector<typename Traits::RangeType> &out;
55 unsigned int first_unused_index;
59 EvalAccess(std::vector<typename Traits::RangeType> &out_)
62 , first_unused_index(0)
67 assert(first_unused_index == out.size());
70 typename Traits::RangeFieldType &
operator[](
unsigned int index)
72 assert(index < out.size());
74 if(first_unused_index <= index)
75 first_unused_index = index+1;
82 template <
typename Traits>
84 std::vector<typename Traits::JacobianType> &out;
87 unsigned int first_unused_index;
93 : out(out_), row(row_)
95 , first_unused_index(0)
100 assert(first_unused_index == out.size());
103 typename Traits::RangeFieldType &
operator[](
unsigned int index)
105 assert(index < out.size());
107 if(first_unused_index <= index)
108 first_unused_index = index+1;
110 return out[index][0][row];
126 template <
typename Traits,
int c>
131 d = Traits::dimDomain - c
139 template <
typename Access>
141 const typename Traits::DomainType &in,
144 const std::array<int, Traits::dimDomain> &derivatives,
147 typename Traits::RangeFieldType prod,
156 for (
int e = bound; e >= 0; --e)
160 int newbound = bound - e;
161 if(e < derivatives[d])
163 eval(in, derivatives, 0, newbound, index, access);
166 for(
int i = e - derivatives[d] + 1; i <= e; ++i)
174 prod *
ipow(in[d], e-derivatives[d]) * coeff,
190 template <
typename Traits>
193 enum { d = Traits::dimDomain-1 };
195 template <
typename Access>
196 static void eval (
const typename Traits::DomainType &in,
197 const std::array<int, Traits::dimDomain> &derivatives,
198 typename Traits::RangeFieldType prod,
199 int bound,
int& index, Access &access)
201 if(bound < derivatives[d])
205 for(
int i = bound - derivatives[d] + 1; i <= bound; ++i)
207 prod *=
ipow(in[d], bound-derivatives[d]) * coeff;
209 access[index] = prod;
230 template<
class D,
class R,
unsigned int d,
unsigned int p,
unsigned diffOrder = p>
236 Dune::FieldMatrix<R,1,d>,diffOrder>
Traits;
246 std::vector<typename Traits::RangeType>& out)
const 248 evaluate<0>(std::array<int, 0>(), in, out);
252 template<
unsigned int k>
253 inline void evaluate (
const std::array<int,k>& directions,
255 std::vector<typename Traits::RangeType>& out)
const 259 std::array<int, d> derivatives;
260 for(
unsigned int i = 0; i < d; ++i) derivatives[i] = 0;
261 for(
unsigned int i = 0; i < k; ++i) ++derivatives[directions[i]];
263 for(
unsigned int lp = 0; lp <= p; ++lp)
271 std::vector<typename Traits::JacobianType>& out)
const 274 std::array<int, d> derivatives;
275 for(
unsigned int i = 0; i < d; ++i)
277 for(
unsigned int i = 0; i < d; ++i)
282 for(
unsigned int lp = 0; lp <= p; ++lp)
297 #endif // DUNE_LOCALFUNCTIONS_MONOMIAL_MONOMIALLOCALBASIS_HH Access output vector of evaluateFunction() and evaluate()
Definition: monomiallocalbasis.hh:52
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
JacobianAccess(std::vector< typename Traits::JacobianType > &out_, unsigned int row_)
Definition: monomiallocalbasis.hh:91
static void eval(const typename Traits::DomainType &in, const std::array< int, Traits::dimDomain > &derivatives, typename Traits::RangeFieldType prod, int bound, int &index, Access &access)
Definition: monomiallocalbasis.hh:196
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: monomiallocalbasis.hh:270
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: monomiallocalbasis.hh:245
static void eval(const typename Traits::DomainType &in, const std::array< int, Traits::dimDomain > &derivatives, typename Traits::RangeFieldType prod, int bound, int &index, Access &access)
Definition: monomiallocalbasis.hh:140
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
Definition: monomiallocalbasis.hh:20
unsigned int size() const
number of shape functions
Definition: monomiallocalbasis.hh:239
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d >, diffOrder > Traits
export type traits for function signature
Definition: monomiallocalbasis.hh:236
~JacobianAccess()
Definition: monomiallocalbasis.hh:99
Access output vector of evaluateJacobian()
Definition: monomiallocalbasis.hh:83
Definition: monomiallocalbasis.hh:21
void evaluate(const std::array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
return given derivative of all components
Definition: monomiallocalbasis.hh:253
T ipow(T base, int exp)
Definition: monomiallocalbasis.hh:37
EvalAccess(std::vector< typename Traits::RangeType > &out_)
Definition: monomiallocalbasis.hh:59
D DomainType
domain type
Definition: localbasis.hh:49
unsigned int order() const
Polynomial order of the shape functions.
Definition: monomiallocalbasis.hh:289
Definition: monomiallocalbasis.hh:127
~EvalAccess()
Definition: monomiallocalbasis.hh:66
Traits::RangeFieldType & operator[](unsigned int index)
Definition: monomiallocalbasis.hh:70
Traits::RangeFieldType & operator[](unsigned int index)
Definition: monomiallocalbasis.hh:103
Constant shape function.
Definition: monomiallocalbasis.hh:231