9 #ifndef TEMPUS_TEST_VANDERPOL_IMEX_EXPLICITMODEL_IMPL_HPP
10 #define TEMPUS_TEST_VANDERPOL_IMEX_EXPLICITMODEL_IMPL_HPP
12 #include "Teuchos_StandardParameterEntryValidators.hpp"
14 #include "Thyra_DefaultSpmdVectorSpace.hpp"
15 #include "Thyra_DetachedVectorView.hpp"
16 #include "Thyra_DetachedMultiVectorView.hpp"
17 #include "Thyra_DefaultSerialDenseLinearOpWithSolveFactory.hpp"
18 #include "Thyra_DefaultMultiVectorLinearOpWithSolve.hpp"
19 #include "Thyra_DefaultLinearOpSource.hpp"
20 #include "Thyra_VectorStdOps.hpp"
21 #include "Thyra_MultiVectorStdOps.hpp"
22 #include "Thyra_DefaultMultiVectorProductVector.hpp"
23 #include "Thyra_DefaultBlockedLinearOp.hpp"
27 namespace Tempus_Test {
29 template <
class Scalar>
33 useProductVector_ = useProductVector;
34 isInitialized_ =
false;
40 acceptModelParams_ =
false;
47 if (useProductVector_ ==
false) {
48 x_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(dim_);
49 f_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(dim_);
55 xSpace_ = Thyra::defaultSpmdVectorSpace<Scalar>(1);
56 for (
int i = 0; i < dim_; ++i) yxSpaceArray.
push_back(xSpace_);
57 x_space_ = Thyra::productVectorSpace<Scalar>(yxSpaceArray);
58 f_space_ = Thyra::productVectorSpace<Scalar>(yxSpaceArray);
62 p_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(np_);
63 g_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(ng_);
64 dxdp_space_ = Thyra::multiVectorProductVectorSpace(x_space_, np_);
66 setParameterList(pList_);
69 template <
class Scalar>
76 template <
class Scalar>
83 template <
class Scalar>
88 "Error, setupInOutArgs_ must be called first!\n");
89 return nominalValues_;
92 template <
class Scalar>
97 if (useProductVector_ ==
false)
98 W_op = Thyra::createMembers(x_space_, dim_);
103 Thyra::createMembers(xSpace_, 1);
105 Thyra::createMembers(xSpace_, 1);
107 Thyra::createMembers(xSpace_, 1);
109 Thyra::createMembers(xSpace_, 1);
110 W_op = Thyra::nonconstBlock2x2(A00, A01, A10, A11);
115 template <
class Scalar>
123 template <
class Scalar>
131 template <
class Scalar>
134 if (isInitialized_) {
142 inArgs.
setSupports(Thyra::ModelEvaluatorBase::IN_ARG_x_dot);
143 inArgs.
setSupports(Thyra::ModelEvaluatorBase::IN_ARG_x);
144 inArgs.
setSupports(Thyra::ModelEvaluatorBase::IN_ARG_t);
145 inArgs.
setSupports(Thyra::ModelEvaluatorBase::IN_ARG_alpha);
146 inArgs.
setSupports(Thyra::ModelEvaluatorBase::IN_ARG_beta);
147 if (acceptModelParams_) {
157 outArgs.
setSupports(Thyra::ModelEvaluatorBase::OUT_ARG_f);
158 outArgs.
setSupports(Thyra::ModelEvaluatorBase::OUT_ARG_W_op);
159 if (acceptModelParams_) {
161 outArgs.
setSupports(Thyra::ModelEvaluatorBase::OUT_ARG_DfDp, 0,
162 Thyra::ModelEvaluatorBase::DERIV_MV_JACOBIAN_FORM);
168 nominalValues_ = inArgs_;
172 nominalValues_.set_t(t0_ic_);
176 x_ic_view[0] = x0_ic_;
177 x_ic_view[1] = x1_ic_;
179 nominalValues_.set_x(x_ic);
181 if (acceptModelParams_) {
185 p_ic_view[0] = epsilon_;
187 nominalValues_.set_p(0, p_ic);
192 x_dot_ic_view[0] = x1_ic_;
193 x_dot_ic_view[1] = ((1.0 - x0_ic_ * x0_ic_) * x1_ic_ - x0_ic_) / epsilon_;
195 nominalValues_.set_x_dot(x_dot_ic);
198 isInitialized_ =
true;
201 template <
class Scalar>
209 if (paramList != Teuchos::null) tmpPL = paramList;
211 this->setMyParamList(tmpPL);
213 bool acceptModelParams = get<bool>(*pl,
"Accept model parameters");
214 bool haveIC = get<bool>(*pl,
"Provide nominal values");
215 bool useDfDpAsTangent = get<bool>(*pl,
"Use DfDp as Tangent");
216 if ((acceptModelParams != acceptModelParams_) || (haveIC != haveIC_)) {
217 isInitialized_ =
false;
219 acceptModelParams_ = acceptModelParams;
221 useDfDpAsTangent_ = useDfDpAsTangent;
222 epsilon_ = get<Scalar>(*pl,
"Coeff epsilon");
223 x0_ic_ = get<Scalar>(*pl,
"IC x0");
224 x1_ic_ = get<Scalar>(*pl,
"IC x1");
225 t0_ic_ = get<Scalar>(*pl,
"IC t0");
229 template <
class Scalar>
236 pl->
set(
"Accept model parameters",
false);
237 pl->
set(
"Provide nominal values",
true);
238 pl->
set(
"Use DfDp as Tangent",
false);
239 Teuchos::setDoubleParameter(
"Coeff epsilon", 1.0e-06,
240 "Coefficient a in model", &*pl);
241 Teuchos::setDoubleParameter(
"IC x0", 2.0,
"Initial Condition for x0", &*pl);
242 Teuchos::setDoubleParameter(
"IC x1", 0.0,
"Initial Condition for x1", &*pl);
243 Teuchos::setDoubleParameter(
"IC t0", 0.0,
"Initial time t0", &*pl);
249 template <
class Scalar>
253 if (!acceptModelParams_) {
254 return Teuchos::null;
259 else if (l == 1 || l == 2)
261 return Teuchos::null;
264 template <
class Scalar>
268 if (!acceptModelParams_) {
269 return Teuchos::null;
275 p_strings->push_back(
"Model Coefficient: epsilon");
277 p_strings->push_back(
"DxDp");
279 p_strings->push_back(
"Dx_dotDp");
283 template <
class Scalar>
291 template <
class Scalar>
298 using Teuchos::rcp_dynamic_cast;
300 "Error, setupInOutArgs_ must be called first!\n");
303 inArgs.
get_x().assert_not_null();
308 Scalar
eps = epsilon_;
310 if (acceptModelParams_) {
312 if (p_in != Teuchos::null) {
316 if (inArgs.
get_p(1) != Teuchos::null) {
317 dxdp_in = rcp_dynamic_cast<
const DMVPV>(inArgs.
get_p(1),
true)
320 if (inArgs.
get_p(2) != Teuchos::null) {
321 dxdotdp_in = rcp_dynamic_cast<
const DMVPV>(inArgs.
get_p(2),
true)
329 if (acceptModelParams_) {
338 f_out_view[0] = x_in_view[1];
339 f_out_view[1] = -x_in_view[0] /
eps;
343 DfDp_out_view(0, 0) = 0.0;
344 DfDp_out_view(1, 0) = x_in_view[0] / (eps *
eps);
347 if (useDfDpAsTangent_ && !
is_null(dxdp_in)) {
349 DfDp_out_view(0, 0) += dxdp(1, 0);
350 DfDp_out_view(1, 0) += -dxdp(0, 0) /
eps;
354 if (useProductVector_ ==
false) {
361 W_view(1, 0) = -beta /
eps;
371 W_block->getNonconstBlock(0, 0),
true);
374 W_block->getNonconstBlock(0, 1),
true);
377 W_block->getNonconstBlock(1, 0),
true);
380 W_block->getNonconstBlock(1, 1),
true);
385 W00_view(0, 0) = 0.0;
386 W01_view(0, 0) = beta;
387 W10_view(0, 0) = -beta /
eps;
388 W11_view(0, 0) = 0.0;
396 x_dot_in = inArgs.
get_x_dot().assert_not_null();
402 f_out_view[0] = x_dot_in_view[0] - x_in_view[1];
403 f_out_view[1] = x_dot_in_view[1] + x_in_view[0] /
eps;
407 DfDp_out_view(0, 0) = 0.0;
408 DfDp_out_view(1, 0) = -x_in_view[0] / (eps *
eps);
411 if (useDfDpAsTangent_ && !
is_null(dxdotdp_in)) {
413 DfDp_out_view(0, 0) += dxdotdp(0, 0);
414 DfDp_out_view(1, 0) += dxdotdp(1, 0);
416 if (useDfDpAsTangent_ && !
is_null(dxdp_in)) {
418 DfDp_out_view(0, 0) += -dxdp(1, 0);
419 DfDp_out_view(1, 0) += dxdp(0, 0) /
eps;
423 if (useProductVector_ ==
false) {
428 W_view(0, 0) = alpha;
429 W_view(0, 1) = -beta;
430 W_view(1, 0) = beta /
eps;
431 W_view(1, 1) = alpha;
440 W_block->getNonconstBlock(0, 0),
true);
443 W_block->getNonconstBlock(0, 1),
true);
446 W_block->getNonconstBlock(1, 0),
true);
449 W_block->getNonconstBlock(1, 1),
true);
454 W00_view(0, 0) = alpha;
455 W01_view(0, 0) = -beta;
456 W10_view(0, 0) = beta /
eps;
457 W11_view(0, 0) = alpha;
465 #endif // TEMPUS_TEST_VANDERPOL_IMEX_EXPLICITMODEL_IMPL_HPP
bool is_null(const boost::shared_ptr< T > &p)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Derivative< Scalar > get_DfDp(int l) const
RCP< const VectorBase< Scalar > > get_x_dot() const
void setModelEvalDescription(const std::string &modelEvalDescription)
void setupInOutArgs_() const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_f_space() const
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
void set_Np_Ng(int Np, int Ng)
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_p_space(int l) const
Evaluation< VectorBase< Scalar > > get_f() const
void setParameterList(Teuchos::RCP< Teuchos::ParameterList > const ¶mList)
#define TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(index, lower_inclusive, upper_exclusive)
void setSupports(EInArgsMembers arg, bool supports=true)
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_W_op() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Thyra::ModelEvaluatorBase::InArgs< Scalar > getNominalValues() const
Thyra::ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
void validateParametersAndSetDefaults(ParameterList const &validParamList, int const depth=1000)
void push_back(const value_type &x)
void setSupports(EOutArgsMembers arg, bool supports=true)
Thyra::ModelEvaluatorBase::OutArgs< Scalar > createOutArgsImpl() const
VanDerPol_IMEX_ExplicitModel(Teuchos::RCP< Teuchos::ParameterList > pList=Teuchos::null, bool useProductVector=false)
void evalModelImpl(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs_bar, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs_bar) const
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_g_space(int j) const
void setModelEvalDescription(const std::string &modelEvalDescription)
RCP< LinearOpBase< Scalar > > get_W_op() const
RCP< const VectorBase< Scalar > > get_x() const
RCP< const VectorBase< Scalar > > get_p(int l) const
RCP< MultiVectorBase< Scalar > > getMultiVector() const