287 virtual void evalModelImpl(
304 typedef ModelEvaluatorDefaultBaseTypes::DefaultDerivLinearOpSupport
305 DefaultDerivLinearOpSupport;
307 typedef ModelEvaluatorDefaultBaseTypes::DefaultDerivMvAdjointSupport
308 DefaultDerivMvAdjointSupport;
310 typedef ModelEvaluatorDefaultBaseTypes::MultiVectorAdjointPair<Scalar>
311 MultiVectorAdjointPair;
327 bool default_W_support_;
334 void lazyInitializeDefaultBase()
const;
336 void assert_l(
const int l)
const;
338 void assert_j(
const int j)
const;
343 static DefaultDerivLinearOpSupport
344 determineDefaultDerivLinearOpSupport(
349 createDefaultLinearOp(
350 const DefaultDerivLinearOpSupport &defaultLinearOpSupport,
356 updateDefaultLinearOpSupport(
358 const DefaultDerivLinearOpSupport &defaultLinearOpSupport
362 getOutArgImplForDefaultLinearOpSupport(
364 const DefaultDerivLinearOpSupport &defaultLinearOpSupport
367 static DefaultDerivMvAdjointSupport
368 determineDefaultDerivMvAdjointSupport(
375 updateDefaultDerivMvAdjointSupport(
377 const DefaultDerivMvAdjointSupport &defaultMvAdjointSupport
391#include "Thyra_ModelEvaluatorHelpers.hpp"
392#include "Thyra_DefaultScaledAdjointLinearOp.hpp"
393#include "Thyra_DetachedMultiVectorView.hpp"
394#include "Teuchos_Assert.hpp"
403template<
class Scalar>
406 lazyInitializeDefaultBase();
407 return prototypeOutArgs_.Np();
411template<
class Scalar>
414 lazyInitializeDefaultBase();
415 return prototypeOutArgs_.Ng();
419template<
class Scalar>
423 lazyInitializeDefaultBase();
424 if (default_W_support_)
425 return this->get_W_factory()->createOp();
430template<
class Scalar>
434 lazyInitializeDefaultBase();
438 const DefaultDerivLinearOpSupport
439 defaultLinearOpSupport = DfDp_default_op_support_[l];
440 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
441 return createDefaultLinearOp(
442 defaultLinearOpSupport,
447 return this->create_DfDp_op_impl(l);
451template<
class Scalar>
455 lazyInitializeDefaultBase();
459 const DefaultDerivLinearOpSupport
460 defaultLinearOpSupport = DgDx_dot_default_op_support_[j];
461 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
462 return createDefaultLinearOp(
463 defaultLinearOpSupport,
464 this->get_g_space(j),
468 return this->create_DgDx_dot_op_impl(j);
472template<
class Scalar>
476 lazyInitializeDefaultBase();
480 const DefaultDerivLinearOpSupport
481 defaultLinearOpSupport = DgDx_default_op_support_[j];
482 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
483 return createDefaultLinearOp(
484 defaultLinearOpSupport,
485 this->get_g_space(j),
489 return this->create_DgDx_op_impl(j);
493template<
class Scalar>
497 lazyInitializeDefaultBase();
502 const DefaultDerivLinearOpSupport
503 defaultLinearOpSupport = DgDp_default_op_support_[j][l];
504 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
505 return createDefaultLinearOp(
506 defaultLinearOpSupport,
507 this->get_g_space(j),
511 return this->create_DgDp_op_impl(j,l);
515template<
class Scalar>
519 lazyInitializeDefaultBase();
520 return prototypeOutArgs_;
524template<
class Scalar>
531 using Teuchos::outArg;
534 lazyInitializeDefaultBase();
536 const int l_Np = outArgs.
Np();
537 const int l_Ng = outArgs.
Ng();
544 assertInArgsEvalObjects(*
this,inArgs);
545 assertOutArgsEvalObjects(*
this,outArgs,&inArgs);
553 MEB::OutArgs<Scalar> outArgsImpl = this->createOutArgsImpl();
558 outArgsImpl.setArgs(outArgs,
true);
561 if (outArgsImpl.supports(MEB::OUT_ARG_f)) {
562 for (
int l = 0; l < l_Np; ++l ) {
563 const DefaultDerivLinearOpSupport defaultLinearOpSupport =
564 DfDp_default_op_support_[l];
565 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
566 outArgsImpl.set_DfDp( l,
567 getOutArgImplForDefaultLinearOpSupport(
568 outArgs.
get_DfDp(l), defaultLinearOpSupport
579 for (
int j = 0; j < l_Ng; ++j ) {
580 const DefaultDerivLinearOpSupport defaultLinearOpSupport =
581 DgDx_dot_default_op_support_[j];
582 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
583 outArgsImpl.set_DgDx_dot( j,
584 getOutArgImplForDefaultLinearOpSupport(
595 for (
int j = 0; j < l_Ng; ++j ) {
596 const DefaultDerivLinearOpSupport defaultLinearOpSupport =
597 DgDx_default_op_support_[j];
598 if (defaultLinearOpSupport.provideDefaultLinearOp()) {
599 outArgsImpl.set_DgDx( j,
600 getOutArgImplForDefaultLinearOpSupport(
601 outArgs.
get_DgDx(j), defaultLinearOpSupport
611 for (
int j = 0; j < l_Ng; ++j ) {
613 DgDp_default_op_support_[j];
615 DgDp_default_mv_support_[j];
616 for (
int l = 0; l < l_Np; ++l ) {
617 const DefaultDerivLinearOpSupport defaultLinearOpSupport =
618 DgDp_default_op_support_j[l];
619 const DefaultDerivMvAdjointSupport defaultMvAdjointSupport =
620 DgDp_default_mv_support_j[l];
621 MEB::Derivative<Scalar> DgDp_j_l;
622 if (!outArgs.
supports(MEB::OUT_ARG_DgDp,j,l).none())
625 defaultLinearOpSupport.provideDefaultLinearOp()
626 && !is_null(DgDp_j_l.getLinearOp())
629 outArgsImpl.set_DgDp( j, l,
630 getOutArgImplForDefaultLinearOpSupport(
631 DgDp_j_l, defaultLinearOpSupport
636 defaultMvAdjointSupport.provideDefaultAdjoint()
637 && !is_null(DgDp_j_l.getMultiVector())
641 DgDp_j_l.getMultiVector();
643 defaultMvAdjointSupport.mvAdjointCopyOrientation()
645 DgDp_j_l.getMultiVectorOrientation()
652 createMembers(DgDp_j_l_mv->domain(), DgDp_j_l_mv->range()->dim());
653 outArgsImpl.set_DgDp( j, l,
654 MEB::Derivative<Scalar>(
656 getOtherDerivativeMultiVectorOrientation(
657 defaultMvAdjointSupport.mvAdjointCopyOrientation()
664 MultiVectorAdjointPair(DgDp_j_l_mv, DgDp_j_l_mv_adj)
681 if ( default_W_support_ && !is_null(W=outArgs.
get_W()) ) {
683 W_factory = this->get_W_factory();
686 uninitializeOp<Scalar>(*W_factory, W.
ptr(), outArg(W_op_const));
688 if (!is_null(W_op_const)) {
698 W_op = this->create_W_op();
700 outArgsImpl.set_W_op(W_op);
709 this->evalModelImpl( inArgs, outArgsImpl );
716 const int numMvAdjointCopies = DgDp_temp_adjoint_copies.
size();
717 for (
int adj_copy_i = 0; adj_copy_i < numMvAdjointCopies; ++adj_copy_i ) {
718 const MultiVectorAdjointPair adjPair =
719 DgDp_temp_adjoint_copies[adj_copy_i];
720 doExplicitMultiVectorAdjoint( *adjPair.mvImplAdjoint, &*adjPair.mvOuter );
726 if ( default_W_support_ && !is_null(W=outArgs.
get_W()) ) {
728 W_factory = this->get_W_factory();
729 W_factory->setOStream(this->getOStream());
730 W_factory->setVerbLevel(this->getVerbLevel());
731 initializeOp<Scalar>(*W_factory, outArgsImpl.get_W_op().
getConst(), W.
ptr());
743template<
class Scalar>
750 isInitialized_ =
false;
751 default_W_support_ =
false;
757 const MEB::InArgs<Scalar> inArgs = this->createInArgs();
758 const MEB::OutArgs<Scalar> outArgsImpl = this->createOutArgsImpl();
765 assertInArgsOutArgsSetup( this->description(), inArgs, outArgsImpl );
772 const int l_Ng = outArgsImpl.Ng();
773 const int l_Np = outArgsImpl.Np();
776 MEB::OutArgsSetup<Scalar> outArgs;
777 outArgs.setModelEvalDescription(this->description());
778 outArgs.set_Np_Ng(l_Np,l_Ng);
779 outArgs.setSupports(outArgsImpl);
782 DfDp_default_op_support_.clear();
783 if (outArgs.supports(MEB::OUT_ARG_f)) {
784 for (
int l = 0; l < l_Np; ++l ) {
785 const MEB::DerivativeSupport DfDp_l_impl_support =
786 outArgsImpl.supports(MEB::OUT_ARG_DfDp,l);
787 const DefaultDerivLinearOpSupport DfDp_l_op_support =
788 determineDefaultDerivLinearOpSupport(DfDp_l_impl_support);
789 DfDp_default_op_support_.push_back(DfDp_l_op_support);
791 MEB::OUT_ARG_DfDp, l,
792 updateDefaultLinearOpSupport(
793 DfDp_l_impl_support, DfDp_l_op_support
800 DgDx_dot_default_op_support_.clear();
801 for (
int j = 0; j < l_Ng; ++j ) {
802 const MEB::DerivativeSupport DgDx_dot_j_impl_support =
803 outArgsImpl.supports(MEB::OUT_ARG_DgDx_dot,j);
804 const DefaultDerivLinearOpSupport DgDx_dot_j_op_support =
805 determineDefaultDerivLinearOpSupport(DgDx_dot_j_impl_support);
806 DgDx_dot_default_op_support_.push_back(DgDx_dot_j_op_support);
808 MEB::OUT_ARG_DgDx_dot, j,
809 updateDefaultLinearOpSupport(
810 DgDx_dot_j_impl_support, DgDx_dot_j_op_support
816 DgDx_default_op_support_.clear();
817 for (
int j = 0; j < l_Ng; ++j ) {
818 const MEB::DerivativeSupport DgDx_j_impl_support =
819 outArgsImpl.supports(MEB::OUT_ARG_DgDx,j);
820 const DefaultDerivLinearOpSupport DgDx_j_op_support =
821 determineDefaultDerivLinearOpSupport(DgDx_j_impl_support);
822 DgDx_default_op_support_.push_back(DgDx_j_op_support);
824 MEB::OUT_ARG_DgDx, j,
825 updateDefaultLinearOpSupport(
826 DgDx_j_impl_support, DgDx_j_op_support
832 DgDp_default_op_support_.clear();
833 DgDp_default_mv_support_.clear();
834 for (
int j = 0; j < l_Ng; ++j ) {
837 for (
int l = 0; l < l_Np; ++l ) {
838 const MEB::DerivativeSupport DgDp_j_l_impl_support =
839 outArgsImpl.supports(MEB::OUT_ARG_DgDp,j,l);
841 const DefaultDerivLinearOpSupport DgDp_j_l_op_support =
842 determineDefaultDerivLinearOpSupport(DgDp_j_l_impl_support);
843 DgDp_default_op_support_[j].push_back(DgDp_j_l_op_support);
845 MEB::OUT_ARG_DgDp, j, l,
846 updateDefaultLinearOpSupport(
847 DgDp_j_l_impl_support, DgDp_j_l_op_support
851 const DefaultDerivMvAdjointSupport DgDp_j_l_mv_support =
852 determineDefaultDerivMvAdjointSupport(
853 DgDp_j_l_impl_support, *this->get_g_space(j), *this->get_p_space(l)
855 DgDp_default_mv_support_[j].push_back(DgDp_j_l_mv_support);
857 MEB::OUT_ARG_DgDp, j, l,
858 updateDefaultDerivMvAdjointSupport(
859 outArgs.supports(MEB::OUT_ARG_DgDp, j, l),
869 default_W_support_ =
false;
870 if ( outArgsImpl.supports(MEB::OUT_ARG_W_op) && !is_null(this->get_W_factory())
871 && !outArgsImpl.supports(MEB::OUT_ARG_W) )
873 default_W_support_ =
true;
874 outArgs.setSupports(MEB::OUT_ARG_W);
875 outArgs.set_W_properties(outArgsImpl.get_W_properties());
882 prototypeOutArgs_ = outArgs;
883 isInitialized_ =
true;
887template<
class Scalar>
890 isInitialized_ =
false;
896template<
class Scalar>
901 MEB::OutArgs<Scalar> outArgs = this->createOutArgsImpl();
903 outArgs.supports(MEB::OUT_ARG_DfDp,l).supports(MEB::DERIV_LINEAR_OP),
905 "Error, The ModelEvaluator subclass "<<this->description()<<
" says that it"
906 " supports the LinearOpBase form of DfDp("<<l<<
") (as determined from its"
907 " OutArgs object created by createOutArgsImpl())"
908 " but this function create_DfDp_op_impl(...) has not been overridden"
909 " to create such an object!"
915template<
class Scalar>
916RCP<LinearOpBase<Scalar> >
917ModelEvaluatorDefaultBase<Scalar>::create_DgDx_dot_op_impl(
int j)
const
919 typedef ModelEvaluatorBase MEB;
920 MEB::OutArgs<Scalar> outArgs = this->createOutArgsImpl();
922 outArgs.supports(MEB::OUT_ARG_DgDx_dot,j).supports(MEB::DERIV_LINEAR_OP),
924 "Error, The ModelEvaluator subclass "<<this->description()<<
" says that it"
925 " supports the LinearOpBase form of DgDx_dot("<<j<<
") (as determined from"
926 " its OutArgs object created by createOutArgsImpl())"
927 " but this function create_DgDx_dot_op_impl(...) has not been overridden"
928 " to create such an object!"
934template<
class Scalar>
935RCP<LinearOpBase<Scalar> >
936ModelEvaluatorDefaultBase<Scalar>::create_DgDx_op_impl(
int j)
const
938 typedef ModelEvaluatorBase MEB;
939 MEB::OutArgs<Scalar> outArgs = this->createOutArgsImpl();
941 outArgs.supports(MEB::OUT_ARG_DgDx,j).supports(MEB::DERIV_LINEAR_OP),
943 "Error, The ModelEvaluator subclass "<<this->description()<<
" says that it"
944 " supports the LinearOpBase form of DgDx("<<j<<
") (as determined from"
945 " its OutArgs object created by createOutArgsImpl())"
946 " but this function create_DgDx_op_impl(...) has not been overridden"
947 " to create such an object!"
953template<
class Scalar>
954RCP<LinearOpBase<Scalar> >
955ModelEvaluatorDefaultBase<Scalar>::create_DgDp_op_impl(
int j,
int l)
const
957 typedef ModelEvaluatorBase MEB;
958 MEB::OutArgs<Scalar> outArgs = this->createOutArgsImpl();
960 outArgs.supports(MEB::OUT_ARG_DgDp,j,l).supports(MEB::DERIV_LINEAR_OP),
962 "Error, The ModelEvaluator subclass "<<this->description()<<
" says that it"
963 " supports the LinearOpBase form of DgDp("<<j<<
","<<l<<
")"
964 " (as determined from its OutArgs object created by createOutArgsImpl())"
965 " but this function create_DgDp_op_impl(...) has not been overridden"
966 " to create such an object!"
972template<
class Scalar>
973RCP<const VectorSpaceBase<Scalar> >
976 return this->get_f_space();
979template<
class Scalar>
983 return this->get_g_space(j);
986#ifdef Thyra_BUILD_HESSIAN_SUPPORT
988template<
class Scalar>
995template<
class Scalar>
996RCP<LinearOpBase<Scalar> >
997ModelEvaluatorDefaultBase<Scalar>::create_hess_f_xp(
int l)
const
1002template<
class Scalar>
1003RCP<LinearOpBase<Scalar> >
1004ModelEvaluatorDefaultBase<Scalar>::create_hess_f_pp(
int l1,
int l2 )
const
1009template<
class Scalar>
1010RCP<LinearOpBase<Scalar> >
1011ModelEvaluatorDefaultBase<Scalar>::create_hess_g_xx(
int j)
const
1016template<
class Scalar>
1017RCP<LinearOpBase<Scalar> >
1018ModelEvaluatorDefaultBase<Scalar>::create_hess_g_xp(
int j,
int l )
const
1023template<
class Scalar>
1024RCP<LinearOpBase<Scalar> >
1025ModelEvaluatorDefaultBase<Scalar>::create_hess_g_pp(
int j,
int l1,
int l2 )
const
1035template<
class Scalar>
1037 :isInitialized_(false), default_W_support_(false)
1044template<
class Scalar>
1047 if (!isInitialized_)
1052template<
class Scalar>
1053void ModelEvaluatorDefaultBase<Scalar>::assert_l(
const int l)
const
1059template<
class Scalar>
1060void ModelEvaluatorDefaultBase<Scalar>::assert_j(
const int j)
const
1069template<
class Scalar>
1070ModelEvaluatorDefaultBaseTypes::DefaultDerivLinearOpSupport
1071ModelEvaluatorDefaultBase<Scalar>::determineDefaultDerivLinearOpSupport(
1072 const ModelEvaluatorBase::DerivativeSupport &derivSupportImpl
1075 typedef ModelEvaluatorBase MEB;
1078 derivSupportImpl.supports(MEB::DERIV_MV_BY_COL)
1080 derivSupportImpl.supports(MEB::DERIV_TRANS_MV_BY_ROW)
1083 !derivSupportImpl.supports(MEB::DERIV_LINEAR_OP)
1086 return DefaultDerivLinearOpSupport(
1087 derivSupportImpl.supports(MEB::DERIV_MV_BY_COL)
1088 ? MEB::DERIV_MV_BY_COL
1089 : MEB::DERIV_TRANS_MV_BY_ROW
1092 return DefaultDerivLinearOpSupport();
1096template<
class Scalar>
1097RCP<LinearOpBase<Scalar> >
1098ModelEvaluatorDefaultBase<Scalar>::createDefaultLinearOp(
1099 const DefaultDerivLinearOpSupport &defaultLinearOpSupport,
1100 const RCP<
const VectorSpaceBase<Scalar> > &fnc_space,
1101 const RCP<
const VectorSpaceBase<Scalar> > &var_space
1104 using Teuchos::rcp_implicit_cast;
1105 typedef LinearOpBase<Scalar> LOB;
1106 typedef ModelEvaluatorBase MEB;
1107 switch(defaultLinearOpSupport.mvImplOrientation()) {
1108 case MEB::DERIV_MV_BY_COL:
1110 return createMembers(fnc_space, var_space->dim());
1111 case MEB::DERIV_TRANS_MV_BY_ROW:
1113 return nonconstAdjoint<Scalar>(
1114 rcp_implicit_cast<LOB>(createMembers(var_space, fnc_space->dim()))
1125template<
class Scalar>
1126ModelEvaluatorBase::DerivativeSupport
1127ModelEvaluatorDefaultBase<Scalar>::updateDefaultLinearOpSupport(
1128 const ModelEvaluatorBase::DerivativeSupport &derivSupportImpl,
1129 const DefaultDerivLinearOpSupport &defaultLinearOpSupport
1132 typedef ModelEvaluatorBase MEB;
1133 MEB::DerivativeSupport derivSupport = derivSupportImpl;
1134 if (defaultLinearOpSupport.provideDefaultLinearOp())
1135 derivSupport.plus(MEB::DERIV_LINEAR_OP);
1136 return derivSupport;
1140template<
class Scalar>
1141ModelEvaluatorBase::Derivative<Scalar>
1142ModelEvaluatorDefaultBase<Scalar>::getOutArgImplForDefaultLinearOpSupport(
1143 const ModelEvaluatorBase::Derivative<Scalar> &deriv,
1144 const DefaultDerivLinearOpSupport &defaultLinearOpSupport
1148 using Teuchos::rcp_dynamic_cast;
1149 typedef ModelEvaluatorBase MEB;
1150 typedef MultiVectorBase<Scalar> MVB;
1151 typedef ScaledAdjointLinearOpBase<Scalar> SALOB;
1156 if (
is_null(deriv.getLinearOp()))
1161 switch(defaultLinearOpSupport.mvImplOrientation()) {
1162 case MEB::DERIV_MV_BY_COL: {
1163 return MEB::Derivative<Scalar>(
1164 rcp_dynamic_cast<MVB>(deriv.getLinearOp(),
true),
1165 MEB::DERIV_MV_BY_COL
1168 case MEB::DERIV_TRANS_MV_BY_ROW: {
1169 return MEB::Derivative<Scalar>(
1170 rcp_dynamic_cast<MVB>(
1171 rcp_dynamic_cast<SALOB>(deriv.getLinearOp(),
true)->getNonconstOrigOp()
1173 MEB::DERIV_TRANS_MV_BY_ROW
1182 return ModelEvaluatorBase::Derivative<Scalar>();
1187template<
class Scalar>
1188ModelEvaluatorDefaultBaseTypes::DefaultDerivMvAdjointSupport
1189ModelEvaluatorDefaultBase<Scalar>::determineDefaultDerivMvAdjointSupport(
1190 const ModelEvaluatorBase::DerivativeSupport &derivSupportImpl,
1191 const VectorSpaceBase<Scalar> &fnc_space,
1192 const VectorSpaceBase<Scalar> &var_space
1195 typedef ModelEvaluatorBase MEB;
1198 const bool implSupportsMv =
1199 ( derivSupportImpl.supports(MEB::DERIV_MV_BY_COL)
1200 || derivSupportImpl.supports(MEB::DERIV_TRANS_MV_BY_ROW) );
1201 const bool implLacksMvOrientSupport =
1202 ( !derivSupportImpl.supports(MEB::DERIV_MV_BY_COL)
1203 || !derivSupportImpl.supports(MEB::DERIV_TRANS_MV_BY_ROW) );
1204 const bool bothSpacesHaveInCoreViews =
1205 ( fnc_space.hasInCoreView() && var_space.hasInCoreView() );
1206 if ( implSupportsMv && implLacksMvOrientSupport && bothSpacesHaveInCoreViews ) {
1207 return DefaultDerivMvAdjointSupport(
1208 derivSupportImpl.supports(MEB::DERIV_MV_BY_COL)
1209 ? MEB::DERIV_TRANS_MV_BY_ROW
1210 : MEB::DERIV_MV_BY_COL
1214 return DefaultDerivMvAdjointSupport();
1218template<
class Scalar>
1219ModelEvaluatorBase::DerivativeSupport
1220ModelEvaluatorDefaultBase<Scalar>::updateDefaultDerivMvAdjointSupport(
1221 const ModelEvaluatorBase::DerivativeSupport &derivSupportImpl,
1222 const DefaultDerivMvAdjointSupport &defaultMvAdjointSupport
1225 typedef ModelEvaluatorBase MEB;
1226 MEB::DerivativeSupport derivSupport = derivSupportImpl;
1227 if (defaultMvAdjointSupport.provideDefaultAdjoint())
1228 derivSupport.plus(defaultMvAdjointSupport.mvAdjointCopyOrientation());
1229 return derivSupport;