3 #ifndef DUNE_RAVIARTTHOMASPREBASIS_HH 4 #define DUNE_RAVIARTTHOMASPREBASIS_HH 9 #include <dune/geometry/type.hh> 15 template <
unsigned int dim,
class Field>
17 template <
unsigned int dim,
class Field>
23 typedef typename MBasisFactory::Object
MBasis;
28 typedef unsigned int Key;
32 template <
class Topology,
class Field >
35 template <
unsigned int dim,
class Field>
37 :
public TopologyFactory< RTPreBasisFactoryTraits< dim, Field > >
43 template <
unsigned int dd,
class FF>
48 template<
class Topology >
52 typename Traits::MBasis *mbasis = Traits::MBasisFactory::template create<Topology>(order+1);
53 typename std::remove_const<Object>::type *tmBasis =
54 new typename std::remove_const<Object>::type(*mbasis);
55 tmBasis->fill(vecMatrix);
59 template <
class Topology,
class Field>
62 static const unsigned int dim = Topology::dimension;
67 MIBasis basis(order+1);
68 FieldVector< MI, dim > x;
69 for(
unsigned int i = 0; i < dim; ++i )
71 std::vector< MI > val( basis.
size() );
75 unsigned int notHomogen = 0;
77 notHomogen = basis.
sizes()[order-1];
78 unsigned int homogen = basis.
sizes()[order]-notHomogen;
79 row_ = (notHomogen*dim+homogen*(dim+1))*dim;
80 row1_ = basis.
sizes()[order]*dim*dim;
81 mat_ =
new Field*[row_];
83 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
85 for (
unsigned int r=0; r<dim; ++r)
87 for (
unsigned int rr=0; rr<dim; ++rr)
89 mat_[row] =
new Field[col_];
90 for (
unsigned int j=0; j<col_; ++j)
100 for (
unsigned int i=0; i<homogen; ++i)
102 for (
unsigned int r=0; r<dim; ++r)
104 mat_[row] =
new Field[col_];
105 for (
unsigned int j=0; j<col_; ++j)
110 MI xval = val[notHomogen+i];
112 for (w=homogen+notHomogen; w<val.size(); ++w)
120 assert(w<val.size());
127 for (
unsigned int i=0; i<rows(); ++i) {
138 template <
class Vector>
139 void row(
const unsigned int row, Vector &vec )
const 141 const unsigned int N = cols();
142 assert( vec.size() == N );
143 for (
unsigned int i=0; i<N; ++i)
152 #endif // DUNE_RAVIARTTHOMASPREBASIS_HH const unsigned int * sizes(unsigned int order) const
Definition: monomialbasis.hh:661
MonomialBasis< Topology, MI > MIBasis
Definition: raviartthomassimplexprebasis.hh:64
Definition: raviartthomassimplexprebasis.hh:33
Definition: monomialbasis.hh:73
static Object * createObject(const Key &order)
Definition: raviartthomassimplexprebasis.hh:49
Definition: basisevaluator.hh:127
Definition: raviartthomassimplexprebasis.hh:18
Traits::Key Key
Definition: raviartthomassimplexprebasis.hh:42
RTPreBasisFactoryTraits< dim, Field > Traits
Definition: raviartthomassimplexprebasis.hh:39
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:695
StandardEvaluator< MBasis > EvalMBasis
Definition: raviartthomassimplexprebasis.hh:24
unsigned int row_
Definition: raviartthomassimplexprebasis.hh:146
unsigned int size() const
Definition: monomialbasis.hh:672
MBasisFactory::Object MBasis
Definition: raviartthomassimplexprebasis.hh:23
unsigned int cols() const
Definition: raviartthomassimplexprebasis.hh:132
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition: field.hh:157
unsigned int Key
Definition: raviartthomassimplexprebasis.hh:28
Traits::Object Object
Definition: raviartthomassimplexprebasis.hh:41
static const unsigned int dimension
Definition: raviartthomassimplexprebasis.hh:20
const Basis Object
Definition: raviartthomassimplexprebasis.hh:27
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
Field ** mat_
Definition: raviartthomassimplexprebasis.hh:147
Definition: monomialbasis.hh:987
MultiIndex< dim, Field > MI
Definition: raviartthomassimplexprebasis.hh:63
Definition: polynomialbasis.hh:265
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition: raviartthomassimplexprebasis.hh:25
RTVecMatrix(unsigned int order)
Definition: raviartthomassimplexprebasis.hh:65
RTPreBasisFactory< dim, Field > Factory
Definition: raviartthomassimplexprebasis.hh:29
MonomialBasisProvider< dim, Field > MBasisFactory
Definition: raviartthomassimplexprebasis.hh:22
Definition: raviartthomassimplexprebasis.hh:44
void row(const unsigned int row, Vector &vec) const
Definition: raviartthomassimplexprebasis.hh:139
Definition: raviartthomassimplexprebasis.hh:16
MonomialBasisProvider< dd, FF > Type
Definition: raviartthomassimplexprebasis.hh:46
Definition: multiindex.hh:23
~RTVecMatrix()
Definition: raviartthomassimplexprebasis.hh:125
unsigned int rows() const
Definition: raviartthomassimplexprebasis.hh:135