3 #ifndef DUNE_MONOMLOCALBASIS_HH
4 #define DUNE_MONOMLOCALBASIS_HH
8 #include <dune/common/fmatrix.hh>
10 #include "../common/localbasis.hh"
18 template<
int d,
int k>
50 template <
typename Traits>
52 std::vector<typename Traits::RangeType> &out;
54 unsigned int first_unused_index;
58 EvalAccess(std::vector<typename Traits::RangeType> &out_)
61 , first_unused_index(0)
66 assert(first_unused_index == out.size());
69 typename Traits::RangeFieldType &
operator[](
unsigned int index)
71 assert(index < out.size());
73 if(first_unused_index <= index)
74 first_unused_index = index+1;
81 template <
typename Traits>
83 std::vector<typename Traits::JacobianType> &out;
86 unsigned int first_unused_index;
92 : out(out_), row(row_)
94 , first_unused_index(0)
99 assert(first_unused_index == out.size());
102 typename Traits::RangeFieldType &
operator[](
unsigned int index)
104 assert(index < out.size());
106 if(first_unused_index <= index)
107 first_unused_index = index+1;
109 return out[index][0][row];
125 template <
typename Traits,
int c>
130 d = Traits::dimDomain - c
138 template <
typename Access>
140 const typename Traits::DomainType &in,
143 const array<int, Traits::dimDomain> &derivatives,
146 typename Traits::RangeFieldType prod,
155 for (
int e = bound; e >= 0; --e)
159 int newbound = bound - e;
160 if(e < derivatives[
d])
162 eval(in, derivatives, 0, newbound, index, access);
165 for(
int i = e - derivatives[d] + 1; i <= e; ++i)
173 prod *
ipow(in[d], e-derivatives[d]) * coeff,
189 template <
typename Traits>
192 enum {
d = Traits::dimDomain-1 };
194 template <
typename Access>
195 static void eval (
const typename Traits::DomainType &in,
196 const array<int, Traits::dimDomain> &derivatives,
197 typename Traits::RangeFieldType prod,
198 int bound,
int& index, Access &access)
200 if(bound < derivatives[
d])
204 for(
int i = bound - derivatives[d] + 1; i <= bound; ++i)
206 prod *=
ipow(in[d], bound-derivatives[d]) * coeff;
208 access[index] = prod;
229 template<
class D,
class R,
unsigned int d,
unsigned int p,
unsigned diffOrder = p>
237 Dune::FieldMatrix<R,1,d>,diffOrder>
Traits;
247 std::vector<typename Traits::RangeType>& out)
const
249 evaluate<0>(array<int, 0>(), in, out);
253 template<
unsigned int k>
254 inline void evaluate (
const array<int,k>& directions,
256 std::vector<typename Traits::RangeType>& out)
const
260 array<int, d> derivatives;
261 for(
unsigned int i = 0; i < d; ++i) derivatives[i] = 0;
262 for(
unsigned int i = 0; i < k; ++i) ++derivatives[directions[i]];
264 for(
unsigned int lp = 0; lp <= p; ++lp)
272 std::vector<typename Traits::JacobianType>& out)
const
275 array<int, d> derivatives;
276 for(
unsigned int i = 0; i < d; ++i)
278 for(
unsigned int i = 0; i < d; ++i)
283 for(
unsigned int lp = 0; lp <= p; ++lp)
298 #endif // DUNE_MONOMLOCALBASIS_HH
Access output vector of evaluateFunction() and evaluate()
Definition: monomlocalbasis.hh:51
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:39
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: monomlocalbasis.hh:246
JacobianAccess(std::vector< typename Traits::JacobianType > &out_, unsigned int row_)
Definition: monomlocalbasis.hh:90
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: monomlocalbasis.hh:237
Constant shape function.
Definition: monomlocalbasis.hh:230
Definition: monomlocalbasis.hh:19
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: monomlocalbasis.hh:271
T ipow(T base, int exp)
Definition: monomlocalbasis.hh:36
The next dimension to try for factors.
Definition: monomlocalbasis.hh:130
~JacobianAccess()
Definition: monomlocalbasis.hh:98
static void eval(const typename Traits::DomainType &in, const array< int, Traits::dimDomain > &derivatives, typename Traits::RangeFieldType prod, int bound, int &index, Access &access)
Definition: monomlocalbasis.hh:195
Access output vector of evaluateJacobian()
Definition: monomlocalbasis.hh:82
EvalAccess(std::vector< typename Traits::RangeType > &out_)
Definition: monomlocalbasis.hh:58
void evaluate(const array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
return given derivative of all components
Definition: monomlocalbasis.hh:254
D DomainType
domain type
Definition: localbasis.hh:51
Definition: monomlocalbasis.hh:20
Definition: monomlocalbasis.hh:126
static void eval(const typename Traits::DomainType &in, const array< int, Traits::dimDomain > &derivatives, typename Traits::RangeFieldType prod, int bound, int &index, Access &access)
Definition: monomlocalbasis.hh:139
~EvalAccess()
Definition: monomlocalbasis.hh:65
unsigned int order() const
Polynomial order of the shape functions.
Definition: monomlocalbasis.hh:290
Traits::RangeFieldType & operator[](unsigned int index)
Definition: monomlocalbasis.hh:69
unsigned int size() const
number of shape functions
Definition: monomlocalbasis.hh:240
Traits::RangeFieldType & operator[](unsigned int index)
Definition: monomlocalbasis.hh:102