3#ifndef DUNE_LOCALFUNCTIONS_MONOMIAL_MONOMIALLOCALBASIS_HH
4#define DUNE_LOCALFUNCTIONS_MONOMIAL_MONOMIALLOCALBASIS_HH
10#include <dune/common/fmatrix.hh>
12#include "../common/localbasis.hh"
20 template<
int d,
int k>
52 template <
typename Traits>
54 std::vector<typename Traits::RangeType> &out;
56 unsigned int first_unused_index;
60 EvalAccess(std::vector<typename Traits::RangeType> &out_)
63 , first_unused_index(0)
68 assert(first_unused_index == out.size());
71 typename Traits::RangeFieldType &
operator[](
unsigned int index)
73 assert(index < out.size());
75 if(first_unused_index <= index)
76 first_unused_index = index+1;
83 template <
typename Traits>
85 std::vector<typename Traits::JacobianType> &out;
88 unsigned int first_unused_index;
94 : out(out_), row(row_)
96 , first_unused_index(0)
101 assert(first_unused_index == out.size());
104 typename Traits::RangeFieldType &
operator[](
unsigned int index)
106 assert(index < out.size());
108 if(first_unused_index <= index)
109 first_unused_index = index+1;
111 return out[index][0][row];
127 template <
typename Traits,
int c>
132 d = Traits::dimDomain - c
140 template <
typename Access>
142 const typename Traits::DomainType &in,
145 const std::array<int, Traits::dimDomain> &derivatives,
148 typename Traits::RangeFieldType prod,
157 for (
int e = bound; e >= 0; --e)
161 int newbound = bound - e;
162 if(e < derivatives[
d])
164 eval(in, derivatives, 0, newbound, index, access);
167 for(
int i = e - derivatives[
d] + 1; i <= e; ++i)
175 prod *
ipow(in[
d], e-derivatives[
d]) * coeff,
191 template <
typename Traits>
194 enum {
d = Traits::dimDomain-1 };
196 template <
typename Access>
197 static void eval (
const typename Traits::DomainType &in,
198 const std::array<int, Traits::dimDomain> &derivatives,
199 typename Traits::RangeFieldType prod,
200 int bound,
int& index, Access &access)
202 if(bound < derivatives[
d])
206 for(
int i = bound - derivatives[
d] + 1; i <= bound; ++i)
208 prod *=
ipow(in[
d], bound-derivatives[
d]) * coeff;
210 access[index] = prod;
231 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>(std::array<int, 0>(), in, out);
259 std::vector<typename Traits::RangeType>& out)
const
261 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
270 std::array<int,1> directions;
271 directions[0] = std::find(
order.begin(),
order.end(), 1)-
order.begin();
272 evaluate<1>(directions, in, out);
277 std::array<int,2> directions;
278 unsigned int counter = 0;
279 auto nonconstOrder =
order;
280 for (
unsigned int i=0; i<d; i++)
282 while (nonconstOrder[i])
284 directions[counter++] = i;
289 evaluate<2>(directions, in, out);
294 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
299 template<
unsigned int k>
300 inline void evaluate (
const std::array<int,k>& directions,
302 std::vector<typename Traits::RangeType>& out)
const
306 std::array<int, d> derivatives;
307 for(
unsigned int i = 0; i < d; ++i) derivatives[i] = 0;
308 for(
unsigned int i = 0; i < k; ++i) ++derivatives[directions[i]];
310 for(
unsigned int lp = 0; lp <= p; ++lp)
318 std::vector<typename Traits::JacobianType>& out)
const
321 std::array<int, d> derivatives;
322 for(
unsigned int i = 0; i < d; ++i)
324 for(
unsigned int i = 0; i < d; ++i)
329 for(
unsigned int lp = 0; lp <= p; ++lp)
Definition brezzidouglasmarini1cube2d.hh:14
T ipow(T base, int exp)
Definition monomiallocalbasis.hh:38
Type traits for LocalBasisVirtualInterface.
Definition localbasis.hh:38
D DomainType
domain type
Definition localbasis.hh:49
Definition monomiallocalbasis.hh:21
@ val
Definition monomiallocalbasis.hh:22
Access output vector of evaluateFunction() and evaluate()
Definition monomiallocalbasis.hh:53
Traits::RangeFieldType & operator[](unsigned int index)
Definition monomiallocalbasis.hh:71
~EvalAccess()
Definition monomiallocalbasis.hh:67
EvalAccess(std::vector< typename Traits::RangeType > &out_)
Definition monomiallocalbasis.hh:60
Access output vector of evaluateJacobian()
Definition monomiallocalbasis.hh:84
~JacobianAccess()
Definition monomiallocalbasis.hh:100
Traits::RangeFieldType & operator[](unsigned int index)
Definition monomiallocalbasis.hh:104
JacobianAccess(std::vector< typename Traits::JacobianType > &out_, unsigned int row_)
Definition monomiallocalbasis.hh:92
Definition monomiallocalbasis.hh:129
@ d
The next dimension to try for factors.
Definition monomiallocalbasis.hh:132
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:141
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:197
Constant shape function.
Definition monomiallocalbasis.hh:233
void partial(const std::array< unsigned int, d > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of any order of all shape functions.
Definition monomiallocalbasis.hh:257
unsigned int size() const
number of shape functions
Definition monomiallocalbasis.hh:240
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:237
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition monomiallocalbasis.hh:246
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition monomiallocalbasis.hh:317
unsigned int order() const
Polynomial order of the shape functions.
Definition monomiallocalbasis.hh:336
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:300