Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Tempus_StepperNewmarkImplicitDForm_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_StepperNewmarkImplicitDForm_decl_hpp
10 #define Tempus_StepperNewmarkImplicitDForm_decl_hpp
11 
12 #include "Tempus_StepperImplicit.hpp"
13 #include "Tempus_WrapperModelEvaluatorSecondOrder.hpp"
14 
15 namespace Tempus {
16 
17 /** \brief Newmark time stepper.
18  *
19  * Here, we implement the Newmark scheme in displacement predictor/corrector
20  * form; see equations (34)-(35) in: A. Mota, W. Klug, M. Ortiz,
21  * "Finite element simulation of firearm injury to the human cranium",
22  * Computational Mechanics 31(1) 115-121 (2003).
23  *
24  * Newmark has two parameters: \f$\beta\f$
25  * and \f$\gamma\f$, both of which need to be in the range \f$[0,1]\f$.
26  * Newmark can be an explicit or implicit method, depending on
27  * the value of the \f$\beta\f$ parameter. If \f$\beta = 0\f$, the method
28  * is explicit; but note that the d-form of the Newmark scheme is not defined
29  * for the explicit case.
30  *
31  * Newmark is second order accurate if \f$\gamma = 0.5\f$; otherwise it
32  * is first order accurate. Some additional properties about the Newmark
33  * Beta scheme can be found
34  * <a href="http://opensees.berkeley.edu/wiki/index.php/Newmark_Method">here</a>.
35  *
36  * The First-Step-As-Last (FSAL) principle is not used with the
37  * Newmark implicit D-Form method.
38  */
39 template <class Scalar>
40 class StepperNewmarkImplicitDForm : virtual public Tempus::StepperImplicit<Scalar> {
41  public:
42 
43  /** \brief Default constructor.
44  *
45  * Requires subsequent setModel(), setSolver() and initialize()
46  * calls before calling takeStep().
47  */
49 
50  /// Constructor
52  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar>>& appModel,
53  const Teuchos::RCP<StepperObserver<Scalar> >& obs,
54  const Teuchos::RCP<Thyra::NonlinearSolverBase<Scalar> >& solver,
55  bool useFSAL,
56  std::string ICConsistency,
57  bool ICConsistencyCheck,
58  bool zeroInitialGuess,
59  std::string schemeName,
60  Scalar beta,
61  Scalar gamma);
62 
63  /// \name Basic stepper methods
64  //@{
65  virtual void
66  setModel(const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar>>& appModel);
67 
68  virtual void setObserver(
69  Teuchos::RCP<StepperObserver<Scalar> > /* obs */ = Teuchos::null){}
70 
71  virtual Teuchos::RCP<StepperObserver<Scalar> > getObserver() const
72  { return Teuchos::null; }
73 
74  /// Set the initial conditions and make them consistent.
75  virtual void setInitialConditions (
76  const Teuchos::RCP<SolutionHistory<Scalar> >& /* solutionHistory */){}
77 
78  /// Take the specified timestep, dt, and return true if successful.
79  virtual void
80  takeStep(const Teuchos::RCP<SolutionHistory<Scalar>>& solutionHistory);
81 
82  /// Get a default (initial) StepperState
83  virtual Teuchos::RCP<Tempus::StepperState<Scalar>>
85  virtual Scalar
86  getOrder() const {
87  if (gamma_ == 0.5)
88  return 2.0;
89  else
90  return 1.0;
91  }
92  virtual Scalar
93  getOrderMin() const {
94  return 1.0;
95  }
96  virtual Scalar
97  getOrderMax() const {
98  return 2.0;
99  }
100  virtual bool isExplicit() const {return false;}
101  virtual bool isImplicit() const {return true;}
102  virtual bool isExplicitImplicit() const
103  {return isExplicit() and isImplicit();}
104  virtual bool isOneStepMethod() const {return true;}
105  virtual bool isMultiStepMethod() const {return !isOneStepMethod();}
106 
107  virtual OrderODE getOrderODE() const {return SECOND_ORDER_ODE;}
108  //@}
109 
110  /// Return W_xDotxDot_coeff = d(xDotDot)/d(x).
111  virtual Scalar getW_xDotDot_coeff (const Scalar dt) const
112  { return Scalar(1.0)/(beta_*dt*dt); }
113  /// Return alpha = d(xDot)/d(x).
114  virtual Scalar getAlpha(const Scalar dt) const { return gamma_/(beta_*dt); }
115  /// Return beta = d(x)/d(x).
116  virtual Scalar getBeta (const Scalar ) const { return Scalar(1.0); }
117 
118  Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const;
119 
120  /// \name Overridden from Teuchos::Describable
121  //@{
122  virtual void describe(Teuchos::FancyOStream& out,
123  const Teuchos::EVerbosityLevel verbLevel) const;
124  //@}
125 
126  virtual bool isValidSetup(Teuchos::FancyOStream & out) const;
127 
128  void
130  Thyra::VectorBase<Scalar>& vPred, const Thyra::VectorBase<Scalar>& v,
131  const Thyra::VectorBase<Scalar>& a, const Scalar dt) const;
132 
133  void
135  Thyra::VectorBase<Scalar>& dPred, const Thyra::VectorBase<Scalar>& d,
136  const Thyra::VectorBase<Scalar>& v, const Thyra::VectorBase<Scalar>& a,
137  const Scalar dt) const;
138 
139  void
141  Thyra::VectorBase<Scalar>& v, const Thyra::VectorBase<Scalar>& vPred,
142  const Thyra::VectorBase<Scalar>& a, const Scalar dt) const;
143 
144  void
146  Thyra::VectorBase<Scalar>& d, const Thyra::VectorBase<Scalar>& dPred,
147  const Thyra::VectorBase<Scalar>& a, const Scalar dt) const;
148 
149  void
151  Thyra::VectorBase<Scalar>& a, const Thyra::VectorBase<Scalar>& dPred,
152  const Thyra::VectorBase<Scalar>& d, const Scalar dt) const;
153 
154 
155  void setSchemeName(std::string schemeName);
156  void setBeta(Scalar beta);
157  void setGamma(Scalar gamma);
158 
159  protected:
160 
161  std::string schemeName_;
162  Scalar beta_;
163  Scalar gamma_;
164 
165  Teuchos::RCP<Teuchos::FancyOStream> out_;
166 
167 };
168 } // namespace Tempus
169 
170 #endif // Tempus_StepperNewmarkImplicitDForm_decl_hpp
void correctVelocity(Thyra::VectorBase< Scalar > &v, const Thyra::VectorBase< Scalar > &vPred, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
void correctDisplacement(Thyra::VectorBase< Scalar > &d, const Thyra::VectorBase< Scalar > &dPred, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
virtual Teuchos::RCP< StepperObserver< Scalar > > getObserver() const
Get Observer.
virtual Teuchos::RCP< Tempus::StepperState< Scalar > > getDefaultStepperState()
Get a default (initial) StepperState.
Stepper integrates second-order ODEs.
Thyra Base interface for implicit time steppers.
virtual void takeStep(const Teuchos::RCP< SolutionHistory< Scalar >> &solutionHistory)
Take the specified timestep, dt, and return true if successful.
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
void predictVelocity(Thyra::VectorBase< Scalar > &vPred, const Thyra::VectorBase< Scalar > &v, const Thyra::VectorBase< Scalar > &a, const Scalar dt) const
virtual void setObserver(Teuchos::RCP< StepperObserver< Scalar > >=Teuchos::null)
Set Observer.
virtual void setInitialConditions(const Teuchos::RCP< SolutionHistory< Scalar > > &)
Set the initial conditions and make them consistent.
StepperObserver class for Stepper class.
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
virtual Scalar getW_xDotDot_coeff(const Scalar dt) const
Return W_xDotxDot_coeff = d(xDotDot)/d(x).
virtual bool isValidSetup(Teuchos::FancyOStream &out) const
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
virtual Scalar getAlpha(const Scalar dt) const
Return alpha = d(xDot)/d(x).
virtual Scalar getBeta(const Scalar) const
Return beta = d(x)/d(x).
virtual void setModel(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar >> &appModel)
void correctAcceleration(Thyra::VectorBase< Scalar > &a, const Thyra::VectorBase< Scalar > &dPred, const Thyra::VectorBase< Scalar > &d, const Scalar dt) const