Rythmos_IntegratorBase.hpp

00001 //@HEADER
00002 // ***********************************************************************
00003 //
00004 //                     Rythmos Package
00005 //                 Copyright (2006) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Todd S. Coffey (tscoffe@sandia.gov)
00025 //
00026 // ***********************************************************************
00027 //@HEADER
00028 
00029 #ifndef Rythmos_INTEGRATOR_BASE_H
00030 #define Rythmos_INTEGRATOR_BASE_H
00031 
00032 #include "Rythmos_InterpolationBufferBase.hpp"
00033 #include "Rythmos_StepperBase.hpp"
00034 #include "Teuchos_as.hpp"
00035 
00036 
00037 namespace Rythmos {
00038 
00039 
00040 namespace Exceptions {
00041 
00042 
00046 class GetFwdPointsFailed : public ::Rythmos::Exceptions::ExceptionBase
00047 {public: GetFwdPointsFailed(const std::string &what):ExceptionBase(what) {}};
00048 
00049 
00050 } // namespace Exceptions
00051 
00052 
00061 template<class Scalar> 
00062 class IntegratorBase : virtual public InterpolationBufferBase<Scalar>
00063 {
00064 public:
00065   
00067   typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
00068 
00070   virtual RCP<IntegratorBase<Scalar> > cloneIntegrator() const
00071     { return Teuchos::null; }
00072 
00074   virtual void setTrailingInterpolationBuffer(
00075     const RCP<InterpolationBufferBase<Scalar> > &trailingInterpBuffer
00076     ) 
00077     { TEST_FOR_EXCEPT(true); }
00078 
00080   virtual bool acceptsTrailingInterpolationBuffer() const 
00081     { return false; }
00082 
00117   virtual void setStepper(
00118     const RCP<StepperBase<Scalar> > &stepper,
00119     const Scalar &finalTime,
00120     const bool landOnFinalTime = true
00121     ) =0;
00122 
00128   virtual Teuchos::RCP<const StepperBase<Scalar> > getStepper() const =0;
00129 
00138   virtual RCP<StepperBase<Scalar> > unSetStepper() =0;
00139 
00182   virtual void getFwdPoints(
00183     const Array<Scalar>& time_vec,
00184     Array<RCP<const Thyra::VectorBase<Scalar> > >* x_vec,
00185     Array<RCP<const Thyra::VectorBase<Scalar> > >* xdot_vec,
00186     Array<ScalarMag>* accuracy_vec
00187     ) =0;
00188 
00197   virtual TimeRange<Scalar> getFwdTimeRange() const =0;
00198 
00199 };
00200 
00201 
00202 // 2007/09/14: rabartl: ToDo: Move these functions into a file
00203 // Rythmos_IntegratorBaseHelpers.hpp.
00204 
00205 
00210 template<class Scalar> 
00211 RCP<const Thyra::VectorBase<Scalar> >
00212 get_fwd_x( IntegratorBase<Scalar>& integrator, const Scalar t )
00213 {
00214   Array<Scalar> time_vec;
00215   time_vec.push_back(t);
00216   Array<RCP<const Thyra::VectorBase<Scalar> > > x_vec;
00217   integrator.getFwdPoints(time_vec,&x_vec,0,0);
00218   return x_vec[0];
00219 }
00220 
00221 
00227 template<class Scalar> 
00228 void get_fwd_x_and_x_dot(
00229   IntegratorBase<Scalar>& integrator,
00230   const Scalar t,
00231   RCP<const Thyra::VectorBase<Scalar> > *x,
00232   RCP<const Thyra::VectorBase<Scalar> > *x_dot
00233   )
00234 {
00235   Array<Scalar> time_vec;
00236   time_vec.push_back(t);
00237   Array<RCP<const Thyra::VectorBase<Scalar> > > x_vec;
00238   Array<RCP<const Thyra::VectorBase<Scalar> > > x_dot_vec;
00239   integrator.getFwdPoints(
00240     time_vec,
00241     x ? &x_vec : 0,
00242     x_dot ? &x_dot_vec : 0,
00243     0
00244     );
00245   if (x) *x = x_vec[0];
00246   if (x_dot) *x_dot = x_dot_vec[0];
00247 }
00248 
00249  
00250 } // namespace Rythmos
00251 
00252 
00253 #endif // Rythmos_INTEGRATOR_BASE_H
00254 

Generated on Tue Oct 20 12:46:08 2009 for Rythmos - Transient Integration for Differential Equations by doxygen 1.4.7