dune-localfunctions  2.6-git
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 
6 #include <fstream>
7 #include <utility>
8 
9 #include <dune/geometry/type.hh>
10 
12 
13 namespace Dune
14 {
15  template <unsigned int dim, class Field>
17  template <unsigned int dim, class Field>
19  {
20  static const unsigned int dimension = dim;
21 
23  typedef typename MBasisFactory::Object MBasis;
26 
27  typedef const Basis Object;
28  typedef unsigned int Key;
30  };
31 
32  template < class Topology, class Field >
33  struct RTVecMatrix;
34 
35  template <unsigned int dim, class Field>
36  struct RTPreBasisFactory
37  : public TopologyFactory< RTPreBasisFactoryTraits< dim, Field > >
38  {
40  static const unsigned int dimension = dim;
41  typedef typename Traits::Object Object;
42  typedef typename Traits::Key Key;
43  template <unsigned int dd, class FF>
45  {
47  };
48  template< class Topology >
49  static Object *createObject ( const Key &order )
50  {
51  RTVecMatrix<Topology,Field> vecMatrix(order);
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);
56  return tmBasis;
57  }
58  };
59  template <class Topology, class Field>
60  struct RTVecMatrix
61  {
62  static const unsigned int dim = Topology::dimension;
65  RTVecMatrix(unsigned int order)
66  {
67  MIBasis basis(order+1);
68  FieldVector< MI, dim > x;
69  for( unsigned int i = 0; i < dim; ++i )
70  x[ i ].set( i, 1 );
71  std::vector< MI > val( basis.size() );
72  basis.evaluate( x, val );
73 
74  col_ = basis.size();
75  unsigned int notHomogen = 0;
76  if (order>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_];
82  int row = 0;
83  for (unsigned int i=0; i<notHomogen+homogen; ++i)
84  {
85  for (unsigned int r=0; r<dim; ++r)
86  {
87  for (unsigned int rr=0; rr<dim; ++rr)
88  {
89  mat_[row] = new Field[col_];
90  for (unsigned int j=0; j<col_; ++j)
91  {
92  mat_[row][j] = 0.;
93  }
94  if (r==rr)
95  mat_[row][i] = 1.;
96  ++row;
97  }
98  }
99  }
100  for (unsigned int i=0; i<homogen; ++i)
101  {
102  for (unsigned int r=0; r<dim; ++r)
103  {
104  mat_[row] = new Field[col_];
105  for (unsigned int j=0; j<col_; ++j)
106  {
107  mat_[row][j] = 0.;
108  }
109  unsigned int w;
110  MI xval = val[notHomogen+i];
111  xval *= x[r];
112  for (w=homogen+notHomogen; w<val.size(); ++w)
113  {
114  if (val[w] == xval)
115  {
116  mat_[row][w] = 1.;
117  break;
118  }
119  }
120  assert(w<val.size());
121  ++row;
122  }
123  }
124  }
126  {
127  for (unsigned int i=0; i<rows(); ++i) {
128  delete [] mat_[i];
129  }
130  delete [] mat_;
131  }
132  unsigned int cols() const {
133  return col_;
134  }
135  unsigned int rows() const {
136  return row_;
137  }
138  template <class Vector>
139  void row( const unsigned int row, Vector &vec ) const
140  {
141  const unsigned int N = cols();
142  assert( vec.size() == N );
143  for (unsigned int i=0; i<N; ++i)
144  field_cast(mat_[row][i],vec[i]);
145  }
146  unsigned int row_,col_,row1_;
147  Field **mat_;
148  };
149 
150 
151 }
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