Rythmos - Transient Integration for Differential Equations Version of the Day
Rythmos_RampingIntegrationControlStrategy_def.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_RAMPING_INTEGRATION_CONTROL_STRATEGY_DEF_HPP
00031 #define RYTHMOS_RAMPING_INTEGRATION_CONTROL_STRATEGY_DEF_HPP
00032 
00033 #include "Rythmos_RampingIntegrationControlStrategy_decl.hpp"
00034 #include "Teuchos_VerboseObjectParameterListHelpers.hpp"
00035 
00036 
00037 namespace Rythmos {
00038 
00039 
00040 template<class Scalar>
00041 RCP<RampingIntegrationControlStrategy<Scalar> >
00042 rampingIntegrationControlStrategy()
00043 {
00044   RCP<RampingIntegrationControlStrategy<Scalar> >
00045     integrationControl =
00046       Teuchos::rcp(new RampingIntegrationControlStrategy<Scalar>());
00047   return integrationControl;
00048 }
00049 
00050 
00051 template<class Scalar>
00052 RCP<RampingIntegrationControlStrategy<Scalar> >
00053 rampingIntegrationControlStrategy( const RCP<ParameterList> &paramList )
00054 {
00055   RCP<RampingIntegrationControlStrategy<Scalar> >
00056     integrationControl =
00057       Teuchos::rcp(new RampingIntegrationControlStrategy<Scalar>());
00058   integrationControl->setParameterList(paramList);
00059   return integrationControl;
00060 }
00061 
00062 
00063 //
00064 // Implementation
00065 //
00066 
00067 
00068 // Static members
00069 
00070 
00071 template<class Scalar>
00072 const std::string
00073 RampingIntegrationControlStrategy<Scalar>::take_variable_steps_name_
00074 = "Take Variable Steps";
00075 
00076 template<class Scalar>
00077 const bool
00078 RampingIntegrationControlStrategy<Scalar>::take_variable_steps_default_
00079 = true;
00080 
00081 
00082 template<class Scalar>
00083 const std::string
00084 RampingIntegrationControlStrategy<Scalar>::num_constant_steps_name_
00085 = "Number of Initial Constant Steps";
00086 
00087 template<class Scalar>
00088 const int
00089 RampingIntegrationControlStrategy<Scalar>::num_constant_steps_default_
00090 = 0;
00091 
00092 
00093 template<class Scalar>
00094 const std::string
00095 RampingIntegrationControlStrategy<Scalar>::num_ramping_steps_name_
00096 = "Number of Ramping Steps";
00097 
00098 template<class Scalar>
00099 const int
00100 RampingIntegrationControlStrategy<Scalar>::num_ramping_steps_default_
00101 = 6;
00102 
00103 
00104 template<class Scalar>
00105 const std::string
00106 RampingIntegrationControlStrategy<Scalar>::initial_dt_name_
00107 = "Initial dt";
00108 
00109 template<class Scalar>
00110 const double
00111 RampingIntegrationControlStrategy<Scalar>::initial_dt_default_
00112 = -1.0;
00113 
00114 
00115 template<class Scalar>
00116 const std::string
00117 RampingIntegrationControlStrategy<Scalar>::min_dt_name_
00118 = "Min dt";
00119 
00120 template<class Scalar>
00121 const double
00122 RampingIntegrationControlStrategy<Scalar>::min_dt_default_
00123 = std::numeric_limits<Scalar>::min();
00124 
00125 
00126 template<class Scalar>
00127 const std::string
00128 RampingIntegrationControlStrategy<Scalar>::max_dt_name_
00129 = "Max dt";
00130 
00131 template<class Scalar>
00132 const double
00133 RampingIntegrationControlStrategy<Scalar>::max_dt_default_
00134 = std::numeric_limits<Scalar>::max();
00135 
00136 
00137 template<class Scalar>
00138 const std::string
00139 RampingIntegrationControlStrategy<Scalar>::ramping_factor_name_
00140 = "Ramping Factor";
00141 
00142 template<class Scalar>
00143 const double
00144 RampingIntegrationControlStrategy<Scalar>::ramping_factor_default_
00145 = 1.0;
00146 
00147 
00148 template<class Scalar>
00149 const std::string
00150 RampingIntegrationControlStrategy<Scalar>::max_step_failures_name_
00151 = "Maximum Number of Step Failures";
00152 
00153 template<class Scalar>
00154 const int
00155 RampingIntegrationControlStrategy<Scalar>::max_step_failures_default_
00156 = 100;
00157 
00158 
00159 
00160 // Constructors/Initializers
00161 
00162 
00163 template<class Scalar>
00164 RampingIntegrationControlStrategy<Scalar>::RampingIntegrationControlStrategy()
00165   : take_variable_steps_(take_variable_steps_default_),
00166     num_constant_steps_(num_constant_steps_default_),
00167     num_ramping_steps_(num_ramping_steps_default_),
00168     initial_dt_(initial_dt_default_),
00169     min_dt_(min_dt_default_),
00170     max_dt_(max_dt_default_),
00171     ramping_factor_(ramping_factor_default_),
00172     num_step_failures_(0),
00173     max_step_failures_(max_step_failures_default_),
00174     current_dt_(-1.0)
00175 {}
00176 
00177 
00178 // Overridden from ParameterListAcceptor
00179 
00180 
00181 template<class Scalar>
00182 void RampingIntegrationControlStrategy<Scalar>::setParameterList(
00183   RCP<ParameterList> const& paramList
00184   )
00185 {
00186   using Teuchos::as;
00187   using Teuchos::get;
00188   typedef Teuchos::ScalarTraits<Scalar> ST;
00189   TEUCHOS_TEST_FOR_EXCEPT(is_null(paramList));
00190   paramList->validateParametersAndSetDefaults(*getValidParameters());
00191   this->setMyParamList(paramList);
00192 
00193   take_variable_steps_ = paramList->get<bool>  (take_variable_steps_name_,
00194                                                 take_variable_steps_default_);
00195   num_constant_steps_  = paramList->get<int>   (num_constant_steps_name_,
00196                                                 num_constant_steps_default_);
00197   num_ramping_steps_   = paramList->get<int>   (num_ramping_steps_name_,
00198                                                 num_ramping_steps_default_);
00199   initial_dt_          = paramList->get<double>(initial_dt_name_,
00200                                                 initial_dt_default_);
00201   min_dt_              = paramList->get<double>(min_dt_name_, min_dt_default_);
00202   max_dt_              = paramList->get<double>(max_dt_name_, max_dt_default_);
00203   ramping_factor_      = paramList->get<double>(ramping_factor_name_,
00204                                                 ramping_factor_default_);
00205   max_step_failures_   = paramList->get<int>   (max_step_failures_name_,
00206                                                 max_step_failures_default_);
00207 
00208   Teuchos::readVerboseObjectSublist(&*paramList,this);
00209 }
00210 
00211 
00212 template<class Scalar>
00213 RCP<const ParameterList>
00214 RampingIntegrationControlStrategy<Scalar>::getValidParameters() const
00215 {
00216   static RCP<const ParameterList> validPL;
00217   if (is_null(validPL) ) {
00218     RCP<ParameterList> pl = Teuchos::parameterList();
00219 
00220     pl->set( take_variable_steps_name_, take_variable_steps_default_,
00221       "Take variable time steps after '" + num_constant_steps_name_ +
00222       "' plus '" + num_ramping_steps_name_ + "' steps.  Variable time "
00223       "stepping allows the Stepper to adjust the time step through a "
00224       "StepControlStrategy after fixed time steps during initial constant "
00225       "steps and ramping steps.  If false, fixed-time steps are taken "
00226       "after ramping.  Fixed time stepping requires the Stepper "
00227       "to take the time step set by this IntegrationControlStrategy.");
00228     pl->set(num_constant_steps_name_, num_constant_steps_default_,
00229       "Number of initial constant steps to take before starting the ramping.");
00230     pl->set(num_ramping_steps_name_, num_ramping_steps_default_,
00231       "Number of ramping steps to take before handing control to "
00232       "variable stepper if '" + take_variable_steps_name_ +
00233       "' is set to true.  Otherwise take fixed-time steps.");
00234     pl->set(initial_dt_name_, initial_dt_default_,
00235       "Initial time step.");
00236     pl->set(min_dt_name_, min_dt_default_,
00237       "Minimum time step.");
00238     pl->set(max_dt_name_, max_dt_default_,
00239       "Maximum time step.");
00240     pl->set(ramping_factor_name_, ramping_factor_default_,
00241       "Time step growth factor used during ramping phase. dt_{n+1} = "
00242       "(ramping factor) * dt_n");
00243     pl->set(max_step_failures_name_, max_step_failures_default_,
00244       "The maximum number of step failures before exiting with error.");
00245     Teuchos::setupVerboseObjectSublist(&*pl);
00246     validPL = pl;
00247   }
00248   return validPL;
00249 }
00250 
00251 
00252 // Overridden from IntegrationControlStrategyBase
00253 
00254 
00255 template<class Scalar>
00256 bool RampingIntegrationControlStrategy<Scalar>::handlesFailedTimeSteps() const
00257 {
00258   return true;
00259 }
00260 
00261 
00262 template<class Scalar>
00263 RCP<IntegrationControlStrategyBase<Scalar> >
00264 RampingIntegrationControlStrategy<Scalar>::cloneIntegrationControlStrategy() const
00265 {
00266   RCP<RampingIntegrationControlStrategy<Scalar> >
00267     integrCtrlStry = rampingIntegrationControlStrategy<Scalar>();
00268   const RCP<const ParameterList> paramList = this->getParameterList();
00269   if (!is_null(paramList))
00270     integrCtrlStry->setParameterList(Teuchos::parameterList(*paramList));
00271 
00272   integrCtrlStry->take_variable_steps_ = this->take_variable_steps_;
00273   integrCtrlStry->num_constant_steps_  = this->num_constant_steps_;
00274   integrCtrlStry->num_ramping_steps_   = this->num_ramping_steps_;
00275   integrCtrlStry->initial_dt_          = this->initial_dt_;
00276   integrCtrlStry->min_dt_              = this->min_dt_;
00277   integrCtrlStry->max_dt_              = this->max_dt_;
00278   integrCtrlStry->ramping_factor_      = this->ramping_factor_;
00279   integrCtrlStry->current_dt_          = this->current_dt_;
00280 
00281   return integrCtrlStry;
00282 }
00283 
00284 
00285 template<class Scalar>
00286 void
00287 RampingIntegrationControlStrategy<Scalar>::resetIntegrationControlStrategy(
00288   const TimeRange<Scalar> &integrationTimeDomain
00289   )
00290 {
00291   typedef Teuchos::ScalarTraits<Scalar> ST;
00292 #ifdef HAVE_RYTHMOS_DEBUG
00293   TEUCHOS_ASSERT(integrationTimeDomain.length() > ST::zero());
00294 #endif
00295   integrationTimeDomain_ = integrationTimeDomain;
00296   if (max_dt_ < ST::zero()) {
00297     max_dt_ = integrationTimeDomain_.length();
00298   }
00299 
00300   current_dt_ = initial_dt_;
00301 }
00302 
00303 
00304 template<class Scalar>
00305 StepControlInfo<Scalar>
00306 RampingIntegrationControlStrategy<Scalar>::getNextStepControlInfo(
00307   const StepperBase<Scalar> &stepper,
00308   const StepControlInfo<Scalar> &stepCtrlInfoLast,
00309   const int timeStepIter
00310   )
00311 {
00312 
00313   typedef Teuchos::ScalarTraits<Scalar> ST;
00314 
00315 #ifdef HAVE_RYTHMOS_DEBUG
00316   TEUCHOS_ASSERT(integrationTimeDomain_.length() > ST::zero());
00317 #endif
00318 
00319   StepControlInfo<Scalar> trialStepCtrlInfo;
00320 
00321   if (timeStepIter < num_constant_steps_)
00322     current_dt_ = initial_dt_;
00323   else if (timeStepIter < num_constant_steps_ + num_ramping_steps_)
00324     current_dt_ *= ramping_factor_;
00325 
00326   current_dt_ = std::min(max_dt_, current_dt_);
00327   current_dt_ = std::max(min_dt_, current_dt_);
00328 
00329   num_step_failures_ = std::min(num_step_failures_-1,0);
00330 
00331   trialStepCtrlInfo.stepSize = current_dt_;
00332   if (take_variable_steps_) {
00333     if (timeStepIter < num_constant_steps_ + num_ramping_steps_)
00334       trialStepCtrlInfo.stepType = STEP_TYPE_FIXED;
00335     else
00336       trialStepCtrlInfo.stepType = STEP_TYPE_VARIABLE;
00337   } else {
00338     trialStepCtrlInfo.stepType = STEP_TYPE_FIXED;
00339   }
00340 
00341   return trialStepCtrlInfo;
00342 }
00343 
00344 
00345 template<class Scalar>
00346 bool RampingIntegrationControlStrategy<Scalar>::resetForFailedTimeStep(
00347   const StepperBase<Scalar> &stepper,
00348   const StepControlInfo<Scalar> &stepCtrlInfoLast,
00349   const int timeStepIter,
00350   const StepControlInfo<Scalar> &stepCtrlInfo
00351   )
00352 {
00353   if (current_dt_ == min_dt_) return false;
00354   num_step_failures_++;
00355   if (num_step_failures_ > max_step_failures_) return false;
00356   current_dt_ = std::max(min_dt_, current_dt_/ramping_factor_);
00357   return true;
00358 }
00359 
00360 
00361 //
00362 // Explicit Instantiation macro
00363 //
00364 // Must be expanded from within the Rythmos namespace!
00365 //
00366 
00367 #define RYTHMOS_RAMPING_INTEGRATION_CONTROL_STRATEGY_INSTANT(SCALAR) \
00368   \
00369   template class RampingIntegrationControlStrategy< SCALAR >; \
00370   \
00371   template RCP<RampingIntegrationControlStrategy< SCALAR > > \
00372   rampingIntegrationControlStrategy(); \
00373   \
00374   template RCP<RampingIntegrationControlStrategy< SCALAR > > \
00375   rampingIntegrationControlStrategy( const RCP<ParameterList> &paramList );
00376 
00377 
00378 } // namespace Rythmos
00379 
00380 
00381 #endif // RYTHMOS_RAMPING_INTEGRATION_CONTROL_STRATEGY_DEF_HPP
 All Classes Functions Variables Typedefs Friends