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_StepperRKBase.hpp"
14 #include "Tempus_StepperExplicit.hpp"
15 #ifndef TEMPUS_HIDE_DEPRECATED_CODE
16  #include "Tempus_StepperRKObserverComposite.hpp"
17 #endif
18 
19 
20 namespace Tempus {
21 
22 /** \brief Explicit Runge-Kutta time stepper.
23  *
24  * For the explicit ODE system,
25  * \f[
26  * \dot{x} = \bar{f}(x,t),
27  * \f]
28  * the general explicit Runge-Kutta method for \f$s\f$-stages can be
29  * written as
30  * \f[
31  * X_{i} = x_{n-1}
32  * + \Delta t\,\sum_{j=1}^{i-1} a_{ij}\,\bar{f}(X_{j},t_{n-1}+c_{j}\Delta t)
33  * \f]
34  * \f[
35  * x_{n} = x_{n-1}
36  * + \Delta t\,\sum_{i=1}^{s}b_{i}\,\bar{f}(X_{i},t_{n-1}+c_{i}\Delta t)
37  * \f]
38  * where \f$X_{i}\f$ are intermediate approximations to the solution
39  * at times, \f$t_{n-1}+c_{i}\Delta t\f$, (stage solutions) which may
40  * be correct to a lower order of accuracy than the solution, \f$x_{n}\f$.
41  * We should note that these lower-order approximations are combined
42  * through \f$b_{i}\f$ so that error terms cancel out and produce a
43  * more accurate solution. Note for explicit RK that \f$a_{ij}=0\f$ for
44  * \f$j \leq i\f$ and does not require any solves.
45  * Note that the stage time derivatives are
46  * \f[
47  * \dot{X}_{i} = \bar{f}(X_{i},t_{n-1}+c_{i}\Delta t),
48  * \f]
49  * and the time derivative by definition is
50  * \f[
51  * \dot{x}_{n} = \bar{f}(x_{n},t_{n}),
52  * \f]
53  *
54  * <b> Algorithm </b>
55  * The single-timestep algorithm for Explicit RK is
56  *
57  * \f{algorithm}{
58  * \renewcommand{\thealgorithm}{}
59  * \caption{Explicit RK with the application-action locations indicated.}
60  * \begin{algorithmic}[1]
61  * \State {\it appAction.execute(solutionHistory, stepper, BEGIN\_STEP)}
62  * \State $X \leftarrow x_{n-1}$ \Comment Set initial guess to last timestep.
63  * \For {$i = 0 \ldots s-1$}
64  * \State {\it appAction.execute(solutionHistory, stepper, BEGIN\_STAGE)}
65  * \State {\it appAction.execute(solutionHistory, stepper, BEFORE\_SOLVE)}
66  * \If { i==0 and useFSAL and (previous step not failed) }
67  * \State tmp = $\dot{X}_0$
68  * \State $\dot{X}_0 = \dot{X}_s$
69  * \State $\dot{X}_s$ = tmp
70  * \State {\bf continue}
71  * \Else
72  * \State $X \leftarrow x_{n-1}
73  * + \Delta t\,\sum_{j=1}^{i-1} a_{ij}\,\dot{X}_j$
74  * \State {\it appAction.execute(solutionHistory, stepper, AFTER\_SOLVE)}
75  * \State {\it appAction.execute(solutionHistory, stepper, BEFORE\_EXPLICIT\_EVAL)}
76  * \State $\dot{X}_i \leftarrow \bar{f}(X_i,t_{n-1}+c_i\Delta t)$
77  * \EndIf
78  * \State {\it appAction.execute(solutionHistory, stepper, END\_STAGE)}
79  * \EndFor
80  * \State $x_n \leftarrow x_{n-1} + \Delta t\,\sum_{i=1}^{s}b_i\,\dot{X}_i$
81  * \State {\it appAction.execute(solutionHistory, stepper, END\_STEP)}
82  * \end{algorithmic}
83  * \f}
84  *
85  * For Explicit RK, FSAL requires \f$c_1 = 0\f$, \f$c_s = 1\f$, and
86  * be stiffly accurate (\f$a_{sj} = b_j\f$). An example of this is
87  * the Bogacki-Shampine 3(2) method.
88  * \f[
89  * \begin{array}{c|cccc} 0 & 0 & & & \\
90  * 1/3 & 1/2 & 0 & & \\
91  * 2/3 & 0 & 3/4 & 0 & \\
92  * 1 & 2/9 & 1/3 & 4/9 & 0 \\ \hline
93  * & 2/9 & 1/3 & 4/9 & 0 \\
94  * & 7/24 & 1/4 & 1/3 & 1/8 \end{array}
95  * \f]
96  */
97 template<class Scalar>
98 class StepperExplicitRK : virtual public Tempus::StepperExplicit<Scalar>,
99  virtual public Tempus::StepperRKBase<Scalar>
100 {
101 
102 public:
103 
104  /// \name Basic stepper methods
105  //@{
106 #ifndef TEMPUS_HIDE_DEPRECATED_CODE
107  virtual void setObserver(
108  Teuchos::RCP<StepperObserver<Scalar> > obs = Teuchos::null);
109 
110  virtual Teuchos::RCP<StepperObserver<Scalar> > getObserver() const
111  { return this->stepperObserver_; }
112 #endif
113  /// Initialize during construction and after changing input parameters.
114  virtual void initialize();
115 
116  /// Set the initial conditions and make them consistent.
117  virtual void setInitialConditions (
118  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory);
119 
120  /// Take the specified timestep, dt, and return true if successful.
121  virtual void takeStep(
122  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory);
123 
124  /// Get a default (initial) StepperState
125  virtual Teuchos::RCP<Tempus::StepperState<Scalar> > getDefaultStepperState();
126  virtual Scalar getInitTimeStep(
127  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory) const;
128 
129  virtual bool isExplicit() const {return true;}
130  virtual bool isImplicit() const {return false;}
131  virtual bool isExplicitImplicit() const
132  {return isExplicit() and isImplicit();}
133  virtual bool isOneStepMethod() const {return true;}
134  virtual bool isMultiStepMethod() const {return !isOneStepMethod();}
135 
136  virtual OrderODE getOrderODE() const {return FIRST_ORDER_ODE;}
137 
138  void getValidParametersBasicERK(Teuchos::RCP<Teuchos::ParameterList> pl) const;
139  virtual std::string getDescription() const = 0;
140  //@}
141 
142  Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const;
143 
144  /// \name Overridden from Teuchos::Describable
145  //@{
146  virtual void describe(Teuchos::FancyOStream & out,
147  const Teuchos::EVerbosityLevel verbLevel) const;
148  //@}
149 
150  virtual bool isValidSetup(Teuchos::FancyOStream & out) const;
151 
152  /// \name Accessors methods
153  //@{
154  /** \brief Use embedded if avialable. */
155  virtual void setUseEmbedded(bool a) { useEmbedded_ = a; }
156  virtual bool getUseEmbedded() const { return useEmbedded_; }
157  virtual bool getUseEmbeddedDefault() const { return false; }
158  //@}
159 
160 
161 protected:
162 
163  /// Default setup for constructor.
164  virtual void setupDefault();
165 
166  /// Setup for constructor.
167 #ifndef TEMPUS_HIDE_DEPRECATED_CODE
168  virtual void setup(
169  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel,
170  const Teuchos::RCP<StepperRKObserverComposite<Scalar> >& obs,
171  bool useFSAL,
172  std::string ICConsistency,
173  bool ICConsistencyCheck,
174  bool useEmbedded);
175 #endif
176  virtual void setup(
177  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel,
178  bool useFSAL,
179  std::string ICConsistency,
180  bool ICConsistencyCheck,
181  bool useEmbedded,
182  const Teuchos::RCP<StepperRKAppAction<Scalar> >& stepperRKAppAction);
183 
184  virtual void setupTableau() = 0;
185 
186 
187  std::vector<Teuchos::RCP<Thyra::VectorBase<Scalar> > > stageXDot_;
188 
189 #ifndef TEMPUS_HIDE_DEPRECATED_CODE
190  Teuchos::RCP<StepperRKObserverComposite<Scalar> > stepperObserver_;
191 #endif
192 
193  // For Embedded RK
195  Teuchos::RCP<Thyra::VectorBase<Scalar> > ee_;
196  Teuchos::RCP<Thyra::VectorBase<Scalar> > abs_u0;
197  Teuchos::RCP<Thyra::VectorBase<Scalar> > abs_u;
198  Teuchos::RCP<Thyra::VectorBase<Scalar> > sc;
199 
200 };
201 
202 } // namespace Tempus
203 
204 #endif // Tempus_StepperExplicitRK_decl_hpp
Explicit Runge-Kutta time stepper.
virtual std::string getDescription() const =0
virtual void setupDefault()
Default setup for constructor.
virtual Scalar getInitTimeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory) const
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
std::vector< Teuchos::RCP< Thyra::VectorBase< Scalar > > > stageXDot_
Teuchos::RCP< Thyra::VectorBase< Scalar > > ee_
virtual bool isValidSetup(Teuchos::FancyOStream &out) const
virtual void setup(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< StepperRKObserverComposite< Scalar > > &obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded)
Setup for constructor.
virtual void setInitialConditions(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Set the initial conditions and make them consistent.
Base class for Runge-Kutta methods, ExplicitRK, DIRK and IMEX.
virtual void setupTableau()=0
Application Action for StepperRKBase.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
virtual Teuchos::RCP< StepperObserver< Scalar > > getObserver() const
Get Observer.
virtual void setObserver(Teuchos::RCP< StepperObserver< Scalar > > obs=Teuchos::null)
Set Observer.
Teuchos::RCP< Thyra::VectorBase< Scalar > > sc
StepperObserver class for Stepper class.
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
void getValidParametersBasicERK(Teuchos::RCP< Teuchos::ParameterList > pl) const
Stepper integrates first-order ODEs.
virtual void takeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Take the specified timestep, dt, and return true if successful.
virtual void initialize()
Initialize during construction and after changing input parameters.
Teuchos::RCP< StepperRKObserverComposite< Scalar > > stepperObserver_
virtual Teuchos::RCP< Tempus::StepperState< Scalar > > getDefaultStepperState()
Get a default (initial) StepperState.
virtual void setUseEmbedded(bool a)
Use embedded if avialable.
Thyra Base interface for implicit time steppers.
Teuchos::RCP< Thyra::VectorBase< Scalar > > abs_u0
Teuchos::RCP< Thyra::VectorBase< Scalar > > abs_u