9 #ifndef Tempus_StepperBDF2_impl_hpp
10 #define Tempus_StepperBDF2_impl_hpp
12 #include "Tempus_WrapperModelEvaluatorBasic.hpp"
14 #include "Tempus_StepperFactory.hpp"
20 template<
class Scalar>
23 this->setStepperName(
"BDF2");
24 this->setStepperType(
"BDF2");
25 this->setUseFSAL(
false);
26 this->setICConsistency(
"None");
27 this->setICConsistencyCheck(
false);
28 this->setZeroInitialGuess(
false);
30 this->setAppAction(Teuchos::null);
31 this->setDefaultSolver();
32 this->setStartUpStepper(
"DIRK 1 Stage Theta Method");
36 template<
class Scalar>
42 std::string ICConsistency,
43 bool ICConsistencyCheck,
44 bool zeroInitialGuess,
47 this->setStepperName(
"BDF2");
48 this->setStepperType(
"BDF2");
49 this->setUseFSAL( useFSAL);
50 this->setICConsistency( ICConsistency);
51 this->setICConsistencyCheck( ICConsistencyCheck);
52 this->setZeroInitialGuess( zeroInitialGuess);
54 this->setAppAction(stepperBDF2AppAction);
55 this->setSolver(solver);
56 this->setStartUpStepper(startUpStepper);
58 if (appModel != Teuchos::null) {
59 this->setModel(appModel);
64 template<
class Scalar>
70 if (startUpStepper_->getModel() == Teuchos::null) {
71 startUpStepper_->setModel(appModel);
72 startUpStepper_->initialize();
75 this->isInitialized_ =
false;
79 template<
class Scalar>
83 if (this->wrapperModel_ != Teuchos::null &&
84 this->wrapperModel_->getAppModel() != Teuchos::null) {
85 startUpStepper_ = sf->createStepper(startupStepperType,
86 this->wrapperModel_->getAppModel());
88 startUpStepper_ = sf->createStepper(startupStepperType);
91 this->isInitialized_ =
false;
95 template<
class Scalar>
99 startUpStepper_ = startUpStepper;
101 if (this->wrapperModel_ != Teuchos::null) {
103 this->wrapperModel_->getAppModel() == Teuchos::null, std::logic_error,
104 "Error - Can not set the startUpStepper to Teuchos::null.\n");
106 if (startUpStepper->getModel() == Teuchos::null &&
107 this->wrapperModel_->getAppModel() != Teuchos::null) {
108 startUpStepper_->setModel(this->wrapperModel_->getAppModel());
109 startUpStepper_->initialize();
113 this->isInitialized_ =
false;
116 template<
class Scalar>
120 if (appAction == Teuchos::null) {
122 stepperBDF2AppAction_ =
126 stepperBDF2AppAction_ = appAction;
132 template<
class Scalar>
139 template<
class Scalar>
145 RCP<SolutionState<Scalar> > initialState = solutionHistory->getCurrentState();
148 if (initialState->getXDot() == Teuchos::null)
149 this->setStepperXDot(initialState->getX()->clone_v());
151 this->setStepperXDot(initialState->getXDot());
157 template<
class Scalar>
161 this->checkInitialized();
165 TEMPUS_FUNC_TIME_MONITOR(
"Tempus::StepperBDF2::takeStep()");
167 int numStates = solutionHistory->getNumStates();
169 RCP<Thyra::VectorBase<Scalar> > xOld;
170 RCP<Thyra::VectorBase<Scalar> > xOldOld;
175 computeStartUp(solutionHistory);
179 "Error in Tempus::StepperBDF2::takeStep(): numStates after \n"
180 <<
"startup stepper must be at least 3, whereas numStates = "
181 << numStates <<
"!\n" <<
"If running with Storage Type = Static, "
182 <<
"make sure Storage Limit > 2.\n");
187 RCP<StepperBDF2<Scalar> > thisStepper = Teuchos::rcpFromRef(*
this);
188 stepperBDF2AppAction_->execute(solutionHistory, thisStepper,
192 RCP<SolutionState<Scalar> > workingState=solutionHistory->getWorkingState();
193 RCP<SolutionState<Scalar> > currentState=solutionHistory->getCurrentState();
195 RCP<Thyra::VectorBase<Scalar> > x = workingState->getX();
196 if (workingState->getXDot() != Teuchos::null)
197 this->setStepperXDot(workingState->getXDot());
198 RCP<Thyra::VectorBase<Scalar> > xDot = this->getStepperXDot();
201 const Scalar time = workingState->getTime();
202 const Scalar dt = workingState->getTimeStep();
203 const Scalar dtOld = currentState->getTimeStep();
205 xOld = solutionHistory->getStateTimeIndexNM1()->getX();
206 xOldOld = solutionHistory->getStateTimeIndexNM2()->getX();
207 order_ = Scalar(2.0);
213 const Scalar alpha = getAlpha(dt, dtOld);
214 const Scalar beta = getBeta (dt);
217 timeDer, dt, alpha, beta));
218 stepperBDF2AppAction_->execute(solutionHistory, thisStepper,
222 this->solveImplicitODE(x, xDot, time, p);
224 stepperBDF2AppAction_->execute(solutionHistory, thisStepper,
228 if (workingState->getXDot() != Teuchos::null)
229 timeDer->compute(x, xDot);
231 workingState->setSolutionStatus(sStatus);
232 workingState->setOrder(getOrder());
233 workingState->computeNorms(currentState);
234 stepperBDF2AppAction_->execute(solutionHistory, thisStepper,
240 template<
class Scalar>
244 auto out = Teuchos::fancyOStream( this->getOStream() );
245 out->setOutputToRootOnly(0);
247 *out <<
"Warning -- Taking a startup step for BDF2 using '"
248 << startUpStepper_->getStepperType()<<
"'!" << std::endl;
251 startUpStepper_->takeStep(solutionHistory);
253 order_ = startUpStepper_->getOrder();
262 template<
class Scalar>
273 template<
class Scalar>
278 auto l_out = Teuchos::fancyOStream( out.
getOStream() );
280 l_out->setOutputToRootOnly(0);
285 *l_out <<
"--- StepperBDF2 ---\n";
286 if (startUpStepper_ != Teuchos::null) {
287 *l_out <<
" startup stepper type = "
288 << startUpStepper_->description() << std::endl;
290 *l_out <<
" startUpStepper_ = "
291 << startUpStepper_ << std::endl;
292 *l_out <<
" startUpStepper_->isInitialized() = "
294 *l_out <<
" stepperBDF2AppAction_ = "
295 << stepperBDF2AppAction_ << std::endl;
296 *l_out <<
"----------------------------" << std::endl;
297 *l_out <<
" order_ = " << order_ << std::endl;
298 *l_out <<
"-------------------" << std::endl;
302 template<
class Scalar>
305 bool isValidSetup =
true;
306 auto l_out = Teuchos::fancyOStream( out.
getOStream() );
307 l_out->setOutputToRootOnly(0);
312 if ( !this->startUpStepper_->isInitialized() ) {
313 isValidSetup =
false;
314 *l_out <<
"The startup stepper is not initialized!\n";
316 if (stepperBDF2AppAction_ == Teuchos::null) {
317 isValidSetup =
false;
318 *l_out <<
"The BDF2 AppAction is not set!\n";
324 template<
class Scalar>
328 auto pl = this->getValidParametersBasicImplicit();
329 pl->set(
"Start Up Stepper Type", startUpStepper_->getStepperType());
336 template<
class Scalar>
343 stepper->setStepperImplicitValues(pl);
345 std::string startUpStepperName =
"DIRK 1 Stage Theta Method";
346 if (pl != Teuchos::null) startUpStepperName =
347 pl->
get<std::string>(
"Start Up Stepper Type", startUpStepperName);
348 stepper->setStartUpStepper(startUpStepperName);
350 if (model != Teuchos::null) {
351 stepper->setModel(model);
352 stepper->initialize();
360 #endif // Tempus_StepperBDF2_impl_hpp
BDF2 (Backward-Difference-Formula-2) time stepper.
virtual void setInitialConditions(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Set the initial conditions and make them consistent.
Time-derivative interface for BDF2.
virtual Teuchos::RCP< Tempus::StepperState< Scalar > > getDefaultStepperState()
Get a default (initial) StepperState.
virtual void setModel(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel)
Set the model.
T & get(const std::string &name, T def_value)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Default modifier for StepperBDF2.
virtual void initialize()
Initialize after construction and changing input parameters.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
virtual void computeStartUp(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Compute the first time step given the supplied startup stepper.
virtual void setModel(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel) override
Set the model.
Thyra Base interface for time steppers.
virtual void setInitialConditions(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory) override
Set the initial conditions and make them consistent.
StepperState is a simple class to hold state information about the stepper.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Thyra Base interface for implicit time steppers.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
virtual void setAppAction(Teuchos::RCP< StepperBDF2AppAction< Scalar > > appAction)
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
virtual bool isValidSetup(Teuchos::FancyOStream &out) const
RCP< std::basic_ostream< char_type, traits_type > > getOStream()
Teuchos::RCP< StepperBDF2< Scalar > > createStepperBDF2(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Application Action for StepperBDF2.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
virtual void takeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Take the specified timestep, dt, and return true if successful.
void setStartUpStepper(std::string startupStepperType)
Set the stepper to use in first step.
virtual void initialize()
Initialize during construction and after changing input parameters.
StepperBDF2()
Default constructor.
std::string toString(const T &t)