5#ifndef DUNE_RAVIARTTHOMASPREBASIS_HH
6#define DUNE_RAVIARTTHOMASPREBASIS_HH
11#include <dune/geometry/type.hh>
17 template < GeometryType::Id geometryId,
class Field >
20 template <
unsigned int dim,
class Field>
24 typedef typename MBasisFactory::Object
MBasis;
29 typedef std::size_t
Key;
31 template <
unsigned int dd,
class FF>
36 template< GeometryType::Id geometryId >
40 MBasis *mbasis = MBasisFactory::template create<geometryId>(order+1);
41 typename std::remove_const<Object>::type *tmBasis =
new typename std::remove_const<Object>::type(*mbasis);
42 tmBasis->fill(vecMatrix);
48 template <GeometryType::Id geometryId,
class Field>
51 static constexpr GeometryType
geometry = geometryId;
78 FieldVector< MI, dim > x;
85 for(
unsigned int i = 0; i <
dim; ++i )
87 std::vector< MI > val( basis.
size() );
95 unsigned int notHomogen = 0;
97 notHomogen = basis.
sizes()[order-1];
100 unsigned int homogen = basis.
sizes()[order]-notHomogen;
128 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
130 for (
unsigned int r=0; r<
dim; ++r)
132 for (
unsigned int rr=0; rr<
dim; ++rr)
136 for (
unsigned int j=0; j<
col_; ++j)
150 for (
unsigned int i=0; i<homogen; ++i)
152 for (
unsigned int r=0; r<
dim; ++r)
156 for (
unsigned int j=0; j<
col_; ++j)
163 MI xval = val[notHomogen+i];
165 for (w=homogen+notHomogen; w<val.size(); ++w)
173 assert(w<val.size());
181 for (
unsigned int i=0; i<
rows(); ++i) {
195 template <
class Vector>
196 void row(
const unsigned int row, Vector &vec )
const
198 const unsigned int N =
cols();
199 assert( vec.size() == N );
200 for (
unsigned int i=0; i<N; ++i)
Definition bdfmcube.hh:18
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition field.hh:159
Definition raviartthomassimplexprebasis.hh:50
static const unsigned int dim
Definition raviartthomassimplexprebasis.hh:52
~RTVecMatrix()
Definition raviartthomassimplexprebasis.hh:179
Field ** mat_
Definition raviartthomassimplexprebasis.hh:204
RTVecMatrix(std::size_t order)
Definition raviartthomassimplexprebasis.hh:55
unsigned int cols() const
Definition raviartthomassimplexprebasis.hh:187
unsigned int row_
Definition raviartthomassimplexprebasis.hh:203
MultiIndex< dim, Field > MI
Definition raviartthomassimplexprebasis.hh:53
void row(const unsigned int row, Vector &vec) const
Definition raviartthomassimplexprebasis.hh:196
MonomialBasis< geometryId, MI > MIBasis
Definition raviartthomassimplexprebasis.hh:54
static constexpr GeometryType geometry
Definition raviartthomassimplexprebasis.hh:51
unsigned int rows() const
Definition raviartthomassimplexprebasis.hh:191
unsigned int col_
Definition raviartthomassimplexprebasis.hh:203
Definition raviartthomassimplexprebasis.hh:22
const Basis Object
Definition raviartthomassimplexprebasis.hh:28
MonomialBasisProvider< dim, Field > MBasisFactory
Definition raviartthomassimplexprebasis.hh:23
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition raviartthomassimplexprebasis.hh:26
MBasisFactory::Object MBasis
Definition raviartthomassimplexprebasis.hh:24
std::size_t Key
Definition raviartthomassimplexprebasis.hh:29
static void release(Object *object)
Definition raviartthomassimplexprebasis.hh:45
StandardEvaluator< MBasis > EvalMBasis
Definition raviartthomassimplexprebasis.hh:25
static Object * create(const Key &order)
Definition raviartthomassimplexprebasis.hh:37
Definition raviartthomassimplexprebasis.hh:33
MonomialBasisProvider< dd, FF > Type
Definition raviartthomassimplexprebasis.hh:34
Definition basisevaluator.hh:131
Definition monomialbasis.hh:440
unsigned int size() const
Definition monomialbasis.hh:476
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition monomialbasis.hh:498
const unsigned int * sizes(unsigned int order) const
Definition monomialbasis.hh:465
Definition monomialbasis.hh:780
Definition multiindex.hh:38
Definition polynomialbasis.hh:348