3#ifndef DUNE_RAVIARTTHOMASPREBASIS_HH
4#define DUNE_RAVIARTTHOMASPREBASIS_HH
9#include <dune/geometry/type.hh>
15 template < GeometryType::Id geometryId,
class Field >
18 template <
unsigned int dim,
class Field>
22 typedef typename MBasisFactory::Object
MBasis;
27 typedef std::size_t
Key;
29 template <
unsigned int dd,
class FF>
34 template< GeometryType::Id geometryId >
38 MBasis *mbasis = MBasisFactory::template create<geometryId>(order+1);
39 typename std::remove_const<Object>::type *tmBasis =
new typename std::remove_const<Object>::type(*mbasis);
40 tmBasis->fill(vecMatrix);
46 template <GeometryType::Id geometryId,
class Field>
49 static constexpr GeometryType
geometry = geometryId;
76 FieldVector< MI, dim > x;
83 for(
unsigned int i = 0; i <
dim; ++i )
85 std::vector< MI > val( basis.
size() );
93 unsigned int notHomogen = 0;
95 notHomogen = basis.
sizes()[order-1];
98 unsigned int homogen = basis.
sizes()[order]-notHomogen;
126 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
128 for (
unsigned int r=0; r<
dim; ++r)
130 for (
unsigned int rr=0; rr<
dim; ++rr)
134 for (
unsigned int j=0; j<
col_; ++j)
148 for (
unsigned int i=0; i<homogen; ++i)
150 for (
unsigned int r=0; r<
dim; ++r)
154 for (
unsigned int j=0; j<
col_; ++j)
161 MI xval = val[notHomogen+i];
163 for (w=homogen+notHomogen; w<val.size(); ++w)
171 assert(w<val.size());
179 for (
unsigned int i=0; i<
rows(); ++i) {
193 template <
class Vector>
194 void row(
const unsigned int row, Vector &vec )
const
196 const unsigned int N =
cols();
197 assert( vec.size() == N );
198 for (
unsigned int i=0; i<N; ++i)
Definition: bdfmcube.hh:16
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition: field.hh:157
Definition: raviartthomassimplexprebasis.hh:48
static const unsigned int dim
Definition: raviartthomassimplexprebasis.hh:50
~RTVecMatrix()
Definition: raviartthomassimplexprebasis.hh:177
Field ** mat_
Definition: raviartthomassimplexprebasis.hh:202
RTVecMatrix(std::size_t order)
Definition: raviartthomassimplexprebasis.hh:53
unsigned int cols() const
Definition: raviartthomassimplexprebasis.hh:185
unsigned int row_
Definition: raviartthomassimplexprebasis.hh:201
MultiIndex< dim, Field > MI
Definition: raviartthomassimplexprebasis.hh:51
void row(const unsigned int row, Vector &vec) const
Definition: raviartthomassimplexprebasis.hh:194
MonomialBasis< geometryId, MI > MIBasis
Definition: raviartthomassimplexprebasis.hh:52
static constexpr GeometryType geometry
Definition: raviartthomassimplexprebasis.hh:49
unsigned int rows() const
Definition: raviartthomassimplexprebasis.hh:189
unsigned int col_
Definition: raviartthomassimplexprebasis.hh:201
Definition: raviartthomassimplexprebasis.hh:20
const Basis Object
Definition: raviartthomassimplexprebasis.hh:26
MonomialBasisProvider< dim, Field > MBasisFactory
Definition: raviartthomassimplexprebasis.hh:21
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition: raviartthomassimplexprebasis.hh:24
MBasisFactory::Object MBasis
Definition: raviartthomassimplexprebasis.hh:22
std::size_t Key
Definition: raviartthomassimplexprebasis.hh:27
static void release(Object *object)
Definition: raviartthomassimplexprebasis.hh:43
StandardEvaluator< MBasis > EvalMBasis
Definition: raviartthomassimplexprebasis.hh:23
static Object * create(const Key &order)
Definition: raviartthomassimplexprebasis.hh:35
Definition: raviartthomassimplexprebasis.hh:31
MonomialBasisProvider< dd, FF > Type
Definition: raviartthomassimplexprebasis.hh:32
Definition: basisevaluator.hh:129
Definition: monomialbasis.hh:438
unsigned int size() const
Definition: monomialbasis.hh:474
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:496
const unsigned int * sizes(unsigned int order) const
Definition: monomialbasis.hh:463
Definition: monomialbasis.hh:789
Definition: multiindex.hh:35
Definition: polynomialbasis.hh:335