Rythmos - Transient Integration for Differential Equations Version of the Day
Rythmos_ForwardEulerStepper_decl.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_FORWARDEULER_STEPPER_DECL_H
00030 #define Rythmos_FORWARDEULER_STEPPER_DECL_H
00031 
00032 #include "Rythmos_StepperBase.hpp"
00033 #include "Rythmos_Types.hpp"
00034 #include "Rythmos_MomentoBase.hpp"
00035 #include "Rythmos_StateSerializerStrategy.hpp"
00036 #include "Thyra_ModelEvaluator.hpp"
00037 #include "Teuchos_ParameterListAcceptorDefaultBase.hpp"
00038 #include "Teuchos_ConstNonconstObjectContainer.hpp"
00039 
00040 namespace Rythmos {
00041 
00046   /*
00047 template<class Scalar>
00048   class ForwardEulerStepperMomento :
00049     virtual public MomentoBase<Scalar>,
00050     virtual public Teuchos::ParameterListAcceptorDefaultBase
00051 {
00052   public:
00053     ForwardEulerStepperMomento() {}
00054     virtual ~ForwardEulerStepperMomento() {}
00055 
00056     void serialize(
00057         const StateSerializerStrategy<Scalar>& stateSerializer,
00058         std::ostream& oStream
00059         ) const
00060     {
00061       using Teuchos::is_null;
00062       TEUCHOS_ASSERT( !is_null(model_) );
00063       if (is_null(solution_vector_)) {
00064         solution_vector_ = Thyra::createMember(model_->get_x_space());
00065       }
00066       if (is_null(residual_vector_)) {
00067         residual_vector_ = Thyra::createMember(model_->get_f_space());
00068       }
00069       if (is_null(solution_vector_old_)) {
00070         solution_vector_old_ = Thyra::createMember(model_->get_x_space());
00071       }
00072       stateSerializer.serializeVectorBase(*solution_vector_,oStream);
00073       stateSerializer.serializeVectorBase(*residual_vector_,oStream);
00074       stateSerializer.serializeVectorBase(*solution_vector_old_,oStream);
00075       stateSerializer.serializeScalar(t_,oStream);
00076       stateSerializer.serializeScalar(t_old_,oStream);
00077       stateSerializer.serializeScalar(dt_,oStream);
00078       stateSerializer.serializeInt(numSteps_,oStream);
00079       stateSerializer.serializeBool(isInitialized_,oStream);
00080       stateSerializer.serializeBool(haveInitialCondition_,oStream);
00081       RCP<ParameterList> pl = parameterList_;
00082       if (Teuchos::is_null(pl)) {
00083         pl = Teuchos::parameterList();
00084       }
00085       stateSerializer.serializeParameterList(*pl,oStream);
00086     }
00087 
00088     void deSerialize(
00089         const StateSerializerStrategy<Scalar>& stateSerializer,
00090         std::istream& iStream
00091         )
00092     {
00093       using Teuchos::outArg;
00094       using Teuchos::is_null;
00095       TEUCHOS_ASSERT( !is_null(model_) );
00096       if (is_null(solution_vector_)) {
00097         solution_vector_ = Thyra::createMember(*model_->get_x_space());
00098       }
00099       if (is_null(residual_vector_)) {
00100         residual_vector_ = Thyra::createMember(*model_->get_f_space());
00101       }
00102       if (is_null(solution_vector_old_)) {
00103         solution_vector_old_ = Thyra::createMember(*model_->get_x_space());
00104       }
00105       stateSerializer.deSerializeVectorBase(outArg(*solution_vector_),iStream);
00106       stateSerializer.deSerializeVectorBase(outArg(*residual_vector_),iStream);
00107       stateSerializer.deSerializeVectorBase(outArg(*solution_vector_old_),iStream);
00108       stateSerializer.deSerializeScalar(outArg(t_),iStream);
00109       stateSerializer.deSerializeScalar(outArg(t_old_),iStream);
00110       stateSerializer.deSerializeScalar(outArg(dt_),iStream);
00111       stateSerializer.deSerializeInt(outArg(numSteps_),iStream);
00112       stateSerializer.deSerializeBool(outArg(isInitialized_),iStream);
00113       stateSerializer.deSerializeBool(outArg(haveInitialCondition_),iStream);
00114       if (is_null(parameterList_)) {
00115         parameterList_ = Teuchos::parameterList();
00116       }
00117       stateSerializer.deSerializeParameterList(outArg(*parameterList_),iStream);
00118     }
00119 
00120     RCP<MomentoBase<Scalar> > clone() const
00121     {
00122       RCP<ForwardEulerStepperMomento<Scalar> > m = rcp(new ForwardEulerStepperMomento<Scalar>());
00123       m->set_solution_vector(solution_vector_);
00124       m->set_residual_vector(residual_vector_);
00125       m->set_solution_vector_old(solution_vector_old_);
00126       m->set_t(t_);
00127       m->set_t_old(t_old_);
00128       m->set_dt(dt_);
00129       m->set_numSteps(numSteps_);
00130       m->set_isInitialized(isInitialized_);
00131       m->set_haveInitialCondition(haveInitialCondition_);
00132       m->set_parameterList(parameterList_);
00133       if (!Teuchos::is_null(this->getMyParamList())) {
00134         m->setParameterList(Teuchos::parameterList(*(this->getMyParamList())));
00135       }
00136       m->setModel(model_);
00137       m->setBasePoint(basePoint_);
00138       // How do I copy the VerboseObject data?  
00139       // 07/10/09 tscoffe:  Its not set up in Teuchos to do this yet
00140       return m;
00141     }
00142 
00143     void set_solution_vector(const RCP<const VectorBase<Scalar> >& solution_vector )
00144     { 
00145       solution_vector_ = Teuchos::null;
00146       if (!Teuchos::is_null(solution_vector)) {
00147         solution_vector_ = solution_vector->clone_v(); 
00148       }
00149     }
00150     RCP<VectorBase<Scalar> > get_solution_vector() const
00151     { return solution_vector_; }
00152 
00153     void set_residual_vector(const RCP<const VectorBase<Scalar> >& residual_vector )
00154     { 
00155       residual_vector_ = Teuchos::null;
00156       if (!Teuchos::is_null(residual_vector)) {
00157         residual_vector_ = residual_vector->clone_v(); 
00158       }
00159     }
00160     RCP<VectorBase<Scalar> > get_residual_vector() const
00161     { return residual_vector_; }
00162 
00163     void set_solution_vector_old(const RCP<const VectorBase<Scalar> >& solution_vector_old )
00164     { 
00165       solution_vector_old_ = Teuchos::null;
00166       if (!Teuchos::is_null(solution_vector_old)) {
00167         solution_vector_old_ = solution_vector_old->clone_v(); 
00168       }
00169     }
00170     RCP<VectorBase<Scalar> > get_solution_vector_old() const
00171     { return solution_vector_old_; }
00172 
00173     void set_t(const Scalar & t)
00174     { t_ = t; }
00175     Scalar get_t() const
00176     { return t_; }
00177 
00178     void set_t_old(const Scalar & t_old)
00179     { t_old_ = t_old; }
00180     Scalar get_t_old() const
00181     { return t_old_; }
00182 
00183     void set_dt(const Scalar & dt)
00184     { dt_ = dt; }
00185     Scalar get_dt() const
00186     { return dt_; }
00187 
00188     void set_numSteps(const int & numSteps)
00189     { numSteps_ = numSteps; }
00190     int get_numSteps() const
00191     { return numSteps_; }
00192 
00193     void set_isInitialized(const bool & isInitialized)
00194     { isInitialized_ = isInitialized; }
00195     bool get_isInitialized() const
00196     { return isInitialized_; }
00197 
00198     void set_haveInitialCondition(const bool & haveInitialCondition)
00199     { haveInitialCondition_ = haveInitialCondition; }
00200     bool get_haveInitialCondition() const
00201     { return haveInitialCondition_; }
00202 
00203     void set_parameterList(const RCP<const ParameterList>& pl)
00204     { 
00205       parameterList_ = Teuchos::null;
00206       if (!Teuchos::is_null(pl)) {
00207         parameterList_ = Teuchos::parameterList(*pl); 
00208       }
00209     }
00210     RCP<ParameterList> get_parameterList() const
00211     { return parameterList_; }
00212 
00213     void setParameterList(const RCP<ParameterList>& paramList)
00214     { this->setMyParamList(paramList); }
00215     RCP<const ParameterList> getValidParameters() const
00216     { return Teuchos::null; }
00217 
00218     void set_model(const RCP<const Thyra::ModelEvaluator<Scalar> >& model)
00219     { model_ = model; }
00220     RCP<const Thyra::ModelEvaluator<Scalar> > get_model() const
00221     { return model_; }
00222 
00223     void set_basePoint(const RCP<const Thyra::ModelEvaluatorBase::InArgs<Scalar> >& basePoint)
00224     { basePoint_ = basePoint; }
00225     RCP<const Thyra::ModelEvaluatorBase::InArgs<Scalar> > get_basePoint() const
00226     { return basePoint_; }
00227 
00228   private:
00229     RCP<Thyra::VectorBase<Scalar> > solution_vector_;
00230     RCP<Thyra::VectorBase<Scalar> > residual_vector_;
00231     RCP<Thyra::VectorBase<Scalar> > solution_vector_old_;
00232     Scalar t_;
00233     Scalar t_old_;
00234     Scalar dt_;
00235     int numSteps_;
00236     bool isInitialized_;
00237     bool haveInitialCondition_;
00238     RCP<ParameterList> parameterList_;
00239      
00240     // Objects that must be set prior to serialization and deSerialization:
00241     RCP<const Thyra::ModelEvaluator<Scalar> > model_;
00242     // Objects that must be set prior to calling ForwardEulerStepper::setMomento: 
00243     RCP<const Thyra::ModelEvaluatorBase::InArgs<Scalar> > basePoint_;
00244 };
00245 */
00246 template<class Scalar>
00247   class ForwardEulerStepperMomento :
00248     virtual public MomentoBase<Scalar>,
00249     virtual public Teuchos::ParameterListAcceptorDefaultBase
00250 {
00251   public:
00252     ForwardEulerStepperMomento();
00253     virtual ~ForwardEulerStepperMomento();
00254 
00255     void serialize(
00256         const StateSerializerStrategy<Scalar>& stateSerializer,
00257         std::ostream& oStream
00258         ) const;
00259 
00260     void deSerialize(
00261         const StateSerializerStrategy<Scalar>& stateSerializer,
00262         std::istream& iStream
00263         );
00264 
00265     RCP<MomentoBase<Scalar> > clone() const;
00266 
00267     void set_solution_vector(const RCP<const VectorBase<Scalar> >& solution_vector );
00268     RCP<VectorBase<Scalar> > get_solution_vector() const;
00269 
00270     void set_residual_vector(const RCP<const VectorBase<Scalar> >& residual_vector );
00271     RCP<VectorBase<Scalar> > get_residual_vector() const;
00272 
00273     void set_solution_vector_old(const RCP<const VectorBase<Scalar> >& solution_vector_old );
00274     RCP<VectorBase<Scalar> > get_solution_vector_old() const;
00275 
00276     void set_t(const Scalar & t);
00277     Scalar get_t() const;
00278 
00279     void set_t_old(const Scalar & t_old);
00280     Scalar get_t_old() const;
00281 
00282     void set_dt(const Scalar & dt);
00283     Scalar get_dt() const;
00284 
00285     void set_numSteps(const int & numSteps);
00286     int get_numSteps() const;
00287 
00288     void set_isInitialized(const bool & isInitialized);
00289     bool get_isInitialized() const;
00290 
00291     void set_haveInitialCondition(const bool & haveInitialCondition);
00292     bool get_haveInitialCondition() const;
00293 
00294     void set_parameterList(const RCP<const ParameterList>& pl);
00295     RCP<ParameterList> get_parameterList() const;
00296 
00297     void setParameterList(const RCP<ParameterList>& paramList);
00298     RCP<const ParameterList> getValidParameters() const;
00299 
00300     void set_model(const RCP<const Thyra::ModelEvaluator<Scalar> >& model);
00301     RCP<const Thyra::ModelEvaluator<Scalar> > get_model() const;
00302 
00303     void set_basePoint(const RCP<const Thyra::ModelEvaluatorBase::InArgs<Scalar> >& basePoint);
00304     RCP<const Thyra::ModelEvaluatorBase::InArgs<Scalar> > get_basePoint() const;
00305 
00306   private:
00307     RCP<Thyra::VectorBase<Scalar> > solution_vector_;
00308     RCP<Thyra::VectorBase<Scalar> > residual_vector_;
00309     RCP<Thyra::VectorBase<Scalar> > solution_vector_old_;
00310     Scalar t_;
00311     Scalar t_old_;
00312     Scalar dt_;
00313     int numSteps_;
00314     bool isInitialized_;
00315     bool haveInitialCondition_;
00316     RCP<ParameterList> parameterList_;
00317      
00318     // Objects that must be set prior to serialization and deSerialization:
00319     RCP<const Thyra::ModelEvaluator<Scalar> > model_;
00320     // Objects that must be set prior to calling ForwardEulerStepper::setMomento: 
00321     RCP<const Thyra::ModelEvaluatorBase::InArgs<Scalar> > basePoint_;
00322 };
00323 
00324 
00326 template<class Scalar>
00327 class ForwardEulerStepper : virtual public StepperBase<Scalar>
00328 {
00329   public:
00330 
00331     typedef Teuchos::ScalarTraits<Scalar> ST;
00332     typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
00333     
00335     ForwardEulerStepper();
00336 
00338     bool supportsCloning() const;
00339 
00341     RCP<StepperBase<Scalar> > cloneStepperAlgorithm() const;
00342 
00344     void setModel(const RCP<const Thyra::ModelEvaluator<Scalar> >& model);
00345 
00347     void setNonconstModel(const RCP<Thyra::ModelEvaluator<Scalar> >& model);
00348 
00350     RCP<const Thyra::ModelEvaluator<Scalar> > getModel() const;
00351 
00353     RCP<Thyra::ModelEvaluator<Scalar> > getNonconstModel();
00354 
00356     void setInitialCondition(
00357       const Thyra::ModelEvaluatorBase::InArgs<Scalar> &initialCondition
00358       );
00359 
00361     Thyra::ModelEvaluatorBase::InArgs<Scalar> getInitialCondition() const;
00362 
00364     RCP<const Thyra::VectorSpaceBase<Scalar> > get_x_space() const;
00365     
00367     ~ForwardEulerStepper();
00368 
00370     Scalar takeStep(Scalar dt, StepSizeType flag);
00371 
00373     const StepStatus<Scalar> getStepStatus() const;
00374 
00376     std::string description() const;
00377 
00379     void describe(
00380       Teuchos::FancyOStream            &out,
00381       const Teuchos::EVerbosityLevel   verbLevel
00382       ) const;
00383     
00386     void addPoints(
00387       const Array<Scalar>& time_vec
00388       ,const Array<RCP<const Thyra::VectorBase<Scalar> > >& x_vec
00389       ,const Array<RCP<const Thyra::VectorBase<Scalar> > >& xdot_vec
00390       );
00391     
00393     void getPoints(
00394       const Array<Scalar>& time_vec
00395       ,Array<RCP<const Thyra::VectorBase<Scalar> > >* x_vec
00396       ,Array<RCP<const Thyra::VectorBase<Scalar> > >* xdot_vec
00397       ,Array<ScalarMag>* accuracy_vec
00398       ) const;
00399 
00401     void setRange(
00402       const TimeRange<Scalar>& range,
00403       const InterpolationBufferBase<Scalar> & IB
00404       );
00405 
00407     TimeRange<Scalar> getTimeRange() const;
00408 
00410     void getNodes(Array<Scalar>* time_vec) const;
00411 
00413     void removeNodes(Array<Scalar>& time_vec);
00414 
00416     int getOrder() const;
00417 
00419 
00420     void setParameterList(RCP<Teuchos::ParameterList> const& paramList);
00421 
00423     RCP<Teuchos::ParameterList> getNonconstParameterList();
00424 
00426     RCP<Teuchos::ParameterList> unsetParameterList();
00427 
00429     RCP<const Teuchos::ParameterList> getValidParameters() const;
00430 
00434     RCP<const MomentoBase<Scalar> > getMomento() const;
00435 
00439     void setMomento( const Ptr<const MomentoBase<Scalar> >& momentoPtr );
00440 
00441   private:
00442 
00443     RCP<const Thyra::ModelEvaluator<Scalar> > model_;
00444     RCP<Thyra::VectorBase<Scalar> > solution_vector_;
00445     RCP<Thyra::VectorBase<Scalar> > residual_vector_;
00446     Scalar t_;
00447     Scalar dt_;
00448     Scalar t_old_;
00449     RCP<Thyra::VectorBase<Scalar> > solution_vector_old_;
00450     Thyra::ModelEvaluatorBase::InArgs<Scalar> basePoint_;
00451     int numSteps_;
00452     bool haveInitialCondition_;
00453 
00454     RCP<Teuchos::ParameterList> parameterList_;
00455     bool isInitialized_;
00456 
00457     // Private member functions:
00458     void defaultInitializAll_();
00459     void initialize_();
00460     void checkConsistentState_();
00461 
00462 };
00463 
00464 // Nonmember constructor
00465 template<class Scalar>
00466 RCP<ForwardEulerStepper<Scalar> > forwardEulerStepper();
00467 
00468 // Nonmember constructor
00469 template<class Scalar>
00470 RCP<ForwardEulerStepper<Scalar> > forwardEulerStepper(const RCP<Thyra::ModelEvaluator<Scalar> >& model);
00471 
00472 } // namespace Rythmos
00473 
00474 #endif //Rythmos_FORWARDEULER_STEPPER_DECL_H
 All Classes Functions Variables Typedefs Friends