Thyra Version of the Day
Loading...
Searching...
No Matches
Thyra_ModelEvaluatorBase_def.hpp
1// @HEADER
2// ***********************************************************************
3//
4// Thyra: Interfaces and Support for Abstract Numerical Algorithms
5// Copyright (2004) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Roscoe A. Bartlett (bartlettra@ornl.gov)
38//
39// ***********************************************************************
40// @HEADER
41
42#ifndef THYRA_MODEL_EVALUATOR_BASE_DEF_HPP
43#define THYRA_MODEL_EVALUATOR_BASE_DEF_HPP
44
45
46#include "Thyra_ModelEvaluatorBase_decl.hpp"
47#include "Thyra_MultiVectorBase.hpp"
48#include "Thyra_VectorBase.hpp"
49#include "Thyra_MultiVectorStdOps.hpp"
50#include "Thyra_VectorStdOps.hpp"
51
52
53namespace Thyra {
54
55
56namespace ModelEvaluatorHelperPack {
57
58
59template<class Scalar>
60inline
61RCP<const Thyra::VectorBase<Scalar> >
62condCloneVec(
63 const RCP<const Thyra::VectorBase<Scalar> > &vec,
64 bool cloneObject
65 )
66{
67 if(cloneObject)
68 return vec->clone_v();
69 return vec;
70}
71
72template<class Scalar>
73inline
74RCP<const Thyra::MultiVectorBase<Scalar> >
75condCloneMultiVec(
76 const RCP<const Thyra::MultiVectorBase<Scalar> > &vec,
77 bool cloneObject
78 )
79{
80 if(cloneObject)
81 return vec->clone_mv();
82 return vec;
83}
84
85inline
86RCP<const Stokhos::ProductEpetraVector >
87condCloneVec_mp(
88 const RCP<const Stokhos::ProductEpetraVector > &vec,
89 bool cloneObject
90 )
91{
92 if(cloneObject)
93 {
94 printf("Warning: clone_v not implemented for ProductEpetraVector: %s %d\n",__FILE__,__LINE__);
95 //return vec->clone_v(); //JF clone_v not implemented for ProductEpetraVector
96 }
97 return vec;
98}
99
100} // namespace ModelEvaluatorHelperPack
101
102
103//
104// ModelEvaluatorBase::InArgs
105//
106
107
108template<class Scalar>
110 :modelEvalDescription_("WARNING! THIS INARGS OBJECT IS UNINITIALIZED!")
111{
114 std::fill_n(&supports_[0],NUM_E_IN_ARGS_MEMBERS,false);
115 t_ = SMT::zero();
116 alpha_ = ST::zero();
117 beta_ = ST::zero();
118 step_size_ = ST::zero();
119 stage_number_ = ST::one();
120}
121
122
123template<class Scalar>
125{ return p_.size(); }
126
127template<class Scalar>
129{ return g_multiplier_.size(); }
130
131template<class Scalar>
133{
135 int(arg)>=NUM_E_IN_ARGS_MEMBERS || int(arg) < 0,std::logic_error
136 ,"model = \'"<<modelEvalDescription_
137 <<"\': Error, arg="<<toString(arg)<<" is invalid!"
138 );
139 return supports_[arg];
140}
141
142template<class Scalar>
144{
145 assert_l(l);
146 return supports_p_mp_[l];
147}
148
149template<class Scalar>
151 const RCP<const VectorBase<Scalar> > &x_dot_dot
152 )
153{ assert_supports(IN_ARG_x_dot_dot); x_dot_dot_ = x_dot_dot; }
154
155template<class Scalar>
158{ assert_supports(IN_ARG_x_dot_dot); return x_dot_dot_; }
159
160template<class Scalar>
162 const RCP<const VectorBase<Scalar> > &x_dot
163 )
164{ assert_supports(IN_ARG_x_dot); x_dot_ = x_dot; }
165
166
167template<class Scalar>
170{ assert_supports(IN_ARG_x_dot); return x_dot_; }
171
173template<class Scalar>
175 const RCP<const VectorBase<Scalar> > &x
177{ assert_supports(IN_ARG_x); x_ = x; }
179
180template<class Scalar>
183{ assert_supports(IN_ARG_x); return x_; }
185
186template<class Scalar>
188 const RCP<const MultiVectorBase<Scalar> > &x_direction
190{ assert_supports(IN_ARG_x); x_direction_ = x_direction; }
192
193template<class Scalar>
195 int l, const RCP<const MultiVectorBase<Scalar> > &p_direction_l
197{ assert_l(l); p_direction_[l] = p_direction_l; }
199
200template<class Scalar>
203{ assert_supports(IN_ARG_x); return x_direction_; }
204
205
206template<class Scalar>
209{ assert_l(l); return p_direction_[l]; }
210
211
212template<class Scalar>
215 )
216{ assert_supports(IN_ARG_x_dot_mp); x_dot_mp_ = x_dot_mp; }
217
218
219template<class Scalar>
222{ assert_supports(IN_ARG_x_dot_mp); return x_dot_mp_; }
223
224
225template<class Scalar>
228 )
229{ assert_supports(IN_ARG_x_mp); x_mp_ = x_mp; }
230
231
232template<class Scalar>
235{ assert_supports(IN_ARG_x_mp); return x_mp_; }
236
237template<class Scalar>
239 const RCP<const VectorBase<Scalar> > &f_multiplier
241{ assert_supports(IN_ARG_x); f_multiplier_ = f_multiplier; }
242
243template<class Scalar>
246{ assert_supports(IN_ARG_x); return f_multiplier_; }
247
248template<class Scalar>
250 int j, const RCP<const VectorBase<Scalar> > &g_multiplier
252{
253 assert_j(j);
254 assert_supports(IN_ARG_x);
255 g_multiplier_[j] = g_multiplier;
256}
258template<class Scalar>
262 assert_j(j);
263 assert_supports(IN_ARG_x);
264 return g_multiplier_[j];
266
267#ifdef HAVE_THYRA_ME_POLYNOMIAL
268
269template<class Scalar>
271 const RCP<const Teuchos::Polynomial< VectorBase<Scalar> > > &x_dot_poly
272 )
273{ assert_supports(IN_ARG_x_dot_poly); x_dot_poly_ = x_dot_poly; }
275
276template<class Scalar>
279{ assert_supports(IN_ARG_x_dot_poly); return x_dot_poly_; }
280
282template<class Scalar>
284 const RCP<const Teuchos::Polynomial< VectorBase<Scalar> > > &x_poly
285 )
286{ assert_supports(IN_ARG_x_poly); x_poly_ = x_poly; }
287
288
289template<class Scalar>
292{ assert_supports(IN_ARG_x_poly); return x_poly_; }
294
295#endif // HAVE_THYRA_ME_POLYNOMIAL
296
297template<class Scalar>
299 int l, const RCP<const VectorBase<Scalar> > &p_l
300 )
301{ assert_l(l); p_[l] = p_l; }
302
303
304template<class Scalar>
307{ assert_l(l); return p_[l]; }
308
309
310template<class Scalar>
313 )
314{ assert_supports(IN_ARG_p_mp, l); p_mp_[l] = p_mp_l; }
315
316template<class Scalar>
319{ assert_supports(IN_ARG_p_mp, l); return p_mp_[l]; }
320
321
322template<class Scalar>
324{ assert_supports(IN_ARG_t); t_ = t; }
325
326
327template<class Scalar>
330{ assert_supports(IN_ARG_t); return t_; }
331
332
333template<class Scalar>
335{ assert_supports(IN_ARG_alpha); alpha_ = alpha; }
336
337
338template<class Scalar>
340{ assert_supports(IN_ARG_alpha); return alpha_; }
341
342
343template<class Scalar>
345{ assert_supports(IN_ARG_beta); beta_ = beta; }
346
347
348template<class Scalar>
350{ assert_supports(IN_ARG_beta); return beta_; }
351
352template<class Scalar>
354{ assert_supports(IN_ARG_W_x_dot_dot_coeff); W_x_dot_dot_coeff_ = W_x_dot_dot_coeff; }
355
356
357template<class Scalar>
359{ assert_supports(IN_ARG_W_x_dot_dot_coeff); return W_x_dot_dot_coeff_; }
360
361template<class Scalar>
363{ assert_supports(IN_ARG_step_size); step_size_ = step_size; }
364
365template<class Scalar>
367{ assert_supports(IN_ARG_step_size); return step_size_; }
368
369template<class Scalar>
371{ assert_supports(IN_ARG_stage_number); return stage_number_; }
372
373
374template<class Scalar>
376{ assert_supports(IN_ARG_stage_number); stage_number_ = stage_number; }
377
378
379template<class Scalar>
381 const InArgs<Scalar>& inArgs, bool ignoreUnsupported, bool cloneObjects
382 )
383{
384 using ModelEvaluatorHelperPack::condCloneVec;
385 using ModelEvaluatorHelperPack::condCloneMultiVec;
386 using ModelEvaluatorHelperPack::condCloneVec_mp;
387 if( inArgs.supports(IN_ARG_x_dot_dot) && nonnull(inArgs.get_x_dot_dot()) ) {
388 if(supports(IN_ARG_x_dot_dot) || !ignoreUnsupported)
389 set_x_dot_dot(condCloneVec(inArgs.get_x_dot_dot(),cloneObjects));
390 }
391 if( inArgs.supports(IN_ARG_x_dot) && nonnull(inArgs.get_x_dot()) ) {
392 if(supports(IN_ARG_x_dot) || !ignoreUnsupported)
393 set_x_dot(condCloneVec(inArgs.get_x_dot(),cloneObjects));
394 }
395 if( inArgs.supports(IN_ARG_x_dot_mp) && nonnull(inArgs.get_x_dot_mp()) ) {
396 if(supports(IN_ARG_x_dot_mp) || !ignoreUnsupported)
397 set_x_dot_mp(condCloneVec_mp(inArgs.get_x_dot_mp(),cloneObjects));
398 }
399 if( inArgs.supports(IN_ARG_x) && nonnull(inArgs.get_x()) ) {
400 if(supports(IN_ARG_x) || !ignoreUnsupported)
401 set_x(condCloneVec(inArgs.get_x(),cloneObjects));
402 }
403 if( inArgs.supports(IN_ARG_x_mp) && nonnull(inArgs.get_x_mp()) ) {
404 if(supports(IN_ARG_x_mp) || !ignoreUnsupported)
405 set_x_mp(condCloneVec_mp(inArgs.get_x_mp(),cloneObjects));
406 }
407 if( inArgs.supports(IN_ARG_x) && nonnull(inArgs.get_x_direction()) ) {
408 if(supports(IN_ARG_x) || !ignoreUnsupported)
409 set_x_direction(condCloneMultiVec(inArgs.get_x_direction(),cloneObjects));
410 }
411 if( inArgs.supports(IN_ARG_x) && nonnull(inArgs.get_f_multiplier()) ) {
412 if(supports(IN_ARG_x) || !ignoreUnsupported)
413 set_f_multiplier(condCloneVec(inArgs.get_f_multiplier(),cloneObjects));
414 }
415 const int min_Ng = TEUCHOS_MIN(this->Ng(),inArgs.Ng());
416 for (int j = 0; j < min_Ng; ++j) {
417 if (nonnull(inArgs.get_g_multiplier(j)))
418 set_g_multiplier(j,condCloneVec(inArgs.get_g_multiplier(j),cloneObjects));
419 }
420#ifdef HAVE_THYRA_ME_POLYNOMIAL
421 if( inArgs.supports(IN_ARG_x_dot_poly) && nonnull(inArgs.get_x_dot_poly()) ) {
422 if(supports(IN_ARG_x_dot_poly) || !ignoreUnsupported) {
424 cloneObjects && "Have not implemented cloning for x_dot_poly yet!" );
425 set_x_dot_poly(inArgs.get_x_dot_poly());
426 }
427 }
428 if( inArgs.supports(IN_ARG_x_poly) && nonnull(inArgs.get_x_poly()) ) {
429 if(supports(IN_ARG_x_poly) || !ignoreUnsupported) {
431 cloneObjects && "Have not implemented cloning for x_poly yet!" );
432 set_x_poly(inArgs.get_x_poly());
433 }
434 }
435#endif // HAVE_THYRA_ME_POLYNOMIAL
436 const int min_Np = TEUCHOS_MIN(this->Np(),inArgs.Np());
437 for (int l = 0; l < min_Np; ++l) {
438 if (nonnull(inArgs.get_p(l)))
439 set_p(l,condCloneVec(inArgs.get_p(l),cloneObjects));
440 }
441 for (int l = 0; l < min_Np; ++l) {
442 if (nonnull(inArgs.get_p_direction(l)))
443 set_p_direction(l,condCloneMultiVec(inArgs.get_p_direction(l),cloneObjects));
444 }
445 for (int l = 0; l < min_Np; ++l) {
446 if (inArgs.supports(IN_ARG_p_mp,l)) {
447 if (nonnull(inArgs.get_p_mp(l)))
448 set_p_mp(l,condCloneVec_mp(inArgs.get_p_mp(l),cloneObjects));
449 }
450 }
451 if (inArgs.supports(IN_ARG_t)) {
452 if(supports(IN_ARG_t) || !ignoreUnsupported)
453 set_t(inArgs.get_t());
454 }
455 if (inArgs.supports(IN_ARG_alpha)) {
456 if(supports(IN_ARG_alpha) || !ignoreUnsupported)
457 set_alpha(inArgs.get_alpha());
458 }
459 if (inArgs.supports(IN_ARG_beta)) {
460 if(supports(IN_ARG_beta) || !ignoreUnsupported)
461 set_beta(inArgs.get_beta());
462 }
464 if(supports(IN_ARG_W_x_dot_dot_coeff) || !ignoreUnsupported)
465 set_W_x_dot_dot_coeff(inArgs.get_W_x_dot_dot_coeff());
466 }
467 if (inArgs.supports(IN_ARG_step_size)) {
468 if(supports(IN_ARG_step_size) || !ignoreUnsupported)
469 set_step_size(inArgs.get_step_size());
470 }
471 if (inArgs.supports(IN_ARG_stage_number)) {
472 if(supports(IN_ARG_stage_number) || !ignoreUnsupported)
473 set_stage_number(inArgs.get_stage_number());
474 }
475 // Extended inArgs
476 if (extended_inargs_.size() > 0)
477 TEUCHOS_TEST_FOR_EXCEPTION(cloneObjects,
478 std::runtime_error,
479 "Extended InArgs does not support cloning!");
480 this->extended_inargs_ = inArgs.extended_inargs_;
481}
482
483
484template<class Scalar>
486 const InArgs<Scalar> &inArgs
487 ) const
488{
489 for ( int inArg_i = 0; inArg_i < NUM_E_IN_ARGS_MEMBERS; ++inArg_i ) {
490 const EInArgsMembers inArg_arg = static_cast<EInArgsMembers>(inArg_i);
491 const std::string inArg_name = toString(inArg_arg);
493 supports(inArg_arg) != inArgs.supports(inArg_arg), std::logic_error,
494 "Error, the input argument "<<inArg_name<<" with support "<<inArgs.supports(inArg_arg)<<"\n"
495 "in the InArgs object for the model:\n\n"
496 " "<<inArgs.modelEvalDescription()<<"\n\n"
497 "is not the same the argument "<<inArg_name<<" with support "<<supports(inArg_arg)<<"\n"
498 "in the InArgs object for the model:\n\n"
499 " "<<modelEvalDescription()<<"\n\n"
500 "and these two InArgs objects are not compatible!"
501 );
502 }
503 TEUCHOS_ASSERT_EQUALITY( this->Np(), inArgs.Np() );
504}
505
506
507template<class Scalar>
509{
510 return modelEvalDescription_;
511}
512
513
514template<class Scalar>
516{
518 std::ostringstream oss;
519 oss
520 << "Thyra::ModelEvaluatorBase::InArgs<"<<ST::name()<<">"
521 << "{"
522 << "model="<<modelEvalDescription_
523 << ",Np="<<Np()
524 << "}";
525 return oss.str();
526}
527
528
529template<class Scalar>
531 Teuchos::FancyOStream &out_arg, const Teuchos::EVerbosityLevel verbLevel
532 ) const
533{
534 using std::endl;
536 using Teuchos::OSTab;
537 using Teuchos::describe;
539 typedef RCP<const VectorBase<Scalar> > CV_ptr;
540
541 if(verbLevel == Teuchos::VERB_NONE)
542 return;
543
545 out = Teuchos::rcp(&out_arg,false);
546 const bool dump_x = includesVerbLevel(verbLevel,Teuchos::VERB_HIGH);
547 const Teuchos::EVerbosityLevel x_verbLevel =
548 dump_x?Teuchos::VERB_EXTREME:verbLevel;
549 const bool print_x_nrm = includesVerbLevel(verbLevel,Teuchos::VERB_LOW);
550 const bool dump_p = includesVerbLevel(verbLevel,Teuchos::VERB_MEDIUM);
551 const Teuchos::EVerbosityLevel p_verbLevel =
552 dump_p?Teuchos::VERB_EXTREME:verbLevel;
553 const bool print_p_nrm = includesVerbLevel(verbLevel,Teuchos::VERB_LOW);
554 OSTab tab(out);
555
556 *out <<"Thyra::ModelEvaluatorBase::InArgs<"<<ST::name()<<">:\n";
557 tab.incrTab();
558
559 *out <<"model = " << modelEvalDescription_ << "\n";
560 *out <<"Np = " << Np() << "\n";
561
562 CV_ptr x_dot_dot;
563 if ( this->supports(IN_ARG_x_dot_dot) && !is_null(x_dot_dot=get_x_dot_dot()) ) {
564 *out << "x_dot_dot = " << Teuchos::describe(*x_dot_dot,x_verbLevel);
565 if (print_x_nrm)
566 *out << "||x_dot_dot|| = " << norm(*x_dot_dot) << endl;
567 }
568
569 CV_ptr x_dot;
570 if ( this->supports(IN_ARG_x_dot) && !is_null(x_dot=get_x_dot()) ) {
571 *out << "x_dot = " << Teuchos::describe(*x_dot,x_verbLevel);
572 if (print_x_nrm)
573 *out << "||x_dot|| = " << norm(*x_dot) << endl;
574 }
575
576 CV_ptr x;
577 if ( this->supports(IN_ARG_x) && !is_null(x=get_x()) ) {
578 *out << "x = " << Teuchos::describe(*x,x_verbLevel);
579 if (print_x_nrm)
580 *out << "||x|| = " << norm(*x) << endl;
581 }
582
583 if (print_x_nrm) {
584 for( int l = 0; l < Np(); ++l ) {
585 CV_ptr p_l;
586 if ( !is_null(p_l = this->get_p(l)) ) {
587 *out << "p("<<l<<") = " << Teuchos::describe(*p_l,p_verbLevel);
588 if (print_p_nrm)
589 *out << "||p("<<l<<")|| = " << norm(*p_l) << endl;
590 }
591 }
592 }
593
594 if (includesVerbLevel(verbLevel,Teuchos::VERB_MEDIUM)) {
595 if (this->supports(IN_ARG_t)) {
596 *out << "t = " << t_ << endl;
597 }
598 if (this->supports(IN_ARG_alpha)) {
599 *out << "alpha = " << alpha_ << endl;
600 }
601 if (this->supports(IN_ARG_beta)) {
602 *out << "beta = " << beta_ << endl;
603 }
604 if (this->supports(IN_ARG_W_x_dot_dot_coeff)) {
605 *out << "W_x_dot_dot_coeff = " << W_x_dot_dot_coeff_ << endl;
606 }
607 if (this->supports(IN_ARG_step_size)) {
608 *out << "step_size = " << step_size_ << endl;
609 }
610 if (this->supports(IN_ARG_stage_number)) {
611 *out << "stage_number = " << stage_number_ << endl;
612 }
613 }
614
615}
616
617
618template<class Scalar>
620 const std::string &modelEvalDescription_in
621 )
622{
623 modelEvalDescription_ = modelEvalDescription_in;
624}
625
626
627template<class Scalar>
629{
630 p_.resize(Np_in);
631 p_direction_.resize(Np_in);
632 p_mp_.resize(Np_in);
633 supports_p_mp_.resize(Np_in);
634 g_multiplier_.resize(Ng_in);
635}
636
637
638template<class Scalar>
640 EInArgsMembers arg, bool supports_in
641 )
642{
644 int(arg)>=NUM_E_IN_ARGS_MEMBERS || int(arg) < 0,std::logic_error
645 ,"model = \'"<<modelEvalDescription_
646 <<"\': Error, arg="<<toString(arg)<<" is invalid!");
647 supports_[arg] = supports_in;
648}
649
650template<class Scalar>
652 EInArgs_p_mp /* arg */, int l, bool supports_in
653 )
654{
655 assert_l(l);
656 supports_p_mp_[l] = supports_in;
657}
658
659
660template<class Scalar>
662 const InArgs<Scalar>& inArgs, const int Np_in
663 )
664{
665 std::copy(
666 &inArgs.supports_[0],
667 &inArgs.supports_[0] + NUM_E_IN_ARGS_MEMBERS, &supports_[0] );
668 this->_set_Np_Ng( Np_in >= 0 ? Np_in : inArgs.Np(), inArgs.Ng() );
669}
670
671
672template<class Scalar>
675 )
676{
677 switch(arg) {
678 case IN_ARG_x: {
679 this->_setSupports(IN_ARG_x_dot_dot,false);
680 this->_setSupports(IN_ARG_x_dot,false);
681 this->_setSupports(IN_ARG_x_dot_poly,false);
682 this->_setSupports(IN_ARG_alpha,false);
683 this->_setSupports(IN_ARG_beta,false);
684 this->_setSupports(IN_ARG_W_x_dot_dot_coeff,false);
685 this->_setSupports(IN_ARG_step_size,false);
686 this->_setSupports(IN_ARG_stage_number,false);
687 break;
688 }
689 default:
691 true ,std::logic_error,
692 "Error, can not handle args other than IN_ARG_x yet!"
693 );
694 }
695 this->_setSupports(arg,false);
696}
697
698
699template<class Scalar>
702 ) const
703{
705 !supports_[arg], std::logic_error
706 ,"Thyra::ModelEvaluatorBase::InArgs<"
707 << Teuchos::ScalarTraits<Scalar>::name() <<">::assert_supports(arg): "
708 "model = \'"<<modelEvalDescription_<<"\': Error, "
709 "The argument arg = " << toString(arg) << " is not supported!"
710 );
711}
712
713template<class Scalar>
715 EInArgs_p_mp /* arg */, int l
716 ) const
717{
718 assert_l(l);
720 !supports_p_mp_[l], std::logic_error
721 ,"Thyra::ModelEvaluatorBase::InArgs<"
722 << Teuchos::ScalarTraits<Scalar>::name() <<">::assert_supports(IN_ARG_p_mp,1): "
723 "model = \'"<<modelEvalDescription_<<"\': Error, "
724 "The argument p_mp(l) with index l = " << l << " is not supported!"
725 );
726}
727
728
729template<class Scalar>
730void ModelEvaluatorBase::InArgs<Scalar>::assert_l(int l) const
731{
733 !( 0 <= l && l < Np() ), std::logic_error
734 ,"Thyra::ModelEvaluatorBase::InArgs<Scalar>::assert_l(l):\n\n"
735 " model = \'"<<modelEvalDescription_<<"\':\n\n"
736 "Error, The parameter l = " << l << " is not in the range [0,"<<Np()<<")!"
737 );
738}
739
740template<class Scalar>
741void ModelEvaluatorBase::InArgs<Scalar>::assert_j(int j) const
742{
744 !( 0 <= j && j < Ng() ), std::logic_error
745 ,"Thyra::ModelEvaluatorBase::InArgs<Scalar>::assert_j(j):\n\n"
746 "model = \'"<<modelEvalDescription_<<"\':\n\n"
747 "Error, The auxiliary function g("<<j<<")"
748 " is not in the range [0,"<<Ng()<<")!"
749 );
750}
751
752//
753// ModelEvaluatorBase::DerivativeMultiVector
754//
755
756
757template<class Scalar>
759{
760 using std::endl;
761 std::ostringstream oss;
762 oss << "DerivativeMultiVector{";
763 if (is_null(getMultiVector())) {
764 oss << "NULL";
765 }
766 else {
767 oss
768 << "multiVec=" << getMultiVector()->description()
769 << ",orientation=" << toString(getOrientation());
770 }
771 oss << "}";
772 return oss.str();
773}
774
775
776template<class Scalar>
779 ) const
780{
781 using std::endl;
782 using Teuchos::describe;
783 Teuchos::OSTab tab1(out);
784 out << "DerivativeMultiVector\n";
785 Teuchos::OSTab tab2(out);
786 out
787 << "multiVec = "
788 << describe(*getMultiVector(),verbLevel)
789 << "orientation = "
790 << toString(getOrientation()) << endl;
791}
792
793
794// 2007/06/12: rabartl: The above description() and describe(...) functions
795// have to be defined here and not in the class DerivativeMultiVector since it
796// relies on the non-member function
797// toString(ModelEvaluatorBase::EDerivativeMultiVectorOrientation) which is
798// defined after the class definition for ModelEvaluatorBase. This was caught
799// by the intel compiler. I am not sure why this worked with gcc.
800
801
802//
803// ModelEvaluatorBase::Derivative
804//
805
806
807template<class Scalar>
808std::string
810{
811 using std::endl;
812 std::ostringstream oss;
813 oss << "Derivative{";
814 if (isEmpty()) {
815 oss << "NULL";
816 }
817 else if (!is_null(getLinearOp())) {
818 oss << "linearOp=" << getLinearOp()->description();
819 }
820 else {
821 oss << "derivMultiVec=" << getDerivativeMultiVector().description();
822 }
823 oss << "}";
824 return oss.str();
825}
826
827
828template<class Scalar>
831 ) const
832{
833 using std::endl;
834 using Teuchos::describe;
835 Teuchos::OSTab tab1(out);
836 out << "Derivative:";
837 if (isEmpty()) {
838 out << " NULL\n";
839 }
840 else if (!is_null(getLinearOp())) {
841 out
842 << endl
843 << "linearOp = " << describe(*getLinearOp(),verbLevel);
844 }
845 else {
846 out
847 << endl
848 << "derivMultiVec = ";
849 getDerivativeMultiVector().describe(out,verbLevel);
850 }
851}
852
853
854//
855// ModelEvaluatorBase::OutArgs
856//
857
858
859template<class Scalar>
861 :modelEvalDescription_("WARNING! THIS OUTARGS OBJECT IS UNINITIALIZED!"),
862 isFailed_(false)
863{
864 std::fill_n(&supports_[0],NUM_E_OUT_ARGS_MEMBERS,false);
865#ifdef Thyra_BUILD_HESSIAN_SUPPORT
866 this->_setHessianSupports(false);
867#endif
868}
869
870
871template<class Scalar>
873{ return DfDp_.size(); }
874
875
876template<class Scalar>
878{ return g_.size(); }
879
880
881template<class Scalar>
884 ) const
885{
887 int(arg)>=NUM_E_OUT_ARGS_MEMBERS || int(arg) < 0,std::logic_error
888 ,"model = \'"<<modelEvalDescription_
889 <<"\': Error, arg="<<toString(arg)<<" is invalid!"
890 );
891 return supports_[arg];
892}
893
894
895template<class Scalar>
898 EOutArgsDfDp /* arg */, int l
899 ) const
900{
901 assert_l(l);
902 return supports_DfDp_[l];
903}
904
906template<class Scalar>
909 EOutArgsDgDx_dot /* arg */, int j
910 ) const
912 assert_j(j);
913 return supports_DgDx_dot_[j];
914}
915
917template<class Scalar>
920 EOutArgsDgDx /* arg */, int j
921 ) const
923 assert_j(j);
924 return supports_DgDx_[j];
926
927
928template<class Scalar>
931 EOutArgsDgDp /* arg */, int j, int l
932 ) const
933{
934 assert_j(j);
935 assert_l(l);
936 return supports_DgDp_[ j*Np() + l ];
937}
938
939
940#ifdef Thyra_BUILD_HESSIAN_SUPPORT
941
942template<class Scalar>
944 EOutArgs_hess_vec_prod_f_xx /* arg */
945 ) const
946{
947 return supports_hess_vec_prod_f_xx_;
948}
949
950template<class Scalar>
952 EOutArgs_hess_vec_prod_f_xp /* arg */, int l
953 ) const
954{
955 assert_l(l);
956 return supports_hess_vec_prod_f_xp_[l];
957}
958
959template<class Scalar>
961 EOutArgs_hess_vec_prod_f_px /* arg */, int l
962 ) const
963{
964 assert_l(l);
965 return supports_hess_vec_prod_f_px_[l];
966}
967
968template<class Scalar>
970 EOutArgs_hess_vec_prod_f_pp /* arg */, int l1, int l2
971 ) const
972{
973 assert_l(l1);
974 assert_l(l2);
975 return supports_hess_vec_prod_f_pp_[ l1*Np() + l2 ];
976}
977
978template<class Scalar>
980 EOutArgs_hess_vec_prod_g_xx /* arg */, int j
981 ) const
982{
983 assert_j(j);
984 return supports_hess_vec_prod_g_xx_[j];
986
987template<class Scalar>
989 EOutArgs_hess_vec_prod_g_xp /* arg */, int j, int l
990 ) const
991{
992 assert_j(j);
993 assert_l(l);
994 return supports_hess_vec_prod_g_xp_[ j*Np() + l ];
995}
996
997template<class Scalar>
999 EOutArgs_hess_vec_prod_g_px /* arg */, int j, int l
1000 ) const
1001{
1002 assert_j(j);
1003 assert_l(l);
1004 return supports_hess_vec_prod_g_px_[ j*Np() + l ];
1005}
1006
1007template<class Scalar>
1009 EOutArgs_hess_vec_prod_g_pp /* arg */, int j, int l1, int l2
1010 ) const
1012 assert_j(j);
1013 assert_l(l1);
1014 assert_l(l2);
1015 return supports_hess_vec_prod_g_pp_[ j*Np()*Np() + l1*Np() + l2 ];
1016}
1018template<class Scalar>
1020 EOutArgs_hess_f_xx /* arg */
1021 ) const
1023 return supports_hess_f_xx_;
1025
1026template<class Scalar>
1028 EOutArgs_hess_f_xp /* arg */, int l
1029 ) const
1030{
1031 assert_l(l);
1032 return supports_hess_f_xp_[l];
1034
1035template<class Scalar>
1037 EOutArgs_hess_f_pp /* arg */, int l1, int l2
1038 ) const
1039{
1040 assert_l(l1);
1041 assert_l(l2);
1042 return supports_hess_f_pp_[ l1*Np() + l2 ];
1044
1045template<class Scalar>
1047 EOutArgs_H_xx /* arg */
1048 ) const
1049{
1050 return supports_H_xx_;
1051}
1053template<class Scalar>
1055 EOutArgs_H_xp /* arg */, int l
1056 ) const
1058 assert_l(l);
1059 return supports_H_xp_[l];
1060}
1061
1062template<class Scalar>
1064 EOutArgs_H_pp /* arg */, int l1, int l2
1065 ) const
1066{
1067 assert_l(l1);
1068 assert_l(l2);
1069 return supports_H_pp_[ l1*Np() + l2 ];
1070}
1071
1072template<class Scalar>
1074 EOutArgs_hess_g_xx /* arg */, int j
1075 ) const
1076{
1077 assert_j(j);
1078 return supports_hess_g_xx_[j];
1079}
1080
1081template<class Scalar>
1083 EOutArgs_hess_g_xp /* arg */, int j, int l
1084 ) const
1085{
1086 assert_j(j);
1087 assert_l(l);
1088 return supports_hess_g_xp_[ j*Np() + l ];
1089}
1090
1091template<class Scalar>
1093 EOutArgs_hess_g_pp /* arg */, int j, int l1, int l2
1094 ) const
1095{
1096 assert_j(j);
1097 assert_l(l1);
1098 assert_l(l2);
1099 return supports_hess_g_pp_[ j*Np()*Np() + l1*Np() + l2 ];
1100}
1101
1102#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
1103
1104
1105template<class Scalar>
1106bool
1108 EOutArgs_g_mp /* arg */, int j
1109 ) const
1110{
1111 assert_j(j);
1112 return supports_g_mp_[j];
1113}
1114
1115template<class Scalar>
1116const ModelEvaluatorBase::DerivativeSupport&
1118 EOutArgsDfDp_mp /* arg */, int l
1119 ) const
1120{
1121 assert_l(l);
1122 return supports_DfDp_mp_[l];
1123}
1124
1125
1126template<class Scalar>
1127const ModelEvaluatorBase::DerivativeSupport&
1129 EOutArgsDgDx_dot_mp /* arg */, int j
1130 ) const
1131{
1132 assert_j(j);
1133 return supports_DgDx_dot_mp_[j];
1134}
1135
1136
1137template<class Scalar>
1138const ModelEvaluatorBase::DerivativeSupport&
1140 EOutArgsDgDx_mp /* arg */, int j
1141 ) const
1142{
1143 assert_j(j);
1144 return supports_DgDx_mp_[j];
1145}
1146
1147
1148template<class Scalar>
1149const ModelEvaluatorBase::DerivativeSupport&
1151 EOutArgsDgDp_mp /* arg */, int j, int l
1152 ) const
1153{
1154 assert_j(j);
1155 assert_l(l);
1156 return supports_DgDp_mp_[ j*Np() + l ];
1157}
1158
1159
1160
1161template<class Scalar>
1164 )
1165{
1166 assert_supports(OUT_ARG_f);
1167 f_ = f;
1168}
1169
1170
1171template<class Scalar>
1174{
1175 assert_supports(OUT_ARG_f);
1176 return f_;
1177}
1178
1180template<class Scalar>
1182 int j, const Evaluation<VectorBase<Scalar> > &g_j
1184{
1185 assert_j(j);
1186 g_[j] = g_j;
1188
1189
1190template<class Scalar>
1193{
1194 assert_j(j);
1195 return g_[j];
1196}
1198
1199template<class Scalar>
1202 )
1203{
1204 assert_supports(OUT_ARG_f_mp);
1205 f_mp_ = f_mp;
1207
1209template<class Scalar>
1212{
1213 assert_supports(OUT_ARG_f_mp);
1214 return f_mp_;
1215}
1216
1217
1218template<class Scalar>
1220 int j, const RCP<Stokhos::ProductEpetraVector> &g_mp_j
1221 )
1222{
1223 assert_supports(OUT_ARG_g_mp,j);
1224 g_mp_[j] = g_mp_j;
1225}
1226
1227
1228template<class Scalar>
1231{
1232 assert_supports(OUT_ARG_g_mp,j);
1233 return g_mp_[j];
1234}
1235
1236
1237template<class Scalar>
1240 )
1241{
1242 assert_supports(OUT_ARG_W);
1243 W_ = W;
1244}
1245
1246
1247template<class Scalar>
1250{
1251 assert_supports(OUT_ARG_W);
1252 return W_;
1253}
1254
1255
1256template<class Scalar>
1260{
1261 assert_supports(OUT_ARG_W_mp);
1262 W_mp_ = W_mp;
1264
1266template<class Scalar>
1269{
1270 assert_supports(OUT_ARG_W_mp);
1271 return W_mp_;
1272}
1273
1274
1275template<class Scalar>
1278 )
1280 assert_supports(OUT_ARG_W_op);
1281 W_op_ = W_op;
1282}
1283
1284
1285template<class Scalar>
1288{
1289 assert_supports(OUT_ARG_W_op);
1290 return W_op_;
1291}
1292
1293
1294template<class Scalar>
1296 const RCP<PreconditionerBase<Scalar> > &W_prec
1297 )
1298{
1299 assert_supports(OUT_ARG_W_prec);
1300 W_prec_ = W_prec;
1301}
1302
1303
1304template<class Scalar>
1307{
1308 assert_supports(OUT_ARG_W_prec);
1309 return W_prec_;
1310}
1311
1312
1313template<class Scalar>
1316{
1317 assert_supports(OUT_ARG_f);
1318 return W_properties_;
1319}
1320
1321
1322template<class Scalar>
1324 int l, const Derivative<Scalar> &DfDp_l
1325 )
1326{
1327 assert_supports(OUT_ARG_DfDp,l,DfDp_l);
1328 DfDp_[l] = DfDp_l;
1329}
1330
1331
1332template<class Scalar>
1335{
1336 assert_supports(OUT_ARG_DfDp,l);
1337 return DfDp_[l];
1338}
1339
1340
1341template<class Scalar>
1344{
1345 assert_supports(OUT_ARG_DfDp,l);
1346 return DfDp_properties_[l];
1347}
1348
1349
1350template<class Scalar>
1352 int l, const MPDerivative &DfDp_mp_l
1353 )
1354{
1355 assert_supports(OUT_ARG_DfDp_mp,l,DfDp_mp_l);
1356 DfDp_mp_[l] = DfDp_mp_l;
1357}
1358
1359
1360template<class Scalar>
1363{
1364 assert_supports(OUT_ARG_DfDp_mp,l);
1365 return DfDp_mp_[l];
1366}
1367
1368
1369template<class Scalar>
1370ModelEvaluatorBase::DerivativeProperties
1371ModelEvaluatorBase::OutArgs<Scalar>::get_DfDp_mp_properties(int l) const
1372{
1373 assert_supports(OUT_ARG_DfDp_mp,l);
1374 return DfDp_mp_properties_[l];
1375}
1376
1377
1378template<class Scalar>
1380 int j, const Derivative<Scalar> &DgDx_dot_j
1381 )
1382{
1383 assert_supports(OUT_ARG_DgDx_dot,j,DgDx_dot_j);
1384 DgDx_dot_[j] = DgDx_dot_j;
1385}
1386
1387
1388template<class Scalar>
1391{
1392 assert_supports(OUT_ARG_DgDx_dot,j);
1393 return DgDx_dot_[j];
1394}
1395
1396
1397template<class Scalar>
1400{
1401 assert_supports(OUT_ARG_DgDx_dot,j);
1402 return DgDx_dot_properties_[j];
1403}
1404
1405
1406template<class Scalar>
1408 int j, const MPDerivative &DgDx_dot_mp_j
1409 )
1410{
1411 assert_supports(OUT_ARG_DgDx_dot_mp,j,DgDx_dot_mp_j);
1412 DgDx_dot_mp_[j] = DgDx_dot_mp_j;
1413}
1414
1415
1416template<class Scalar>
1419{
1420 assert_supports(OUT_ARG_DgDx_dot_mp,j);
1421 return DgDx_dot_mp_[j];
1422}
1423
1424
1425template<class Scalar>
1426ModelEvaluatorBase::DerivativeProperties
1427ModelEvaluatorBase::OutArgs<Scalar>::get_DgDx_dot_mp_properties(int j) const
1428{
1429 assert_supports(OUT_ARG_DgDx_dot_mp,j);
1430 return DgDx_dot_mp_properties_[j];
1431}
1432
1433
1434template<class Scalar>
1436 int j, const Derivative<Scalar> &DgDx_j
1437 )
1438{
1439 assert_supports(OUT_ARG_DgDx,j,DgDx_j);
1440 DgDx_[j] = DgDx_j;
1441}
1442
1443
1444template<class Scalar>
1447{
1448 assert_supports(OUT_ARG_DgDx,j);
1449 return DgDx_[j];
1450}
1451
1452
1453template<class Scalar>
1456{
1457 assert_supports(OUT_ARG_DgDx,j);
1458 return DgDx_properties_[j];
1459}
1460
1461
1462template<class Scalar>
1464 int j, const MPDerivative &DgDx_mp_j
1465 )
1466{
1467 assert_supports(OUT_ARG_DgDx_mp,j,DgDx_mp_j);
1468 DgDx_mp_[j] = DgDx_mp_j;
1469}
1470
1471
1472template<class Scalar>
1475{
1476 assert_supports(OUT_ARG_DgDx_mp,j);
1477 return DgDx_mp_[j];
1478}
1479
1480
1481template<class Scalar>
1482ModelEvaluatorBase::DerivativeProperties
1483ModelEvaluatorBase::OutArgs<Scalar>::get_DgDx_mp_properties(int j) const
1484{
1485 assert_supports(OUT_ARG_DgDx_mp,j);
1486 return DgDx_mp_properties_[j];
1487}
1488
1489
1490template<class Scalar>
1492 int j, int l, const Derivative<Scalar> &DgDp_j_l
1493 )
1494{
1495 assert_supports(OUT_ARG_DgDp,j,l,DgDp_j_l);
1496 DgDp_[ j*Np() + l ] = DgDp_j_l;
1497}
1498
1499
1500template<class Scalar>
1503{
1504 assert_supports(OUT_ARG_DgDp,j,l);
1505 return DgDp_[ j*Np() + l ];
1506}
1507
1508
1509template<class Scalar>
1512{
1513 assert_supports(OUT_ARG_DgDp,j,l);
1514 return DgDp_properties_[ j*Np() + l ];
1515}
1516
1517
1518template<class Scalar>
1520 int j, int l, const MPDerivative &DgDp_mp_j_l
1521 )
1522{
1523 assert_supports(OUT_ARG_DgDp_mp,j,l,DgDp_mp_j_l);
1524 DgDp_mp_[ j*Np() + l ] = DgDp_mp_j_l;
1525}
1526
1527
1528template<class Scalar>
1531{
1532 assert_supports(OUT_ARG_DgDp_mp,j,l);
1533 return DgDp_mp_[ j*Np() + l ];
1534}
1535
1536
1537template<class Scalar>
1538ModelEvaluatorBase::DerivativeProperties
1539ModelEvaluatorBase::OutArgs<Scalar>::get_DgDp_mp_properties(int j, int l) const
1540{
1541 assert_supports(OUT_ARG_DgDp_mp,j,l);
1542 return DgDp_mp_properties_[ j*Np() + l ];
1543}
1544
1545
1546#ifdef Thyra_BUILD_HESSIAN_SUPPORT
1547
1548template<class Scalar>
1549void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_f_xx(
1550 const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_f_xx
1551 )
1552{
1553 assert_supports(OUT_ARG_hess_vec_prod_f_xx);
1554 hess_vec_prod_f_xx_ = hess_vec_prod_f_xx;
1555}
1556
1557template<class Scalar>
1558void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_f_xp(
1559 int l, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_f_xp_l
1560 )
1561{
1562 assert_supports(OUT_ARG_hess_vec_prod_f_xp,l);
1563 hess_vec_prod_f_xp_[l] = hess_vec_prod_f_xp_l;
1564}
1565
1566template<class Scalar>
1567void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_f_px(
1568 int l, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_f_px_l
1569 )
1570{
1571 assert_supports(OUT_ARG_hess_vec_prod_f_px,l);
1572 hess_vec_prod_f_px_[l] = hess_vec_prod_f_px_l;
1573}
1574
1575template<class Scalar>
1576void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_f_pp(
1577 int l1, int l2, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_f_pp_l1_l2
1578 )
1579{
1580 assert_supports(OUT_ARG_hess_vec_prod_f_pp,l1,l2);
1581 hess_vec_prod_f_pp_[ l1*Np() + l2 ] = hess_vec_prod_f_pp_l1_l2;
1582}
1583
1584template<class Scalar>
1585RCP<MultiVectorBase<Scalar> >
1586ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_f_xx() const
1587{
1588 assert_supports(OUT_ARG_hess_vec_prod_f_xx);
1589 return hess_vec_prod_f_xx_;
1590}
1591
1592template<class Scalar>
1593RCP<MultiVectorBase<Scalar> >
1594ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_f_xp(int l) const
1595{
1596 assert_supports(OUT_ARG_hess_vec_prod_f_xp,l);
1597 return hess_vec_prod_f_xp_[l];
1598}
1599
1600template<class Scalar>
1601RCP<MultiVectorBase<Scalar> >
1602ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_f_px(int l) const
1603{
1604 assert_supports(OUT_ARG_hess_vec_prod_f_px,l);
1605 return hess_vec_prod_f_px_[l];
1606}
1607
1608template<class Scalar>
1609RCP<MultiVectorBase<Scalar> >
1610ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_f_pp(int l1, int l2) const
1611{
1612 assert_supports(OUT_ARG_hess_vec_prod_f_pp,l1,l2);
1613 return hess_vec_prod_f_pp_[ l1*Np() + l2 ];
1614}
1615
1616
1617template<class Scalar>
1618void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_g_xx(
1619 int j, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_g_xx_j
1620 )
1621{
1622 assert_supports(OUT_ARG_hess_vec_prod_g_xx,j);
1623 hess_vec_prod_g_xx_[j] = hess_vec_prod_g_xx_j;
1624}
1625
1626template<class Scalar>
1627void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_g_xp(
1628 int j, int l, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_g_xp_j_l
1629 )
1630{
1631 assert_supports(OUT_ARG_hess_vec_prod_g_xp,j,l);
1632 hess_vec_prod_g_xp_[ j*Np() + l ] = hess_vec_prod_g_xp_j_l;
1633}
1634
1635template<class Scalar>
1636void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_g_px(
1637 int j, int l, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_g_px_j_l
1638 )
1639{
1640 assert_supports(OUT_ARG_hess_vec_prod_g_px,j,l);
1641 hess_vec_prod_g_px_[ j*Np() + l ] = hess_vec_prod_g_px_j_l;
1642}
1643
1644template<class Scalar>
1645void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_vec_prod_g_pp(
1646 int j, int l1, int l2, const RCP<MultiVectorBase<Scalar> > &hess_vec_prod_g_pp_j_l1_l2
1647 )
1648{
1649 assert_supports(OUT_ARG_hess_vec_prod_g_pp,j,l1,l2);
1650 hess_vec_prod_g_pp_[ j*Np()*Np() + l1*Np() + l2 ] = hess_vec_prod_g_pp_j_l1_l2;
1651}
1652
1653template<class Scalar>
1654RCP<MultiVectorBase<Scalar> >
1655ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_g_xx(int j) const
1656{
1657 assert_supports(OUT_ARG_hess_vec_prod_g_xx,j);
1658 return hess_vec_prod_g_xx_[j];
1659}
1660
1661template<class Scalar>
1662RCP<MultiVectorBase<Scalar> >
1663ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_g_xp(int j, int l) const
1664{
1665 assert_supports(OUT_ARG_hess_vec_prod_g_xp,j,l);
1666 return hess_vec_prod_g_xp_[ j*Np() + l ];
1667}
1668
1669template<class Scalar>
1670RCP<MultiVectorBase<Scalar> >
1671ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_g_px(int j, int l) const
1672{
1673 assert_supports(OUT_ARG_hess_vec_prod_g_px,j,l);
1674 return hess_vec_prod_g_px_[ j*Np() + l ];
1675}
1676
1677template<class Scalar>
1678RCP<MultiVectorBase<Scalar> >
1679ModelEvaluatorBase::OutArgs<Scalar>::get_hess_vec_prod_g_pp(int j, int l1, int l2) const
1680{
1681 assert_supports(OUT_ARG_hess_vec_prod_g_pp,j,l1,l2);
1682 return hess_vec_prod_g_pp_[ j*Np()*Np() + l1*Np() + l2 ];
1683}
1684
1685template<class Scalar>
1686void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_f_xx(
1687 const RCP<LinearOpBase<Scalar> > &hess_f_xx
1688 )
1689{
1690 assert_supports(OUT_ARG_hess_f_xx);
1691 hess_f_xx_ = hess_f_xx;
1692}
1693
1694template<class Scalar>
1695void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_f_xp(
1696 int l, const RCP<LinearOpBase<Scalar> > &hess_f_xp_l
1697 )
1698{
1699 assert_supports(OUT_ARG_hess_f_xp,l);
1700 hess_f_xp_[l] = hess_f_xp_l;
1701}
1702
1703template<class Scalar>
1704void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_f_pp(
1705 int l1, int l2, const RCP<LinearOpBase<Scalar> > &hess_f_pp_l1_l2
1706 )
1707{
1708 assert_supports(OUT_ARG_hess_f_pp,l1,l2);
1709 hess_f_pp_[ l1*Np() + l2 ] = hess_f_pp_l1_l2;
1710}
1711
1712template<class Scalar>
1713RCP<LinearOpBase<Scalar> >
1714ModelEvaluatorBase::OutArgs<Scalar>::get_hess_f_xx() const
1715{
1716 assert_supports(OUT_ARG_hess_f_xx);
1717 return hess_f_xx_;
1718}
1719
1720template<class Scalar>
1721RCP<LinearOpBase<Scalar> >
1722ModelEvaluatorBase::OutArgs<Scalar>::get_hess_f_xp(int l) const
1723{
1724 assert_supports(OUT_ARG_hess_f_xp,l);
1725 return hess_f_xp_[l];
1726}
1727
1728template<class Scalar>
1729RCP<LinearOpBase<Scalar> >
1730ModelEvaluatorBase::OutArgs<Scalar>::get_hess_f_pp(int l1, int l2) const
1731{
1732 assert_supports(OUT_ARG_hess_f_pp,l1,l2);
1733 return hess_f_pp_[ l1*Np() + l2 ];
1734}
1735
1736template<class Scalar>
1737void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_g_xx(
1738 int j, const RCP<LinearOpBase<Scalar> > &hess_g_xx_j
1739 )
1740{
1741 assert_supports(OUT_ARG_hess_g_xx,j);
1742 hess_g_xx_[j] = hess_g_xx_j;
1743}
1744
1745template<class Scalar>
1746void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_g_xp(
1747 int j, int l, const RCP<LinearOpBase<Scalar> > &hess_g_xp_j_l
1748 )
1749{
1750 assert_supports(OUT_ARG_hess_g_xp,j,l);
1751 hess_g_xp_[ j*Np() + l ] = hess_g_xp_j_l;
1752}
1753
1754template<class Scalar>
1755void ModelEvaluatorBase::OutArgs<Scalar>::set_hess_g_pp(
1756 int j, int l1, int l2, const RCP<LinearOpBase<Scalar> > &hess_g_pp_j_l1_l2
1757 )
1758{
1759 assert_supports(OUT_ARG_hess_g_pp,j,l1,l2);
1760 hess_g_pp_[ j*Np()*Np() + l1*Np() + l2 ] = hess_g_pp_j_l1_l2;
1761}
1762
1763template<class Scalar>
1764RCP<LinearOpBase<Scalar> >
1765ModelEvaluatorBase::OutArgs<Scalar>::get_hess_g_xx(int j) const
1766{
1767 assert_supports(OUT_ARG_hess_g_xx,j);
1768 return hess_g_xx_[j];
1769}
1770
1771template<class Scalar>
1772RCP<LinearOpBase<Scalar> >
1773ModelEvaluatorBase::OutArgs<Scalar>::get_hess_g_xp(int j, int l) const
1774{
1775 assert_supports(OUT_ARG_hess_g_xp,j,l);
1776 return hess_g_xp_[ j*Np() + l ];
1777}
1778
1779template<class Scalar>
1780RCP<LinearOpBase<Scalar> >
1781ModelEvaluatorBase::OutArgs<Scalar>::get_hess_g_pp(int j, int l1, int l2) const
1782{
1783 assert_supports(OUT_ARG_hess_g_pp,j,l1,l2);
1784 return hess_g_pp_[ j*Np()*Np() + l1*Np() + l2 ];
1785}
1786
1787template<class Scalar>
1788void ModelEvaluatorBase::OutArgs<Scalar>::set_H_xx(
1789 const RCP<LinearOpBase<Scalar> > &H_xx
1790 )
1791{
1792 assert_supports(OUT_ARG_H_xx);
1793 H_xx_ = H_xx;
1794}
1795
1796template<class Scalar>
1797void ModelEvaluatorBase::OutArgs<Scalar>::set_H_xp(
1798 int l, const RCP<LinearOpBase<Scalar> > &H_xp_l
1799 )
1800{
1801 assert_supports(OUT_ARG_H_xp,l);
1802 H_xp_[l] = H_xp_l;
1803}
1804
1805template<class Scalar>
1806void ModelEvaluatorBase::OutArgs<Scalar>::set_H_pp(
1807 int l1, int l2, const RCP<LinearOpBase<Scalar> > &H_pp_l1_l2
1808 )
1809{
1810 assert_supports(OUT_ARG_H_pp,l1,l2);
1811 H_pp_[ l1*Np() + l2 ] = H_pp_l1_l2;
1812}
1813
1814template<class Scalar>
1815RCP<LinearOpBase<Scalar> >
1816ModelEvaluatorBase::OutArgs<Scalar>::get_H_xx() const
1817{
1818 assert_supports(OUT_ARG_H_xx);
1819 return H_xx_;
1820}
1821
1822template<class Scalar>
1823RCP<LinearOpBase<Scalar> >
1824ModelEvaluatorBase::OutArgs<Scalar>::get_H_xp(int l) const
1825{
1826 assert_supports(OUT_ARG_H_xp,l);
1827 return H_xp_[l];
1828}
1829
1830template<class Scalar>
1831RCP<LinearOpBase<Scalar> >
1832ModelEvaluatorBase::OutArgs<Scalar>::get_H_pp(int l1, int l2) const
1833{
1834 assert_supports(OUT_ARG_H_pp,l1,l2);
1835 return H_pp_[ l1*Np() + l2 ];
1836}
1837
1838#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
1839
1840
1841#ifdef HAVE_THYRA_ME_POLYNOMIAL
1842
1843
1844template<class Scalar>
1845void ModelEvaluatorBase::OutArgs<Scalar>::set_f_poly(
1846 const RCP<Teuchos::Polynomial< VectorBase<Scalar> > > &f_poly
1847 )
1848{
1849 f_poly_ = f_poly;
1850}
1851
1852
1853template<class Scalar>
1854RCP<Teuchos::Polynomial< VectorBase<Scalar> > >
1855ModelEvaluatorBase::OutArgs<Scalar>::get_f_poly() const
1856{
1857 return f_poly_;
1858}
1859
1860
1861#endif // HAVE_THYRA_ME_POLYNOMIAL
1862
1863
1864template<class Scalar>
1866 const OutArgs<Scalar>& inputOutArgs, bool ignoreUnsupported
1867 )
1868{
1869 typedef ModelEvaluatorBase MEB;
1870 const int min_Np = TEUCHOS_MIN(this->Np(),inputOutArgs.Np());
1871 const int min_Ng = TEUCHOS_MIN(this->Ng(),inputOutArgs.Ng());
1872 // f
1873 if ( inputOutArgs.supports(OUT_ARG_f) && nonnull(inputOutArgs.get_f()) ) {
1874 if ( supports(OUT_ARG_f) || !ignoreUnsupported )
1875 set_f(inputOutArgs.get_f());
1876 }
1877 if ( inputOutArgs.supports(OUT_ARG_f_mp) && nonnull(inputOutArgs.get_f_mp()) ) {
1878 if ( supports(OUT_ARG_f_mp) || !ignoreUnsupported )
1879 set_f_mp(inputOutArgs.get_f_mp());
1880 }
1881#ifdef HAVE_THYRA_ME_POLYNOMIAL
1882 // f_poly
1883 if ( inputOutArgs.supports(OUT_ARG_f_poly) && nonnull(inputOutArgs.get_f_poly()) ) {
1884 if ( supports(OUT_ARG_f_poly) || !ignoreUnsupported )
1885 set_f_poly(inputOutArgs.get_f_poly());
1886 }
1887#endif // HAVE_THYRA_ME_POLYNOMIAL
1888 // g(j)
1889 for ( int j = 0; j < min_Ng; ++j ) {
1890 if ( nonnull(inputOutArgs.get_g(j)) )
1891 set_g(j,inputOutArgs.get_g(j));
1892 }
1893 for ( int j = 0; j < min_Ng; ++j ) {
1894 if ( inputOutArgs.supports(OUT_ARG_g_mp,j) && nonnull(inputOutArgs.get_g_mp(j)) ) {
1895 if ( supports(OUT_ARG_g_mp,j) || !ignoreUnsupported )
1896 set_g_mp(j,inputOutArgs.get_g_mp(j));
1897 }
1898 }
1899 // W
1900 if( inputOutArgs.supports(OUT_ARG_W) && nonnull(inputOutArgs.get_W()) ) {
1901 if ( supports(OUT_ARG_W) || !ignoreUnsupported )
1902 set_W(inputOutArgs.get_W());
1903 }
1904 if( inputOutArgs.supports(OUT_ARG_W_mp) && nonnull(inputOutArgs.get_W_mp()) ) {
1905 if ( supports(OUT_ARG_W_mp) || !ignoreUnsupported )
1906 set_W_mp(inputOutArgs.get_W_mp());
1907 }
1908 // W_op
1909 if( inputOutArgs.supports(OUT_ARG_W_op) && nonnull(inputOutArgs.get_W_op()) ) {
1910 if ( supports(OUT_ARG_W_op) || !ignoreUnsupported )
1911 set_W_op(inputOutArgs.get_W_op());
1912 }
1913 // W_prec
1914 if( inputOutArgs.supports(OUT_ARG_W_prec) && nonnull(inputOutArgs.get_W_prec()) ) {
1915 if ( supports(OUT_ARG_W_prec) || !ignoreUnsupported )
1916 set_W_prec(inputOutArgs.get_W_prec());
1917 }
1918 // DfDp(l)
1919 for ( int l = 0; l < min_Np; ++l ) {
1920 MEB::Derivative<Scalar> DfDp_l;
1921 if ( !inputOutArgs.supports(OUT_ARG_DfDp,l).none()
1922 && !(DfDp_l=inputOutArgs.get_DfDp(l)).isEmpty() )
1923 {
1924 if ( DfDp_l.isSupportedBy(supports(OUT_ARG_DfDp,l)) || !ignoreUnsupported )
1925 set_DfDp(l,DfDp_l);
1926 }
1927 }
1928 for ( int l = 0; l < min_Np; ++l ) {
1929 MEB::MPDerivative DfDp_mp_l;
1930 if ( !inputOutArgs.supports(OUT_ARG_DfDp_mp,l).none()
1931 && !(DfDp_mp_l=inputOutArgs.get_DfDp_mp(l)).isEmpty() )
1932 {
1933 if ( DfDp_mp_l.isSupportedBy(supports(OUT_ARG_DfDp_mp,l)) || !ignoreUnsupported )
1934 set_DfDp_mp(l,DfDp_mp_l);
1935 }
1936 }
1937 // DgDx_dot(j) and DgDx(j)
1938 for ( int j = 0; j < min_Ng; ++j ) {
1939 // DgDx_dot(j)
1940 MEB::Derivative<Scalar> DgDx_dot_j;
1941 if ( !inputOutArgs.supports(OUT_ARG_DgDx_dot,j).none()
1942 && !(DgDx_dot_j=inputOutArgs.get_DgDx_dot(j)).isEmpty() )
1943 {
1944 if( DgDx_dot_j.isSupportedBy(supports(OUT_ARG_DgDx_dot,j)) || !ignoreUnsupported )
1945 set_DgDx_dot(j,DgDx_dot_j);
1946 }
1947 // DgDx(j)
1948 MEB::Derivative<Scalar> DgDx_j;
1949 if ( !inputOutArgs.supports(OUT_ARG_DgDx,j).none()
1950 && !(DgDx_j=inputOutArgs.get_DgDx(j)).isEmpty() ) {
1951 if ( DgDx_j.isSupportedBy(supports(OUT_ARG_DgDx,j)) || !ignoreUnsupported )
1952 set_DgDx(j,DgDx_j);
1953 }
1954 }
1955 for ( int j = 0; j < min_Ng; ++j ) {
1956 // DgDx_dot(j)
1957 MEB::MPDerivative DgDx_dot_mp_j;
1958 if ( !inputOutArgs.supports(OUT_ARG_DgDx_dot_mp,j).none()
1959 && !(DgDx_dot_mp_j=inputOutArgs.get_DgDx_dot_mp(j)).isEmpty() )
1960 {
1961 if( DgDx_dot_mp_j.isSupportedBy(supports(OUT_ARG_DgDx_dot_mp,j)) || !ignoreUnsupported )
1962 set_DgDx_dot_mp(j,DgDx_dot_mp_j);
1963 }
1964 // DgDx(j)
1965 MEB::MPDerivative DgDx_mp_j;
1966 if ( !inputOutArgs.supports(OUT_ARG_DgDx_mp,j).none()
1967 && !(DgDx_mp_j=inputOutArgs.get_DgDx_mp(j)).isEmpty() ) {
1968 if ( DgDx_mp_j.isSupportedBy(supports(OUT_ARG_DgDx_mp,j)) || !ignoreUnsupported )
1969 set_DgDx_mp(j,DgDx_mp_j);
1970 }
1971 }
1972 // DgDp(j,l)
1973 for ( int l = 0; l < min_Np; ++l ) {
1974 for ( int j = 0; j < min_Ng; ++j ) {
1975 MEB::Derivative<Scalar> DgDp_j_l;
1976 if ( !inputOutArgs.supports(OUT_ARG_DgDp,j,l).none()
1977 && !(DgDp_j_l=inputOutArgs.get_DgDp(j,l)).isEmpty() )
1978 {
1979 if ( DgDp_j_l.isSupportedBy(supports(OUT_ARG_DgDp,j,l)) || !ignoreUnsupported )
1980 set_DgDp(j,l,DgDp_j_l);
1981 }
1982 }
1983 }
1984 for ( int l = 0; l < min_Np; ++l ) {
1985 for ( int j = 0; j < min_Ng; ++j ) {
1986 MEB::MPDerivative DgDp_mp_j_l;
1987 if ( !inputOutArgs.supports(OUT_ARG_DgDp_mp,j,l).none()
1988 && !(DgDp_mp_j_l=inputOutArgs.get_DgDp_mp(j,l)).isEmpty() )
1989 {
1990 if ( DgDp_mp_j_l.isSupportedBy(supports(OUT_ARG_DgDp_mp,j,l)) || !ignoreUnsupported )
1991 set_DgDp_mp(j,l,DgDp_mp_j_l);
1992 }
1993 }
1994 }
1995
1996#ifdef Thyra_BUILD_HESSIAN_SUPPORT
1997
1998 // hess_vec_prod_f
1999 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_f_xx) && nonnull(inputOutArgs.get_hess_vec_prod_f_xx()) ) {
2000 if ( supports(OUT_ARG_hess_vec_prod_f_xx) || !ignoreUnsupported )
2001 set_hess_vec_prod_f_xx(inputOutArgs.get_hess_vec_prod_f_xx());
2002 }
2003 for ( int l1 = 0; l1 < min_Np; ++l1 ) {
2004 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_f_xp,l1) && nonnull(inputOutArgs.get_hess_vec_prod_f_xp(l1)) ) {
2005 if ( supports(OUT_ARG_hess_vec_prod_f_xp,l1) || !ignoreUnsupported )
2006 set_hess_vec_prod_f_xp(l1,inputOutArgs.get_hess_vec_prod_f_xp(l1));
2007 }
2008 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_f_px,l1) && nonnull(inputOutArgs.get_hess_vec_prod_f_px(l1)) ) {
2009 if ( supports(OUT_ARG_hess_vec_prod_f_px,l1) || !ignoreUnsupported )
2010 set_hess_vec_prod_f_px(l1,inputOutArgs.get_hess_vec_prod_f_px(l1));
2011 }
2012 for ( int l2 = 0; l2 < min_Np; ++l2 ) {
2013 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_f_pp,l1,l2) && nonnull(inputOutArgs.get_hess_vec_prod_f_pp(l1,l2)) ) {
2014 if ( supports(OUT_ARG_hess_vec_prod_f_pp,l1,l2) || !ignoreUnsupported )
2015 set_hess_vec_prod_f_pp(l1,l2,inputOutArgs.get_hess_vec_prod_f_pp(l1,l2));
2016 }
2017 }
2018 }
2019
2020 // hess_vec_prod_g
2021 for ( int j = 0; j < min_Ng; ++j ) {
2022 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_g_xx,j) && nonnull(inputOutArgs.get_hess_vec_prod_g_xx(j)) ) {
2023 if ( supports(OUT_ARG_hess_vec_prod_g_xx,j) || !ignoreUnsupported )
2024 set_hess_vec_prod_g_xx(j,inputOutArgs.get_hess_vec_prod_g_xx(j));
2025 }
2026 for ( int l1 = 0; l1 < min_Np; ++l1 ) {
2027 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_g_xp,j,l1) && nonnull(inputOutArgs.get_hess_vec_prod_g_xp(j,l1)) ) {
2028 if ( supports(OUT_ARG_hess_vec_prod_g_xp,j,l1) || !ignoreUnsupported )
2029 set_hess_vec_prod_g_xp(j,l1,inputOutArgs.get_hess_vec_prod_g_xp(j,l1));
2030 }
2031 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_g_px,j,l1) && nonnull(inputOutArgs.get_hess_vec_prod_g_px(j,l1)) ) {
2032 if ( supports(OUT_ARG_hess_vec_prod_g_px,j,l1) || !ignoreUnsupported )
2033 set_hess_vec_prod_g_px(j,l1,inputOutArgs.get_hess_vec_prod_g_px(j,l1));
2034 }
2035 for ( int l2 = 0; l2 < min_Np; ++l2 ) {
2036 if( inputOutArgs.supports(OUT_ARG_hess_vec_prod_g_pp,j,l1,l2) && nonnull(inputOutArgs.get_hess_vec_prod_g_pp(j,l1,l2)) ) {
2037 if ( supports(OUT_ARG_hess_vec_prod_g_pp,j,l1,l2) || !ignoreUnsupported )
2038 set_hess_vec_prod_g_pp(j,l1,l2,inputOutArgs.get_hess_vec_prod_g_pp(j,l1,l2));
2039 }
2040 }
2041 }
2042 }
2043
2044 // hess_f
2045 if( inputOutArgs.supports(OUT_ARG_hess_f_xx) && nonnull(inputOutArgs.get_hess_f_xx()) ) {
2046 if ( supports(OUT_ARG_hess_f_xx) || !ignoreUnsupported )
2047 set_hess_f_xx(inputOutArgs.get_hess_f_xx());
2048 }
2049 for ( int l1 = 0; l1 < min_Np; ++l1 ) {
2050 if( inputOutArgs.supports(OUT_ARG_hess_f_xp,l1) && nonnull(inputOutArgs.get_hess_f_xp(l1)) ) {
2051 if ( supports(OUT_ARG_hess_f_xp,l1) || !ignoreUnsupported )
2052 set_hess_f_xp(l1,inputOutArgs.get_hess_f_xp(l1));
2053 }
2054 for ( int l2 = 0; l2 < min_Np; ++l2 ) {
2055 if( inputOutArgs.supports(OUT_ARG_hess_f_pp,l1,l2) && nonnull(inputOutArgs.get_hess_f_pp(l1,l2)) ) {
2056 if ( supports(OUT_ARG_hess_f_pp,l1,l2) || !ignoreUnsupported )
2057 set_hess_f_pp(l1,l2,inputOutArgs.get_hess_f_pp(l1,l2));
2058 }
2059 }
2060 }
2061
2062 // hess_g
2063 for ( int j = 0; j < min_Ng; ++j ) {
2064 if( inputOutArgs.supports(OUT_ARG_hess_g_xx,j) && nonnull(inputOutArgs.get_hess_g_xx(j)) ) {
2065 if ( supports(OUT_ARG_hess_g_xx,j) || !ignoreUnsupported )
2066 set_hess_g_xx(j,inputOutArgs.get_hess_g_xx(j));
2067 }
2068 for ( int l1 = 0; l1 < min_Np; ++l1 ) {
2069 if( inputOutArgs.supports(OUT_ARG_hess_g_xp,j,l1) && nonnull(inputOutArgs.get_hess_g_xp(j,l1)) ) {
2070 if ( supports(OUT_ARG_hess_g_xp,j,l1) || !ignoreUnsupported )
2071 set_hess_g_xp(j,l1,inputOutArgs.get_hess_g_xp(j,l1));
2072 }
2073 for ( int l2 = 0; l2 < min_Np; ++l2 ) {
2074 if( inputOutArgs.supports(OUT_ARG_hess_g_pp,j,l1,l2) && nonnull(inputOutArgs.get_hess_g_pp(j,l1,l2)) ) {
2075 if ( supports(OUT_ARG_hess_g_pp,j,l1,l2) || !ignoreUnsupported )
2076 set_hess_g_pp(j,l1,l2,inputOutArgs.get_hess_g_pp(j,l1,l2));
2077 }
2078 }
2079 }
2080 }
2081
2082 // H
2083 if( inputOutArgs.supports(OUT_ARG_H_xx) && nonnull(inputOutArgs.get_H_xx()) ) {
2084 if ( supports(OUT_ARG_H_xx) || !ignoreUnsupported )
2085 set_H_xx(inputOutArgs.get_H_xx());
2086 }
2087 for ( int l1 = 0; l1 < min_Np; ++l1 ) {
2088 if( inputOutArgs.supports(OUT_ARG_H_xp,l1) && nonnull(inputOutArgs.get_H_xp(l1)) ) {
2089 if ( supports(OUT_ARG_H_xp,l1) || !ignoreUnsupported )
2090 set_H_xp(l1,inputOutArgs.get_H_xp(l1));
2091 }
2092 for ( int l2 = 0; l2 < min_Np; ++l2 ) {
2093 if( inputOutArgs.supports(OUT_ARG_H_pp,l1,l2) && nonnull(inputOutArgs.get_H_pp(l1,l2)) ) {
2094 if ( supports(OUT_ARG_H_pp,l1,l2) || !ignoreUnsupported )
2095 set_H_pp(l1,l2,inputOutArgs.get_H_pp(l1,l2));
2096 }
2097 }
2098 }
2099
2100#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2101
2102 // Extended outArgs
2103 this->extended_outargs_ = inputOutArgs.extended_outargs_;
2104
2105 // ToDo: Add more args as needed!
2106}
2107
2108
2109template<class Scalar>
2111{
2113 isFailed_ = true;
2114 if( this->supports(OUT_ARG_f) && nonnull(this->get_f()) ) {
2115 assign(this->get_f().ptr(),ST::nan());
2116 }
2117 for( int j = 0; j < this->Ng(); ++j ) {
2118 if (nonnull(this->get_g(j)))
2119 assign(this->get_g(j).ptr(),ST::nan());
2120 }
2121 // ToDo: Set other objects to NaN as well!
2122}
2123
2124
2125template<class Scalar>
2127{
2128 return isFailed_;
2129}
2130
2131
2132template<class Scalar>
2134{
2135 if (!is_null(f_))
2136 return false;
2137 if (!is_null(W_))
2138 return false;
2139 if (!is_null(W_op_))
2140 return false;
2141 for ( int l = 0; l < Np(); ++l ) {
2142 if (!DfDp_[l].isEmpty())
2143 return false;
2144 }
2145#ifdef HAVE_THYRA_ME_POLYNOMIAL
2146 if (!is_null(f_poly_))
2147 return false;
2148#endif // HAVE_THYRA_ME_POLYNOMIAL
2149 for ( int j = 0; j < Ng(); ++j ) {
2150 if (!is_null(g_[j]))
2151 return false;
2152 if (!DgDx_dot_[j].isEmpty())
2153 return false;
2154 if (!DgDx_[j].isEmpty())
2155 return false;
2156 for ( int l = 0; l < Np(); ++l ) {
2157 if (!DgDp_[j*Np()+l].isEmpty())
2158 return false;
2159 }
2160 }
2161 return true;
2162}
2163
2164
2165template<class Scalar>
2167 const OutArgs<Scalar> &outArgs
2168 ) const
2169{
2170
2171 for ( int outArg_i = 0; outArg_i < NUM_E_OUT_ARGS_MEMBERS; ++outArg_i ) {
2172 const EOutArgsMembers outArg_arg = static_cast<EOutArgsMembers>(outArg_i);
2173 const std::string outArg_name = toString(outArg_arg);
2175 supports(outArg_arg) != outArgs.supports(outArg_arg), std::logic_error,
2176 "Error, the output argument "<<outArg_name<<" with support "<<outArgs.supports(outArg_arg)<<"\n"
2177 "in the OutArgs object for the model:\n\n"
2178 " "<<outArgs.modelEvalDescription()<<"\n\n"
2179 "is not the same the argument "<<outArg_name<<" with support "<<supports(outArg_arg)<<"\n"
2180 "in the OutArgs object for the model:\n\n"
2181 " "<<modelEvalDescription()<<"\n\n"
2182 "and these two OutArgs objects are not compatible!"
2183 );
2184 }
2185
2186 const int l_Np = this->Np();
2187 const int l_Ng = this->Ng();
2188 TEUCHOS_ASSERT_EQUALITY( l_Np, outArgs.Np() );
2189 TEUCHOS_ASSERT_EQUALITY( l_Ng, outArgs.Ng() );
2190
2191 if (supports(OUT_ARG_f)) {
2192 for ( int l = 0; l < l_Np; ++l ) {
2194 !supports(OUT_ARG_DfDp,l).isSameSupport(outArgs.supports(OUT_ARG_DfDp,l)),
2195 std::logic_error,
2196 "Error, the support for DfDp("<<l<<") is not the same for the models\n\n"
2197 " "<<outArgs.modelEvalDescription()<<"\n\n"
2198 "and:\n\n"
2199 " "<<modelEvalDescription()<<"\n\n"
2200 "and these two OutArgs objects are not compatible!"
2201 );
2202 }
2203 }
2204
2205 for ( int j = 0; j < l_Ng; ++j ) {
2207 !supports(OUT_ARG_DgDx_dot,j).isSameSupport(outArgs.supports(OUT_ARG_DgDx_dot,j)),
2208 std::logic_error,
2209 "Error, the support for DgDx_dot("<<j<<") is not the same for the models\n\n"
2210 " "<<outArgs.modelEvalDescription()<<"\n\n"
2211 "and:\n\n"
2212 " "<<modelEvalDescription()<<"\n\n"
2213 "and these two OutArgs objects are not compatible!"
2214 );
2216 !supports(OUT_ARG_DgDx,j).isSameSupport(outArgs.supports(OUT_ARG_DgDx,j)),
2217 std::logic_error,
2218 "Error, the support for DgDx("<<j<<") is not the same for the models\n\n"
2219 " "<<outArgs.modelEvalDescription()<<"\n\n"
2220 "and:\n\n"
2221 " "<<modelEvalDescription()<<"\n\n"
2222 "and these two OutArgs objects are not compatible!"
2223 );
2224 for ( int l = 0; l < l_Np; ++l ) {
2226 !supports(OUT_ARG_DgDp,j,l).isSameSupport(outArgs.supports(OUT_ARG_DgDp,j,l)),
2227 std::logic_error,
2228 "Error, the support for DgDp("<<j<<","<<l<<") is not the same for the models\n\n"
2229 " "<<outArgs.modelEvalDescription()<<"\n\n"
2230 "and:\n\n"
2231 " "<<modelEvalDescription()<<"\n\n"
2232 "and these two OutArgs objects are not compatible!"
2233 );
2234 }
2235 }
2236}
2237
2238
2239template<class Scalar>
2241{
2242 return modelEvalDescription_;
2243}
2244
2245
2246template<class Scalar>
2248{
2250 std::ostringstream oss;
2251 oss
2252 << "Thyra::ModelEvaluatorBase::OutArgs<"<<ST::name()<<">"
2253 << "{"
2254 << "model="<<modelEvalDescription_
2255 << ",Np="<<Np()
2256 << ",Ng="<<Ng()
2257 << "}";
2258 return oss.str();
2259}
2260
2261
2262template<class Scalar>
2264 Teuchos::FancyOStream &out_arg, const Teuchos::EVerbosityLevel verbLevel
2265 ) const
2266{
2267
2268 using Teuchos::OSTab;
2269 using Teuchos::describe;
2271 typedef RCP<const VectorBase<Scalar> > CV_ptr;
2272 typedef RCP<const LinearOpBase<Scalar> > CLO_ptr;
2273 typedef RCP<const LinearOpWithSolveBase<Scalar> > CLOWS_ptr;
2274 typedef ModelEvaluatorBase MEB;
2275 typedef MEB::Derivative<Scalar> Deriv;
2276
2277 if( verbLevel == Teuchos::VERB_NONE && verbLevel == Teuchos::VERB_DEFAULT )
2278 return;
2279
2281 out = Teuchos::rcp(&out_arg,false);
2282 OSTab tab(out);
2283
2284 *out <<"Thyra::ModelEvaluatorBase::OutArgs<"<<ST::name()<<">:\n";
2285 tab.incrTab();
2286
2287 *out <<"model = " << modelEvalDescription_ << "\n";
2288 *out <<"Np = " << Np() << "\n";
2289 *out <<"Ng = " << Ng() << "\n";
2290
2291 CV_ptr f;
2292 if (this->supports(OUT_ARG_f) && !is_null(f=get_f()) ) {
2293 *out << "f = " << Teuchos::describe(*f,verbLevel);
2294 }
2295
2296 for( int j = 0; j < Ng(); ++j ) {
2297 CV_ptr g_j;
2298 if (!is_null(g_j=this->get_g(j)))
2299 *out << "g("<<j<<") = " << Teuchos::describe(*g_j,verbLevel);
2300 }
2301
2302 CLOWS_ptr W;
2303 if ( this->supports(OUT_ARG_W) && !is_null(W=get_W()) ) {
2304 *out << "W = " << Teuchos::describe(*W,verbLevel);
2305 }
2306
2307 CLO_ptr W_op;
2308 if ( this->supports(OUT_ARG_W_op) && !is_null(W_op=get_W_op()) ) {
2309 *out << "W_op = " << Teuchos::describe(*W_op,verbLevel);
2310 }
2311
2312 for( int l = 0; l < Np(); ++l ) {
2313 Deriv DfDp_l;
2314 if (
2315 !this->supports(OUT_ARG_DfDp,l).none()
2316 && !(DfDp_l=get_DfDp(l)).isEmpty()
2317 )
2318 {
2319 *out << "DfDp("<<l<<") = ";
2320 DfDp_l.describe(*out,verbLevel);
2321 }
2322 }
2323
2324 for( int j = 0; j < Ng(); ++j ) {
2325
2326 Deriv DgDx_dot_j;
2327 if (
2328 !this->supports(OUT_ARG_DgDx_dot,j).none()
2329 && !(DgDx_dot_j=get_DgDx_dot(j)).isEmpty()
2330 )
2331 {
2332 *out << "DgDx_dot("<<j<<") = ";
2333 DgDx_dot_j.describe(*out,verbLevel);
2334 }
2335
2336 Deriv DgDx_j;
2337 if (
2338 !this->supports(OUT_ARG_DgDx,j).none()
2339 && !(DgDx_j=get_DgDx(j)).isEmpty()
2340 )
2341 {
2342 *out << "DgDx("<<j<<") = ";
2343 DgDx_j.describe(*out,verbLevel);
2344 }
2345
2346 for( int l = 0; l < Np(); ++l ) {
2347
2348 Deriv DgDp_j_l;
2349 if (
2350 !this->supports(OUT_ARG_DgDp,j,l).none()
2351 && !(DgDp_j_l=get_DgDp(j,l)).isEmpty()
2352 )
2353 {
2354 *out << "DgDp("<<j<<","<<l<<") = ";
2355 DgDp_j_l.describe(*out,verbLevel);
2356 }
2357 }
2358
2359 }
2360
2361 // ToDo: Add output for more objects?
2362
2363}
2364
2365
2366// protected
2367
2368
2369template<class Scalar>
2371 const std::string &modelEvalDescription_in
2372 )
2373{
2374 modelEvalDescription_ = modelEvalDescription_in;
2375}
2376
2377template<class Scalar>
2379{
2380 if(Np_in) {
2381 supports_DfDp_.resize(Np_in);
2382 DfDp_.resize(Np_in); std::fill_n(DfDp_.begin(),Np_in,Derivative<Scalar>());
2383 DfDp_properties_.resize(Np_in); std::fill_n(DfDp_properties_.begin(),Np_in,DerivativeProperties());
2384
2385 supports_DfDp_mp_.resize(Np_in);
2386 DfDp_mp_.resize(Np_in); std::fill_n(DfDp_mp_.begin(),Np_in,MPDerivative());
2387 DfDp_mp_properties_.resize(Np_in); std::fill_n(DfDp_mp_properties_.begin(),Np_in,DerivativeProperties());
2388 }
2389 if(Ng_in) {
2390 g_.resize(Ng_in); std::fill_n(g_.begin(),Ng_in,Teuchos::null);
2391 supports_DgDx_dot_.resize(Ng_in);
2392 DgDx_dot_.resize(Ng_in); std::fill_n(DgDx_dot_.begin(),Ng_in,Derivative<Scalar>());
2393 DgDx_dot_properties_.resize(Ng_in); std::fill_n(DgDx_dot_properties_.begin(),Ng_in,DerivativeProperties());
2394 supports_DgDx_.resize(Ng_in);
2395 DgDx_.resize(Ng_in); std::fill_n(DgDx_.begin(),Ng_in,Derivative<Scalar>());
2396 DgDx_properties_.resize(Ng_in); std::fill_n(DgDx_properties_.begin(),Ng_in,DerivativeProperties());
2397
2398#ifdef Thyra_BUILD_HESSIAN_SUPPORT
2399 supports_hess_vec_prod_g_xx_.resize(Ng_in);
2400 hess_vec_prod_g_xx_.resize(Ng_in); std::fill_n(hess_vec_prod_g_xx_.begin(),Ng_in,RCP<MultiVectorBase<Scalar> >());
2401
2402 supports_hess_g_xx_.resize(Ng_in);
2403 hess_g_xx_.resize(Ng_in); std::fill_n(hess_g_xx_.begin(),Ng_in,RCP<LinearOpBase<Scalar> >());
2404#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2405
2406 g_mp_.resize(Ng_in); std::fill_n(g_mp_.begin(),Ng_in,Teuchos::null);
2407 supports_g_mp_.resize(Ng_in);
2408 supports_DgDx_dot_mp_.resize(Ng_in);
2409 DgDx_dot_mp_.resize(Ng_in); std::fill_n(DgDx_dot_mp_.begin(),Ng_in,MPDerivative());
2410 DgDx_dot_mp_properties_.resize(Ng_in); std::fill_n(DgDx_dot_mp_properties_.begin(),Ng_in,DerivativeProperties());
2411 supports_DgDx_mp_.resize(Ng_in);
2412 DgDx_mp_.resize(Ng_in); std::fill_n(DgDx_mp_.begin(),Ng_in,MPDerivative());
2413 DgDx_mp_properties_.resize(Ng_in); std::fill_n(DgDx_mp_properties_.begin(),Ng_in,DerivativeProperties());
2414 }
2415#ifdef Thyra_BUILD_HESSIAN_SUPPORT
2416 if(Np_in) {
2417 const int Np = Np_in;
2418 const int NpNp = Np_in*Np_in;
2419
2420 supports_hess_vec_prod_f_xp_.resize(Np);
2421 hess_vec_prod_f_xp_.resize(Np); std::fill_n(hess_vec_prod_f_xp_.begin(),Np,RCP<MultiVectorBase<Scalar> >());
2422
2423 supports_hess_vec_prod_f_px_.resize(Np);
2424 hess_vec_prod_f_px_.resize(Np); std::fill_n(hess_vec_prod_f_px_.begin(),Np,RCP<MultiVectorBase<Scalar> >());
2425
2426 supports_hess_vec_prod_f_pp_.resize(NpNp);
2427 hess_vec_prod_f_pp_.resize(NpNp); std::fill_n(hess_vec_prod_f_pp_.begin(),NpNp,RCP<MultiVectorBase<Scalar> >());
2428
2429 supports_hess_f_xp_.resize(Np);
2430 hess_f_xp_.resize(Np); std::fill_n(hess_f_xp_.begin(),Np,RCP<LinearOpBase<Scalar> >());
2431
2432 supports_hess_f_pp_.resize(NpNp);
2433 hess_f_pp_.resize(NpNp); std::fill_n(hess_f_pp_.begin(),NpNp,RCP<LinearOpBase<Scalar> >());
2434
2435 supports_H_xp_.resize(Np);
2436 H_xp_.resize(Np); std::fill_n(H_xp_.begin(),Np,RCP<LinearOpBase<Scalar> >());
2437
2438 supports_H_pp_.resize(NpNp);
2439 H_pp_.resize(NpNp); std::fill_n(H_pp_.begin(),NpNp,RCP<LinearOpBase<Scalar> >());
2440 }
2441#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2442 if(Np_in && Ng_in) {
2443 const int NpNg = Np_in*Ng_in;
2444 const int NpNpNg = Np_in*Np_in*Ng_in;
2445 supports_DgDp_.resize(NpNg);
2446 DgDp_.resize(NpNg); std::fill_n(DgDp_.begin(),NpNg,Derivative<Scalar>());
2447 DgDp_properties_.resize(NpNg); std::fill_n(DgDp_properties_.begin(),NpNg,DerivativeProperties());
2448
2449#ifdef Thyra_BUILD_HESSIAN_SUPPORT
2450 supports_hess_vec_prod_g_xp_.resize(NpNg);
2451 hess_vec_prod_g_xp_.resize(NpNg); std::fill_n(hess_vec_prod_g_xp_.begin(),NpNg,RCP<MultiVectorBase<Scalar> >());
2452
2453 supports_hess_vec_prod_g_px_.resize(NpNg);
2454 hess_vec_prod_g_px_.resize(NpNg); std::fill_n(hess_vec_prod_g_px_.begin(),NpNg,RCP<MultiVectorBase<Scalar> >());
2455
2456 supports_hess_vec_prod_g_pp_.resize(NpNpNg);
2457 hess_vec_prod_g_pp_.resize(NpNpNg); std::fill_n(hess_vec_prod_g_pp_.begin(),NpNpNg,RCP<MultiVectorBase<Scalar> >());
2458
2459 supports_hess_g_xp_.resize(NpNg);
2460 hess_g_xp_.resize(NpNg); std::fill_n(hess_g_xp_.begin(),NpNg,RCP<LinearOpBase<Scalar> >());
2461
2462 supports_hess_g_pp_.resize(NpNpNg);
2463 hess_g_pp_.resize(NpNpNg); std::fill_n(hess_g_pp_.begin(),NpNpNg,RCP<LinearOpBase<Scalar> >());
2464#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2465
2466 supports_DgDp_mp_.resize(NpNg);
2467 DgDp_mp_.resize(NpNg); std::fill_n(DgDp_mp_.begin(),NpNg,MPDerivative());
2468 DgDp_mp_properties_.resize(NpNg); std::fill_n(DgDp_mp_properties_.begin(),NpNg,DerivativeProperties());
2469 }
2470}
2471
2472
2473template<class Scalar>
2475 EOutArgsMembers arg, bool supports_in )
2476{
2478 int(arg)>=NUM_E_OUT_ARGS_MEMBERS || int(arg) < 0,std::logic_error
2479 ,"model = \'"<<modelEvalDescription_
2480 <<"\': Error, arg="<<toString(arg)<<" is invalid!"
2481 );
2482 supports_[arg] = supports_in;
2483}
2484
2485
2486template<class Scalar>
2488 EOutArgsDfDp /* arg */, int l, const DerivativeSupport& supports_in
2489 )
2490{
2491 assert_supports(OUT_ARG_f);
2492 assert_l(l);
2493 supports_DfDp_[l] = supports_in;
2494}
2495
2496
2497template<class Scalar>
2499 EOutArgsDgDx_dot /* arg */, int j, const DerivativeSupport& supports_in
2500 )
2501{
2502 assert_j(j);
2503 supports_DgDx_dot_[j] = supports_in;
2504}
2505
2506
2507template<class Scalar>
2509 EOutArgsDgDx /* arg */, int j, const DerivativeSupport& supports_in
2510 )
2511{
2512 assert_j(j);
2513 supports_DgDx_[j] = supports_in;
2514}
2515
2516
2517#ifdef Thyra_BUILD_HESSIAN_SUPPORT
2518
2519template<class Scalar>
2521 EOutArgs_hess_vec_prod_f_xx /* arg */, bool supports_in
2522 )
2523{
2524 supports_hess_vec_prod_f_xx_ = supports_in;
2525}
2526
2527template<class Scalar>
2529 EOutArgs_hess_vec_prod_f_xp /* arg */, int l, bool supports_in
2530 )
2531{
2532 assert_l(l);
2533 supports_hess_vec_prod_f_xp_[l] = supports_in;
2534}
2535
2536template<class Scalar>
2538 EOutArgs_hess_vec_prod_f_px /* arg */, int l, bool supports_in
2539 )
2540{
2541 assert_l(l);
2542 supports_hess_vec_prod_f_px_[l] = supports_in;
2543}
2544
2545template<class Scalar>
2547 EOutArgs_hess_vec_prod_f_pp /* arg */, int l1, int l2, bool supports_in
2548 )
2549{
2550 assert_l(l1);
2551 assert_l(l2);
2552 supports_hess_vec_prod_f_pp_[ l1*Np()+ l2 ] = supports_in;
2553}
2554
2555#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2556
2557
2558template<class Scalar>
2560 EOutArgsDgDp /* arg */, int j, int l, const DerivativeSupport& supports_in
2561 )
2562{
2563 assert_j(j);
2564 assert_l(l);
2565 supports_DgDp_[ j*Np()+ l ] = supports_in;
2566}
2567
2568
2569#ifdef Thyra_BUILD_HESSIAN_SUPPORT
2570
2571template<class Scalar>
2573 EOutArgs_hess_vec_prod_g_xx /* arg */, int j, bool supports_in
2574 )
2575{
2576 assert_j(j);
2577 supports_hess_vec_prod_g_xx_[j] = supports_in;
2578}
2579
2580template<class Scalar>
2582 EOutArgs_hess_vec_prod_g_xp /* arg */, int j, int l, bool supports_in
2583 )
2584{
2585 assert_j(j);
2586 assert_l(l);
2587 supports_hess_vec_prod_g_xp_[ j*Np()+ l ] = supports_in;
2588}
2589
2590template<class Scalar>
2592 EOutArgs_hess_vec_prod_g_px /* arg */, int j, int l, bool supports_in
2593 )
2594{
2595 assert_j(j);
2596 assert_l(l);
2597 supports_hess_vec_prod_g_px_[ j*Np()+ l ] = supports_in;
2598}
2599
2600template<class Scalar>
2602 EOutArgs_hess_vec_prod_g_pp /* arg */, int j, int l1, int l2, bool supports_in
2603 )
2604{
2605 assert_j(j);
2606 assert_l(l1);
2607 assert_l(l2);
2608 supports_hess_vec_prod_g_pp_[ j*Np()*Np()+ l1*Np() + l2 ] = supports_in;
2609}
2610
2611
2612template<class Scalar>
2614 EOutArgs_hess_f_xx /* arg */, bool supports_in
2615 )
2616{
2617 supports_hess_f_xx_ = supports_in;
2618}
2619
2620template<class Scalar>
2622 EOutArgs_hess_f_xp /* arg */, int l, bool supports_in
2623 )
2624{
2625 assert_l(l);
2626 supports_hess_f_xp_[l] = supports_in;
2627}
2628
2629template<class Scalar>
2631 EOutArgs_hess_f_pp /* arg */, int l1, int l2, bool supports_in
2632 )
2633{
2634 assert_l(l1);
2635 assert_l(l2);
2636 supports_hess_f_pp_[ l1*Np()+ l2 ] = supports_in;
2637}
2638
2639template<class Scalar>
2641 EOutArgs_hess_g_xx /* arg */, int j, bool supports_in
2642 )
2643{
2644 assert_j(j);
2645 supports_hess_g_xx_[j] = supports_in;
2646}
2647
2648template<class Scalar>
2650 EOutArgs_hess_g_xp /* arg */, int j, int l, bool supports_in
2651 )
2652{
2653 assert_j(j);
2654 assert_l(l);
2655 supports_hess_g_xp_[ j*Np()+ l ] = supports_in;
2656}
2657
2658template<class Scalar>
2660 EOutArgs_hess_g_pp /* arg */, int j, int l1, int l2, bool supports_in
2661 )
2662{
2663 assert_j(j);
2664 assert_l(l1);
2665 assert_l(l2);
2666 supports_hess_g_pp_[ j*Np()*Np()+ l1*Np() + l2 ] = supports_in;
2667}
2668
2669template<class Scalar>
2671 EOutArgs_H_xx /* arg */, bool supports_in
2672 )
2673{
2674 supports_H_xx_ = supports_in;
2675}
2676
2677template<class Scalar>
2679 EOutArgs_H_xp /* arg */, int l, bool supports_in
2680 )
2681{
2682 assert_l(l);
2683 supports_H_xp_[l] = supports_in;
2684}
2685
2686template<class Scalar>
2688 EOutArgs_H_pp /* arg */, int l1, int l2, bool supports_in
2689 )
2690{
2691 assert_l(l1);
2692 assert_l(l2);
2693 supports_H_pp_[ l1*Np()+ l2 ] = supports_in;
2694}
2695
2696#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2697
2698
2699template<class Scalar>
2701 EOutArgs_g_mp /* arg */, int j, bool supports_in
2702 )
2703{
2704 //assert_supports(OUT_ARG_g_mp,j);
2705 assert_j(j);
2706 supports_g_mp_[j] = supports_in;
2707}
2708
2709
2710template<class Scalar>
2712 EOutArgsDfDp_mp /* arg */, int l, const DerivativeSupport& supports_in
2713 )
2714{
2715 assert_supports(OUT_ARG_f_mp); //JF this is not in epetraext ME
2716 assert_l(l);
2717 supports_DfDp_mp_[l] = supports_in;
2718}
2719
2720
2721template<class Scalar>
2723 EOutArgsDgDx_dot_mp /* arg */, int j, const DerivativeSupport& supports_in
2724 )
2725{
2726 assert_j(j);
2727 supports_DgDx_dot_mp_[j] = supports_in;
2728}
2729
2730
2731template<class Scalar>
2733 EOutArgsDgDx_mp /* arg */, int j, const DerivativeSupport& supports_in
2734 )
2735{
2736 assert_j(j);
2737 supports_DgDx_mp_[j] = supports_in;
2738}
2739
2740
2741template<class Scalar>
2743 EOutArgsDgDp_mp /* arg */, int j, int l, const DerivativeSupport& supports_in
2744 )
2745{
2746 assert_j(j);
2747 assert_l(l);
2748 supports_DgDp_mp_[ j*Np()+ l ] = supports_in;
2749}
2750
2751
2752template<class Scalar>
2754 const DerivativeProperties &properties
2755 )
2756{
2757 W_properties_ = properties;
2758}
2759
2760
2761template<class Scalar>
2763 int l, const DerivativeProperties &properties
2764 )
2765{
2766 assert_supports(OUT_ARG_DfDp,l);
2767 DfDp_properties_[l] = properties;
2768}
2769
2770
2771template<class Scalar>
2773 int j, const DerivativeProperties &properties
2774 )
2775{
2776 assert_supports(OUT_ARG_DgDx_dot,j);
2777 DgDx_dot_properties_[j] = properties;
2778}
2779
2780
2781template<class Scalar>
2783 int j, const DerivativeProperties &properties
2784 )
2785{
2786 assert_supports(OUT_ARG_DgDx,j);
2787 DgDx_properties_[j] = properties;
2788}
2789
2790
2791template<class Scalar>
2793 int j, int l, const DerivativeProperties &properties
2794 )
2795{
2796 assert_supports(OUT_ARG_DgDp,j,l);
2797 DgDp_properties_[ j*Np()+ l ] = properties;
2798}
2799
2800
2801template<class Scalar>
2803 int l, const DerivativeProperties &properties
2804 )
2805{
2806 assert_supports(OUT_ARG_DfDp_mp,l);
2807 DfDp_mp_properties_[l] = properties;
2808}
2809
2810
2811template<class Scalar>
2813 int j, const DerivativeProperties &properties
2814 )
2815{
2816 assert_supports(OUT_ARG_DgDx_dot_mp,j);
2817 DgDx_dot_mp_properties_[j] = properties;
2818}
2819
2820
2821template<class Scalar>
2822void ModelEvaluatorBase::OutArgs<Scalar>::_set_DgDx_mp_properties(
2823 int j, const DerivativeProperties &properties
2824 )
2825{
2826 assert_supports(OUT_ARG_DgDx_mp,j);
2827 DgDx_mp_properties_[j] = properties;
2828}
2829
2830
2831template<class Scalar>
2832void ModelEvaluatorBase::OutArgs<Scalar>::_set_DgDp_mp_properties(
2833 int j, int l, const DerivativeProperties &properties
2834 )
2835{
2836 assert_supports(OUT_ARG_DgDp_mp,j,l);
2837 DgDp_mp_properties_[ j*Np()+ l ] = properties;
2838}
2839
2840
2841template<class Scalar>
2843 const OutArgs<Scalar>& inputOutArgs
2844 )
2845{
2846 typedef ModelEvaluatorBase MEB;
2847 const int l_Np = TEUCHOS_MIN(this->Np(),inputOutArgs.Np());
2848 const int l_Ng = TEUCHOS_MIN(this->Ng(),inputOutArgs.Ng());
2849 std::copy(
2850 &inputOutArgs.supports_[0],
2851 &inputOutArgs.supports_[0] + NUM_E_OUT_ARGS_MEMBERS, &supports_[0] );
2852 for( int l = 0; l < l_Np; ++l ) {
2853 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DfDp,l);
2854 if (!ds.none()) {
2855 this->_setSupports(MEB::OUT_ARG_DfDp,l,ds);
2856 this->_set_DfDp_properties(l,inputOutArgs.get_DfDp_properties(l));
2857 }
2858 }
2859 for( int l = 0; l < l_Np; ++l ) {
2860 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DfDp_mp,l);
2861 if (!ds.none()) {
2862 this->_setSupports(MEB::OUT_ARG_DfDp_mp,l,ds);
2863 this->_set_DfDp_mp_properties(l,inputOutArgs.get_DfDp_mp_properties(l));
2864 }
2865 }
2866 for( int j = 0; j < l_Ng; ++j ) {
2867 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DgDx_dot,j);
2868 this->_setSupports(MEB::OUT_ARG_DgDx_dot,j,ds);
2869 if(!ds.none()) this->_set_DgDx_dot_properties(j,inputOutArgs.get_DgDx_dot_properties(j));
2870 }
2871 for( int j = 0; j < l_Ng; ++j ) {
2872 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DgDx_dot_mp,j);
2873 this->_setSupports(MEB::OUT_ARG_DgDx_dot_mp,j,ds);
2874 if(!ds.none()) this->_set_DgDx_dot_mp_properties(j,inputOutArgs.get_DgDx_dot_mp_properties(j));
2875 }
2876 for( int j = 0; j < l_Ng; ++j ) {
2877 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DgDx,j);
2878 this->_setSupports(MEB::OUT_ARG_DgDx,j,ds);
2879 if(!ds.none()) this->_set_DgDx_properties(j,inputOutArgs.get_DgDx_properties(j));
2880 }
2881 for( int j = 0; j < l_Ng; ++j ) {
2882 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DgDx_mp,j);
2883 this->_setSupports(MEB::OUT_ARG_DgDx_mp,j,ds);
2884 if(!ds.none()) this->_set_DgDx_mp_properties(j,inputOutArgs.get_DgDx_mp_properties(j));
2885 }
2886 for( int j = 0; j < l_Ng; ++j ) for( int l = 0; l < l_Np; ++l ) {
2887 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DgDp,j,l);
2888 this->_setSupports(MEB::OUT_ARG_DgDp,j,l,ds);
2889 if(!ds.none()) this->_set_DgDp_properties(j,l,inputOutArgs.get_DgDp_properties(j,l));
2890 }
2891 for( int j = 0; j < l_Ng; ++j ) for( int l = 0; l < l_Np; ++l ) {
2892 DerivativeSupport ds = inputOutArgs.supports(MEB::OUT_ARG_DgDp_mp,j,l);
2893 this->_setSupports(MEB::OUT_ARG_DgDp_mp,j,l,ds);
2894 if(!ds.none()) this->_set_DgDp_mp_properties(j,l,inputOutArgs.get_DgDp_mp_properties(j,l));
2895 }
2896 if(this->supports(OUT_ARG_W) || this->supports(OUT_ARG_W_op))
2897 this->_set_W_properties(inputOutArgs.get_W_properties());
2898 if(this->supports(OUT_ARG_W_mp))
2899 this->_set_W_properties(inputOutArgs.get_W_properties()); //JF should this be W_mp_properties?
2900
2901#ifdef Thyra_BUILD_HESSIAN_SUPPORT
2902 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_xx,inputOutArgs.supports(OUT_ARG_hess_vec_prod_f_xx));
2903 for( int l1 = 0; l1 < l_Np; ++l1 ) {
2904 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_xp,l1,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_f_xp,l1));
2905 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_px,l1,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_f_px,l1));
2906 for( int l2 = 0; l2 < l_Np; ++l2 ) {
2907 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_pp,l1,l2,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_f_pp,l1,l2));
2908 }
2909 }
2910 for( int j = 0; j < l_Ng; ++j ) {
2911 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_xx,j,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_g_xx,j));
2912 for( int l1 = 0; l1 < l_Np; ++l1 ) {
2913 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_xp,j,l1,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_g_xp,j,l1));
2914 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_px,j,l1,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_g_px,j,l1));
2915 for( int l2 = 0; l2 < l_Np; ++l2 ) {
2916 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_pp,j,l1,l2,inputOutArgs.supports(MEB::OUT_ARG_hess_vec_prod_g_pp,j,l1,l2));
2917 }
2918 }
2919 }
2920 this->_setSupports(MEB::OUT_ARG_hess_f_xx,inputOutArgs.supports(OUT_ARG_hess_f_xx));
2921 for( int l1 = 0; l1 < l_Np; ++l1 ) {
2922 this->_setSupports(MEB::OUT_ARG_hess_f_xp,l1,inputOutArgs.supports(MEB::OUT_ARG_hess_f_xp,l1));
2923 for( int l2 = 0; l2 < l_Np; ++l2 ) {
2924 this->_setSupports(MEB::OUT_ARG_hess_f_pp,l1,l2,inputOutArgs.supports(MEB::OUT_ARG_hess_f_pp,l1,l2));
2925 }
2926 }
2927 for( int j = 0; j < l_Ng; ++j ) {
2928 this->_setSupports(MEB::OUT_ARG_hess_g_xx,j,inputOutArgs.supports(MEB::OUT_ARG_hess_g_xx,j));
2929 for( int l1 = 0; l1 < l_Np; ++l1 ) {
2930 this->_setSupports(MEB::OUT_ARG_hess_g_xp,j,l1,inputOutArgs.supports(MEB::OUT_ARG_hess_g_xp,j,l1));
2931 for( int l2 = 0; l2 < l_Np; ++l2 ) {
2932 this->_setSupports(MEB::OUT_ARG_hess_g_pp,j,l1,l2,inputOutArgs.supports(MEB::OUT_ARG_hess_g_pp,j,l1,l2));
2933 }
2934 }
2935 }
2936 this->_setSupports(MEB::OUT_ARG_H_xx,inputOutArgs.supports(OUT_ARG_H_xx));
2937 for( int l1 = 0; l1 < l_Np; ++l1 ) {
2938 this->_setSupports(MEB::OUT_ARG_H_xp,l1,inputOutArgs.supports(MEB::OUT_ARG_H_xp,l1));
2939 for( int l2 = 0; l2 < l_Np; ++l2 ) {
2940 this->_setSupports(MEB::OUT_ARG_H_pp,l1,l2,inputOutArgs.supports(MEB::OUT_ARG_H_pp,l1,l2));
2941 }
2942 }
2943#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
2944
2945 extended_outargs_ = inputOutArgs.extended_outargs_;
2946}
2947
2948
2949template<class Scalar>
2951 EInArgsMembers arg
2952 )
2953{
2954 switch(arg) {
2955 case IN_ARG_x: {
2956 const int l_Ng = this->Ng();
2957 for( int j = 0; j < l_Ng; ++j ) {
2958 this->_setSupports(OUT_ARG_DgDx_dot,j,DerivativeSupport());
2959 this->_setSupports(OUT_ARG_DgDx,j,DerivativeSupport());
2960 }
2961 break;
2962 }
2963 case IN_ARG_x_mp: {
2964 const int l_Ng = this->Ng();
2965 for( int j = 0; j < l_Ng; ++j ) {
2966 this->_setSupports(OUT_ARG_DgDx_dot_mp,j,DerivativeSupport());
2967 this->_setSupports(OUT_ARG_DgDx_mp,j,DerivativeSupport());
2968 }
2969 break;
2970 }
2971 default:
2973 true ,std::logic_error,
2974 "Error, can not handle args other than IN_ARG_x yet!"
2975 );
2976 }
2977}
2978
2979
2980template<class Scalar>
2982 EOutArgsMembers arg
2983 )
2984{
2985 switch(arg) {
2986 case OUT_ARG_f: {
2987 this->_setSupports(OUT_ARG_W,false);
2988 this->_setSupports(OUT_ARG_W_op,false);
2989 this->_setSupports(OUT_ARG_f_poly,false);
2990 const int l_Np = this->Np();
2991 for( int l = 0; l < l_Np; ++l )
2992 this->_setSupports(OUT_ARG_DfDp,l,DerivativeSupport());
2993 break;
2994 }
2995 case OUT_ARG_f_mp: {
2996 this->_setSupports(OUT_ARG_W_mp,false);
2997 this->_setSupports(OUT_ARG_W_op,false);
2998 this->_setSupports(OUT_ARG_f_poly,false);
2999 const int l_Np = this->Np();
3000 for( int l = 0; l < l_Np; ++l )
3001 this->_setSupports(OUT_ARG_DfDp_mp,l,DerivativeSupport());
3002 break;
3003 }
3004 default:
3006 true ,std::logic_error,
3007 "Error, can not handle args other than OUT_ARG_f yet!"
3008 );
3009 }
3010 this->_setSupports(arg,false);
3011}
3012
3013#ifdef Thyra_BUILD_HESSIAN_SUPPORT
3014template<class Scalar>
3016 const bool supports
3017 )
3018{
3019 typedef ModelEvaluatorBase MEB;
3020 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_xx,supports);
3021 this->_setSupports(MEB::OUT_ARG_hess_f_xx,supports);
3022 this->_setSupports(MEB::OUT_ARG_H_xx,supports);
3023 for (int l1=0; l1<this->Np(); ++l1)
3024 {
3025 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_xp,l1,supports);
3026 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_px,l1,supports);
3027 this->_setSupports(MEB::OUT_ARG_hess_f_xp,l1,supports);
3028 this->_setSupports(MEB::OUT_ARG_H_xp,l1,supports);
3029 for (int l2=0; l2<this->Np(); ++l2)
3030 {
3031 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_f_pp,l1,l2,supports);
3032 this->_setSupports(MEB::OUT_ARG_hess_f_pp,l1,l2,supports);
3033 this->_setSupports(MEB::OUT_ARG_H_pp,l1,l2,supports);
3034 }
3035 }
3036
3037 for (int j=0; j<this->Ng(); ++j)
3038 {
3039 this->_setSupports(MEB::OUT_ARG_hess_g_xx,j,supports);
3040 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_xx,j,supports);
3041 for (int l1=0; l1<this->Np(); ++l1)
3042 {
3043 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_xp,j,l1,supports);
3044 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_px,j,l1,supports);
3045 this->_setSupports(MEB::OUT_ARG_hess_g_xp,j,l1,supports);
3046 for (int l2=0; l2<this->Np(); ++l2)
3047 {
3048 this->_setSupports(MEB::OUT_ARG_hess_vec_prod_g_pp,j,l1,l2,supports);
3049 this->_setSupports(MEB::OUT_ARG_hess_g_pp,j,l1,l2,supports);
3050 }
3051 }
3052 }
3053}
3054#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
3055
3056// private
3057
3058
3059template<class Scalar>
3060void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(EOutArgsMembers arg) const
3061{
3063 !this->supports(arg), std::logic_error
3064 ,"Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(arg):\n\n"
3065 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3066 "Error, The argument arg = " << toString(arg) << " is not supported!"
3067 );
3068}
3069
3070
3071template<class Scalar>
3072void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3073 EOutArgsDfDp arg, int l, const Derivative<Scalar> &deriv
3074 ) const
3075{
3076 const DerivativeSupport derivSupport = this->supports(arg,l);
3078 !deriv.isSupportedBy(derivSupport), std::logic_error,
3079 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DfDp,l):\n\n"
3080 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3081 "Error, The argument DfDp("<<l<<") = " << deriv.description() << "\n"
3082 "is not supported!\n\n"
3083 "The supported types include " << derivSupport.description() << "!"
3084 );
3085}
3086
3087
3088template<class Scalar>
3089void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3090 EOutArgsDgDx_dot arg, int j, const Derivative<Scalar> &deriv
3091 ) const
3092{
3093 const DerivativeSupport derivSupport = this->supports(arg,j);
3095 !deriv.isSupportedBy(derivSupport), std::logic_error,
3096 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DgDx_dot,j):\n\n"
3097 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3098 "Error, The argument DgDx_dot("<<j<<") = " << deriv.description() << "\n"
3099 "is not supported!\n\n"
3100 "The supported types include " << derivSupport.description() << "!"
3101 );
3102}
3103
3104
3105template<class Scalar>
3106void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3107 EOutArgsDgDx arg, int j, const Derivative<Scalar> &deriv
3108 ) const
3109{
3110 const DerivativeSupport derivSupport = this->supports(arg,j);
3112 !deriv.isSupportedBy(derivSupport), std::logic_error,
3113 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DgDx,j):\n\n"
3114 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3115 "Error, The argument DgDx("<<j<<") = " << deriv.description() << "\n"
3116 "is not supported!\n\n"
3117 "The supported types include " << derivSupport.description() << "!"
3118 );
3119}
3120
3121
3122template<class Scalar>
3123void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3124 EOutArgsDgDp arg, int j, int l, const Derivative<Scalar> &deriv
3125 ) const
3126{
3127 const DerivativeSupport derivSupport = this->supports(arg,j,l);
3129 !deriv.isSupportedBy(derivSupport), std::logic_error,
3130 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DgDp,j,l):\n\n"
3131 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3132 "Error, The argument DgDp("<<j<<","<<l<<") = " << deriv.description() << "\n"
3133 "is not supported!\n\n"
3134 "The supported types include " << derivSupport.description() << "!"
3135 );
3136}
3137
3138
3139#ifdef Thyra_BUILD_HESSIAN_SUPPORT
3140
3141template<class Scalar>
3142void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3143 EOutArgs_hess_vec_prod_f_xx arg
3144 ) const
3145{
3146 const bool support = this->supports(arg);
3148 !support, std::logic_error,
3149 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_f_xx):\n\n"
3150 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3151 "Error, The argument hess_vec_prod_f_xx() is not supported!"
3152 );
3153}
3154
3155template<class Scalar>
3156void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3157 EOutArgs_hess_vec_prod_f_xp arg, int l
3158 ) const
3159{
3160 const bool support = this->supports(arg,l);
3162 !support, std::logic_error,
3163 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_f_xp,l):\n\n"
3164 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3165 "Error, The argument hess_vec_prod_f_xp("<<l<<") is not supported!"
3166 );
3167}
3168
3169template<class Scalar>
3170void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3171 EOutArgs_hess_vec_prod_f_px arg, int l
3172 ) const
3173{
3174 const bool support = this->supports(arg,l);
3176 !support, std::logic_error,
3177 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_f_px,l):\n\n"
3178 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3179 "Error, The argument hess_vec_prod_f_px("<<l<<") is not supported!"
3180 );
3181}
3182
3183template<class Scalar>
3184void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3185 EOutArgs_hess_vec_prod_f_pp arg, int l1, int l2
3186 ) const
3187{
3188 const bool support = this->supports(arg,l1,l2);
3190 !support, std::logic_error,
3191 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_f_pp,l1,l2):\n\n"
3192 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3193 "Error, The argument hess_vec_prod_f_pp("<<l1<<","<<l2<<") is not supported!"
3194 );
3195}
3196
3197template<class Scalar>
3198void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3199 EOutArgs_hess_vec_prod_g_xx arg, int j
3200 ) const
3201{
3202 const bool support = this->supports(arg,j);
3204 !support, std::logic_error,
3205 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_g_xx,j):\n\n"
3206 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3207 "Error, The argument hess_vec_prod_g_xx("<<j<<") is not supported!"
3208 );
3209}
3210
3211template<class Scalar>
3212void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3213 EOutArgs_hess_vec_prod_g_xp arg, int j, int l
3214 ) const
3215{
3216 const bool support = this->supports(arg,j,l);
3218 !support, std::logic_error,
3219 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_g_xp,j,l):\n\n"
3220 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3221 "Error, The argument hess_vec_prod_g_xp("<<j<<","<<l<<") is not supported!"
3222 );
3223}
3224
3225template<class Scalar>
3226void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3227 EOutArgs_hess_vec_prod_g_px arg, int j, int l
3228 ) const
3229{
3230 const bool support = this->supports(arg,j,l);
3232 !support, std::logic_error,
3233 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_g_px,j,l):\n\n"
3234 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3235 "Error, The argument hess_vec_prod_g_px("<<j<<","<<l<<") is not supported!"
3236 );
3237}
3238
3239template<class Scalar>
3240void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3241 EOutArgs_hess_vec_prod_g_pp arg, int j, int l1, int l2
3242 ) const
3243{
3244 const bool support = this->supports(arg,j,l1,l2);
3246 !support, std::logic_error,
3247 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_vec_prod_g_pp,j,l1,l2):\n\n"
3248 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3249 "Error, The argument hess_vec_prod_g_pp("<<j<<","<<l1<<","<<l2<<") is not supported!"
3250 );
3251}
3252
3253template<class Scalar>
3254void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3255 EOutArgs_hess_f_xx arg
3256 ) const
3257{
3258 const bool support = this->supports(arg);
3260 !support, std::logic_error,
3261 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_f_xx):\n\n"
3262 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3263 "Error, The argument hess_f_xx() is not supported!"
3264 );
3265}
3266
3267template<class Scalar>
3268void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3269 EOutArgs_hess_f_xp arg, int l
3270 ) const
3271{
3272 const bool support = this->supports(arg,l);
3274 !support, std::logic_error,
3275 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_f_xp,l):\n\n"
3276 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3277 "Error, The argument hess_f_xp() is not supported!"
3278 );
3279}
3280
3281template<class Scalar>
3282void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3283 EOutArgs_hess_f_pp arg, int l1, int l2
3284 ) const
3285{
3286 const bool support = this->supports(arg,l1,l2);
3288 !support, std::logic_error,
3289 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_f_pp,l1,l2):\n\n"
3290 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3291 "Error, The argument hess_f_pp() is not supported!"
3292 );
3293}
3294
3295template<class Scalar>
3296void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3297 EOutArgs_hess_g_xx arg, int j
3298 ) const
3299{
3300 const bool support = this->supports(arg,j);
3302 !support, std::logic_error,
3303 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_g_xx,j):\n\n"
3304 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3305 "Error, The argument hess_g_xx() is not supported!"
3306 );
3307}
3308
3309template<class Scalar>
3310void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3311 EOutArgs_hess_g_xp arg, int j, int l
3312 ) const
3313{
3314 const bool support = this->supports(arg,j,l);
3316 !support, std::logic_error,
3317 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_g_xp,j,l):\n\n"
3318 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3319 "Error, The argument hess_g_xp() is not supported!"
3320 );
3321}
3322
3323template<class Scalar>
3324void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3325 EOutArgs_hess_g_pp arg, int j, int l1, int l2
3326 ) const
3327{
3328 const bool support = this->supports(arg,j,l1,l2);
3330 !support, std::logic_error,
3331 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_hess_g_pp,j,l1,l2):\n\n"
3332 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3333 "Error, The argument hess_g_pp() is not supported!"
3334 );
3335}
3336
3337template<class Scalar>
3338void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3339 EOutArgs_H_xx arg
3340 ) const
3341{
3342 const bool support = this->supports(arg);
3344 !support, std::logic_error,
3345 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_H_xx):\n\n"
3346 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3347 "Error, The argument H_xx() is not supported!"
3348 );
3349}
3350
3351template<class Scalar>
3352void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3353 EOutArgs_H_xp arg, int l
3354 ) const
3355{
3356 const bool support = this->supports(arg,l);
3358 !support, std::logic_error,
3359 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_H_xp,l):\n\n"
3360 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3361 "Error, The argument H_xp() is not supported!"
3362 );
3363}
3364
3365template<class Scalar>
3366void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3367 EOutArgs_H_pp arg, int l1, int l2
3368 ) const
3369{
3370 const bool support = this->supports(arg,l1,l2);
3372 !support, std::logic_error,
3373 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_H_pp,l1,l2):\n\n"
3374 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3375 "Error, The argument H_pp() is not supported!"
3376 );
3377}
3378
3379#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
3380
3381
3382template<class Scalar>
3383void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3384 EOutArgs_g_mp /* arg */, int j
3385 ) const
3386{
3387 assert_j(j);
3389 !supports_g_mp_[j], std::logic_error,
3390 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_g_mp,j):\n\n"
3391 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3392 "Error, The argument g_mp("<<j<<") \n"
3393 "is not supported!\n\n"
3394 );
3395}
3396
3397
3398template<class Scalar>
3399void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3400 EOutArgsDfDp_mp arg, int l, const MPDerivative &deriv
3401 ) const
3402{
3403 const DerivativeSupport derivSupport = this->supports(arg,l);
3405 !deriv.isSupportedBy(derivSupport), std::logic_error,
3406 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DfDp_mp,l):\n\n"
3407 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3408 "Error, The argument DfDp_mp("<<l<<") = " << deriv.description() << "\n"
3409 "is not supported!\n\n"
3410 "The supported types include " << derivSupport.description() << "!"
3411 );
3412}
3413
3414
3415template<class Scalar>
3416void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3417 EOutArgsDgDx_dot_mp arg, int j, const MPDerivative &deriv
3418 ) const
3419{
3420 const DerivativeSupport derivSupport = this->supports(arg,j);
3422 !deriv.isSupportedBy(derivSupport), std::logic_error,
3423 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DgDx_dot_mp,j):\n\n"
3424 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3425 "Error, The argument DgDx_dot_mp("<<j<<") = " << deriv.description() << "\n"
3426 "is not supported!\n\n"
3427 "The supported types include " << derivSupport.description() << "!"
3428 );
3429}
3430
3431
3432template<class Scalar>
3433void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3434 EOutArgsDgDx_mp arg, int j, const MPDerivative &deriv
3435 ) const
3436{
3437 const DerivativeSupport derivSupport = this->supports(arg,j);
3439 !deriv.isSupportedBy(derivSupport), std::logic_error,
3440 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DgDx_mp,j):\n\n"
3441 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3442 "Error, The argument DgDx_mp("<<j<<") = " << deriv.description() << "\n"
3443 "is not supported!\n\n"
3444 "The supported types include " << derivSupport.description() << "!"
3445 );
3446}
3447
3448
3449template<class Scalar>
3450void ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(
3451 EOutArgsDgDp_mp arg, int j, int l, const MPDerivative &deriv
3452 ) const
3453{
3454 const DerivativeSupport derivSupport = this->supports(arg,j,l);
3456 !deriv.isSupportedBy(derivSupport), std::logic_error,
3457 "Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_supports(OUT_ARG_DgDp_mp,j,l):\n\n"
3458 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3459 "Error, The argument DgDp_mp("<<j<<","<<l<<") = " << deriv.description() << "\n"
3460 "is not supported!\n\n"
3461 "The supported types include " << derivSupport.description() << "!"
3462 );
3463}
3464
3465
3466template<class Scalar>
3467void ModelEvaluatorBase::OutArgs<Scalar>::assert_l(int l) const
3468{
3470 !( 0 <= l && l < Np() ), std::logic_error
3471 ,"Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_l(l):\n\n"
3472 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3473 "Error, The parameter subvector p("<<l<<")"
3474 " is not in the range [0,"<<Np()<<")!"
3475 );
3476}
3477
3478
3479template<class Scalar>
3480void ModelEvaluatorBase::OutArgs<Scalar>::assert_j(int j) const
3481{
3483 !( 0 <= j && j < Ng() ), std::logic_error
3484 ,"Thyra::ModelEvaluatorBase::OutArgs<Scalar>::assert_j(j):\n\n"
3485 "model = \'"<<modelEvalDescription_<<"\':\n\n"
3486 "Error, The auxiliary function g("<<j<<")"
3487 " is not in the range [0,"<<Ng()<<")!"
3488 );
3489}
3490
3491
3492//
3493// ModelEvaluatorBase::InArgsSetup
3494//
3495
3496
3497template<class Scalar>
3500
3501
3502template<class Scalar>
3506
3507
3508template<class Scalar>
3510 const std::string &modelEvalDescription_in )
3511{
3512 this->_setModelEvalDescription(modelEvalDescription_in);
3513}
3514
3515
3516template<class Scalar>
3518{ this->_set_Np_Ng(Np_in, 0); }
3519
3520template<class Scalar>
3522{ this->_set_Np_Ng(Np_in, Ng_in); }
3523
3524template<class Scalar>
3526{ this->_setSupports(arg,supports_in); }
3527
3528
3529template<class Scalar>
3530void ModelEvaluatorBase::InArgsSetup<Scalar>::setSupports( EInArgs_p_mp arg, int l, bool supports_in )
3531{ this->_setSupports(arg,l,supports_in); }
3532
3533
3534template<class Scalar>
3536 const InArgs<Scalar>& inArgs, const int Np_in
3537 )
3538{
3539 this->_setSupports(inArgs, Np_in);
3540}
3541
3542
3543template<class Scalar>
3545 EInArgsMembers arg
3546 )
3547{
3548 this->_setUnsupportsAndRelated(arg);
3549}
3550
3551
3552//
3553// ModelEvaluatorBase::OutArgsSetup
3554//
3555
3556
3557template<class Scalar>
3560
3561
3562template<class Scalar>
3564 const OutArgs<Scalar>& inputOutArgs
3565 )
3566 :OutArgs<Scalar>(inputOutArgs)
3567{}
3568
3569
3570template<class Scalar>
3572 const std::string &modelEvalDescription_in
3573 )
3574{ this->_setModelEvalDescription(modelEvalDescription_in); }
3575
3576
3577template<class Scalar>
3579{ this->_set_Np_Ng(Np_in, Ng_in); }
3580
3581
3582template<class Scalar>
3584 EOutArgsMembers arg, bool supports_in
3585 )
3586{ this->_setSupports(arg,supports_in); }
3587
3588
3589template<class Scalar>
3591 EOutArgsDfDp arg, int l, const DerivativeSupport& supports_in
3592 )
3593{ this->_setSupports(arg,l,supports_in); }
3594
3595
3596template<class Scalar>
3598 EOutArgsDgDx_dot arg, int j, const DerivativeSupport& supports_in
3599 )
3600{ this->_setSupports(arg,j,supports_in); }
3601
3602
3603template<class Scalar>
3605 EOutArgsDgDx arg, int j, const DerivativeSupport& supports_in
3606 )
3607{ this->_setSupports(arg,j,supports_in); }
3608
3609
3610#ifdef Thyra_BUILD_HESSIAN_SUPPORT
3611
3612template<class Scalar>
3614 EOutArgs_hess_vec_prod_f_xx arg, bool supports_in
3615 )
3616{ this->_setSupports(arg,supports_in); }
3617
3618template<class Scalar>
3620 EOutArgs_hess_vec_prod_f_xp arg, int l, bool supports_in
3621 )
3622{ this->_setSupports(arg,l,supports_in); }
3623
3624template<class Scalar>
3626 EOutArgs_hess_vec_prod_f_px arg, int l, bool supports_in
3627 )
3628{ this->_setSupports(arg,l,supports_in); }
3629
3630template<class Scalar>
3632 EOutArgs_hess_vec_prod_f_pp arg, int l1, int l2, bool supports_in
3633 )
3634{ this->_setSupports(arg,l1,l2,supports_in); }
3635
3636#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
3637
3638
3639template<class Scalar>
3641 EOutArgsDgDp arg, int j, int l, const DerivativeSupport& supports_in
3642 )
3643{ this->_setSupports(arg,j,l,supports_in); }
3644
3645
3646#ifdef Thyra_BUILD_HESSIAN_SUPPORT
3647
3648template<class Scalar>
3650 EOutArgs_hess_vec_prod_g_xx arg, int j, bool supports_in
3651 )
3652{ this->_setSupports(arg,j,supports_in); }
3653
3654template<class Scalar>
3656 EOutArgs_hess_vec_prod_g_xp arg, int j, int l, bool supports_in
3657 )
3658{ this->_setSupports(arg,j,l,supports_in); }
3659
3660template<class Scalar>
3662 EOutArgs_hess_vec_prod_g_px arg, int j, int l, bool supports_in
3663 )
3664{ this->_setSupports(arg,j,l,supports_in); }
3665
3666template<class Scalar>
3668 EOutArgs_hess_vec_prod_g_pp arg, int j, int l1, int l2, bool supports_in
3669 )
3670{ this->_setSupports(arg,j,l1,l2,supports_in); }
3671
3672template<class Scalar>
3674 EOutArgs_hess_f_xx arg, bool supports_in
3675 )
3676{ this->_setSupports(arg,supports_in); }
3677
3678template<class Scalar>
3680 EOutArgs_hess_f_xp arg, int l, bool supports_in
3681 )
3682{ this->_setSupports(arg,l,supports_in); }
3683
3684template<class Scalar>
3686 EOutArgs_hess_f_pp arg, int l1, int l2, bool supports_in
3687 )
3688{ this->_setSupports(arg,l1,l2,supports_in); }
3689
3690template<class Scalar>
3692 EOutArgs_hess_g_xx arg, int j, bool supports_in
3693 )
3694{ this->_setSupports(arg,j,supports_in); }
3695
3696template<class Scalar>
3698 EOutArgs_hess_g_xp arg, int j, int l, bool supports_in
3699 )
3700{ this->_setSupports(arg,j,l,supports_in); }
3701
3702template<class Scalar>
3704 EOutArgs_hess_g_pp arg, int j, int l1, int l2, bool supports_in
3705 )
3706{ this->_setSupports(arg,j,l1,l2,supports_in); }
3707
3708template<class Scalar>
3710 EOutArgs_H_xx arg, bool supports_in
3711 )
3712{ this->_setSupports(arg,supports_in); }
3713
3714template<class Scalar>
3716 EOutArgs_H_xp arg, int l, bool supports_in
3717 )
3718{ this->_setSupports(arg,l,supports_in); }
3719
3720template<class Scalar>
3722 EOutArgs_H_pp arg, int l1, int l2, bool supports_in
3723 )
3724{ this->_setSupports(arg,l1,l2,supports_in); }
3725
3726#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
3727
3728
3729template<class Scalar>
3731 EOutArgs_g_mp arg, int j, bool supports_in
3732 )
3733{ this->_setSupports(arg,j,supports_in); }
3734
3735
3736template<class Scalar>
3738 EOutArgsDfDp_mp arg, int l, const DerivativeSupport& supports_in
3739 )
3740{ this->_setSupports(arg,l,supports_in); }
3741
3742
3743template<class Scalar>
3745 EOutArgsDgDx_dot_mp arg, int j, const DerivativeSupport& supports_in
3746 )
3747{ this->_setSupports(arg,j,supports_in); }
3748
3749
3750template<class Scalar>
3752 EOutArgsDgDx_mp arg, int j, const DerivativeSupport& supports_in
3753 )
3754{ this->_setSupports(arg,j,supports_in); }
3755
3756
3757template<class Scalar>
3759 EOutArgsDgDp_mp arg, int j, int l, const DerivativeSupport& supports_in
3760 )
3761{ this->_setSupports(arg,j,l,supports_in); }
3762
3763
3764template<class Scalar>
3766 const DerivativeProperties &properties
3767 )
3768{ this->_set_W_properties(properties); }
3769
3770
3771template<class Scalar>
3773 int l, const DerivativeProperties &properties
3774 )
3775{ this->_set_DfDp_properties(l,properties); }
3776
3777
3778template<class Scalar>
3780 int j, const DerivativeProperties &properties
3781 )
3782{ this->_set_DgDx_dot_properties(j,properties); }
3783
3784
3785template<class Scalar>
3787 int j, const DerivativeProperties &properties
3788 )
3789{ this->_set_DgDx_properties(j,properties); }
3790
3791
3792template<class Scalar>
3794 int j, int l, const DerivativeProperties &properties
3795 )
3796{ this->_set_DgDp_properties(j,l,properties); }
3797
3798
3799template<class Scalar>
3801 int l, const DerivativeProperties &properties
3802 )
3803{ this->_set_DfDp_mp_properties(l,properties); }
3804
3805
3806template<class Scalar>
3808 int j, const DerivativeProperties &properties
3809 )
3810{ this->_set_DgDx_dot_mp_properties(j,properties); }
3811
3812
3813template<class Scalar>
3814void ModelEvaluatorBase::OutArgsSetup<Scalar>::set_DgDx_mp_properties(
3815 int j, const DerivativeProperties &properties
3816 )
3817{ this->_set_DgDx_mp_properties(j,properties); }
3818
3819
3820template<class Scalar>
3821void ModelEvaluatorBase::OutArgsSetup<Scalar>::set_DgDp_mp_properties(
3822 int j, int l, const DerivativeProperties &properties
3823 )
3824{ this->_set_DgDp_mp_properties(j,l,properties); }
3825
3826
3827template<class Scalar>
3829 const OutArgs<Scalar>& inputOutArgs
3830 )
3831{ this->_setSupports(inputOutArgs); }
3832
3833
3834template<class Scalar>
3836 EInArgsMembers arg
3837 )
3838{ this->_setUnsupportsAndRelated(arg); }
3839
3840
3841template<class Scalar>
3843 EOutArgsMembers arg
3844 )
3845{ this->_setUnsupportsAndRelated(arg); }
3846
3847
3848#ifdef Thyra_BUILD_HESSIAN_SUPPORT
3849template<class Scalar>
3851 const bool supports
3852 )
3853{ this->_setHessianSupports(supports); }
3854#endif // ifdef Thyra_BUILD_HESSIAN_SUPPORT
3855
3856} // namespace Thyra
3857
3858
3859
3860//
3861// Explicit instantiation macro
3862//
3863// Must be expanded from within the Thyra namespace!
3864//
3865
3866
3867#define THYRA_MODEL_EVALUATOR_BASE_INSTANT(SCALAR) \
3868 \
3869 template class ModelEvaluatorBase::InArgs<SCALAR >; \
3870 \
3871 template std::string \
3872 ModelEvaluatorBase::Derivative<SCALAR >::description() const; \
3873 \
3874 template \
3875 void ModelEvaluatorBase::Derivative<SCALAR >::describe( \
3876 Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel \
3877 ) const; \
3878 \
3879 template class ModelEvaluatorBase::OutArgs<SCALAR >; \
3880 \
3881 template class ModelEvaluatorBase::InArgsSetup<SCALAR >; \
3882 \
3883 template class ModelEvaluatorBase::OutArgsSetup<SCALAR >;
3884
3885
3886#endif // THYRA_MODEL_EVALUATOR_BASE_DEF_HPP
virtual void describe(FancyOStream &out, const EVerbosityLevel verbLevel=verbLevel_default) const
Base class for all linear operators.
Base class for all linear operators that can support a high-level solve operation.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Determines the forms of a general derivative that are supported.
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Type to embed evaluation accuracy with an RCP-managed object.
void setModelEvalDescription(const std::string &modelEvalDescription)
void setSupports(EInArgsMembers arg, bool supports=true)
Concrete aggregate class for all input arguments computable by a ModelEvaluator subclass object.
void set_stage_number(Scalar stage_number)
Precondition: supports(IN_ARG_stage_number)==true.
Scalar get_W_x_dot_dot_coeff() const
Precondition: supports(IN_ARG_W_x_dot_dot_coeff)==true.
void set_p(int l, const RCP< const VectorBase< Scalar > > &p_l)
Set p(l) where 0 <= l && l < this->Np().
int Ng() const
Return the number of axillary response functions g(j)(...) supported (Ng >= 0).
Scalar get_alpha() const
Precondition: supports(IN_ARG_alpha)==true.
RCP< const VectorBase< Scalar > > get_f_multiplier() const
Precondition: supports(IN_ARG_x)==true.
RCP< const MultiVectorBase< Scalar > > get_x_direction() const
Precondition: supports(IN_ARG_x)==true.
RCP< const VectorBase< Scalar > > get_p(int l) const
Get p(l) where 0 <= l && l < this->Np().
void set_W_x_dot_dot_coeff(Scalar W_x_dot_dot_coeff)
Precondition: supports(IN_ARG_W_x_dot_dot_coeff)==true.
void setArgs(const InArgs< Scalar > &inArgs, bool ignoreUnsupported=false, bool cloneObjects=false)
Set non-null arguments (does not overwrite non-NULLs with NULLs) .
void set_x_mp(const RCP< const Stokhos::ProductEpetraVector > &x_mp)
Precondition: supports(IN_ARG_x_mp)==true.
void set_p_direction(int l, const RCP< const MultiVectorBase< Scalar > > &p_direction_l)
Precondition: supports(IN_ARG_x)==true.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Create a more detailed description along about this object and the ModelEvaluator that created it.
RCP< const Stokhos::ProductEpetraVector > get_x_mp() const
Precondition: supports(IN_ARG_x_mp)==true.
RCP< const MultiVectorBase< Scalar > > get_p_direction(int l) const
Get p(l) where 0 <= l && l < this->Np().
void set_x_dot_mp(const RCP< const Stokhos::ProductEpetraVector > &x_dot_mp)
Precondition: supports(IN_ARG_x_dot_mp)==true.
void set_f_multiplier(const RCP< const VectorBase< Scalar > > &f_multiplier)
Precondition: supports(IN_ARG_x)==true.
void assertSameSupport(const InArgs< Scalar > &inArgs) const
Assert that two InArgs objects have the same support.
void set_x(const RCP< const VectorBase< Scalar > > &x)
Precondition: supports(IN_ARG_x)==true.
Scalar get_beta() const
Precondition: supports(IN_ARG_beta)==true.
void set_beta(Scalar beta)
Precondition: supports(IN_ARG_beta)==true.
void set_step_size(Scalar step_size)
Precondition: supports(IN_ARG_step_size)==true.
void set_x_direction(const RCP< const MultiVectorBase< Scalar > > &x_direction)
Precondition: supports(IN_ARG_x)==true.
RCP< const VectorBase< Scalar > > get_x() const
Precondition: supports(IN_ARG_x)==true.
bool supports() const
Determines if an extended input argument of type ObjectType is supported.
Teuchos::ScalarTraits< Scalar >::magnitudeType ScalarMag
.
ScalarMag get_t() const
.Precondition: supports(IN_ARG_t)==true
void set_x_dot(const RCP< const VectorBase< Scalar > > &x_dot)
Precondition: supports(IN_ARG_x_dot)==true.
Scalar get_stage_number() const
Precondition: supports(IN_ARG_stage_number)==true.
RCP< const VectorBase< Scalar > > get_g_multiplier(int j) const
Precondition: supports(IN_ARG_x)==true.
void set_alpha(Scalar alpha)
Precondition: supports(IN_ARG_alpha)==true.
void set_t(ScalarMag t)
Precondition: supports(IN_ARG_t)==true.
RCP< const VectorBase< Scalar > > get_x_dot() const
Precondition: supports(IN_ARG_x_dot)==true.
bool supports(EInArgsMembers arg) const
Determines if an input argument is supported or not.
void set_g_multiplier(int j, const RCP< const VectorBase< Scalar > > &g_multiplier)
Precondition: supports(IN_ARG_x)==true.
void set_x_dot_dot(const RCP< const VectorBase< Scalar > > &x_dot_dot)
Precondition: supports(IN_ARG_x_dot_dot)==true.
void _setSupports(EInArgsMembers arg, bool supports)
Scalar get_step_size() const
Precondition: supports(IN_ARG_step_size)==true.
void _setModelEvalDescription(const std::string &modelEvalDescription)
RCP< const Stokhos::ProductEpetraVector > get_x_dot_mp() const
Precondition: supports(IN_ARG_x_dotmp)==true.
int Np() const
Return the number of parameter subvectors p(l) supported (Np >= 0).
RCP< const VectorBase< Scalar > > get_x_dot_dot() const
Precondition: supports(IN_ARG_x_dot_dot)==true.
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
Protected subclass of OutArgs that only ModelEvaluator subclasses can access to set up the selection ...
void setSupports(EOutArgsMembers arg, bool supports=true)
void set_DgDx_dot_properties(int j, const DerivativeProperties &properties)
void set_DfDp_properties(int l, const DerivativeProperties &properties)
void set_W_properties(const DerivativeProperties &properties)
void set_DgDx_properties(int j, const DerivativeProperties &properties)
void set_DgDp_properties(int j, int l, const DerivativeProperties &properties)
void setModelEvalDescription(const std::string &modelEvalDescription)
Concrete aggregate class for all output arguments computable by a ModelEvaluator subclass object.
void _set_DfDp_properties(int l, const DerivativeProperties &properties)
Derivative< Scalar > get_DfDp(int l) const
Precondition: supports(OUT_ARG_DfDp,l)==true.
void _setSupports(EOutArgsMembers arg, bool supports)
void _set_DgDx_dot_properties(int j, const DerivativeProperties &properties)
RCP< LinearOpWithSolveBase< Scalar > > get_W() const
Precondition: supports(OUT_ARG_W)==true.
DerivativeProperties get_DfDp_properties(int l) const
Return the know properties of DfDp(l) (precondition: supports(OUT_ARG_DfDp,l)==true).
int Ng() const
Return the number of axillary response functions g(j)(...) supported (Ng >= 0).
void setArgs(const OutArgs< Scalar > &outArgs, bool ignoreUnsupported=false)
Set all arguments fron outArgs into *this.
void set_f_mp(const RCP< Stokhos::ProductEpetraVector > &f_mp)
Precondition: supports(OUT_ARG_f_mp)==true.
void set_f(const Evaluation< VectorBase< Scalar > > &f)
Precondition: supports(OUT_ARG_f)==true.
void _set_DgDx_properties(int j, const DerivativeProperties &properties)
bool supports() const
Determines if an extended output argument of type ObjectType is supported.
Derivative< Scalar > get_DgDp(int j, int l) const
Precondition: supports(OUT_ARG_DgDp,j,l)==true.
void assertSameSupport(const OutArgs< Scalar > &outArgs) const
Assert that two OutArgs objects have the same support.
Evaluation< VectorBase< Scalar > > get_f() const
Precondition: supports(OUT_ARG_f)==true.
void set_W(const RCP< LinearOpWithSolveBase< Scalar > > &W)
Precondition: supports(OUT_ARG_W)==true.
void set_DfDp(int l, const Derivative< Scalar > &DfDp_l)
Precondition: supports(OUT_ARG_DfDp,l)==true.
DerivativeProperties get_W_properties() const
Return the known properties of W (precondition: supports(OUT_ARG_f)==true).
void set_W_prec(const RCP< PreconditionerBase< Scalar > > &W_prec)
Precondition: supports(OUT_ARG_W_op)==true.
RCP< Stokhos::ProductEpetraOperator > get_W_mp() const
Precondition: supports(OUT_ARG_W_mp)==true.
void set_W_mp(const RCP< Stokhos::ProductEpetraOperator > &W_mp)
Precondition: supports(OUT_ARG_W_mp)==true.
void _set_DgDp_properties(int j, int l, const DerivativeProperties &properties)
RCP< LinearOpBase< Scalar > > get_W_op() const
Precondition: supports(OUT_ARG_W_op)==true.
void set_DgDx(int j, const Derivative< Scalar > &DgDx_j)
Precondition: supports(OUT_ARG_DgDx,j)==true.
void set_DgDx_dot(int j, const Derivative< Scalar > &DgDx_dot_j)
Precondition: supports(OUT_ARG_DgDx_dot,j)==true.
DerivativeProperties get_DgDx_dot_properties(int j) const
Return the know properties of DgDx_dot(j) (precondition: supports(OUT_ARG_DgDx_dot,...
Derivative< Scalar > get_DgDx_dot(int j) const
Precondition: supports(OUT_ARG_DgDx_dot,j)==true.
DerivativeProperties get_DgDx_properties(int j) const
Return the know properties of DgDx(j) (precondition: supports(OUT_ARG_DgDx,j)==true).
Derivative< Scalar > get_DgDx(int j) const
Precondition: supports(OUT_ARG_DgDx,j)==true.
void set_W_op(const RCP< LinearOpBase< Scalar > > &W_op)
Precondition: supports(OUT_ARG_W_op)==true.
DerivativeProperties get_DgDp_properties(int j, int l) const
Return the know properties of DgDp(j,l) (precondition: supports(OUT_ARG_DgDp,j,l)==true).
bool supports(EOutArgsMembers arg) const
Determine if an input argument is supported or not.
RCP< PreconditionerBase< Scalar > > get_W_prec() const
Precondition: supports(OUT_ARG_W_op)==true.
RCP< Stokhos::ProductEpetraVector > get_g_mp(int j) const
Precondition: supports(OUT_ARG_g_mp)==true..
void set_DgDp(int j, int l, const Derivative< Scalar > &DgDp_j_l)
Precondition: supports(OUT_ARG_DgDp,j,l)==true.
int Np() const
Return the number of parameter subvectors p(l) supported (Np >= 0).
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Create a more detailed description along about this object and the ModelEvaluator that created it.
void _setModelEvalDescription(const std::string &modelEvalDescription)
void _set_W_properties(const DerivativeProperties &properties)
void setFailed() const
Set that the evaluation as a whole failed.
RCP< Stokhos::ProductEpetraVector > get_f_mp() const
Precondition: supports(OUT_ARG_f_mp)==true.
void set_g(int j, const Evaluation< VectorBase< Scalar > > &g_j)
Precondition: supports(OUT_ARG_g)==true.
void set_g_mp(int j, const RCP< Stokhos::ProductEpetraVector > &g_mp_j)
Precondition: supports(OUT_ARG_g_mp)==true.
bool isFailed() const
Return if the evaluation failed or not.
Evaluation< VectorBase< Scalar > > get_g(int j) const
Precondition: supports(OUT_ARG_g)==true..
Base subclass for ModelEvaluator that defines some basic types.
Interface for a collection of column vectors called a multi-vector.
Simple interface class to access a precreated preconditioner as one or more linear operators objects ...
Abstract interface for finite-dimensional dense vectors.
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
T_To & dyn_cast(T_From &from)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
TEUCHOSCORE_LIB_DLL_EXPORT bool includesVerbLevel(const EVerbosityLevel verbLevel, const EVerbosityLevel requestedVerbLevel, const bool isDefaultLevel=false)
Simple public strict containing properties of a derivative object.