Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Tempus_StepperHHTAlpha_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ****************************************************************************
3 // Tempus: Copyright (2017) Sandia Corporation
4 //
5 // Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
6 // ****************************************************************************
7 // @HEADER
8 
9 #ifndef Tempus_StepperHHTAlpha_decl_hpp
10 #define Tempus_StepperHHTAlpha_decl_hpp
11 
12 #include "Tempus_StepperImplicit.hpp"
13 #include "Tempus_WrapperModelEvaluatorSecondOrder.hpp"
14 
15 namespace Tempus {
16 
17 
18 /** \brief HHT-Alpha time stepper.
19  *
20  * Here, we implement the HHT-Alpha scheme in predictor/corrector form;
21  * see equations (10) and (13)-(19) in: G.M. Hulbert, J. Chung,
22  * "Explicit time integration algorithms for structural dynamics with
23  * optimal numerical dissipation", Comput. Methods Appl. Mech. Engrg.
24  * 137 175-188 (1996).
25  *
26  * There are four parameters in the scheme: \f$\alpha_m\f$, \f$\alpha_f\f$,
27  * \f$\beta\f$ and \f$\gamma\f$, all of which must be in the range \f$[0,1]\f$.
28  * When \f$\alpha_m=\alpha_f = 0\f$, the scheme reduces to the Newmark Beta
29  * scheme (see Tempus::StepperNewmark for details). Like the Newmark Beta
30  * scheme, the HHT-Alpha scheme can be either first or second order accurate,
31  * and either explicit or implicit.
32  *
33  * Although the general form of the scheme has been implemented in Tempus,
34  * it has only been verified for the case when \f$\alpha_m=\alpha_f = 0\f$
35  * (corresponding to the Newmark Beta) scheme, so other values for these
36  * parameters are not allowed at the present time. Also, note that, like
37  * the Newmark Beta stepper, the linear solve for the explicit version of
38  * this scheme has not been optimized (the mass matrix is not lumped).
39  *
40  * The First-Step-As-Last (FSAL) principle is not used with the
41  * HHT-Alpha method.
42  */
43 template<class Scalar>
44 class StepperHHTAlpha : virtual public Tempus::StepperImplicit<Scalar>
45 {
46 public:
47 
48  /** \brief Default constructor.
49  *
50  * - Constructs with a default ParameterList.
51  * - Can reset ParameterList with setParameterList().
52  * - Requires subsequent setModel() and initialize() calls before calling
53  * takeStep().
54  */
56 
57  /// Constructor
59  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel,
60  Teuchos::RCP<Teuchos::ParameterList> pList = Teuchos::null);
61 
62  /// \name Basic stepper methods
63  //@{
64  virtual void setModel(
65  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& appModel);
66 
67  virtual void setObserver(
68  Teuchos::RCP<StepperObserver<Scalar> > /* obs */ = Teuchos::null){}
69 
70  /// Initialize during construction and after changing input parameters.
71  virtual void initialize();
72 
73  /// Set the initial conditions and make them consistent.
74  virtual void setInitialConditions (
75  const Teuchos::RCP<SolutionHistory<Scalar> >& /* solutionHistory */){}
76 
77  /// Take the specified timestep, dt, and return true if successful.
78  virtual void takeStep(
79  const Teuchos::RCP<SolutionHistory<Scalar> >& solutionHistory);
80 
81  /// Get a default (initial) StepperState
82  virtual Teuchos::RCP<Tempus::StepperState<Scalar> > getDefaultStepperState();
83  virtual Scalar getOrder() const {
84  if (gamma_ == 0.5) return 2.0;
85  else return 1.0;
86  }
87  virtual Scalar getOrderMin() const {return 1.0;}
88  virtual Scalar getOrderMax() const {return 2.0;}
89 
90  virtual bool isExplicit() const {return false;}
91  virtual bool isImplicit() const {return true;}
92  virtual bool isExplicitImplicit() const
93  {return isExplicit() and isImplicit();}
94  virtual bool isOneStepMethod() const {return true;}
95  virtual bool isMultiStepMethod() const {return !isOneStepMethod();}
96 
97  virtual OrderODE getOrderODE() const {return SECOND_ORDER_ODE;}
98  //@}
99 
100  /// Return W_xDotxDot_coeff = d(xDotDot)/d(x).
101  virtual Scalar getW_xDotDot_coeff (const Scalar dt) const
102  { return Scalar(1.0)/(beta_*dt*dt); }
103  /// Return alpha = d(xDot)/d(x).
104  virtual Scalar getAlpha(const Scalar dt) const { return gamma_/(beta_*dt); }
105  /// Return beta = d(x)/d(x).
106  virtual Scalar getBeta (const Scalar ) const { return Scalar(1.0); }
107 
108  /// \name ParameterList methods
109  //@{
110  void setParameterList(const Teuchos::RCP<Teuchos::ParameterList> & pl);
111  Teuchos::RCP<Teuchos::ParameterList> getNonconstParameterList();
112  Teuchos::RCP<Teuchos::ParameterList> unsetParameterList();
113  Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const;
114  Teuchos::RCP<Teuchos::ParameterList> getDefaultParameters() const;
115  //@}
116 
117  /// \name Overridden from Teuchos::Describable
118  //@{
119  virtual std::string description() const;
120  virtual void describe(Teuchos::FancyOStream & out,
121  const Teuchos::EVerbosityLevel verbLevel) const;
122  //@}
123 
124  void predictVelocity(Thyra::VectorBase<Scalar>& vPred,
125  const Thyra::VectorBase<Scalar>& v,
126  const Thyra::VectorBase<Scalar>& a,
127  const Scalar dt) const;
128 
129  void predictDisplacement(Thyra::VectorBase<Scalar>& dPred,
130  const Thyra::VectorBase<Scalar>& d,
131  const Thyra::VectorBase<Scalar>& v,
132  const Thyra::VectorBase<Scalar>& a,
133  const Scalar dt) const;
134 
135  void predictVelocity_alpha_f(Thyra::VectorBase<Scalar>& vPred,
136  const Thyra::VectorBase<Scalar>& v) const;
137 
138  void predictDisplacement_alpha_f(Thyra::VectorBase<Scalar>& dPred,
139  const Thyra::VectorBase<Scalar>& d) const;
140 
141  void correctAcceleration(Thyra::VectorBase<Scalar>& a_n_plus1,
142  const Thyra::VectorBase<Scalar>& a_n) const;
143 
144  void correctVelocity(Thyra::VectorBase<Scalar>& v,
145  const Thyra::VectorBase<Scalar>& vPred,
146  const Thyra::VectorBase<Scalar>& a,
147  const Scalar dt) const;
148 
149  void correctDisplacement(Thyra::VectorBase<Scalar>& d,
150  const Thyra::VectorBase<Scalar>& dPred,
151  const Thyra::VectorBase<Scalar>& a,
152  const Scalar dt) const;
153 
154 private:
155 
156  Scalar beta_;
157  Scalar gamma_;
158  Scalar alpha_f_;
159  Scalar alpha_m_;
160 
161  Teuchos::RCP<Teuchos::FancyOStream> out_;
162 
163 };
164 } // namespace Tempus
165 
166 #endif // Tempus_StepperHHTAlpha_decl_hpp
void correctDisplacement(Thyra::VectorBase< Scalar > &d, const Thyra::VectorBase< Scalar > &dPred, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
virtual Scalar getW_xDotDot_coeff(const Scalar dt) const
Return W_xDotxDot_coeff = d(xDotDot)/d(x).
virtual bool isExplicitImplicit() const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
void correctVelocity(Thyra::VectorBase< Scalar > &v, const Thyra::VectorBase< Scalar > &vPred, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
virtual OrderODE getOrderODE() const
virtual Scalar getOrderMin() const
void correctAcceleration(Thyra::VectorBase< Scalar > &a_n_plus1, const Thyra::VectorBase< Scalar > &a_n) const
virtual Teuchos::RCP< Tempus::StepperState< Scalar > > getDefaultStepperState()
Get a default (initial) StepperState.
StepperHHTAlpha()
Default constructor.
virtual void setObserver(Teuchos::RCP< StepperObserver< Scalar > >=Teuchos::null)
Set Observer.
Stepper integrates second-order ODEs.
virtual Scalar getAlpha(const Scalar dt) const
Return alpha = d(xDot)/d(x).
Thyra Base interface for implicit time steppers.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &pl)
void predictVelocity_alpha_f(Thyra::VectorBase< Scalar > &vPred, const Thyra::VectorBase< Scalar > &v) const
StepperObserver class for Stepper class.
Teuchos::RCP< SolutionHistory< Scalar > > solutionHistory(Teuchos::RCP< Teuchos::ParameterList > pList=Teuchos::null)
Nonmember constructor.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
void predictDisplacement_alpha_f(Thyra::VectorBase< Scalar > &dPred, const Thyra::VectorBase< Scalar > &d) const
virtual Scalar getBeta(const Scalar) const
Return beta = d(x)/d(x).
void predictDisplacement(Thyra::VectorBase< Scalar > &dPred, const Thyra::VectorBase< Scalar > &d, const Thyra::VectorBase< Scalar > &v, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
virtual Scalar getOrderMax() const
virtual std::string description() const
Teuchos::RCP< Teuchos::FancyOStream > out_
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
void predictVelocity(Thyra::VectorBase< Scalar > &vPred, const Thyra::VectorBase< Scalar > &v, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
virtual void setInitialConditions(const Teuchos::RCP< SolutionHistory< Scalar > > &)
Set the initial conditions and make them consistent.
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
virtual bool isMultiStepMethod() const
virtual void initialize()
Initialize during construction and after changing input parameters.
Teuchos::RCP< Teuchos::ParameterList > getDefaultParameters() const
virtual void setModel(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel)
virtual void takeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Take the specified timestep, dt, and return true if successful.