dune-localfunctions  2.7.1
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 < class Topology, class Field >
16  struct RTVecMatrix;
17 
18  template <unsigned int dim, class Field>
20  {
22  typedef typename MBasisFactory::Object MBasis;
25 
26  typedef const Basis Object;
27  typedef std::size_t Key;
28 
29  template <unsigned int dd, class FF>
31  {
33  };
34  template< class Topology >
35  static Object *create ( const Key &order )
36  {
37  RTVecMatrix<Topology,Field> vecMatrix(order);
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);
42  return tmBasis;
43  }
44  static void release( Object *object ) { delete object; }
45  };
46  template <class Topology, class Field>
47  struct RTVecMatrix
48  {
49  static const unsigned int dim = Topology::dimension;
52  RTVecMatrix(std::size_t order)
53  {
54  MIBasis basis(order+1);
55  FieldVector< MI, dim > x;
56  for( unsigned int i = 0; i < dim; ++i )
57  x[ i ].set( i, 1 );
58  std::vector< MI > val( basis.size() );
59  basis.evaluate( x, val );
60 
61  col_ = basis.size();
62  unsigned int notHomogen = 0;
63  if (order>0)
64  notHomogen = basis.sizes()[order-1];
65  unsigned int homogen = basis.sizes()[order]-notHomogen;
66  row_ = (notHomogen*dim+homogen*(dim+1))*dim;
67  row1_ = basis.sizes()[order]*dim*dim;
68  mat_ = new Field*[row_];
69  int row = 0;
70  for (unsigned int i=0; i<notHomogen+homogen; ++i)
71  {
72  for (unsigned int r=0; r<dim; ++r)
73  {
74  for (unsigned int rr=0; rr<dim; ++rr)
75  {
76  mat_[row] = new Field[col_];
77  for (unsigned int j=0; j<col_; ++j)
78  {
79  mat_[row][j] = 0.;
80  }
81  if (r==rr)
82  mat_[row][i] = 1.;
83  ++row;
84  }
85  }
86  }
87  for (unsigned int i=0; i<homogen; ++i)
88  {
89  for (unsigned int r=0; r<dim; ++r)
90  {
91  mat_[row] = new Field[col_];
92  for (unsigned int j=0; j<col_; ++j)
93  {
94  mat_[row][j] = 0.;
95  }
96  unsigned int w;
97  MI xval = val[notHomogen+i];
98  xval *= x[r];
99  for (w=homogen+notHomogen; w<val.size(); ++w)
100  {
101  if (val[w] == xval)
102  {
103  mat_[row][w] = 1.;
104  break;
105  }
106  }
107  assert(w<val.size());
108  ++row;
109  }
110  }
111  }
113  {
114  for (unsigned int i=0; i<rows(); ++i) {
115  delete [] mat_[i];
116  }
117  delete [] mat_;
118  }
119  unsigned int cols() const {
120  return col_;
121  }
122  unsigned int rows() const {
123  return row_;
124  }
125  template <class Vector>
126  void row( const unsigned int row, Vector &vec ) const
127  {
128  const unsigned int N = cols();
129  assert( vec.size() == N );
130  for (unsigned int i=0; i<N; ++i)
131  field_cast(mat_[row][i],vec[i]);
132  }
133  unsigned int row_,col_,row1_;
134  Field **mat_;
135  };
136 
137 
138 }
139 #endif // DUNE_RAVIARTTHOMASPREBASIS_HH
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