Rythmos - Transient Integration for Differential Equations Version of the Day
Loading...
Searching...
No Matches
Rythmos_DefaultIntegrator_decl.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_DEFAULT_INTEGRATOR_DECL_HPP
30#define RYTHMOS_DEFAULT_INTEGRATOR_DECL_HPP
31
32
33#include "Rythmos_IntegrationControlStrategyAcceptingIntegratorBase.hpp"
34#include "Rythmos_InterpolationBufferAppenderAcceptingIntegratorBase.hpp"
35#include "Rythmos_TrailingInterpolationBufferAcceptingIntegratorBase.hpp"
36#include "Rythmos_IntegrationObserverBase.hpp"
37#include "Rythmos_StepControlInfo.hpp"
38#include "Teuchos_ParameterListAcceptorDefaultBase.hpp"
39
40
41
42
43namespace Rythmos {
44
45
49template<class Scalar>
54 virtual public Teuchos::ParameterListAcceptorDefaultBase
55{
56public:
57
59 typedef typename ScalarTraits<Scalar>::magnitudeType ScalarMag;
60
63
66
69 const RCP<IntegrationObserverBase<Scalar> > &integrationObserver
70 );
71
74
77 const RCP<InterpolationBufferAppenderBase<Scalar> > &interpBufferAppender
78 );
79
81 RCP<const InterpolationBufferAppenderBase<Scalar> >
83
85 RCP<InterpolationBufferAppenderBase<Scalar> >
87
89 RCP<InterpolationBufferAppenderBase<Scalar> >
91
93
96
99 const RCP<IntegrationControlStrategyBase<Scalar> > &integrationControlStrategy
100 );
101
103 RCP<IntegrationControlStrategyBase<Scalar> >
105
107 RCP<const IntegrationControlStrategyBase<Scalar> >
109
111
114
116 void setParameterList(RCP<ParameterList> const& paramList);
117
119 RCP<const ParameterList> getValidParameters() const;
120
122
125
127 RCP<IntegratorBase<Scalar> > cloneIntegrator() const;
128
130 void setStepper(
131 const RCP<StepperBase<Scalar> > &stepper,
132 const Scalar &finalTime,
133 const bool landOnFinalTime = true
134 );
135
137 RCP<StepperBase<Scalar> > unSetStepper();
138
140 RCP<const StepperBase<Scalar> > getStepper() const;
141
143 RCP<StepperBase<Scalar> > getNonconstStepper() const;
144
147
150 const RCP<InterpolationBufferBase<Scalar> > &trailingInterpBuffer
151 );
152
154 RCP<InterpolationBufferBase<Scalar> >
156
158 RCP<const InterpolationBufferBase<Scalar> >
160
162 RCP<InterpolationBufferBase<Scalar> >
164
166
168 void getFwdPoints(
169 const Array<Scalar>& time_vec,
170 Array<RCP<const Thyra::VectorBase<Scalar> > >* x_vec,
171 Array<RCP<const Thyra::VectorBase<Scalar> > >* xdot_vec,
172 Array<ScalarMag>* accuracy_vec
173 );
174
176 TimeRange<Scalar> getFwdTimeRange() const;
177
179
182
184 RCP<const Thyra::VectorSpaceBase<Scalar> > get_x_space() const;
185
187 void addPoints(
188 const Array<Scalar>& time_vec,
189 const Array<RCP<const Thyra::VectorBase<Scalar> > >& x_vec,
190 const Array<RCP<const Thyra::VectorBase<Scalar> > >& xdot_vec
191 );
192
194 void getPoints(
195 const Array<Scalar>& time_vec,
196 Array<RCP<const Thyra::VectorBase<Scalar> > >* x_vec,
197 Array<RCP<const Thyra::VectorBase<Scalar> > >* xdot_vec,
198 Array<ScalarMag>* accuracy_vec
199 ) const;
200
202 TimeRange<Scalar> getTimeRange() const;
203
205 void getNodes(Array<Scalar>* time_vec) const;
206
208 void removeNodes(Array<Scalar>& time_vec);
209
211 int getOrder() const;
212
214
215private:
216
217 // ////////////////////////
218 // Private data members
219
220 RCP<IntegrationControlStrategyBase<Scalar> > integrationControlStrategy_;
221 RCP<IntegrationObserverBase<Scalar> > integrationObserver_;
222
223 RCP<InterpolationBufferBase<Scalar> > trailingInterpBuffer_;
224 RCP<InterpolationBufferAppenderBase<Scalar> > interpBufferAppender_;
225
226 RCP<StepperBase<Scalar> > stepper_;
227 TimeRange<Scalar> integrationTimeDomain_;
228 bool landOnFinalTime_;
229
230 int maxNumTimeSteps_;
231
232 int currTimeStepIndex_;
233 StepControlInfo<Scalar> stepCtrlInfoLast_;
234
235 static const std::string maxNumTimeSteps_name_;
236 static const int maxNumTimeSteps_default_;
237
238 // /////////////////////////
239 // Private member functions
240
241 void finalizeSetup();
242
243 bool advanceStepperToTime( const Scalar& t );
244
245};
246
247
252template<class Scalar>
253RCP<DefaultIntegrator<Scalar> >
255
256
261template<class Scalar>
262RCP<DefaultIntegrator<Scalar> >
264 const RCP<IntegrationControlStrategyBase<Scalar> > &integrationControlStrategy,
265 const RCP<IntegrationObserverBase<Scalar> > &integrationObserver
266 );
267
268
273template<class Scalar>
274RCP<DefaultIntegrator<Scalar> >
276 const RCP<IntegrationControlStrategyBase<Scalar> > &integrationControlStrategy
277 );
278
279
284template<class Scalar>
285RCP<DefaultIntegrator<Scalar> >
287 const RCP<IntegrationObserverBase<Scalar> > &integrationObserver
288 );
289
290// 2007/08/30: rabartl: Above, note that I had to name the nonmember
291// constructors taking an single RCP argument different names from each other
292// in order to get around the classic ambiguity problem with implicit
293// conversions of smart pointers.
294
295
296} // namespace Rythmos
297
298
299#endif //RYTHMOS_DEFAULT_INTEGRATOR_DECL_HPP
A concrete subclass for IntegratorBase that allows a good deal of customization.
RCP< DefaultIntegrator< Scalar > > defaultIntegrator(const RCP< IntegrationControlStrategyBase< Scalar > > &integrationControlStrategy, const RCP< IntegrationObserverBase< Scalar > > &integrationObserver)
void setStepper(const RCP< StepperBase< Scalar > > &stepper, const Scalar &finalTime, const bool landOnFinalTime=true)
RCP< const IntegrationControlStrategyBase< Scalar > > getIntegrationControlStrategy() const
RCP< InterpolationBufferBase< Scalar > > getNonconstTrailingInterpolationBuffer()
RCP< StepperBase< Scalar > > getNonconstStepper() const
RCP< InterpolationBufferAppenderBase< Scalar > > getNonconstInterpolationBufferAppender()
RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const
void removeNodes(Array< Scalar > &time_vec)
void setTrailingInterpolationBuffer(const RCP< InterpolationBufferBase< Scalar > > &trailingInterpBuffer)
RCP< const StepperBase< Scalar > > getStepper() const
RCP< DefaultIntegrator< Scalar > > defaultIntegrator()
void addPoints(const Array< Scalar > &time_vec, const Array< RCP< const Thyra::VectorBase< Scalar > > > &x_vec, const Array< RCP< const Thyra::VectorBase< Scalar > > > &xdot_vec)
RCP< StepperBase< Scalar > > unSetStepper()
RCP< IntegratorBase< Scalar > > cloneIntegrator() const
void getFwdPoints(const Array< Scalar > &time_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *x_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *xdot_vec, Array< ScalarMag > *accuracy_vec)
RCP< DefaultIntegrator< Scalar > > observedDefaultIntegrator(const RCP< IntegrationObserverBase< Scalar > > &integrationObserver)
void setIntegrationObserver(const RCP< IntegrationObserverBase< Scalar > > &integrationObserver)
void getNodes(Array< Scalar > *time_vec) const
void setInterpolationBufferAppender(const RCP< InterpolationBufferAppenderBase< Scalar > > &interpBufferAppender)
RCP< const ParameterList > getValidParameters() const
RCP< const InterpolationBufferAppenderBase< Scalar > > getInterpolationBufferAppender()
void setIntegrationControlStrategy(const RCP< IntegrationControlStrategyBase< Scalar > > &integrationControlStrategy)
void setParameterList(RCP< ParameterList > const &paramList)
void getPoints(const Array< Scalar > &time_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *x_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *xdot_vec, Array< ScalarMag > *accuracy_vec) const
RCP< InterpolationBufferBase< Scalar > > unSetTrailingInterpolationBuffer()
RCP< IntegrationControlStrategyBase< Scalar > > getNonconstIntegrationControlStrategy()
RCP< const InterpolationBufferBase< Scalar > > getTrailingInterpolationBuffer() const
TimeRange< Scalar > getFwdTimeRange() const
RCP< InterpolationBufferAppenderBase< Scalar > > unSetInterpolationBufferAppender()
RCP< DefaultIntegrator< Scalar > > controlledDefaultIntegrator(const RCP< IntegrationControlStrategyBase< Scalar > > &integrationControlStrategy)
ScalarTraits< Scalar >::magnitudeType ScalarMag
Mix-in interface for integrator objects that accept an integration control strategy object to be used...
Mix-in interface for integrator objects that accept an interpolationBufferAppender object to be used ...
Mix-in interface for integrator objects that accept a trailing interpolation buffer object to be used...