dune-fem  2.4.1-rc
butchertable.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_BUTCHERTABLE_HH
2 #define DUNE_FEM_SOLVER_RUNGEKUTTA_BUTCHERTABLE_HH
3 
4 //- dune-common includes
5 #include <dune/common/dynmatrix.hh>
6 #include <dune/common/dynvector.hh>
7 
8 namespace DuneODE
9 {
10 
11  // SimpleButcherTable
12  // ------------------
13 
14  template< class Field >
16  {
18 
19  public:
20  typedef Field FieldType;
21 
22  SimpleButcherTable ( int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c )
23  : stages_( stages ), order_( order ),
24  a_( a ), b_( b ), c_( c )
25  {}
26 
27  Dune::DynamicMatrix< FieldType > A () const { return makeMatrix( stages_, stages_, a_ ); }
28  Dune::DynamicVector< FieldType > b () const { return makeVector( stages_, b_ ); }
29  Dune::DynamicVector< FieldType > c () const { return makeVector( stages_, c_ ); }
30 
31  int order () const { return order_; }
32  int stages () const { return stages_; }
33 
34  protected:
35  static Dune::DynamicMatrix< FieldType > makeMatrix ( int m, int n, const FieldType *data )
36  {
37  Dune::DynamicMatrix< FieldType > A( m, n );
38  for( int i = 0; i < m; ++i )
39  std::copy( data + i*n, data + (i+1)*n, A[ i ].begin() );
40  return A;
41  }
42 
43  static Dune::DynamicVector< FieldType > makeVector ( int n, const FieldType *data )
44  {
45  Dune::DynamicVector< FieldType > v( n );
46  std::copy( data, data + n, v.begin() );
47  return v;
48  }
49 
51  const FieldType *a_, *b_, *c_;
52  };
53 
54 
55 
56  // implicit butcher tables
57  // -----------------------
58 
63 
64 
65 
66  // semiimplicit butcher tables
67  // ---------------------------
68 
76 
77  // ROW butcher tables
78  // -----------------------
79 
80  template< class Field >
82  {
85 
86  public:
87  typedef Field FieldType;
88 
89  ROWSimpleButcherTable ( int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c, const FieldType *a2 )
90  : Base(stages,order,a,b,c)
91  , a2_(a2)
92  {}
93 
94  Dune::DynamicMatrix< FieldType > B () const { return Base::makeMatrix( stages_, stages_, a2_ ); }
95 
96  private:
97  using Base::stages_;
98  const FieldType *a2_;
99  };
102 } // namespace DuneODE
103 
104 #endif // #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_BUTCHERTABLE_HH
static Dune::DynamicMatrix< FieldType > makeMatrix(int m, int n, const FieldType *data)
Definition: butchertable.hh:35
Definition: butchertable.hh:15
int order() const
Definition: butchertable.hh:31
static Dune::DynamicVector< FieldType > makeVector(int n, const FieldType *data)
Definition: butchertable.hh:43
SimpleButcherTable< double > implicitEulerButcherTable()
Definition: butchertable.cc:58
SimpleButcherTable< double > semiImplicitEulerButcherTable(bool expl)
Definition: butchertable.cc:93
Definition: multistep.hh:16
SimpleButcherTable< double > semiImplicitSSP222ButcherTable(bool expl)
Definition: butchertable.cc:189
SimpleButcherTable< double > gauss2ButcherTable()
Definition: butchertable.cc:73
ROWSimpleButcherTable< double > row3ButcherTable()
Definition: butchertable.cc:362
const FieldType * c_
Definition: butchertable.hh:51
ROWSimpleButcherTable(int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c, const FieldType *a2)
Definition: butchertable.hh:89
SimpleButcherTable< double > semiImplicitARK46ButcherTable(bool expl)
SimpleButcherTable< double > semiImplicit33ButcherTable(bool expl)
Definition: butchertable.cc:159
Dune::DynamicMatrix< FieldType > B() const
Definition: butchertable.hh:94
Dune::DynamicMatrix< FieldType > A() const
Definition: butchertable.hh:27
SimpleButcherTable< double > implicit3ButcherTable()
Definition: butchertable.cc:44
Field FieldType
Definition: butchertable.hh:87
int stages() const
Definition: butchertable.hh:32
Definition: butchertable.hh:81
SimpleButcherTable< double > semiImplicitARK34ButcherTable(bool expl)
Field FieldType
Definition: butchertable.hh:20
SimpleButcherTable< double > implicit34ButcherTable()
Definition: butchertable.cc:24
SimpleButcherTable(int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c)
Definition: butchertable.hh:22
const FieldType * b_
Definition: butchertable.hh:51
int order_
Definition: butchertable.hh:50
const FieldType * a_
Definition: butchertable.hh:51
Dune::DynamicVector< FieldType > c() const
Definition: butchertable.hh:29
int stages_
Definition: butchertable.hh:50
Dune::DynamicVector< FieldType > b() const
Definition: butchertable.hh:28
SimpleButcherTable< double > semiImplicit23ButcherTable(bool expl)
Definition: butchertable.cc:124
SimpleButcherTable< double > semiImplicitIERK45ButcherTable(bool expl)
Definition: butchertable.cc:326
ROWSimpleButcherTable< double > row2ButcherTable()
Definition: butchertable.cc:341