Rythmos - Transient Integration for Differential Equations Version of the Day
Rythmos_StepperAsModelEvaluator.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 
00030 #ifndef RYTHMOS_STEPPER_AS_MODEL_EVALUATOR_HPP
00031 #define RYTHMOS_STEPPER_AS_MODEL_EVALUATOR_HPP
00032 
00033 
00034 #include "Thyra_ResponseOnlyModelEvaluatorBase.hpp"
00035 #include "Thyra_ModelEvaluatorDelegatorBase.hpp"
00036 
00037 #include "Rythmos_StepperBase.hpp"
00038 #include "Rythmos_IntegratorBase.hpp"
00039 #include "Teuchos_StandardMemberCompositionMacros.hpp"
00040 #include "Teuchos_Assert.hpp"
00041 
00042 
00043 namespace Rythmos {
00044 
00045 
00052 template<class Scalar>
00053 class StepperAsModelEvaluator
00054   : virtual public Thyra::ResponseOnlyModelEvaluatorBase<Scalar>
00055 {
00056 public:
00057 
00060 
00062   StepperAsModelEvaluator();
00063 
00065   void initialize(
00066     const RCP<StepperBase<Scalar> > &stepper,
00067     const RCP<IntegratorBase<Scalar> > &integrator,
00068     const Thyra::ModelEvaluatorBase::InArgs<Scalar> &initialCondition
00069     );
00070 
00072 
00075 
00077   int Np() const;
00079   int Ng() const;
00081   RCP<const Thyra::VectorSpaceBase<Scalar> >
00082   get_p_space(int l) const;
00084   RCP<const Thyra::VectorSpaceBase<Scalar> >
00085   get_g_space(int j) const;
00087   Thyra::ModelEvaluatorBase::InArgs<Scalar> createInArgs() const;
00088 
00090 
00091 private:
00092 
00095 
00097   Thyra::ModelEvaluatorBase::OutArgs<Scalar> createOutArgsImpl() const;
00099   void evalModelImpl(
00100     const Thyra::ModelEvaluatorBase::InArgs<Scalar>& inArgs,
00101     const Thyra::ModelEvaluatorBase::OutArgs<Scalar>& outArgs
00102     ) const;
00103 
00105 
00106 private:
00107 
00108   // //////////////////////
00109   // Private types
00110 
00111   typedef Array<RCP<const Thyra::VectorSpaceBase<Scalar> > > SpaceArray_t;
00112   
00113 
00114   // //////////////////////
00115   // Private data members
00116 
00117   RCP<StepperBase<Scalar> > stepper_;
00118   RCP<IntegratorBase<Scalar> > integrator_;
00119   Thyra::ModelEvaluatorBase::InArgs<Scalar> initialCondition_;
00120 
00121   int Np_;
00122   int Ng_;
00123 
00124   SpaceArray_t g_space_;
00125   SpaceArray_t p_space_;
00126 
00127   mutable  Thyra::ModelEvaluatorBase::InArgs<Scalar> currentInitialCondition_;
00128 
00129 };
00130 
00131 
00133 template<class Scalar>
00134 RCP<StepperAsModelEvaluator<Scalar> >
00135 stepperAsModelEvaluator(
00136   const RCP<StepperBase<Scalar> > &stepper,
00137   const RCP<IntegratorBase<Scalar> > &integrator,
00138   const Thyra::ModelEvaluatorBase::InArgs<Scalar> &initialCondition
00139   )
00140 {
00141   using Teuchos::rcp;
00142   RCP<StepperAsModelEvaluator<Scalar> >
00143     stepperAsModelEvaluator = rcp(new StepperAsModelEvaluator<Scalar>());
00144   stepperAsModelEvaluator->initialize(stepper,integrator,initialCondition);
00145   return stepperAsModelEvaluator;
00146 }
00147 
00148 
00149 // ////////////////////////
00150 // Definitions
00151 
00152 
00153 // Constructors, Initialization, Misc.
00154 
00155 
00156 template<class Scalar>
00157 StepperAsModelEvaluator<Scalar>::StepperAsModelEvaluator()
00158   : Np_(0), Ng_(0)
00159 {}
00160 
00161 
00162 template<class Scalar>
00163 void StepperAsModelEvaluator<Scalar>::initialize(
00164   const RCP<StepperBase<Scalar> > &stepper,
00165   const RCP<IntegratorBase<Scalar> > &integrator,
00166   const Thyra::ModelEvaluatorBase::InArgs<Scalar> &initialCondition
00167   )
00168 {
00169 
00170 #ifdef RYTHMOS_DEBUG
00171   TEST_FOR_EXCEPT(is_null(stepper));
00172   TEST_FOR_EXCEPT(is_null(stepper->getModel()));
00173   TEST_FOR_EXCEPT(is_null(integrator));
00174 #endif
00175   stepper_ = stepper;
00176   integrator_ = integrator;
00177   initialCondition_ = initialCondition;
00178   currentInitialCondition_ = initialCondition;
00179 
00180   const RCP<const Thyra::ModelEvaluator<Scalar> >
00181     stepperModel = stepper_->getModel();
00182 
00183   Np_ = stepperModel->Np();
00184   p_space_.clear();
00185   for ( int l = 0; l < Np_; ++l ) {
00186     p_space_.push_back(stepperModel->get_p_space(l));
00187   }
00188 
00189   Ng_ = 1;
00190   g_space_.clear();
00191   g_space_.push_back(stepper_->getModel()->get_x_space());
00192 
00193 }
00194 
00195 
00196 // Public functions overridden from ModelEvaulator
00197 
00198 
00199 template<class Scalar>
00200 int StepperAsModelEvaluator<Scalar>::Np() const
00201 {
00202   return Np_;
00203 }
00204 
00205 
00206 template<class Scalar>
00207 int StepperAsModelEvaluator<Scalar>::Ng() const
00208 {
00209   return Ng_;
00210 }
00211 
00212 
00213 template<class Scalar>
00214 RCP<const Thyra::VectorSpaceBase<Scalar> >
00215 StepperAsModelEvaluator<Scalar>::get_p_space(int l) const
00216 {
00217 #ifdef RYTHMOS_DEBUG
00218   TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE( l, 0, Np_ );
00219 #endif
00220   return p_space_[l];
00221 }
00222 
00223 
00224 template<class Scalar>
00225 RCP<const Thyra::VectorSpaceBase<Scalar> >
00226 StepperAsModelEvaluator<Scalar>::get_g_space(int j) const
00227 {
00228 #ifdef RYTHMOS_DEBUG
00229   TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE( j, 0, Ng_ );
00230 #endif
00231   return g_space_[j];
00232 }
00233 
00234 
00235 template<class Scalar>
00236 Thyra::ModelEvaluatorBase::InArgs<Scalar>
00237 StepperAsModelEvaluator<Scalar>::createInArgs() const
00238 {
00239   typedef Thyra::ModelEvaluatorBase MEB;
00240   MEB::InArgsSetup<Scalar> inArgs;
00241   inArgs.setModelEvalDescription(this->description());
00242   inArgs.set_Np(Np_);
00243   inArgs.setSupports(MEB::IN_ARG_t);
00244   return inArgs;
00245 }
00246 
00247 
00248 // Private functions overridden from ModelEvaulatorDefaultBase
00249 
00250 
00251 template<class Scalar>
00252 Thyra::ModelEvaluatorBase::OutArgs<Scalar>
00253 StepperAsModelEvaluator<Scalar>::createOutArgsImpl() const
00254 {
00255   typedef Thyra::ModelEvaluatorBase MEB;
00256   MEB::OutArgsSetup<Scalar> outArgs;
00257   outArgs.setModelEvalDescription(this->description());
00258   outArgs.set_Np_Ng(Np_,Ng_);
00259   return outArgs;
00260 }
00261 
00262 
00263 template<class Scalar>
00264 void StepperAsModelEvaluator<Scalar>::evalModelImpl(
00265   const Thyra::ModelEvaluatorBase::InArgs<Scalar>& inArgs,
00266   const Thyra::ModelEvaluatorBase::OutArgs<Scalar>& outArgs
00267   ) const
00268 {
00269 
00270   using Teuchos::as;
00271   using Teuchos::describe;
00272   typedef Teuchos::ScalarTraits<Scalar> ST;
00273   typedef Teuchos::VerboseObjectTempState<InterpolationBufferBase<Scalar> > VOTSSB;
00274 
00275   THYRA_MODEL_EVALUATOR_DECORATOR_EVAL_MODEL_GEN_BEGIN(
00276     "Rythmos::StepperAsModelEvaluator", inArgs, outArgs, Teuchos::null
00277     );
00278   VOTSSB integrator_outputTempState(integrator_,out,incrVerbLevel(verbLevel,-1));
00279   //VOTSSB stepper_outputTempState(stepper_,out,incrVerbLevel(verbLevel,-1));
00280 
00281   // InArgs
00282 
00283   const Scalar finalTime = inArgs.get_t();
00284   for ( int l = 0; l < Np_; ++l ) {
00285     currentInitialCondition_.set_p(l,inArgs.get_p(l));
00286   }
00287 
00288   // OutArgs
00289 
00290   RCP<Thyra::VectorBase<Scalar> >
00291     g_out = outArgs.get_g(0);
00292 
00293   TEST_FOR_EXCEPT_MSG(
00294     is_null(g_out), "You must ask for g(0) when you call this function!"
00295     );
00296 
00297 #ifdef RYTHMOS_DEBUG
00298 
00299   THYRA_ASSERT_VEC_SPACES(
00300     "StepperAsModelEvaluator<Scalar>::evalModel(...)",
00301     *g_out->space(), *stepper_->get_x_space() );
00302 
00303 #endif
00304   
00305   // Set up the integrator
00306 
00307   stepper_->setInitialCondition(currentInitialCondition_);
00308   integrator_->setStepper(stepper_,finalTime);
00309 
00310   // Compute the desired response
00311 
00312   if (!is_null(g_out)) {
00313 
00314     // Get x and xdot at the end time
00315     Array<Scalar> time_vec = Teuchos::tuple<Scalar>(finalTime);
00316     Array<RCP<const Thyra::VectorBase<Scalar> > > x_vec;
00317     integrator_->getFwdPoints( time_vec, &x_vec, 0, 0 );
00318     
00319     Thyra::V_V( g_out.ptr(), *x_vec[0] );
00320 
00321   }
00322 
00323   THYRA_MODEL_EVALUATOR_DECORATOR_EVAL_MODEL_END();
00324 
00325 }
00326 
00327 
00328 } // namespace Rythmos
00329 
00330 
00331 #endif // RYTHMOS_STEPPER_AS_MODEL_EVALUATOR_HPP
 All Classes Functions Variables Typedefs Friends