Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
Tempus_StepperOptimizationInterface.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_Stepper_Optimization_Interface_hpp
10#define Tempus_Stepper_Optimization_Interface_hpp
11
12//Teuchos
13#include "Teuchos_Array.hpp"
14#include "Teuchos_RCP.hpp"
15
16// Thyra
17#include "Thyra_VectorBase.hpp"
18#include "Thyra_LinearOpBase.hpp"
19#include "Thyra_LinearOpWithSolveBase.hpp"
20#include "Tempus_config.hpp"
21
22namespace Tempus {
23
39template<class Scalar>
41{
42public:
43
45
47
49 virtual int stencilLength() const = 0;
50
52 virtual void computeStepResidual(
54 const Teuchos::Array< Teuchos::RCP<const Thyra::VectorBase<Scalar> > >& x,
55 const Teuchos::Array<Scalar>& t,
57 const int param_index) const = 0;
58
60
64 virtual void computeStepJacobian(
66 const Teuchos::Array< Teuchos::RCP<const Thyra::VectorBase<Scalar> > >& x,
67 const Teuchos::Array<Scalar>& t,
69 const int param_index,
70 const int deriv_index) const = 0;
71
75 const Teuchos::Array< Teuchos::RCP<const Thyra::VectorBase<Scalar> > >& x,
76 const Teuchos::Array<Scalar>& t,
78 const int param_index) const = 0;
79
81
84 virtual void computeStepSolver(
86 const Teuchos::Array< Teuchos::RCP<const Thyra::VectorBase<Scalar> > >& x,
87 const Teuchos::Array<Scalar>& t,
89 const int param_index) const = 0;
90};
91
92} // namespace Tempus
93
94#endif // Tempus_Stepper_hpp
Stepper interface to support full-space optimization.
virtual int stencilLength() const =0
Return the number of solution vectors in the time step stencil.
virtual void computeStepResidual(Thyra::VectorBase< Scalar > &residual, const Teuchos::Array< Teuchos::RCP< const Thyra::VectorBase< Scalar > > > &x, const Teuchos::Array< Scalar > &t, const Thyra::VectorBase< Scalar > &p, const int param_index) const =0
Compute time step residual.
virtual void computeStepParamDeriv(Thyra::LinearOpBase< Scalar > &deriv, const Teuchos::Array< Teuchos::RCP< const Thyra::VectorBase< Scalar > > > &x, const Teuchos::Array< Scalar > &t, const Thyra::VectorBase< Scalar > &p, const int param_index) const =0
Compute time step derivative w.r.t. model parameters.
virtual void computeStepJacobian(Thyra::LinearOpBase< Scalar > &jacobian, const Teuchos::Array< Teuchos::RCP< const Thyra::VectorBase< Scalar > > > &x, const Teuchos::Array< Scalar > &t, const Thyra::VectorBase< Scalar > &p, const int param_index, const int deriv_index) const =0
Compute time step Jacobian.
virtual void computeStepSolver(Thyra::LinearOpWithSolveBase< Scalar > &jacobian_solver, const Teuchos::Array< Teuchos::RCP< const Thyra::VectorBase< Scalar > > > &x, const Teuchos::Array< Scalar > &t, const Thyra::VectorBase< Scalar > &p, const int param_index) const =0
Compute time step Jacobian solver.