Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
Tempus_IntegratorObserverBasic_impl.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_IntegratorObserverBasic_impl_hpp
10#define Tempus_IntegratorObserverBasic_impl_hpp
11
12#include "Tempus_Stepper.hpp"
13
14namespace Tempus {
15
16template<class Scalar>
18
19template<class Scalar>
21
22template<class Scalar>
24observeStartIntegrator(const Integrator<Scalar>& integrator){
25
26 std::time_t begin = std::time(nullptr);
27 const Teuchos::RCP<Teuchos::FancyOStream> out = integrator.getOStream();
28 out->setOutputToRootOnly(0);
29 Teuchos::OSTab ostab(out,0,"ScreenOutput");
30 *out << "\nTempus - IntegratorBasic\n"
31 << std::asctime(std::localtime(&begin)) << "\n"
32 << " Stepper = " << integrator.getStepper()->description() << "\n"
33 << " Simulation Time Range [" << integrator.getTimeStepControl()->getInitTime()
34 << ", " << integrator.getTimeStepControl()->getFinalTime() << "]\n"
35 << " Simulation Index Range [" << integrator.getTimeStepControl()->getInitIndex()
36 << ", " << integrator.getTimeStepControl()->getFinalIndex() << "]\n"
37 << "============================================================================\n"
38 << " Step Time dt Abs Error Rel Error Order nFail dCompTime"
39 << std::endl;
40}
41
42template<class Scalar>
44observeStartTimeStep(const Integrator<Scalar>& /* integrator */){}
45
46template<class Scalar>
48observeNextTimeStep(const Integrator<Scalar>& /* integrator */){}
49
50template<class Scalar>
52observeBeforeTakeStep(const Integrator<Scalar>& /* integrator */){}
53
54template<class Scalar>
56observeAfterTakeStep(const Integrator<Scalar>& /* integrator */){}
57
58template<class Scalar>
60observeAfterCheckTimeStep(const Integrator<Scalar>& /* integrator */){}
61
62template<class Scalar>
64observeEndTimeStep(const Integrator<Scalar>& integrator){
65
66 using Teuchos::RCP;
67 auto cs = integrator.getSolutionHistory()->getCurrentState();
68
69 if ((cs->getOutputScreen() == true) ||
70 (cs->getOutput() == true) ||
71 (cs->getTime() == integrator.getTimeStepControl()->getFinalTime())) {
72
73 const Scalar steppertime = integrator.getStepperTimer()->totalElapsedTime();
74 // reset the stepper timer
75 integrator.getStepperTimer()->reset();
76
77 const Teuchos::RCP<Teuchos::FancyOStream> out = integrator.getOStream();
78 out->setOutputToRootOnly(0);
79 Teuchos::OSTab ostab(out, 0, "ScreenOutput");
80 *out<<std::scientific
81 <<std::setw( 6)<<std::setprecision(3)<<cs->getIndex()
82 <<std::setw(11)<<std::setprecision(3)<<cs->getTime()
83 <<std::setw(11)<<std::setprecision(3)<<cs->getTimeStep()
84 <<std::setw(11)<<std::setprecision(3)<<cs->getErrorAbs()
85 <<std::setw(11)<<std::setprecision(3)<<cs->getErrorRel()
86 <<std::fixed <<std::setw( 7)<<std::setprecision(1)<<cs->getOrder()
87 <<std::scientific<<std::setw( 7)<<std::setprecision(3)<<cs->getNFailures()
88 <<std::setw(11)<<std::setprecision(3)<<steppertime
89 <<std::endl;
90 }
91
92}
93
94template<class Scalar>
96observeEndIntegrator(const Integrator<Scalar>& integrator){
97
98 std::string exitStatus;
99 //const Scalar runtime = integrator.getIntegratorTimer()->totalElapsedTime();
100 if (integrator.getSolutionHistory()->getCurrentState()->getSolutionStatus() ==
101 Status::FAILED || integrator.getStatus() == Status::FAILED) {
102 exitStatus = "Time integration FAILURE!";
103 } else {
104 exitStatus = "Time integration complete.";
105 }
106 std::time_t end = std::time(nullptr);
107 const Scalar runtime = integrator.getIntegratorTimer()->totalElapsedTime();
108 const Teuchos::RCP<Teuchos::FancyOStream> out = integrator.getOStream();
109 out->setOutputToRootOnly(0);
110 Teuchos::OSTab ostab(out,0,"ScreenOutput");
111 *out << "============================================================================\n"
112 << " Total runtime = " << runtime << " sec = "
113 << runtime/60.0 << " min\n"
114 << std::asctime(std::localtime(&end))
115 << "\nNumber of Accepted Steps = " << integrator.getSolutionHistory()->getCurrentState()->getIndex()
116 << "\nNumber of Failures = " << integrator.getSolutionHistory()->getCurrentState()->getMetaData()->getNRunningFailures()
117 << "\n"
118 << exitStatus << "\n"
119 << std::endl;
120}
121
122} // namespace Tempus
123#endif // Tempus_IntegratorObserverBasic_impl_hpp
virtual void observeEndIntegrator(const Integrator< Scalar > &integrator) override
Observe the end of the time integrator.
virtual void observeEndTimeStep(const Integrator< Scalar > &integrator) override
Observe the end of the time step loop.
virtual void observeStartIntegrator(const Integrator< Scalar > &integrator) override
Observe the beginning of the time integrator.
virtual void observeBeforeTakeStep(const Integrator< Scalar > &integrator) override
Observe before Stepper takes step.
virtual void observeNextTimeStep(const Integrator< Scalar > &integrator) override
Observe after the next time step size is selected.
virtual void observeAfterTakeStep(const Integrator< Scalar > &integrator) override
Observe after Stepper takes step.
virtual void observeStartTimeStep(const Integrator< Scalar > &integrator) override
Observe the beginning of the time step loop.
virtual void observeAfterCheckTimeStep(const Integrator< Scalar > &integrator) override
Observe after checking time step. Observer can still fail the time step here.