Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Tempus_StepperExplicitRK_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ****************************************************************************
3 // Tempus: Copyright (2017) Sandia Corporation
4 //
5 // Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
6 // ****************************************************************************
7 // @HEADER
8 
9 #ifndef Tempus_StepperExplicitRK_decl_hpp
10 #define Tempus_StepperExplicitRK_decl_hpp
11 
12 #include "Tempus_config.hpp"
13 #include "Tempus_StepperExplicit.hpp"
16 
17 
18 namespace Tempus {
19 
20 /** \brief Explicit Runge-Kutta time stepper.
21  *
22  * For the explicit ODE system,
23  * \f[
24  * \dot{x} = \bar{f}(x,t),
25  * \f]
26  * the general explicit Runge-Kutta method for \f$s\f$-stages can be
27  * written as
28  * \f[
29  * X_{i} = x_{n-1}
30  * + \Delta t\,\sum_{j=1}^{i-1} a_{ij}\,\bar{f}(X_{j},t_{n-1}+c_{j}\Delta t)
31  * \f]
32  * \f[
33  * x_{n} = x_{n-1}
34  * + \Delta t\,\sum_{i=1}^{s}b_{i}\,\bar{f}(X_{i},t_{n-1}+c_{i}\Delta t)
35  * \f]
36  * where \f$X_{i}\f$ are intermediate approximations to the solution
37  * at times, \f$t_{n-1}+c_{i}\Delta t\f$, (stage solutions) which may
38  * be correct to a lower order of accuracy than the solution, \f$x_{n}\f$.
39  * We should note that these lower-order approximations are combined
40  * through \f$b_{i}\f$ so that error terms cancel out and produce a
41  * more accurate solution. Note for explicit RK that \f$a_{ij}=0\f$ for
42  * \f$j \leq i\f$ and does not require any solves.
43  * Note that the stage time derivatives are
44  * \f[
45  * \dot{X}_{i} = \bar{f}(X_{i},t_{n-1}+c_{i}\Delta t),
46  * \f]
47  * and the time derivative by definition is
48  * \f[
49  * \dot{x}_{n} = \bar{f}(x_{n},t_{n}),
50  * \f]
51  *
52  * <b> Algorithm </b>
53  * The single-timestep algorithm for Explicit RK is simply,
54  * - for \f$i = 1 \ldots s\f$ do
55  * - \f$X_i \leftarrow x_{n-1}
56  * + \Delta t\,\sum_{j=1}^{i-1} a_{ij}\,\dot{X}_j\f$
57  * - Evaluate \f$\bar{f}(X_{i},t_{n-1}+c_{i}\Delta t)\f$
58  * - \f$\dot{X}_i \leftarrow \bar{f}(X_i,t_{n-1}+c_i\Delta t)\f$
59  * - end for
60  * - \f$x_n \leftarrow x_{n-1} + \Delta t\,\sum_{i=1}^{s}b_i\,\dot{X}_i\f$
61  *
62  * When using the First-Step-As-Last (FSAL) priniciple, where one can
63  * reuse the last function evaulation as the first evaluation of the next
64  * time step, the algorithm is only slightly more complicated.
65  * - for \f$i = 1 \ldots s\f$ do
66  * - if ( i==1 && useFSAL && (previous step not failed) )
67  * - tmp = \f$\dot{X}_1\f$
68  * - \f$\dot{X}_1 = \dot{X}_s\f$
69  * - \f$\dot{X}_s\f$ = tmp
70  * - else
71  * - \f$X_i \leftarrow x_{n-1}
72  * + \Delta t\,\sum_{j=1}^{i-1} a_{ij}\,\dot{X}_j\f$
73  * - Evaluate \f$\bar{f}(X_{i},t_{n-1}+c_{i}\Delta t)\f$
74  * - \f$\dot{X}_i \leftarrow \bar{f}(X_i,t_{n-1}+c_i\Delta t)\f$
75  * - end for
76  * - \f$x_n \leftarrow x_{n-1} + \Delta t\,\sum_{i=1}^{s}b_i\,\dot{X}_i\f$
77  *
78  * For Explicit RK, FSAL requires \f$c_1 = 0\f$, \f$c_s = 1\f$, and
79  * be stiffly accurate (\f$a_{sj} = b_j\f$). An example of this is
80  * the Bogacki-Shampine 3(2) method.
81  * \f[
82  * \begin{array}{c|cccc} 0 & 0 & & & \\
83  * 1/3 & 1/2 & 0 & & \\
84  * 2/3 & 0 & 3/4 & 0 & \\
85  * 1 & 2/9 & 1/3 & 4/9 & 0 \\ \hline
86  * & 2/9 & 1/3 & 4/9 & 0
87  * & 7/24 & 1/4 & 1/3 & 1/8 \end{array}
88  * \f]
89  */
90 template<class Scalar>
91 class StepperExplicitRK : virtual public Tempus::StepperExplicit<Scalar>
92 {
93 public:
94 
95  /** \brief Default constructor.
96  *
97  * - Constructs with a default ParameterList.
98  * - Can reset ParameterList with setParameterList().
99  * - Requires subsequent setModel() and initialize() calls before calling
100  * takeStep().
101  */
103 
104  /// Constructor to specialize Stepper parameters.
106  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel,
107  Teuchos::RCP<Teuchos::ParameterList> pList);
108 
109  /// Constructor to use default Stepper parameters.
111  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel,
112  std::string stepperType = "RK Explicit 4 Stage");
113 
114  /// Constructor for StepperFactory.
116  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel,
117  std::string stepperType, Teuchos::RCP<Teuchos::ParameterList> pList);
118 
119  /// \name Basic stepper methods
120  //@{
121  virtual void setObserver(
122  Teuchos::RCP<StepperObserver<Scalar> > obs = Teuchos::null);
123 
124  virtual void setTableau(std::string stepperType);
125 
126  virtual void setTableau(
127  Teuchos::RCP<Teuchos::ParameterList> pList = Teuchos::null);
128 
129  virtual void setTableau(
130  Teuchos::RCP<const RKButcherTableau<Scalar> > ERK_ButcherTableau);
131 
132  /// Initialize during construction and after changing input parameters.
133  virtual void initialize();
134 
135  /// Set the initial conditions and make them consistent.
136  virtual void setInitialConditions (
137  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory);
138 
139  /// Take the specified timestep, dt, and return true if successful.
140  virtual void takeStep(
141  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory);
142 
143  virtual std::string getStepperType() const
144  { return this->stepperPL_->template get<std::string>("Stepper Type"); }
145 
146  /// Get a default (initial) StepperState
147  virtual Teuchos::RCP<Tempus::StepperState<Scalar> > getDefaultStepperState();
148  virtual Scalar getOrder() const {return ERK_ButcherTableau_->order();}
149  virtual Scalar getOrderMin() const {return ERK_ButcherTableau_->orderMin();}
150  virtual Scalar getOrderMax() const {return ERK_ButcherTableau_->orderMax();}
151  virtual Scalar getInitTimeStep(
152  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory) const;
153 
154  virtual bool isExplicit() const {return true;}
155  virtual bool isImplicit() const {return false;}
156  virtual bool isExplicitImplicit() const
157  {return isExplicit() and isImplicit();}
158  virtual bool isOneStepMethod() const {return true;}
159  virtual bool isMultiStepMethod() const {return !isOneStepMethod();}
160 
161  virtual OrderODE getOrderODE() const {return FIRST_ORDER_ODE;}
162  //@}
163 
164  /// \name ParameterList methods
165  //@{
166  void setParameterList(const Teuchos::RCP<Teuchos::ParameterList> & pl);
167  Teuchos::RCP<Teuchos::ParameterList> getNonconstParameterList();
168  Teuchos::RCP<Teuchos::ParameterList> unsetParameterList();
169  Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const;
170  Teuchos::RCP<Teuchos::ParameterList> getDefaultParameters() const;
171  //@}
172 
173  /// \name Overridden from Teuchos::Describable
174  //@{
175  virtual std::string description() const;
176  virtual void describe(Teuchos::FancyOStream & out,
177  const Teuchos::EVerbosityLevel verbLevel) const;
178  //@}
179 
180 protected:
181 
182  Teuchos::RCP<const RKButcherTableau<Scalar> > ERK_ButcherTableau_;
183 
184  std::vector<Teuchos::RCP<Thyra::VectorBase<Scalar> > > stageXDot_;
185  Teuchos::RCP<Thyra::VectorBase<Scalar> > stageX_;
186 
187  Teuchos::RCP<StepperExplicitRKObserver<Scalar> > stepperExplicitRKObserver_;
188 
189  // For Embedded RK
190  Teuchos::RCP<Thyra::VectorBase<Scalar> > ee_;
191  Teuchos::RCP<Thyra::VectorBase<Scalar> > abs_u0;
192  Teuchos::RCP<Thyra::VectorBase<Scalar> > abs_u;
193  Teuchos::RCP<Thyra::VectorBase<Scalar> > sc;
194 
195 };
196 
197 } // namespace Tempus
198 
199 #endif // Tempus_StepperExplicitRK_decl_hpp
Explicit Runge-Kutta time stepper.
Teuchos::RCP< Teuchos::ParameterList > getDefaultParameters() const
Teuchos::RCP< StepperExplicitRKObserver< Scalar > > stepperExplicitRKObserver_
virtual Scalar getInitTimeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory) const
Teuchos::RCP< const RKButcherTableau< Scalar > > ERK_ButcherTableau_
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
virtual std::string description() const
std::vector< Teuchos::RCP< Thyra::VectorBase< Scalar > > > stageXDot_
Teuchos::RCP< Thyra::VectorBase< Scalar > > ee_
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
virtual void setInitialConditions(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Set the initial conditions and make them consistent.
Teuchos::RCP< Teuchos::ParameterList > stepperPL_
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
virtual void setObserver(Teuchos::RCP< StepperObserver< Scalar > > obs=Teuchos::null)
Set Observer.
Teuchos::RCP< Thyra::VectorBase< Scalar > > sc
StepperObserver class for Stepper class.
Teuchos::RCP< SolutionHistory< Scalar > > solutionHistory(Teuchos::RCP< Teuchos::ParameterList > pList=Teuchos::null)
Nonmember constructor.
Teuchos::RCP< Thyra::VectorBase< Scalar > > stageX_
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
Stepper integrates first-order ODEs.
virtual void takeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Take the specified timestep, dt, and return true if successful.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &pl)
virtual void initialize()
Initialize during construction and after changing input parameters.
virtual void setTableau(std::string stepperType)
virtual Teuchos::RCP< Tempus::StepperState< Scalar > > getDefaultStepperState()
Get a default (initial) StepperState.
Thyra Base interface for implicit time steppers.
Teuchos::RCP< Thyra::VectorBase< Scalar > > abs_u0
virtual std::string getStepperType() const
Teuchos::RCP< Thyra::VectorBase< Scalar > > abs_u