dune-localfunctions  2.3.1-rc1
raviartthomassimplexprebasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_RAVIARTTHOMASPREBASIS_HH
4 #define DUNE_RAVIARTTHOMASPREBASIS_HH
5 #include <fstream>
6 #include <utility>
7 
8 #include <dune/geometry/genericgeometry/topologytypes.hh>
9 
11 
12 namespace Dune
13 {
14  template <unsigned int dim, class Field>
16  template <unsigned int dim, class Field>
18  {
19  static const unsigned int dimension = dim;
20 
22  typedef typename MBasisFactory::Object MBasis;
25 
26  typedef const Basis Object;
27  typedef unsigned int Key;
29  };
30 
31  template < class Topology, class Field >
32  struct RTVecMatrix;
33 
34  template <unsigned int dim, class Field>
35  struct RTPreBasisFactory
36  : public TopologyFactory< RTPreBasisFactoryTraits< dim, Field > >
37  {
39  static const unsigned int dimension = dim;
40  typedef typename Traits::Object Object;
41  typedef typename Traits::Key Key;
42  template <unsigned int dd, class FF>
44  {
46  };
47  template< class Topology >
48  static Object *createObject ( const Key &order )
49  {
50  RTVecMatrix<Topology,Field> vecMatrix(order);
51  typename Traits::MBasis *mbasis = Traits::MBasisFactory::template create<Topology>(order+1);
52  typename remove_const<Object>::type *tmBasis =
53  new typename remove_const<Object>::type(*mbasis);
54  tmBasis->fill(vecMatrix);
55  return tmBasis;
56  }
57  };
58  template <class Topology, class Field>
59  struct RTVecMatrix
60  {
61  static const unsigned int dim = Topology::dimension;
64  RTVecMatrix(unsigned int order)
65  {
66  MIBasis basis(order+1);
67  FieldVector< MI, dim > x;
68  for( unsigned int i = 0; i < dim; ++i )
69  x[ i ].set( i, 1 );
70  std::vector< MI > val( basis.size() );
71  basis.evaluate( x, val );
72 
73  col_ = basis.size();
74  unsigned int notHomogen = 0;
75  if (order>0)
76  notHomogen = basis.sizes()[order-1];
77  unsigned int homogen = basis.sizes()[order]-notHomogen;
78  row_ = (notHomogen*dim+homogen*(dim+1))*dim;
79  row1_ = basis.sizes()[order]*dim*dim;
80  mat_ = new Field*[row_];
81  int row = 0;
82  for (unsigned int i=0; i<notHomogen+homogen; ++i)
83  {
84  for (unsigned int r=0; r<dim; ++r)
85  {
86  for (unsigned int rr=0; rr<dim; ++rr)
87  {
88  mat_[row] = new Field[col_];
89  for (unsigned int j=0; j<col_; ++j)
90  {
91  mat_[row][j] = 0.;
92  }
93  if (r==rr)
94  mat_[row][i] = 1.;
95  ++row;
96  }
97  }
98  }
99  for (unsigned int i=0; i<homogen; ++i)
100  {
101  for (unsigned int r=0; r<dim; ++r)
102  {
103  mat_[row] = new Field[col_];
104  for (unsigned int j=0; j<col_; ++j)
105  {
106  mat_[row][j] = 0.;
107  }
108  unsigned int w;
109  MI xval = val[notHomogen+i];
110  xval *= x[r];
111  for (w=homogen+notHomogen; w<val.size(); ++w)
112  {
113  if (val[w] == xval)
114  {
115  mat_[row][w] = 1.;
116  break;
117  }
118  }
119  assert(w<val.size());
120  ++row;
121  }
122  }
123  }
125  {
126  for (unsigned int i=0; i<rows(); ++i) {
127  delete [] mat_[i];
128  }
129  delete [] mat_;
130  }
131  unsigned int cols() const {
132  return col_;
133  }
134  unsigned int rows() const {
135  return row_;
136  }
137  template <class Vector>
138  void row( const unsigned int row, Vector &vec ) const
139  {
140  const unsigned int N = cols();
141  assert( vec.size() == N );
142  for (unsigned int i=0; i<N; ++i)
143  field_cast(mat_[row][i],vec[i]);
144  }
145  unsigned int row_,col_,row1_;
146  Field **mat_;
147  };
148 
149 
150 }
151 #endif // DUNE_RAVIARTTHOMASPREBASIS_HH
RTPreBasisFactoryTraits< dim, Field > Traits
Definition: raviartthomassimplexprebasis.hh:38
Definition: basisevaluator.hh:128
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition: raviartthomassimplexprebasis.hh:24
const Basis Object
Definition: raviartthomassimplexprebasis.hh:26
StandardEvaluator< MBasis > EvalMBasis
Definition: raviartthomassimplexprebasis.hh:23
Traits::Object Object
Definition: raviartthomassimplexprebasis.hh:40
unsigned int Key
Definition: raviartthomassimplexprebasis.hh:27
RTPreBasisFactory< dim, Field > Factory
Definition: raviartthomassimplexprebasis.hh:28
Definition: raviartthomassimplexprebasis.hh:17
MonomialBasisProvider< dim, Field > MBasisFactory
Definition: raviartthomassimplexprebasis.hh:21
RTVecMatrix(unsigned int order)
Definition: raviartthomassimplexprebasis.hh:64
MBasisFactory::Object MBasis
Definition: raviartthomassimplexprebasis.hh:22
static const unsigned int dimension
Definition: raviartthomassimplexprebasis.hh:39
static Object * createObject(const Key &order)
Definition: raviartthomassimplexprebasis.hh:48
Definition: multiindex.hh:23
Definition: raviartthomassimplexprebasis.hh:43
Definition: monomialbasis.hh:981
Definition: raviartthomassimplexprebasis.hh:32
Field ** mat_
Definition: raviartthomassimplexprebasis.hh:146
Definition: raviartthomassimplexprebasis.hh:15
void row(const unsigned int row, Vector &vec) const
Definition: raviartthomassimplexprebasis.hh:138
const unsigned int size() const
Definition: monomialbasis.hh:666
Definition: polynomialbasis.hh:251
static const unsigned int dimension
Definition: raviartthomassimplexprebasis.hh:19
~RTVecMatrix()
Definition: raviartthomassimplexprebasis.hh:124
unsigned int col_
Definition: raviartthomassimplexprebasis.hh:145
unsigned int row_
Definition: raviartthomassimplexprebasis.hh:145
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition: field.hh:157
unsigned int row1_
Definition: raviartthomassimplexprebasis.hh:145
Traits::Key Key
Definition: raviartthomassimplexprebasis.hh:41
MonomialBasis< Topology, MI > MIBasis
Definition: raviartthomassimplexprebasis.hh:63
unsigned int cols() const
Definition: raviartthomassimplexprebasis.hh:131
Definition: monomialbasis.hh:73
static const unsigned int dim
Definition: raviartthomassimplexprebasis.hh:61
unsigned int rows() const
Definition: raviartthomassimplexprebasis.hh:134
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:689
MonomialBasisProvider< dd, FF > Type
Definition: raviartthomassimplexprebasis.hh:45
const unsigned int * sizes(unsigned int order) const
Definition: monomialbasis.hh:655
MultiIndex< dim, Field > MI
Definition: raviartthomassimplexprebasis.hh:62