dune-fem  2.8-git
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 #include <algorithm>
5 
6 #include <dune/common/dynmatrix.hh>
7 #include <dune/common/dynvector.hh>
8 
9 namespace DuneODE
10 {
11 
12  // SimpleButcherTable
13  // ------------------
14 
15  template< class Field >
17  {
19 
20  public:
21  typedef Field FieldType;
22 
23  SimpleButcherTable ( int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c )
24  : stages_( stages ), order_( order ),
25  a_( a ), b_( b ), c_( c )
26  {}
27 
28  Dune::DynamicMatrix< FieldType > A () const { return makeMatrix( stages_, stages_, a_ ); }
29  Dune::DynamicVector< FieldType > b () const { return makeVector( stages_, b_ ); }
30  Dune::DynamicVector< FieldType > c () const { return makeVector( stages_, c_ ); }
31 
32  int order () const { return order_; }
33  int stages () const { return stages_; }
34 
35  protected:
36  static Dune::DynamicMatrix< FieldType > makeMatrix ( int m, int n, const FieldType *data )
37  {
38  Dune::DynamicMatrix< FieldType > A( m, n );
39  for( int i = 0; i < m; ++i )
40  std::copy( data + i*n, data + (i+1)*n, A[ i ].begin() );
41  return A;
42  }
43 
44  static Dune::DynamicVector< FieldType > makeVector ( int n, const FieldType *data )
45  {
46  Dune::DynamicVector< FieldType > v( n );
47  std::copy( data, data + n, v.begin() );
48  return v;
49  }
50 
52  const FieldType *a_, *b_, *c_;
53  };
54 
55 
56 
57  // explicit butcher tables
58  // -----------------------
64 
65  // implicit butcher tables
66  // -----------------------
67 
72 
73 
74 
75  // semiimplicit butcher tables
76  // ---------------------------
77 
85 
86  // ROW butcher tables
87  // -----------------------
88 
89  template< class Field >
91  {
94 
95  public:
96  typedef Field FieldType;
97 
98  ROWSimpleButcherTable ( int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c, const FieldType *a2 )
99  : Base(stages,order,a,b,c)
100  , a2_(a2)
101  {}
102 
103  Dune::DynamicMatrix< FieldType > B () const { return Base::makeMatrix( stages_, stages_, a2_ ); }
104 
105  private:
106  using Base::stages_;
107  const FieldType *a2_;
108  };
109  ROWSimpleButcherTable< double > row2ButcherTable ();
110  ROWSimpleButcherTable< double > row3ButcherTable ();
111 } // namespace DuneODE
112 
113 #endif // #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_BUTCHERTABLE_HH
Definition: multistep.hh:17
SimpleButcherTable< double > implicit34ButcherTable()
Definition: butchertable.cc:121
SimpleButcherTable< double > expl6ButcherTable()
Definition: butchertable.cc:76
SimpleButcherTable< double > semiImplicit33ButcherTable(bool expl)
Definition: butchertable.cc:256
SimpleButcherTable< double > implicitEulerButcherTable()
Definition: butchertable.cc:155
SimpleButcherTable< double > semiImplicitEulerButcherTable(bool expl)
Definition: butchertable.cc:190
SimpleButcherTable< double > explicitEulerButcherTable()
Definition: butchertable.cc:16
ROWSimpleButcherTable< double > row3ButcherTable()
Definition: butchertable.cc:371
SimpleButcherTable< double > implicit3ButcherTable()
Definition: butchertable.cc:141
SimpleButcherTable< double > tvd3ButcherTable()
Definition: butchertable.cc:43
SimpleButcherTable< double > semiImplicitARK34ButcherTable(bool expl)
SimpleButcherTable< double > tvd2ButcherTable()
Definition: butchertable.cc:29
SimpleButcherTable< double > semiImplicit23ButcherTable(bool expl)
Definition: butchertable.cc:221
SimpleButcherTable< double > semiImplicitIERK45ButcherTable(bool expl)
Definition: butchertable.cc:335
SimpleButcherTable< double > semiImplicitSSP222ButcherTable(bool expl)
Definition: butchertable.cc:286
ROWSimpleButcherTable< double > row2ButcherTable()
Definition: butchertable.cc:350
SimpleButcherTable< double > rk4ButcherTable()
Definition: butchertable.cc:58
SimpleButcherTable< double > gauss2ButcherTable()
Definition: butchertable.cc:170
SimpleButcherTable< double > semiImplicitARK46ButcherTable(bool expl)
Definition: butchertable.hh:17
Dune::DynamicMatrix< FieldType > A() const
Definition: butchertable.hh:28
static Dune::DynamicVector< FieldType > makeVector(int n, const FieldType *data)
Definition: butchertable.hh:44
const FieldType * b_
Definition: butchertable.hh:52
int order_
Definition: butchertable.hh:51
static Dune::DynamicMatrix< FieldType > makeMatrix(int m, int n, const FieldType *data)
Definition: butchertable.hh:36
const FieldType * a_
Definition: butchertable.hh:52
SimpleButcherTable(int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c)
Definition: butchertable.hh:23
const FieldType * c_
Definition: butchertable.hh:52
Field FieldType
Definition: butchertable.hh:21
int order() const
Definition: butchertable.hh:32
Dune::DynamicVector< FieldType > b() const
Definition: butchertable.hh:29
int stages_
Definition: butchertable.hh:51
int stages() const
Definition: butchertable.hh:33
Dune::DynamicVector< FieldType > c() const
Definition: butchertable.hh:30
Definition: butchertable.hh:91
Dune::DynamicMatrix< FieldType > B() const
Definition: butchertable.hh:103
ROWSimpleButcherTable(int stages, int order, const FieldType *a, const FieldType *b, const FieldType *c, const FieldType *a2)
Definition: butchertable.hh:98
Field FieldType
Definition: butchertable.hh:96