Rythmos - Transient Integration for Differential Equations Version of the Day
Loading...
Searching...
No Matches
Rythmos_LoggingIntegrationObserver.hpp
1//@HEADER
2// ***********************************************************************
3//
4// Rythmos Package
5// Copyright (2006) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// This library is free software; you can redistribute it and/or modify
11// it under the terms of the GNU Lesser General Public License as
12// published by the Free Software Foundation; either version 2.1 of the
13// License, or (at your option) any later version.
14//
15// This library is distributed in the hope that it will be useful, but
16// WITHOUT ANY WARRANTY; without even the implied warranty of
17// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18// Lesser General Public License for more details.
19//
20// You should have received a copy of the GNU Lesser General Public
21// License along with this library; if not, write to the Free Software
22// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23// USA
24// Questions? Contact Todd S. Coffey (tscoffe@sandia.gov)
25//
26// ***********************************************************************
27//@HEADER
28
29#ifndef Rythmos_LOGGING_INTEGRATION_OBSERVER_HPP
30#define Rythmos_LOGGING_INTEGRATION_OBSERVER_HPP
31
32#include "Rythmos_IntegrationObserverBase.hpp"
33#include "Teuchos_RCP.hpp"
34#include <map>
35#include <list>
36#include <string>
37
38namespace Rythmos {
39
40
43template<class Scalar>
45{
46public:
47
49
50 void resetLogCounters();
51
52 Teuchos::RCP<const std::map<std::string,int> > getCounters();
53
54 Teuchos::RCP<const std::list<std::string> > getOrder();
55
58
59 RCP<IntegrationObserverBase<Scalar> > cloneIntegrationObserver() const;
60
61 void
62 resetIntegrationObserver(const TimeRange<Scalar> &integrationTimeDomain);
63
64 void observeStartTimeIntegration(const StepperBase<Scalar> &stepper);
65
66 void observeEndTimeIntegration(const StepperBase<Scalar> &stepper);
67
69 const StepperBase<Scalar> &stepper,
70 const StepControlInfo<Scalar> &stepCtrlInfo,
71 const int timeStepIter
72 );
73
75 const StepperBase<Scalar> &stepper,
76 const StepControlInfo<Scalar> &stepCtrlInfo,
77 const int timeStepIter
78 );
79
81 const StepperBase<Scalar> &stepper,
82 const StepControlInfo<Scalar> &stepCtrlInfo,
83 const int timeStepIter
84 );
85
87
93
94 const std::string nameCloneIntegrationObserver_;
95 const std::string nameResetIntegrationObserver_;
96 const std::string nameObserveStartTimeIntegration_;
97 const std::string nameObserveEndTimeIntegration_;
98 const std::string nameObserveStartTimeStep_;
99 const std::string nameObserveCompletedTimeStep_;
100 const std::string nameObserveFailedTimeStep_;
101
103
104
105private:
106
112 void logCall(const std::string call) const;
113
114private:
115
116 Teuchos::RCP< std::map<std::string,int> > counters_;
117 Teuchos::RCP< std::list<std::string> > order_;
118
119};
120
121
126template<class Scalar>
127Teuchos::RCP<LoggingIntegrationObserver<Scalar> >
129{
130 const Teuchos::RCP<LoggingIntegrationObserver<Scalar> > lio =
131 Teuchos::rcp(new LoggingIntegrationObserver<Scalar>);
132
133 return lio;
134}
135
136
137// //////////////////////////////////////////////////////
138// Implementations
139
140template<typename Scalar>
142 nameCloneIntegrationObserver_("cloneIntegrationObserver"),
143 nameResetIntegrationObserver_("resetIntegrationObserver"),
144 nameObserveStartTimeIntegration_("observeStartTimeIntegration"),
145 nameObserveEndTimeIntegration_("observeEndTimeIntegration"),
146 nameObserveStartTimeStep_("observeStartTimeStep"),
147 nameObserveCompletedTimeStep_("observeCompletedTimeStep"),
148 nameObserveFailedTimeStep_("observeFailedTimeStep")
149{
150 counters_ = Teuchos::rcp(new std::map<std::string,int>);
151 order_ = Teuchos::rcp(new std::list<std::string>);
152 this->resetLogCounters();
153}
154
155template<typename Scalar>
156void LoggingIntegrationObserver<Scalar>::
157resetLogCounters()
158{
159 (*counters_)[nameCloneIntegrationObserver_] = 0;
160 (*counters_)[nameResetIntegrationObserver_] = 0;
161 (*counters_)[nameObserveStartTimeIntegration_] = 0;
162 (*counters_)[nameObserveEndTimeIntegration_] = 0;
163 (*counters_)[nameObserveStartTimeStep_] = 0;
164 (*counters_)[nameObserveCompletedTimeStep_] = 0;
165 (*counters_)[nameObserveFailedTimeStep_] = 0;
166 order_->clear();
167}
168
169template<typename Scalar>
170RCP<IntegrationObserverBase<Scalar> >
172{
173 logCall(nameCloneIntegrationObserver_);
174 Teuchos::RCP<IntegrationObserverBase<Scalar> > observer =
175 Teuchos::rcp(new LoggingIntegrationObserver<Scalar>(*this));
176 return observer;
177}
178
179template<typename Scalar>
180void
182resetIntegrationObserver(const TimeRange<Scalar> &integrationTimeDomain)
183{
184 logCall(nameResetIntegrationObserver_);
185}
186
187template<typename Scalar>
189observeStartTimeIntegration(const StepperBase<Scalar> &stepper)
190{
191 logCall(nameObserveStartTimeIntegration_);
192}
193
194template<typename Scalar>
196observeEndTimeIntegration(const StepperBase<Scalar> &stepper)
197{
198 logCall(nameObserveEndTimeIntegration_);
199}
200
201template<typename Scalar>
203 const StepperBase<Scalar> &stepper,
204 const StepControlInfo<Scalar> &stepCtrlInfo,
205 const int timeStepIter
206 )
207{
208 logCall(nameObserveStartTimeStep_);
209}
210
211template<typename Scalar>
213 const StepperBase<Scalar> &stepper,
214 const StepControlInfo<Scalar> &stepCtrlInfo,
215 const int timeStepIter
216 )
217{
218 logCall(nameObserveCompletedTimeStep_);
219}
220
221template<typename Scalar>
223 const StepperBase<Scalar> &stepper,
224 const StepControlInfo<Scalar> &stepCtrlInfo,
225 const int timeStepIter
226 )
227{
228 logCall(nameObserveFailedTimeStep_);
229}
230
231template<typename Scalar>
232RCP<const std::map<std::string,int> > LoggingIntegrationObserver<Scalar>::getCounters()
233{
234 return counters_;
235}
236
237template<typename Scalar>
238RCP<const std::list<std::string> > LoggingIntegrationObserver<Scalar>::getOrder()
239{
240 return order_;
241}
242
243
244template<typename Scalar>
245void LoggingIntegrationObserver<Scalar>::logCall(const std::string call) const
246{
247 (*counters_)[call] += 1;
248 order_->push_back(call);
249}
250
251
252} // namespace Rythmos
253
254
255#endif //Rythmos_LOGGING_INTEGRATION_OBSERVER_HPP
Base class for strategy objects that observe and time integration by observing the stepper object.
Logging IntegrationOberserver that counts calls to observer functions and lists their order.
Teuchos::RCP< LoggingIntegrationObserver< Scalar > > createLoggingIntegrationObserver()
Nonmember constructor.
void observeEndTimeIntegration(const StepperBase< Scalar > &stepper)
Observe the end of a time integration loop.
void observeCompletedTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
Observe a successfully completed integration step.
void observeFailedTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
Observer a failed integration step.
void observeStartTimeIntegration(const StepperBase< Scalar > &stepper)
Observe the beginning of a time integration loop.
void resetIntegrationObserver(const TimeRange< Scalar > &integrationTimeDomain)
Reset the observer to prepair it to observe another integration.
RCP< IntegrationObserverBase< Scalar > > cloneIntegrationObserver() const
Clone this integration observer if supported .
void observeStartTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
Observer the beginning of an integration step.