109 RCP<const Thyra::LinearOpWithSolveFactoryBase<Scalar> > W_factory = this->get_W_factory();
110 RCP<Thyra::LinearOpBase<Scalar> > matrix = this->create_W_op();
117 RCP<Thyra::MultiVectorBase<Scalar> > multivec = Teuchos::rcp_dynamic_cast<Thyra::MultiVectorBase<Scalar> >(matrix,
true);
119 RCP<Thyra::VectorBase<Scalar> > vec = Thyra::createMember(x_space_);
121 Thyra::DetachedVectorView<Scalar> vec_view( *vec );
125 V_V(multivec->col(0).ptr(),*vec);
127 Thyra::DetachedVectorView<Scalar> vec_view( *vec );
131 V_V(multivec->col(1).ptr(),*vec);
134 RCP<Thyra::LinearOpWithSolveBase<Scalar> > W =
135 Thyra::linearOpWithSolve<Scalar>(*W_factory, matrix);
189 const Thyra::ModelEvaluatorBase::InArgs<Scalar> &inArgs,
190 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs
194 TEUCHOS_TEST_FOR_EXCEPTION( !isInitialized_, std::logic_error,
195 "Error, setupInOutArgs_ must be called first!\n");
197 const RCP<const Thyra::VectorBase<Scalar> > x_in =
198 inArgs.get_x().assert_not_null();
199 Thyra::ConstDetachedVectorView<Scalar> x_in_view( *x_in );
202 Scalar epsilon = epsilon_;
203 if (acceptModelParams_) {
204 const RCP<const Thyra::VectorBase<Scalar> > p_in =
205 inArgs.get_p(0).assert_not_null();
206 Thyra::ConstDetachedVectorView<Scalar> p_in_view( *p_in );
207 epsilon = p_in_view[0];
210 Scalar beta = inArgs.get_beta();
212 const RCP<Thyra::VectorBase<Scalar> > f_out = outArgs.get_f();
213 const RCP<Thyra::LinearOpBase<Scalar> > W_out = outArgs.get_W_op();
214 RCP<Thyra::MultiVectorBase<Scalar> > DfDp_out;
215 if (acceptModelParams_) {
216 Thyra::ModelEvaluatorBase::Derivative<Scalar> DfDp = outArgs.get_DfDp(0);
217 DfDp_out = DfDp.getMultiVector();
220 if (inArgs.get_x_dot().is_null()) {
223 if (!is_null(f_out)) {
224 Thyra::DetachedVectorView<Scalar> f_out_view( *f_out );
225 f_out_view[0] = x_in_view[1];
227 ((1.0-x_in_view[0]*x_in_view[0])*x_in_view[1]-x_in_view[0])/epsilon;
229 if (!is_null(W_out)) {
230 RCP<Thyra::MultiVectorBase<Scalar> > matrix =
231 Teuchos::rcp_dynamic_cast<Thyra::MultiVectorBase<Scalar> >(W_out,
true);
232 Thyra::DetachedMultiVectorView<Scalar> matrix_view( *matrix );
233 matrix_view(0,0) = 0.0;
234 matrix_view(0,1) = +beta;
236 -beta*(2.0*x_in_view[0]*x_in_view[1]+1.0)/epsilon;
238 beta*(1.0 - x_in_view[0]*x_in_view[0])/epsilon;
241 if (!is_null(DfDp_out)) {
242 Thyra::DetachedMultiVectorView<Scalar> DfDp_out_view( *DfDp_out );
243 DfDp_out_view(0,0) = 0.0;
245 -((1.0-x_in_view[0]*x_in_view[0])*x_in_view[1]-x_in_view[0])
251 RCP<const Thyra::VectorBase<Scalar> > x_dot_in;
252 x_dot_in = inArgs.get_x_dot().assert_not_null();
253 Scalar alpha = inArgs.get_alpha();
254 if (!is_null(f_out)) {
255 Thyra::DetachedVectorView<Scalar> f_out_view( *f_out );
256 Thyra::ConstDetachedVectorView<Scalar> x_dot_in_view( *x_dot_in );
257 f_out_view[0] = x_dot_in_view[0] - x_in_view[1];
258 f_out_view[1] = x_dot_in_view[1]
259 - ((1.0-x_in_view[0]*x_in_view[0])*x_in_view[1]-x_in_view[0])/epsilon;
261 if (!is_null(W_out)) {
262 RCP<Thyra::MultiVectorBase<Scalar> > matrix =
263 Teuchos::rcp_dynamic_cast<Thyra::MultiVectorBase<Scalar> >(W_out,
true);
264 Thyra::DetachedMultiVectorView<Scalar> matrix_view( *matrix );
265 matrix_view(0,0) = alpha;
266 matrix_view(0,1) = -beta;
268 beta*(2.0*x_in_view[0]*x_in_view[1]+1.0)/epsilon;
269 matrix_view(1,1) = alpha
270 - beta*(1.0 - x_in_view[0]*x_in_view[0])/epsilon;
273 if (!is_null(DfDp_out)) {
274 Thyra::DetachedMultiVectorView<Scalar> DfDp_out_view( *DfDp_out );
275 DfDp_out_view(0,0) = 0.0;
277 ((1.0-x_in_view[0]*x_in_view[0])*x_in_view[1]-x_in_view[0])
326 if (isInitialized_) {
332 Thyra::ModelEvaluatorBase::InArgsSetup<Scalar> inArgs;
333 inArgs.setModelEvalDescription(this->description());
334 inArgs.setSupports( Thyra::ModelEvaluatorBase::IN_ARG_t );
335 inArgs.setSupports( Thyra::ModelEvaluatorBase::IN_ARG_x );
336 inArgs.setSupports( Thyra::ModelEvaluatorBase::IN_ARG_beta );
337 inArgs.setSupports( Thyra::ModelEvaluatorBase::IN_ARG_x_dot );
338 inArgs.setSupports( Thyra::ModelEvaluatorBase::IN_ARG_alpha );
339 if (acceptModelParams_) {
347 Thyra::ModelEvaluatorBase::OutArgsSetup<Scalar> outArgs;
348 outArgs.setModelEvalDescription(this->description());
349 outArgs.setSupports( Thyra::ModelEvaluatorBase::OUT_ARG_f );
350 outArgs.setSupports( Thyra::ModelEvaluatorBase::OUT_ARG_W_op );
351 if (acceptModelParams_) {
352 outArgs.set_Np_Ng(Np_,Ng_);
353 outArgs.setSupports( Thyra::ModelEvaluatorBase::OUT_ARG_DfDp,0,
354 Thyra::ModelEvaluatorBase::DERIV_MV_JACOBIAN_FORM );
360 nominalValues_ = inArgs_;
364 nominalValues_.set_t(t0_ic_);
365 const RCP<Thyra::VectorBase<Scalar> > x_ic = createMember(x_space_);
367 Thyra::DetachedVectorView<Scalar> x_ic_view( *x_ic );
368 x_ic_view[0] = x0_ic_;
369 x_ic_view[1] = x1_ic_;
371 nominalValues_.set_x(x_ic);
372 if (acceptModelParams_) {
373 const RCP<Thyra::VectorBase<Scalar> > p_ic = createMember(p_space_);
375 Thyra::DetachedVectorView<Scalar> p_ic_view( *p_ic );
376 p_ic_view[0] = epsilon_;
378 nominalValues_.set_p(0,p_ic);
380 const RCP<Thyra::VectorBase<Scalar> > x_dot_ic = createMember(x_space_);
382 Thyra::DetachedVectorView<Scalar> x_dot_ic_view( *x_dot_ic );
383 x_dot_ic_view[0] = x1_ic_;
384 x_dot_ic_view[1] = ((1.0-x0_ic_*x0_ic_)*x1_ic_-x0_ic_)/epsilon_;
386 nominalValues_.set_x_dot(x_dot_ic);
389 isInitialized_ =
true;
399 using Teuchos::ParameterList;
400 Teuchos::RCP<ParameterList> tmpPL = Teuchos::rcp(
new ParameterList(
"VanDerPolModel"));
401 if (paramList != Teuchos::null) tmpPL = paramList;
402 tmpPL->validateParametersAndSetDefaults(*this->getValidParameters());
403 this->setMyParamList(tmpPL);
404 Teuchos::RCP<ParameterList> pl = this->getMyNonconstParamList();
405 bool acceptModelParams = get<bool>(*pl,
"Accept model parameters");
406 bool haveIC = get<bool>(*pl,
"Provide nominal values");
407 if ( (acceptModelParams != acceptModelParams_) ||
410 isInitialized_ =
false;
412 acceptModelParams_ = acceptModelParams;
414 epsilon_ = get<Scalar>(*pl,
"Coeff epsilon");
415 x0_ic_ = get<Scalar>(*pl,
"IC x0");
416 x1_ic_ = get<Scalar>(*pl,
"IC x1");
417 t0_ic_ = get<Scalar>(*pl,
"IC t0");