3 #ifndef DUNE_RAVIARTTHOMASPREBASIS_HH
4 #define DUNE_RAVIARTTHOMASPREBASIS_HH
9 #include <dune/geometry/type.hh>
15 template <
class Topology,
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<
class Topology >
38 MBasis *mbasis = MBasisFactory::template create<Topology>(order+1);
39 typename std::remove_const<Object>::type *tmBasis =
40 new typename std::remove_const<Object>::type(*mbasis);
41 tmBasis->fill(vecMatrix);
46 template <
class Topology,
class Field>
49 static const unsigned int dim = Topology::dimension;
55 FieldVector< MI, dim > x;
56 for(
unsigned int i = 0; i <
dim; ++i )
58 std::vector< MI > val( basis.
size() );
62 unsigned int notHomogen = 0;
64 notHomogen = basis.
sizes()[order-1];
65 unsigned int homogen = basis.
sizes()[order]-notHomogen;
70 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
72 for (
unsigned int r=0; r<
dim; ++r)
74 for (
unsigned int rr=0; rr<
dim; ++rr)
77 for (
unsigned int j=0; j<
col_; ++j)
87 for (
unsigned int i=0; i<homogen; ++i)
89 for (
unsigned int r=0; r<
dim; ++r)
92 for (
unsigned int j=0; j<
col_; ++j)
97 MI xval = val[notHomogen+i];
99 for (w=homogen+notHomogen; w<val.size(); ++w)
107 assert(w<val.size());
114 for (
unsigned int i=0; i<
rows(); ++i) {
125 template <
class Vector>
126 void row(
const unsigned int row, Vector &vec )
const
128 const unsigned int N =
cols();
129 assert( vec.size() == N );
130 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:49
void row(const unsigned int row, Vector &vec) const
Definition: raviartthomassimplexprebasis.hh:126
MultiIndex< dim, Field > MI
Definition: raviartthomassimplexprebasis.hh:50
~RTVecMatrix()
Definition: raviartthomassimplexprebasis.hh:112
unsigned int cols() const
Definition: raviartthomassimplexprebasis.hh:119
unsigned int row_
Definition: raviartthomassimplexprebasis.hh:133
unsigned int rows() const
Definition: raviartthomassimplexprebasis.hh:122
unsigned int col_
Definition: raviartthomassimplexprebasis.hh:133
RTVecMatrix(std::size_t order)
Definition: raviartthomassimplexprebasis.hh:52
MonomialBasis< Topology, MI > MIBasis
Definition: raviartthomassimplexprebasis.hh:51
unsigned int row1_
Definition: raviartthomassimplexprebasis.hh:133
Field ** mat_
Definition: raviartthomassimplexprebasis.hh:134
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
static Object * create(const Key &order)
Definition: raviartthomassimplexprebasis.hh:35
MBasisFactory::Object MBasis
Definition: raviartthomassimplexprebasis.hh:22
std::size_t Key
Definition: raviartthomassimplexprebasis.hh:27
static void release(Object *object)
Definition: raviartthomassimplexprebasis.hh:44
StandardEvaluator< MBasis > EvalMBasis
Definition: raviartthomassimplexprebasis.hh:23
Definition: raviartthomassimplexprebasis.hh:31
MonomialBasisProvider< dd, FF > Type
Definition: raviartthomassimplexprebasis.hh:32
Definition: basisevaluator.hh:129
Definition: monomialbasis.hh:506
unsigned int size() const
Definition: monomialbasis.hh:541
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:564
const unsigned int * sizes(unsigned int order) const
Definition: monomialbasis.hh:530
Definition: monomialbasis.hh:845
Definition: multiindex.hh:35
Definition: polynomialbasis.hh:335