9 #ifndef TEMPUS_TEST_DAHLQUIST_TEST_MODEL_IMPL_HPP
10 #define TEMPUS_TEST_DAHLQUIST_TEST_MODEL_IMPL_HPP
12 #include "Thyra_DefaultSpmdVectorSpace.hpp"
13 #include "Thyra_DetachedVectorView.hpp"
14 #include "Thyra_DetachedMultiVectorView.hpp"
15 #include "Thyra_DefaultSerialDenseLinearOpWithSolveFactory.hpp"
16 #include "Thyra_DefaultMultiVectorLinearOpWithSolve.hpp"
17 #include "Thyra_DefaultLinearOpSource.hpp"
18 #include "Thyra_VectorStdOps.hpp"
19 #include "Thyra_MultiVectorStdOps.hpp"
20 #include "Thyra_DefaultMultiVectorProductVector.hpp"
24 namespace Tempus_Test {
26 template <
class Scalar>
29 constructDahlquistTestModel(-1.0,
false);
32 template <
class Scalar>
35 constructDahlquistTestModel(lambda, includeXDot);
38 template <
class Scalar>
43 includeXDot_ = includeXDot;
44 isInitialized_ =
false;
46 xDotIC_ = Scalar(lambda);
53 acceptModelParams_ =
false;
56 x_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(dim_);
57 f_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(dim_);
63 MEB::InArgsSetup<Scalar> inArgs;
64 inArgs.setModelEvalDescription(this->description());
65 inArgs.setSupports(MEB::IN_ARG_t);
66 inArgs.setSupports(MEB::IN_ARG_x);
67 inArgs.setSupports(MEB::IN_ARG_x_dot);
69 inArgs.setSupports(MEB::IN_ARG_beta);
70 inArgs.setSupports(MEB::IN_ARG_alpha);
71 if (acceptModelParams_) inArgs.set_Np(Np_);
78 MEB::OutArgsSetup<Scalar> outArgs;
79 outArgs.setModelEvalDescription(this->description());
80 outArgs.setSupports(MEB::OUT_ARG_f);
82 outArgs.setSupports(MEB::OUT_ARG_W_op);
83 if (acceptModelParams_) {
84 outArgs.set_Np_Ng(Np_, Ng_);
85 outArgs.setSupports(MEB::OUT_ARG_DfDp, 0, MEB::DERIV_MV_JACOBIAN_FORM);
86 outArgs.setSupports(MEB::OUT_ARG_DgDp, 0, 0, MEB::DERIV_MV_JACOBIAN_FORM);
87 outArgs.setSupports(MEB::OUT_ARG_DgDx, 0, MEB::DERIV_MV_GRADIENT_FORM);
93 nominalValues_ = inArgs_;
95 nominalValues_.set_t(Scalar(0.0));
101 nominalValues_.set_x(x_ic);
108 x_dot_ic_view[0] = xDotIC_;
110 nominalValues_.set_x_dot(x_dot_ic);
113 isInitialized_ =
true;
116 template <
class Scalar>
122 inArgs.
set_t(exact_t);
128 exact_x_view[0] = exp(lambda_ * exact_t);
130 inArgs.
set_x(exact_x);
135 createMember(x_space_);
138 exact_x_dot_view[0] = lambda_ * exp(lambda_ * exact_t);
146 template <
class Scalar>
153 template <
class Scalar>
160 template <
class Scalar>
164 return nominalValues_;
167 template <
class Scalar>
174 template <
class Scalar>
181 template <
class Scalar>
187 using Teuchos::rcp_dynamic_cast;
191 "Error, setupInOutArgs_ must be called first!\n");
202 f_out_view[0] = lambda_ * x_in_view[0];
206 true, std::logic_error,
207 "Error -- Dahlquist Test Model requires f_out!\n");
215 matrix_view(0, 0) = lambda_;
221 x_dot_in = inArgs.
get_x_dot().assert_not_null();
228 f_out_view[0] = x_dot_in_view[0] - lambda_ * x_in_view[0];
236 matrix_view(0, 0) = alpha - beta * lambda_;
242 template <
class Scalar>
248 this->get_W_factory();
258 vec_view[0] = lambda_;
260 V_V(multivec->col(0).
ptr(), *vec);
264 Thyra::linearOpWithSolve<Scalar>(*W_factory, matrix);
268 template <
class Scalar>
274 Thyra::createMembers(x_space_, dim_);
278 template <
class Scalar>
283 Thyra::defaultSerialDenseLinearOpWithSolveFactory<Scalar>();
287 template <
class Scalar>
291 if (!acceptModelParams_) {
292 return Teuchos::null;
297 else if (l == 1 || l == 2)
299 return Teuchos::null;
302 template <
class Scalar>
306 if (!acceptModelParams_) {
307 return Teuchos::null;
313 p_strings->push_back(
"Model Coefficient: a");
318 p_strings->push_back(
"DxDp");
320 p_strings->push_back(
"Dx_dotDp");
324 template <
class Scalar>
333 #endif // TEMPUS_TEST_DAHLQUIST_TEST_MODEL_IMPL_HPP
Thyra::ModelEvaluatorBase::InArgs< Scalar > getNominalValues() const
void evalModelImpl(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs_bar, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs_bar) const
bool is_null(const boost::shared_ptr< T > &p)
RCP< const VectorBase< Scalar > > get_x_dot() const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_W_op() const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Evaluation< VectorBase< Scalar > > get_f() const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const
#define TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(index, lower_inclusive, upper_exclusive)
Thyra::ModelEvaluatorBase::OutArgs< Scalar > createOutArgsImpl() const
Thyra::ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
void set_x(const RCP< const VectorBase< Scalar > > &x)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Thyra::ModelEvaluatorBase::InArgs< Scalar > getExactSolution(double t) const
Exact solution.
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_f_space() const
Teuchos::RCP< const Thyra::LinearOpWithSolveFactoryBase< Scalar > > get_W_factory() const
void set_x_dot(const RCP< const VectorBase< Scalar > > &x_dot)
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_g_space(int j) const
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_p_space(int l) const
Teuchos::RCP< Thyra::LinearOpWithSolveBase< Scalar > > create_W() const
void constructDahlquistTestModel(Scalar lambda, bool includeXDot)
RCP< LinearOpBase< Scalar > > get_W_op() const
RCP< const VectorBase< Scalar > > get_x() const