Rythmos - Transient Integration for Differential Equations  Version of the Day
 All Classes Functions Variables Typedefs Pages
Rythmos_IntegrationObserverBase.hpp
1 
2 #ifndef RYTHMOS_INTEGRATION_OBSERVER_BASE_HPP
3 #define RYTHMOS_INTEGRATION_OBSERVER_BASE_HPP
4 
5 #include "Rythmos_Types.hpp"
6 #include "Teuchos_Describable.hpp"
7 #include "Teuchos_VerboseObject.hpp"
8 
9 namespace Rythmos {
10 
11 
12 // Forwards
13 template<class Scalar> class TimeRange;
14 template<class Scalar> class StepControlInfo;
15 template<class Scalar> class StepperBase;
16 
17 
23 template<class Scalar>
25  : virtual public Teuchos::Describable,
26  virtual public Teuchos::VerboseObject<IntegrationObserverBase<Scalar> >
27 {
28 public:
29 
35  virtual RCP<IntegrationObserverBase<Scalar> >
36  cloneIntegrationObserver() const = 0;
37 
49  virtual void resetIntegrationObserver(
50  const TimeRange<Scalar> &integrationTimeDomain
51  // ToDo: Pass in the initial condition to the observer
52  ) = 0;
53 
68  virtual void observeStartTimeIntegration(
69  const StepperBase<Scalar> &stepper);
70 
85  virtual void observeEndTimeIntegration(
86  const StepperBase<Scalar> &stepper);
87 
88  // ToDo: add observeStartTimeStep(stepper, ...)
113  virtual void observeStartTimeStep(
114  const StepperBase<Scalar> &stepper,
115  const StepControlInfo<Scalar> &stepCtrlInfo,
116  const int timeStepIter
117  );
118 
150  virtual void observeCompletedTimeStep(
151  const StepperBase<Scalar> &stepper,
152  const StepControlInfo<Scalar> &stepCtrlInfo,
153  const int timeStepIter
154  ) = 0;
155 
191  virtual void observeFailedTimeStep(
192  const StepperBase<Scalar> &stepper,
193  const StepControlInfo<Scalar> &stepCtrlInfo,
194  const int timeStepIter
195  );
196 
197 };
198 
199 
205 template<class Scalar>
206 bool isInitialTimeStep(
207  const TimeRange<Scalar> &currentTimeRange,
208  const TimeRange<Scalar> &fullTimeRange
209  )
210 {
211  typedef Teuchos::ScalarTraits<Scalar> ST;
212  return compareTimeValues(currentTimeRange.lower(), fullTimeRange.lower()) == ST::zero();
213 }
214 
215 
221 template<class Scalar>
222 bool isFinalTimeStep(
223  const TimeRange<Scalar> &currentTimeRange,
224  const TimeRange<Scalar> &fullTimeRange
225  )
226 {
227  typedef Teuchos::ScalarTraits<Scalar> ST;
228  return compareTimeValues(currentTimeRange.upper(), fullTimeRange.upper()) >= ST::zero();
229 }
230 
231 
232 // /////////////////////////////////////////////////////////////
233 /* Default implementations for backwards compatibility
234 
235  Some of the observer methods were added after rythmos was released.
236  Even though all methods should be pure virtual, we provide a
237  default implementation here for the recently added methods to
238  maintain backwards compatibility. These should be removed in the
239  future.
240 */
241 
242 template<class Scalar>
245 {
246 
247 }
248 
249 template<class Scalar>
252 {
253 
254 }
255 
256 template<class Scalar>
259  const StepperBase<Scalar> &/* stepper */,
260  const StepControlInfo<Scalar> &/* stepCtrlInfo */,
261  const int /* timeStepIter */
262  )
263 {
264 
265 }
266 
267 template<class Scalar>
270  const StepperBase<Scalar> &/* stepper */,
271  const StepControlInfo<Scalar> &/* stepCtrlInfo */,
272  const int /* timeStepIter */
273  )
274 {
275 
276 }
277 
278 
279 } // namespace Rythmos
280 
281 
282 #endif // RYTHMOS_INTEGRATION_OBSERVER_BASE_HPP
Base class for defining stepper functionality.
virtual void observeEndTimeIntegration(const StepperBase< Scalar > &stepper)
Observe the end of a time integration loop.
Simple struct to aggregate integration/stepper control information.
virtual void observeFailedTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
Observer a failed integration step.
Base class for strategy objects that observe and time integration by observing the stepper object...
virtual void observeStartTimeIntegration(const StepperBase< Scalar > &stepper)
Observe the beginning of a time integration loop.
virtual void resetIntegrationObserver(const TimeRange< Scalar > &integrationTimeDomain)=0
Reset the observer to prepair it to observe another integration.
virtual void observeStartTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
Observer the beginning of an integration step.
virtual RCP< IntegrationObserverBase< Scalar > > cloneIntegrationObserver() const =0
Clone this integration observer if supported .
virtual void observeCompletedTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)=0
Observe a successfully completed integration step.