Rythmos - Transient Integration for Differential Equations Version of the Day
Loading...
Searching...
No Matches
Rythmos_CompositeIntegrationObserver.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_COMPOSITE_INTEGRATION_OBSERVER_HPP
30#define RYTHMOS_COMPOSITE_INTEGRATION_OBSERVER_HPP
31
32
33#include "Rythmos_IntegrationObserverBase.hpp"
34#include "Teuchos_as.hpp"
35
36
37namespace Rythmos {
38
39
44template<class Scalar>
46 : public IntegrationObserverBase<Scalar>
47{
48public:
49
52
55
57 void addObserver(
58 const RCP<IntegrationObserverBase<Scalar> > &observer
59 );
60
61 // ToDo: Add functions to add observers
62
64
67
69 virtual RCP<IntegrationObserverBase<Scalar> >
71
73 virtual void resetIntegrationObserver(
74 const TimeRange<Scalar> &integrationTimeDomain
75 );
76
77 void observeStartTimeIntegration(const StepperBase<Scalar> &stepper);
78
79 void observeEndTimeIntegration(const StepperBase<Scalar> &stepper);
80
82 const StepperBase<Scalar> &stepper,
83 const StepControlInfo<Scalar> &stepCtrlInfo,
84 const int timeStepIter
85 );
86
88 virtual void observeCompletedTimeStep(
89 const StepperBase<Scalar> &stepper,
90 const StepControlInfo<Scalar> &stepCtrlInfo,
91 const int timeStepIter
92 );
93
95 virtual void observeFailedTimeStep(
96 const StepperBase<Scalar> &stepper,
97 const StepControlInfo<Scalar> &stepCtrlInfo,
98 const int timeStepIter
99 );
100
102
103private:
104
105 Array<RCP<IntegrationObserverBase<Scalar> > > observers_;
106
107};
108
109
114template<class Scalar>
115RCP<CompositeIntegrationObserver<Scalar> > createCompositeIntegrationObserver()
116{
117 RCP<CompositeIntegrationObserver<Scalar> >
118 frsco(new CompositeIntegrationObserver<Scalar>());
119 return frsco;
120}
121
122
123//
124// Implementations
125//
126
127
128// Constructors/Initializers/Accessors
129
130
131template<class Scalar>
134
135
136template<class Scalar>
138 const RCP<IntegrationObserverBase<Scalar> > &observer
139 )
140{
141 observers_.push_back(observer);
142}
143
144
145// Overridden from IntegrationObserverBase
146
147
148template<class Scalar>
149RCP<IntegrationObserverBase<Scalar> >
151{
152 using Teuchos::as;
153 RCP<CompositeIntegrationObserver<Scalar> >
154 compositeObserver = createCompositeIntegrationObserver<Scalar>();
155 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
156 compositeObserver->addObserver(observers_[i]->cloneIntegrationObserver());
157 }
158 return compositeObserver;
159}
160
161
162template<class Scalar>
164 const TimeRange<Scalar> &integrationTimeDomain
165 )
166{
167 using Teuchos::as;
168 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
169 observers_[i]->resetIntegrationObserver(integrationTimeDomain);
170 }
171}
172
173template<class Scalar>
175 const StepperBase<Scalar> &stepper
176 )
177{
178 using Teuchos::as;
179
180 const RCP<FancyOStream> out = this->getOStream();
181 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
182
183 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
184 RCP<IntegrationObserverBase<Scalar> > observer = observers_[i];
185 observer->setOStream(out);
186 observer->setVerbLevel(verbLevel);
187 observer->observeStartTimeIntegration(stepper);
188 }
189}
190
191template<class Scalar>
193 const StepperBase<Scalar> &stepper
194 )
195{
196 using Teuchos::as;
197
198 const RCP<FancyOStream> out = this->getOStream();
199 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
200
201 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
202 RCP<IntegrationObserverBase<Scalar> > observer = observers_[i];
203 observer->setOStream(out);
204 observer->setVerbLevel(verbLevel);
205 observer->observeEndTimeIntegration(stepper);
206 }
207}
208
209template<class Scalar>
211 const StepperBase<Scalar> &stepper,
212 const StepControlInfo<Scalar> &stepCtrlInfo,
213 const int timeStepIter
214 )
215{
216 using Teuchos::as;
217
218 const RCP<FancyOStream> out = this->getOStream();
219 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
220
221 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
222 RCP<IntegrationObserverBase<Scalar> > observer = observers_[i];
223 observer->setOStream(out);
224 observer->setVerbLevel(verbLevel);
225 observer->observeStartTimeStep(stepper,stepCtrlInfo,timeStepIter);
226 }
227}
228
229template<class Scalar>
231 const StepperBase<Scalar> &stepper,
232 const StepControlInfo<Scalar> &stepCtrlInfo,
233 const int timeStepIter
234 )
235{
236 using Teuchos::as;
237
238 const RCP<FancyOStream> out = this->getOStream();
239 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
240
241 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
242 RCP<IntegrationObserverBase<Scalar> > observer = observers_[i];
243 observer->setOStream(out);
244 observer->setVerbLevel(verbLevel);
245 observer->observeCompletedTimeStep(stepper,stepCtrlInfo,timeStepIter);
246 }
247}
248
249
250template<class Scalar>
252 const StepperBase<Scalar> &stepper,
253 const StepControlInfo<Scalar> &stepCtrlInfo,
254 const int timeStepIter
255 )
256{
257 using Teuchos::as;
258
259 const RCP<FancyOStream> out = this->getOStream();
260 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
261
262 for (int i = 0; i < as<int>(observers_.size()); ++i ) {
263 RCP<IntegrationObserverBase<Scalar> > observer = observers_[i];
264 observer->setOStream(out);
265 observer->setVerbLevel(verbLevel);
266 observer->observeFailedTimeStep(stepper,stepCtrlInfo,timeStepIter);
267 }
268}
269
270
271} // namespace Rythmos
272
273
274#endif //RYTHMOS_COMPOSITE_INTEGRATOR_OBSERVER_HPP
virtual void observeCompletedTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
virtual RCP< IntegrationObserverBase< Scalar > > cloneIntegrationObserver() const
RCP< CompositeIntegrationObserver< Scalar > > createCompositeIntegrationObserver()
Non-member constructor.
virtual void resetIntegrationObserver(const TimeRange< Scalar > &integrationTimeDomain)
void observeStartTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
Observer the beginning of an integration step.
void observeEndTimeIntegration(const StepperBase< Scalar > &stepper)
Observe the end of a time integration loop.
void observeStartTimeIntegration(const StepperBase< Scalar > &stepper)
Observe the beginning of a time integration loop.
virtual void observeFailedTimeStep(const StepperBase< Scalar > &stepper, const StepControlInfo< Scalar > &stepCtrlInfo, const int timeStepIter)
void addObserver(const RCP< IntegrationObserverBase< Scalar > > &observer)
Base class for strategy objects that observe and time integration by observing the stepper object.