Rythmos - Transient Integration for Differential Equations  Version of the Day
Rythmos_ImplicitBDFStepperRampingStepControl_def.hpp
1 //@HEADER
2 // ***********************************************************************
3 //
4 // Rythmos Package
5 // Copyright (2006) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact Todd S. Coffey (tscoffe@sandia.gov)
25 //
26 // ***********************************************************************
27 //@HEADER
28 
29 #ifndef Rythmos_IMPLICITBDF_STEPPER_RAMPING_STEP_CONTROL_DEF_H
30 #define Rythmos_IMPLICITBDF_STEPPER_RAMPING_STEP_CONTROL_DEF_H
31 
32 #include "Rythmos_ImplicitBDFStepper.hpp"
33 #include "Rythmos_ImplicitBDFStepperErrWtVecCalc.hpp"
34 #include "Teuchos_StandardParameterEntryValidators.hpp"
35 
36 namespace Rythmos {
37 
38 template<class Scalar>
39 ImplicitBDFStepperRampingStepControl<Scalar>::
40 ImplicitBDFStepperRampingStepControl() :
41  stepControlState_(UNINITIALIZED)
42 {
43 
44 }
45 
46 template<class Scalar>
47 void ImplicitBDFStepperRampingStepControl<Scalar>::setStepControlState_(
48  StepControlStrategyState newState)
49 {
50  if (stepControlState_ == UNINITIALIZED) {
51  TEUCHOS_TEST_FOR_EXCEPT(newState != BEFORE_FIRST_STEP);
52  } else if (stepControlState_ == BEFORE_FIRST_STEP) {
53  TEUCHOS_TEST_FOR_EXCEPT(newState != MID_STEP);
54  } else if (stepControlState_ == MID_STEP) {
55  TEUCHOS_TEST_FOR_EXCEPT(newState != AFTER_CORRECTION);
56  } else if (stepControlState_ == AFTER_CORRECTION) {
57  TEUCHOS_TEST_FOR_EXCEPT(newState != READY_FOR_NEXT_STEP);
58  } else if (stepControlState_ == READY_FOR_NEXT_STEP) {
59  TEUCHOS_TEST_FOR_EXCEPT(newState != MID_STEP);
60  }
61  stepControlState_ = newState;
62 }
63 
64 template<class Scalar>
65 StepControlStrategyState
67 {
68  return(stepControlState_);
69 }
70 
71 template<class Scalar>
73 {
74  TEUCHOS_TEST_FOR_EXCEPT(!((stepControlState_ == BEFORE_FIRST_STEP) ||
75  (stepControlState_ == READY_FOR_NEXT_STEP)));
76 
77  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
78  "updateCoeffs_() is not implemented!");
79 }
80 
81 template<class Scalar>
83  const StepperBase<Scalar>& stepper)
84 {
85  // Initialize can be called from the stepper when setInitialCondition
86  // is called.
87  using Teuchos::as;
88  typedef Teuchos::ScalarTraits<Scalar> ST;
89  using Thyra::createMember;
90 
91  // Set initial time:
92  TimeRange<Scalar> stepperRange = stepper.getTimeRange();
93  TEUCHOS_TEST_FOR_EXCEPTION(
94  !stepperRange.isValid(),
95  std::logic_error,
96  "Error, Stepper does not have valid time range for initialization "
97  "of ImplicitBDFStepperRampingStepControl!\n");
98 
99  if (is_null(parameterList_)) {
100  RCP<Teuchos::ParameterList> emptyParameterList =
101  Teuchos::rcp(new Teuchos::ParameterList);
102  this->setParameterList(emptyParameterList);
103  }
104 
105  if (is_null(errWtVecCalc_)) {
106  RCP<ImplicitBDFStepperErrWtVecCalc<Scalar> > IBDFErrWtVecCalc =
107  rcp(new ImplicitBDFStepperErrWtVecCalc<Scalar>());
108  errWtVecCalc_ = IBDFErrWtVecCalc;
109  }
110 
111  stepControlState_ = UNINITIALIZED;
112 
113  requestedStepSize_ = Scalar(-1.0);
114  currentStepSize_ = initialStepSize_;
115  currentOrder_ = 1;
116  nextStepSize_ = initialStepSize_;
117  nextOrder_ = 1;
118  numberOfSteps_ = 0;
119  totalNumberOfFailedSteps_ = 0;
120  countOfConstantStepsAfterFailure_ = 0;
121 
122  if (is_null(delta_)) {
123  delta_ = createMember(stepper.get_x_space());
124  }
125  if (is_null(errWtVec_)) {
126  errWtVec_ = createMember(stepper.get_x_space());
127  }
128  V_S(delta_.ptr(),ST::zero());
129 
130  if ( doOutput_(Teuchos::VERB_HIGH) ) {
131  RCP<Teuchos::FancyOStream> out = this->getOStream();
132  Teuchos::OSTab ostab(out,1,"initialize");
133  *out << "currentOrder_ = " << currentOrder_ << std::endl;
134  *out << "numberOfSteps_ = " << numberOfSteps_ << std::endl;
135  }
136 
137  setStepControlState_(BEFORE_FIRST_STEP);
138 
139 }
140 
141 template<class Scalar>
143  const StepperBase<Scalar>& stepper,
144  const Scalar& stepSize,
145  const StepSizeType& stepSizeType)
146 {
147  typedef Teuchos::ScalarTraits<Scalar> ST;
148 
149  TEUCHOS_TEST_FOR_EXCEPT(!((stepControlState_ == UNINITIALIZED) ||
150  (stepControlState_ == BEFORE_FIRST_STEP) ||
151  (stepControlState_ == READY_FOR_NEXT_STEP) ||
152  (stepControlState_ == MID_STEP)));
153 
154  TEUCHOS_TEST_FOR_EXCEPTION(
155  ((stepSizeType == STEP_TYPE_FIXED) && (stepSize == ST::zero())),
156  std::logic_error,
157  "Error, step size type == STEP_TYPE_FIXED, "
158  "but requested step size == 0!\n");
159 
160  bool didInitialization = false;
161  if (stepControlState_ == UNINITIALIZED) {
162  initialize(stepper);
163  didInitialization = true;
164  }
165 
166  // errWtVecSet_ is called during initialize
167  if (!didInitialization) {
168  const ImplicitBDFStepper<Scalar>& implicitBDFStepper =
169  Teuchos::dyn_cast<const ImplicitBDFStepper<Scalar> >(stepper);
170  const Thyra::VectorBase<Scalar>& xHistory =
171  implicitBDFStepper.getxHistory(0);
172  errWtVecCalc_->errWtVecSet(&*errWtVec_,xHistory,relErrTol_,absErrTol_);
173  }
174 
175  requestedStepSize_ = stepSize;
176  stepSizeType_ = stepSizeType;
177 }
178 
179 template<class Scalar>
181  const StepperBase<Scalar>& stepper, Scalar* stepSize,
182  StepSizeType* stepSizeType, int* order)
183 {
184  TEUCHOS_TEST_FOR_EXCEPT(!((stepControlState_ == BEFORE_FIRST_STEP) ||
185  (stepControlState_ == MID_STEP) ||
186  (stepControlState_ == READY_FOR_NEXT_STEP) )
187  );
188 
189  if (stepControlState_ == BEFORE_FIRST_STEP) {
190  nextStepSize_ = initialStepSize_;
191  nextOrder_ = 1;
192  }
193 
194  // Now starting a step - rotate next values into current values
195  if (stepSizeType_ == STEP_TYPE_FIXED)
196  currentStepSize_ = requestedStepSize_;
197  else
198  currentStepSize_ = nextStepSize_;
199 
200  currentOrder_ = nextOrder_;
201 
202  // Limit the step size to the requested step size
203  currentStepSize_ = std::min(requestedStepSize_, currentStepSize_);
204 
205  *stepSize = currentStepSize_;
206  *stepSizeType = stepSizeType_;
207  *order = currentOrder_;
208 
209  if (stepControlState_ != MID_STEP) {
210  setStepControlState_(MID_STEP);
211  }
212 
213  // Output
214  if (doOutput_(Teuchos::VERB_MEDIUM)){
215  Teuchos::FancyOStream& out = *this->getOStream();
216  Teuchos::OSTab ostab1(out,2,"** nextStepSize_ **");
217  out << "Values returned to stepper:" << std::endl;
218  Teuchos::OSTab ostab2(out,2,"** nextStepSize_ **");
219  out << "currentStepSize_ = " << currentStepSize_ << std::endl;
220  out << "currentOrder_ = " << currentOrder_ << std::endl;
221  out << "requestedStepSize_ = " << requestedStepSize_ << std::endl;
222  }
223 
224 }
225 
226 template<class Scalar>
228  const StepperBase<Scalar>& stepper
229  ,const RCP<const Thyra::VectorBase<Scalar> >& soln
230  ,const RCP<const Thyra::VectorBase<Scalar> >& ee
231  ,int solveStatus)
232 {
233  TEUCHOS_TEST_FOR_EXCEPT(stepControlState_ != MID_STEP);
234 
235  TEUCHOS_TEST_FOR_EXCEPTION(is_null(ee), std::logic_error,
236  "Error, ee == Teuchos::null!\n");
237 
238  ee_ = ee;
239 
240  newtonConvergenceStatus_ = solveStatus;
241 
242  if ( doOutput_(Teuchos::VERB_MEDIUM) && newtonConvergenceStatus_ < 0) {
243  RCP<Teuchos::FancyOStream> out = this->getOStream();
244  Teuchos::OSTab ostab(out,1,"setCorrection");
245  *out << "\nImplicitBDFStepperRampingStepControl::setCorrection(): "
246  << "Nonlinear Solver Failed!\n";
247  }
248 
249  setStepControlState_(AFTER_CORRECTION);
250 }
251 
252 template<class Scalar>
254  const StepperBase<Scalar>& stepper, Scalar* LETValue)
255 {
256  using Teuchos::as;
257  typedef Teuchos::ScalarTraits<Scalar> ST;
258 
259  TEUCHOS_TEST_FOR_EXCEPT(stepControlState_ != AFTER_CORRECTION);
260 
261 
262  if ( doOutput_(Teuchos::VERB_HIGH) ) {
263  RCP<Teuchos::FancyOStream> out = this->getOStream();
264  Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
265  Teuchos::OSTab ostab(out,1,"acceptStep");
266  *out << "ee_ = " << std::endl;
267  ee_->describe(*out,verbLevel);
268  *out << "errWtVec_ = " << std::endl;
269  errWtVec_->describe(*out,verbLevel);
270  }
271 
272  Scalar enorm = wRMSNorm_(*errWtVec_,*ee_);
273 
274  Scalar LET = ck_ * enorm;
275 
276  if (LETValue) {
277  *LETValue = LET;
278  *LETValue = Scalar(0.0);
279  }
280 
281  if (newtonConvergenceStatus_ < 0)
282  return false;
283 
284  bool return_status = false;
285 
286  if (LET < ST::one() || !useLETToDetermineConvergence_)
287  return_status = true;
288 
289  if ( doOutput_(Teuchos::VERB_HIGH) ) {
290  RCP<Teuchos::FancyOStream> out = this->getOStream();
291  Teuchos::OSTab ostab(out,1,"acceptStep");
292  *out << "return_status = " << return_status << std::endl;
293  *out << "Local Truncation Error Check: (ck*enorm) < 1: (" << LET
294  << ") <?= 1" << std::endl;
295  if ( doOutput_(Teuchos::VERB_EXTREME) ) {
296  *out << "ck_ = " << ck_ << std::endl;
297  *out << "enorm = " << enorm << std::endl;
298  }
299  }
300 
301  return(return_status);
302 }
303 
304 template<class Scalar>
306  const StepperBase<Scalar>& stepper)
307 {
308  TEUCHOS_TEST_FOR_EXCEPT(stepControlState_ != AFTER_CORRECTION);
309  using Teuchos::as;
310  // typedef Teuchos::ScalarTraits<Scalar> ST; // unused
311 
312  if ( doOutput_(Teuchos::VERB_HIGH) ) {
313  RCP<Teuchos::FancyOStream> out = this->getOStream();
314 
315  Teuchos::OSTab ostab1(out,2,"completeStep_");
316  *out << "\n** Begin completeStep() **" << std::endl;
317  Teuchos::OSTab ostab2(out,2,"** Begin completeStep_ **");
318  *out << "numberOfSteps_ = " << numberOfSteps_ << std::endl;
319  *out << "numConstantSteps_ = " << numConstantSteps_ << std::endl;
320  *out << "currentStepSize_ = " << currentStepSize_ << std::endl;
321  *out << "nextStepSize_ = " << nextStepSize_ << std::endl;
322  *out << "currentOrder_ = " << currentOrder_ << std::endl;
323  *out << "nextOrder_ = " << nextOrder_ << std::endl;
324  *out << "stepSizeIncreaseFactor_ = " << stepSizeIncreaseFactor_ <<std::endl;
325  *out << "countOfConstantStepsAfterFailure_ = "
326  << countOfConstantStepsAfterFailure_ << std::endl;
327  }
328 
329  numberOfSteps_ ++;
330 
331  if (countOfConstantStepsAfterFailure_ > 0) {
332  // We track the number of consecutive time step failures so that
333  // if we have a bunch of nonlinear failures, lets keep the time
334  // step constant for a while before we start to ramp again. This
335  // keeps us from oscillating between ramping and cutting step
336  // sizes and wasting resources.
337 
338  nextStepSize_ = currentStepSize_;
339  nextOrder_ = currentOrder_;
340 
341  // Decrement failure counter
342  countOfConstantStepsAfterFailure_ =
343  std::max( (countOfConstantStepsAfterFailure_ - 1), 0);
344 
345  if ( doOutput_(Teuchos::VERB_HIGH) ) {
346  RCP<Teuchos::FancyOStream> out = this->getOStream();
347  Teuchos::OSTab ostab(out,1,"completeStep_");
348  *out << "\nNext Step Size held constant due to previous failed steps!\n";
349  *out << "countOfConstantStepsAfterFailure_ = "
350  << countOfConstantStepsAfterFailure_ << std::endl;
351  }
352  }
353  else {
354 
355  // Phase 1: Constant step size at 1st order
356  if (numberOfSteps_ < numConstantSteps_) {
357  if (currentStepSize_ < initialStepSize_)
358  nextStepSize_ = std::min(initialStepSize_,
359  currentStepSize_ * stepSizeIncreaseFactor_);
360  nextOrder_ = 1;
361  }
362  // Phase 2: Constant step size, ramping the order
363  else if (currentOrder_ < maxOrder_) {
364  if (currentStepSize_ < initialStepSize_)
365  nextStepSize_ = std::min(initialStepSize_,
366  currentStepSize_ * stepSizeIncreaseFactor_);
367  else
368  nextStepSize_ = currentStepSize_;
369 
370  nextOrder_ = currentOrder_ + 1;
371  }
372  // Phase 3: Ramping dt to max step size, highest order
373  else if ( (numberOfSteps_ >= numConstantSteps_) &&
374  (currentOrder_ == maxOrder_) ) {
375  nextStepSize_ = std::min(maxStepSize_,
376  currentStepSize_ * stepSizeIncreaseFactor_);
377  nextOrder_ = maxOrder_;
378  }
379  else {
380  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
381  "RampingStepControlStrategy logic is broken. Please contact "
382  "developers. Aborting run!");
383  }
384 
385  if (restrictStepSizeByNumberOfNonlinearIterations_) {
386  const Rythmos::ImplicitBDFStepper<Scalar>* ibdfStepper =
387  dynamic_cast<const Rythmos::ImplicitBDFStepper<Scalar>* >(&stepper);
388  TEUCHOS_ASSERT(ibdfStepper != NULL);
389  TEUCHOS_ASSERT(nonnull(ibdfStepper->getNonlinearSolveStatus().extraParameters));
390  int numberOfNonlinearIterations = ibdfStepper->getNonlinearSolveStatus().extraParameters->template get<int>("Number of Iterations");
391  if (numberOfNonlinearIterations >= numberOfNonlinearIterationsForStepSizeRestriction_) {
392  nextStepSize_ = currentStepSize_;
393  }
394  }
395 
396 
397  } // if (countOfConstantStepsAfterFailure_ > 0)
398 
399  setStepControlState_(READY_FOR_NEXT_STEP);
400 
401  if ( doOutput_(Teuchos::VERB_HIGH) ) {
402  RCP<Teuchos::FancyOStream> out = this->getOStream();
403  Teuchos::OSTab ostab1(out,2,"** completeStep_ **");
404  *out << "** End of completeStep() **" << std::endl;
405  Teuchos::OSTab ostab2(out,2,"** End completeStep_ **");
406  *out << "numberOfSteps_ = " << numberOfSteps_ << std::endl;
407  *out << "numConstantSteps_ = " << numConstantSteps_ << std::endl;
408  *out << "currentStepSize_ = " << currentStepSize_ << std::endl;
409  *out << "nextStepSize_ = " << nextStepSize_ << std::endl;
410  *out << "currentOrder_ = " << currentOrder_ << std::endl;
411  *out << "nextOrder_ = " << nextOrder_ << std::endl;
412  *out << "stepSizeIncreaseFactor_ = " << stepSizeIncreaseFactor_ <<std::endl;
413  *out << "countOfConstantStepsAfterFailure_ = "
414  << countOfConstantStepsAfterFailure_ << std::endl;
415  }
416 }
417 
418 template<class Scalar>
419 AttemptedStepStatusFlag
421  const StepperBase<Scalar>& stepper)
422 {
423  TEUCHOS_TEST_FOR_EXCEPT(stepControlState_ != AFTER_CORRECTION);
424 
425  using Teuchos::as;
426 
427  ++totalNumberOfFailedSteps_;
428  ++countOfConstantStepsAfterFailure_;
429 
430  // If time step size is already at the min time step, then quit
431  if (currentStepSize_ <= minStepSize_)
432  return (REP_ERR_FAIL);
433 
434  // Otherwise, cut the time step and keep order the same
435  nextStepSize_ = std::max(minStepSize_,
436  (currentStepSize_ * stepSizeDecreaseFactor_) );
437 
438  setStepControlState_(READY_FOR_NEXT_STEP);
439 
440  return(PREDICT_AGAIN);
441 }
442 
443 template<class Scalar>
445  Teuchos::FancyOStream &out,
446  const Teuchos::EVerbosityLevel verbLevel
447  ) const
448 {
449 
450  using Teuchos::as;
451 
452  if ( (as<int>(verbLevel) == as<int>(Teuchos::VERB_DEFAULT) ) ||
453  (as<int>(verbLevel) >= as<int>(Teuchos::VERB_LOW) )
454  ) {
455  out << this->description() << "::describe" << std::endl;
456  }
457  else if (as<int>(verbLevel) >= as<int>(Teuchos::VERB_LOW)) {
458  out << "currentStepSize_ = " << currentStepSize_ << std::endl;
459  out << "currentOrder_ = " << currentOrder_ << std::endl;
460  }
461  else if (as<int>(verbLevel) >= as<int>(Teuchos::VERB_MEDIUM)) {
462  }
463  else if (as<int>(verbLevel) >= as<int>(Teuchos::VERB_HIGH)) {
464  out << "ee_ = ";
465  if (ee_ == Teuchos::null) {
466  out << "Teuchos::null" << std::endl;
467  } else {
468  ee_->describe(out,verbLevel);
469  }
470  out << "delta_ = ";
471  if (delta_ == Teuchos::null) {
472  out << "Teuchos::null" << std::endl;
473  } else {
474  delta_->describe(out,verbLevel);
475  }
476  out << "errWtVec_ = ";
477  if (errWtVec_ == Teuchos::null) {
478  out << "Teuchos::null" << std::endl;
479  } else {
480  errWtVec_->describe(out,verbLevel);
481  }
482  }
483 }
484 
485 template<class Scalar>
487  RCP<Teuchos::ParameterList> const& paramList
488  )
489 {
490 
491  using Teuchos::as;
492  // typedef Teuchos::ScalarTraits<Scalar> ST; // unused
493 
494  TEUCHOS_TEST_FOR_EXCEPT(paramList == Teuchos::null);
495 
496  parameterList_ = Teuchos::parameterList(*paramList);
497 
498  parameterList_->validateParametersAndSetDefaults(*this->getValidParameters());
499 
500  Teuchos::ParameterList& p = *parameterList_;
501 
502  numConstantSteps_ = p.get<int>("Number of Constant First Order Steps");
503  initialStepSize_ = p.get<Scalar>("Initial Step Size");
504  minStepSize_ = p.get<Scalar>("Min Step Size");
505  maxStepSize_ = p.get<Scalar>("Max Step Size");
506  stepSizeIncreaseFactor_ = p.get<Scalar>("Step Size Increase Factor");
507  stepSizeDecreaseFactor_ = p.get<Scalar>("Step Size Decrease Factor");
508 
509  minOrder_ = p.get<int>("Min Order");
510  TEUCHOS_TEST_FOR_EXCEPTION(
511  !((1 <= minOrder_) && (minOrder_ <= 5)), std::logic_error,
512  "Error, minOrder_ = " << minOrder_ << " is not in range [1,5]!\n"
513  );
514  maxOrder_ = p.get<int>("Max Order");
515  TEUCHOS_TEST_FOR_EXCEPTION(
516  !((1 <= maxOrder_) && (maxOrder_ <= 5)), std::logic_error,
517  "Error, maxOrder_ = " << maxOrder_ << " is not in range [1,5]!\n"
518  );
519 
520  absErrTol_ = p.get<Scalar>("Absolute Error Tolerance");
521  relErrTol_ = p.get<Scalar>("Relative Error Tolerance");
522 
523  {
524  std::string let_acceptance =
525  p.get<std::string>("Use LET To Determine Step Acceptance");
526  useLETToDetermineConvergence_ = (let_acceptance == "TRUE");
527 
528  // Currently the using LET for step acceptance is not supported
529  // since we can't calculate the LETValue. Once this is
530  // implemented, delete the line below.
531  TEUCHOS_TEST_FOR_EXCEPTION(useLETToDetermineConvergence_, std::logic_error,
532  "Error - the flag \"Use LET To Determine Step Acceptance\" is set to "
533  "\"TRUE\" but the local error computation is currently not supported. "
534  "Please set this flag to \"FALSE\" for now.");
535  }
536 
537  if (p.get<std::string>("Restrict Step Size Increase by Number of Nonlinear Iterations") == "TRUE")
538  restrictStepSizeByNumberOfNonlinearIterations_ = true;
539  else if (p.get<std::string>("Restrict Step Size Increase by Number of Nonlinear Iterations") == "FALSE")
540  restrictStepSizeByNumberOfNonlinearIterations_ = false;
541 
542  numberOfNonlinearIterationsForStepSizeRestriction_ =
543  p.get<int>("Number of Nonlinear Iterations for Step Size Restriction");
544 
545  if ( doOutput_(Teuchos::VERB_HIGH) ) {
546  RCP<Teuchos::FancyOStream> out = this->getOStream();
547  Teuchos::OSTab ostab(out,1,"setParameterList");
548  out->precision(15);
549  *out << "minOrder_ = " << minOrder_ << std::endl;
550  *out << "maxOrder_ = " << maxOrder_ << std::endl;
551  *out << "relErrTol_ = " << relErrTol_ << std::endl;
552  *out << "absErrTol_ = " << absErrTol_ << std::endl;
553  *out << "stepSizeType = " << stepSizeType_ << std::endl;
554  *out << "stopTime_ = " << stopTime_ << std::endl;
555  }
556 
557 }
558 
559 template<class Scalar>
560 RCP<const Teuchos::ParameterList>
562 {
563  using Teuchos::RCP;
564  using Teuchos::rcp;
565  using Teuchos::ParameterList;
566 
567  static RCP<ParameterList> p;
568 
569  if (is_null(p)) {
570 
571  p = rcp(new ParameterList);
572 
573  p->set<int>("Number of Constant First Order Steps", 10,
574  "Number of constant steps to take before handing control to "
575  "variable stepper.");
576  p->set<Scalar>("Initial Step Size", Scalar(1.0e-3),
577  "Initial time step size and target step size to take during the "
578  "initial constant step phase (could be reduced due to step failures).");
579  p->set<Scalar>("Min Step Size", Scalar(1.0e-7), "Minimum time step size.");
580  p->set<Scalar>("Max Step Size", Scalar(1.0), "Maximum time step size.");
581  p->set<Scalar>("Step Size Increase Factor", Scalar(1.2),
582  "Time step growth factor used after a successful time step. dt_{n+1} = "
583  "(increase factor) * dt_n");
584  p->set<Scalar>("Step Size Decrease Factor", Scalar(0.5),
585  "Time step reduction factor used for a failed time step. dt_{n+1} = "
586  "(decrease factor) * dt_n");
587  p->set<int>("Min Order", 1, "Minimum order to run at.");
588  p->set<int>("Max Order", 5, "Maximum order to run at.");
589  p->set<Scalar>("Absolute Error Tolerance", Scalar(1.0e-5),
590  "abstol value used in WRMS calculation.");
591  p->set<Scalar>("Relative Error Tolerance", Scalar(1.0e-3),
592  "reltol value used in WRMS calculation.");
593  Teuchos::setStringToIntegralParameter<int>(
594  "Use LET To Determine Step Acceptance",
595  "FALSE",
596  "If set to TRUE, then acceptance of step dependes on LET in addition "
597  "to Nonlinear solver converging.",
598  Teuchos::tuple<std::string>("TRUE","FALSE"),
599  p.get());
600  Teuchos::setStringToIntegralParameter<int>(
601  "Restrict Step Size Increase by Number of Nonlinear Iterations",
602  "FALSE",
603  "If set to TRUE, then the step size will not be allowed to increase "
604  "if the number of nonlinear iterations was greater than or equal to the "
605  "specified value.",
606  Teuchos::tuple<std::string>("TRUE","FALSE"),
607  p.get());
608  p->set<int>("Number of Nonlinear Iterations for Step Size Restriction",
609  2,
610  "If \" Restrct Step Size Increase by Number of Nonlinear Iterations\" is "
611  "true, the step size will not be allowed to increase if the number of nonlinear "
612  "iterations was greater than or equal to the specified value.");
613  }
614 
615  return (p);
616 }
617 
618 template<class Scalar>
619 RCP<Teuchos::ParameterList>
621 {
622  RCP<Teuchos::ParameterList> temp_param_list = parameterList_;
623  parameterList_ = Teuchos::null;
624  return(temp_param_list);
625 }
626 
627 template<class Scalar>
628 RCP<Teuchos::ParameterList>
630 {
631  return(parameterList_);
632 }
633 
634 template<class Scalar>
636  const StepperBase<Scalar>& stepper)
637 {
638  if (stepControlState_ == UNINITIALIZED) {
639  initialize(stepper);
640  }
641  const ImplicitBDFStepper<Scalar>& bdfstepper =
642  Teuchos::dyn_cast<const ImplicitBDFStepper<Scalar> >(stepper);
643  int desiredOrder = bdfstepper.getOrder();
644  TEUCHOS_TEST_FOR_EXCEPT(!((1 <= desiredOrder) &&
645  (desiredOrder <= maxOrder_)));
646  if (stepControlState_ == BEFORE_FIRST_STEP) {
647  TEUCHOS_TEST_FOR_EXCEPTION(
648  desiredOrder > 1,
649  std::logic_error,
650  "Error, this ImplicitBDF stepper has not taken a step yet, so it "
651  "cannot take a step of order " << desiredOrder << " > 1!\n");
652  }
653  TEUCHOS_TEST_FOR_EXCEPT(!(desiredOrder <= currentOrder_+1));
654  currentOrder_ = desiredOrder;
655 
656  if ( doOutput_(Teuchos::VERB_EXTREME) ) {
657  RCP<Teuchos::FancyOStream> out = this->getOStream();
658  Teuchos::OSTab ostab(out,1,"setStepControlData");
659  *out << "currentOrder_ = " << currentOrder_ << std::endl;
660  }
661 }
662 
663 template<class Scalar>
665 {
666  return true;
667 }
668 
669 
670 template<class Scalar>
671 RCP<StepControlStrategyBase<Scalar> >
673 {
674 
675  RCP<ImplicitBDFStepperRampingStepControl<Scalar> > stepControl =
677 
678  if (!is_null(parameterList_)) {
679  stepControl->setParameterList(parameterList_);
680  }
681 
682  return stepControl;
683 }
684 
685 template<class Scalar>
687  const RCP<ErrWtVecCalcBase<Scalar> >& errWtVecCalc)
688 {
689  TEUCHOS_TEST_FOR_EXCEPT(is_null(errWtVecCalc));
690  errWtVecCalc_ = errWtVecCalc;
691 }
692 
693 template<class Scalar>
694 RCP<const ErrWtVecCalcBase<Scalar> >
696 {
697  return(errWtVecCalc_);
698 }
699 
700 template<class Scalar>
702  const Thyra::VectorBase<Scalar>& weight,
703  const Thyra::VectorBase<Scalar>& vector) const
704 {
705  return(norm_2(weight,vector));
706 }
707 
708 template<class Scalar>
710 {
711  TEUCHOS_TEST_FOR_EXCEPTION(
712  stepControlState_ == UNINITIALIZED, std::logic_error,
713  "Error, attempting to call getMinOrder before intiialization!\n"
714  );
715  return(minOrder_);
716 }
717 
718 template<class Scalar>
720 {
721  TEUCHOS_TEST_FOR_EXCEPTION(
722  stepControlState_ == UNINITIALIZED, std::logic_error,
723  "Error, attempting to call getMaxOrder before initialization!\n"
724  );
725  return(maxOrder_);
726 }
727 
728 template<class Scalar>
730  Teuchos::EVerbosityLevel verbLevel)
731 {
732  Teuchos::EVerbosityLevel currentObjectVerbLevel = this->getVerbLevel();
733 
734  if ( Teuchos::as<int>(currentObjectVerbLevel) >= Teuchos::as<int>(verbLevel) )
735  return true;
736 
737  return false;
738 }
739 
740 template<class Scalar>
742 {
743  return numberOfSteps_;
744 }
745 
746 template<class Scalar>
748 {
749  return totalNumberOfFailedSteps_;
750 }
751 
752 template<class Scalar>
754 {
755  return currentStepSize_;
756 }
757 
758 template<class Scalar>
760 {
761  return currentOrder_;
762 }
763 
764 
765 //
766 // Explicit Instantiation macro
767 //
768 // Must be expanded from within the Rythmos namespace!
769 //
770 
771 #define RYTHMOS_IMPLICITBDF_STEPPER_RAMPING_STEPCONTROL_INSTANT(SCALAR) \
772  template class ImplicitBDFStepperRampingStepControl< SCALAR >;
773 
774 
775 } // namespace Rythmos
776 #endif // Rythmos_IMPLICITBDF_STEPPER_RAMPING_STEP_CONTROL_DEF_H
777 
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Base class for defining stepper functionality.
void setParameterList(RCP< Teuchos::ParameterList > const &paramList)
const Thyra::SolveStatus< Scalar > & getNonlinearSolveStatus() const
Returns the Thyra::SolveStatus object from the last nonlinear solve.
bool acceptStep(const StepperBase< Scalar > &stepper, Scalar *LETValue)
void setRequestedStepSize(const StepperBase< Scalar > &stepper, const Scalar &stepSize, const StepSizeType &stepSizeType)
virtual RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const =0
Return the space for x and x_dot.
virtual TimeRange< Scalar > getTimeRange() const =0
Return the range of time values where interpolation calls can be performed.
void setErrWtVecCalc(const RCP< ErrWtVecCalcBase< Scalar > > &errWtVecCalc)
void setCorrection(const StepperBase< Scalar > &stepper, const RCP< const Thyra::VectorBase< Scalar > > &soln, const RCP< const Thyra::VectorBase< Scalar > > &ee, int solveStatus)
RCP< StepControlStrategyBase< Scalar > > cloneStepControlStrategyAlgorithm() const
void nextStepSize(const StepperBase< Scalar > &stepper, Scalar *stepSize, StepSizeType *stepSizeType, int *order)
const Thyra::VectorBase< Scalar > & getxHistory(int index) const
AttemptedStepStatusFlag rejectStep(const StepperBase< Scalar > &stepper)