43 #ifndef __Panzer_ModelEvaluator_impl_hpp__ 44 #define __Panzer_ModelEvaluator_impl_hpp__ 46 #include "Teuchos_DefaultComm.hpp" 47 #include "Teuchos_ArrayRCP.hpp" 49 #include "PanzerDiscFE_config.hpp" 67 #include "Thyra_TpetraThyraWrappers.hpp" 68 #include "Thyra_SpmdVectorBase.hpp" 69 #include "Thyra_DefaultSpmdVector.hpp" 70 #include "Thyra_DefaultSpmdVectorSpace.hpp" 71 #include "Thyra_DefaultMultiVectorProductVectorSpace.hpp" 72 #include "Thyra_DefaultMultiVectorProductVector.hpp" 73 #include "Thyra_MultiVectorStdOps.hpp" 74 #include "Thyra_VectorStdOps.hpp" 76 #include "Tpetra_CrsMatrix.hpp" 80 template<
typename Scalar>
85 const std::vector<Teuchos::RCP<Teuchos::Array<std::string> > >& p_names,
86 const std::vector<Teuchos::RCP<Teuchos::Array<double> > >& p_values,
88 const Teuchos::RCP<panzer::GlobalData>& global_data,
89 bool build_transient_support,
92 , num_me_parameters_(0)
94 , fd_perturb_size_(1e-7)
95 , require_in_args_refresh_(true)
96 , require_out_args_refresh_(true)
97 , responseLibrary_(rLibrary)
98 , global_data_(global_data)
99 , build_transient_support_(build_transient_support)
101 , solverFactory_(solverFactory)
102 , oneTimeDirichletBeta_on_(false)
103 , oneTimeDirichletBeta_(0.0)
104 , build_volume_field_managers_(true)
105 , build_bc_field_managers_(true)
106 , active_evaluation_types_(Sacado::mpl::size<
panzer::
Traits::EvalTypes>::value, true)
110 using Teuchos::rcp_dynamic_cast;
111 using Teuchos::tuple;
113 using Thyra::createMember;
115 TEUCHOS_ASSERT(
lof_!=Teuchos::null);
118 ae_tm_.buildObjects(builder);
127 x_space_ = tof->getThyraDomainSpace();
128 f_space_ = tof->getThyraRangeSpace();
133 for(std::size_t i=0;i<p_names.size();i++)
141 template<
typename Scalar>
145 const Teuchos::RCP<panzer::GlobalData>& global_data,
146 bool build_transient_support,
double t_init)
148 , num_me_parameters_(0)
150 , fd_perturb_size_(1e-7)
151 , require_in_args_refresh_(true)
152 , require_out_args_refresh_(true)
153 , global_data_(global_data)
154 , build_transient_support_(build_transient_support)
156 , solverFactory_(solverFactory)
157 , oneTimeDirichletBeta_on_(false)
158 , oneTimeDirichletBeta_(0.0)
159 , build_volume_field_managers_(true)
160 , build_bc_field_managers_(true)
161 , active_evaluation_types_(Sacado::mpl::size<
panzer::
Traits::EvalTypes>::value, true)
164 using Teuchos::rcp_dynamic_cast;
166 TEUCHOS_ASSERT(
lof_!=Teuchos::null);
175 x_space_ = tof->getThyraDomainSpace();
176 f_space_ = tof->getThyraRangeSpace();
186 template<
typename Scalar>
190 TEUCHOS_ASSERT(
false);
195 template<
typename Scalar>
196 Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar> >
203 template<
typename Scalar>
204 Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar> >
210 template<
typename Scalar>
211 Teuchos::RCP<const Teuchos::Array<std::string> >
214 TEUCHOS_TEST_FOR_EXCEPTION(!(i>=0 && i<num_me_parameters_),std::runtime_error,
215 "panzer::ModelEvaluator::get_p_names: Requested parameter index out of range.");
217 if (i < Teuchos::as<int>(parameters_.size()))
218 return parameters_[i]->names;
219 else if (i < Teuchos::as<int>(parameters_.size()+tangent_space_.size())) {
220 Teuchos::RCP<Teuchos::Array<std::string> > names = rcp(
new Teuchos::Array<std::string>);
221 int param_index = i-parameters_.size();
222 std::ostringstream ss;
223 ss <<
"TANGENT VECTOR: " << param_index;
224 names->push_back(ss.str());
227 else if (build_transient_support_ && i < Teuchos::as<int>(parameters_.size()+2*tangent_space_.size())) {
228 Teuchos::RCP<Teuchos::Array<std::string> > names = rcp(
new Teuchos::Array<std::string>);
229 int param_index = i-parameters_.size()-tangent_space_.size();
230 std::ostringstream ss;
231 ss <<
"TIME DERIVATIVE TANGENT VECTOR: " << param_index;
232 names->push_back(ss.str());
236 return Teuchos::null;
239 template<
typename Scalar>
240 Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar> >
243 TEUCHOS_TEST_FOR_EXCEPTION(!(i>=0 && i<num_me_parameters_),std::runtime_error,
244 "panzer::ModelEvaluator::get_p_space: Requested parameter index out of range.");
246 if (i < Teuchos::as<int>(parameters_.size()))
247 return parameters_[i]->space;
248 else if (i < Teuchos::as<int>(parameters_.size()+tangent_space_.size()))
249 return tangent_space_[i-parameters_.size()];
250 else if (build_transient_support_ && i < Teuchos::as<int>(parameters_.size()+2*tangent_space_.size()))
251 return tangent_space_[i-parameters_.size()-tangent_space_.size()];
253 return Teuchos::null;
256 template<
typename Scalar>
257 Teuchos::ArrayView<const std::string>
260 TEUCHOS_TEST_FOR_EXCEPTION(!(i>=0 && i<Teuchos::as<int>(responses_.size())),std::runtime_error,
261 "panzer::ModelEvaluator::get_g_names: Requested response index out of range.");
263 return Teuchos::ArrayView<const std::string>(&(responses_[i]->name),1);
266 template<
typename Scalar>
270 TEUCHOS_ASSERT(i>=0 &&
271 static_cast<typename std::vector<Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar>
> >::size_type>(i)<responses_.size());
273 return responses_[i]->name;
276 template<
typename Scalar>
277 Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar> >
280 TEUCHOS_ASSERT(i>=0 &&
281 static_cast<typename std::vector<Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar>
> >::size_type>(i)<responses_.size());
283 return responses_[i]->space;
286 template<
typename Scalar>
287 Thyra::ModelEvaluatorBase::InArgs<Scalar>
290 return getNominalValues();
293 template<
typename Scalar>
294 Thyra::ModelEvaluatorBase::InArgs<Scalar>
298 using Teuchos::rcp_dynamic_cast;
300 if(require_in_args_refresh_) {
301 typedef Thyra::ModelEvaluatorBase MEB;
308 MEB::InArgsSetup<Scalar> nomInArgs;
309 nomInArgs = nominalValues_;
310 nomInArgs.setSupports(nominalValues_);
313 nomInArgs.set_Np(num_me_parameters_);
314 for(std::size_t p=0;p<parameters_.size();p++) {
316 nomInArgs.set_p(p,parameters_[p]->initial_value);
322 nominalValues_ = nomInArgs;
326 require_in_args_refresh_ =
false;
328 return nominalValues_;
331 template<
typename Scalar>
335 typedef Thyra::ModelEvaluatorBase MEB;
341 MEB::InArgsSetup<Scalar> nomInArgs;
342 nomInArgs.setModelEvalDescription(this->description());
343 nomInArgs.setSupports(MEB::IN_ARG_x);
344 Teuchos::RCP<Thyra::VectorBase<Scalar> > x_nom = Thyra::createMember(x_space_);
345 Thyra::assign(x_nom.ptr(),0.0);
346 nomInArgs.set_x(x_nom);
347 if(build_transient_support_) {
348 nomInArgs.setSupports(MEB::IN_ARG_x_dot,
true);
349 nomInArgs.setSupports(MEB::IN_ARG_t,
true);
350 nomInArgs.setSupports(MEB::IN_ARG_alpha,
true);
351 nomInArgs.setSupports(MEB::IN_ARG_beta,
true);
352 nomInArgs.setSupports(MEB::IN_ARG_step_size,
true);
353 nomInArgs.setSupports(MEB::IN_ARG_stage_number,
true);
355 Teuchos::RCP<Thyra::VectorBase<Scalar> > x_dot_nom = Thyra::createMember(x_space_);
356 Thyra::assign(x_dot_nom.ptr(),0.0);
357 nomInArgs.set_x_dot(x_dot_nom);
358 nomInArgs.set_t(t_init_);
359 nomInArgs.set_alpha(0.0);
360 nomInArgs.set_beta(0.0);
362 nomInArgs.set_step_size(0.0);
363 nomInArgs.set_stage_number(1.0);
367 nomInArgs.set_Np(num_me_parameters_);
368 std::size_t v_index = 0;
369 for(std::size_t p=0;p<parameters_.size();p++) {
370 nomInArgs.set_p(p,parameters_[p]->initial_value);
371 if (!parameters_[p]->is_distributed) {
372 Teuchos::RCP<Thyra::VectorBase<Scalar> > v_nom_x = Thyra::createMember(*tangent_space_[v_index]);
373 Thyra::assign(v_nom_x.ptr(),0.0);
374 nomInArgs.set_p(v_index+parameters_.size(),v_nom_x);
375 if (build_transient_support_) {
376 Teuchos::RCP<Thyra::VectorBase<Scalar> > v_nom_xdot = Thyra::createMember(*tangent_space_[v_index]);
377 Thyra::assign(v_nom_xdot.ptr(),0.0);
378 nomInArgs.set_p(v_index+parameters_.size()+tangent_space_.size(),v_nom_xdot);
384 nominalValues_ = nomInArgs;
387 template <
typename Scalar>
391 build_volume_field_managers_ = value;
394 template <
typename Scalar>
398 build_bc_field_managers_ = value;
401 template <
typename Scalar>
403 setupModel(
const Teuchos::RCP<panzer::WorksetContainer> & wc,
404 const std::vector<Teuchos::RCP<panzer::PhysicsBlock> >& physicsBlocks,
405 const std::vector<panzer::BC> & bcs,
410 const Teuchos::ParameterList& closure_models,
411 const Teuchos::ParameterList& user_data,
412 bool writeGraph,
const std::string & graphPrefix,
413 const Teuchos::ParameterList& me_params)
417 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::ModelEvaluator::setupModel()",setupModel);
423 Teuchos::RCP<panzer::FieldManagerBuilder> fmb;
425 PANZER_FUNC_TIME_MONITOR_DIFF(
"allocate FieldManagerBuilder",allocFMB);
427 fmb->setActiveEvaluationTypes(active_evaluation_types_);
430 PANZER_FUNC_TIME_MONITOR_DIFF(
"fmb->setWorksetContainer()",setupWorksets);
431 fmb->setWorksetContainer(wc);
433 if (build_volume_field_managers_) {
434 PANZER_FUNC_TIME_MONITOR_DIFF(
"fmb->setupVolumeFieldManagers()",setupVolumeFieldManagers);
435 fmb->setupVolumeFieldManagers(physicsBlocks,volume_cm_factory,closure_models,*lof_,user_data);
437 if (build_bc_field_managers_) {
438 PANZER_FUNC_TIME_MONITOR_DIFF(
"fmb->setupBCFieldManagers()",setupBCFieldManagers);
439 fmb->setupBCFieldManagers(bcs,physicsBlocks,eqset_factory,bc_cm_factory,bc_factory,closure_models,*lof_,user_data);
444 if (build_volume_field_managers_)
445 fmb->writeVolumeGraphvizDependencyFiles(graphPrefix, physicsBlocks);
446 if (build_bc_field_managers_)
447 fmb->writeBCGraphvizDependencyFiles(graphPrefix+
"BC_");
451 PANZER_FUNC_TIME_MONITOR_DIFF(
"AssemblyEngine_TemplateBuilder::buildObjects()",AETM_BuildObjects);
453 ae_tm_.buildObjects(builder);
461 PANZER_FUNC_TIME_MONITOR_DIFF(
"build response library",buildResponses);
463 responseLibrary_->initialize(wc,lof_->getRangeGlobalIndexer(),lof_);
465 buildResponses(physicsBlocks,eqset_factory,volume_cm_factory,closure_models,user_data,writeGraph,graphPrefix+
"Responses_");
466 buildDistroParamDfDp_RL(wc,physicsBlocks,bcs,eqset_factory,bc_factory,volume_cm_factory,closure_models,user_data,writeGraph,graphPrefix+
"Response_DfDp_");
467 buildDistroParamDgDp_RL(wc,physicsBlocks,bcs,eqset_factory,bc_factory,volume_cm_factory,closure_models,user_data,writeGraph,graphPrefix+
"Response_DgDp_");
470 fd_perturb_size_ = 1.0e-7;
471 if (me_params.isParameter(
"FD Forward Sensitivities"))
472 do_fd_dfdp_ = me_params.get<
bool>(
"FD Forward Sensitivities");
473 if (me_params.isParameter(
"FD Perturbation Size"))
474 fd_perturb_size_ = me_params.get<
double>(
"FD Perturbation Size");
478 template <
typename Scalar>
485 using Teuchos::rcp_dynamic_cast;
486 using Teuchos::rcp_const_cast;
487 typedef Thyra::ModelEvaluatorBase MEB;
490 if(Teuchos::is_null(ghostedContainer_)) {
491 ghostedContainer_ = lof_->buildGhostedLinearObjContainer();
498 bool is_transient =
false;
499 if (inArgs.supports(MEB::IN_ARG_x_dot ))
500 is_transient = !Teuchos::is_null(inArgs.get_x_dot());
502 if(Teuchos::is_null(xContainer_))
503 xContainer_ = lof_->buildReadOnlyDomainContainer();
504 if(Teuchos::is_null(xdotContainer_) && is_transient)
505 xdotContainer_ = lof_->buildReadOnlyDomainContainer();
507 const RCP<const Thyra::VectorBase<Scalar> > x = inArgs.get_x();
508 RCP<const Thyra::VectorBase<Scalar> > x_dot;
511 TEUCHOS_TEST_FOR_EXCEPTION(is_transient && !build_transient_support_, std::runtime_error,
512 "ModelEvaluator was not built with transient support enabled!");
514 ae_inargs.
container_ = lof_->buildLinearObjContainer();
516 ae_inargs.
alpha = 0.0;
517 ae_inargs.
beta = 1.0;
519 if (build_transient_support_) {
520 x_dot = inArgs.get_x_dot();
521 ae_inargs.
alpha = inArgs.get_alpha();
522 ae_inargs.
beta = inArgs.get_beta();
523 ae_inargs.
time = inArgs.get_t();
525 ae_inargs.
step_size= inArgs.get_step_size();
534 int num_param_vecs = parameters_.size();
535 for (
int i=0; i<num_param_vecs; i++) {
537 RCP<const Thyra::VectorBase<Scalar> > paramVec = inArgs.get_p(i);
538 if ( paramVec!=Teuchos::null && !parameters_[i]->is_distributed) {
541 Teuchos::ArrayRCP<const Scalar> p_data;
542 rcp_dynamic_cast<
const Thyra::SpmdVectorBase<Scalar> >(paramVec,
true)->getLocalData(Teuchos::ptrFromRef(p_data));
544 for (
unsigned int j=0; j < parameters_[i]->scalar_value.size(); j++) {
545 parameters_[i]->scalar_value[j].baseValue = p_data[j];
546 parameters_[i]->scalar_value[j].family->setRealValueForAllTypes(parameters_[i]->scalar_value[j].baseValue);
549 else if ( paramVec!=Teuchos::null && parameters_[i]->is_distributed) {
552 std::string key = (*parameters_[i]->names)[0];
553 RCP<GlobalEvaluationData> ged = distrParamGlobalEvaluationData_.getDataObject(key);
555 TEUCHOS_ASSERT(ged!=Teuchos::null);
560 if(loc_pair_ged!=Teuchos::null) {
562 RCP<ThyraObjContainer<Scalar> > th_ged = rcp_dynamic_cast<
ThyraObjContainer<Scalar> >(loc_pair_ged->getGlobalLOC(),
true);
566 TEUCHOS_ASSERT(ro_ged!=Teuchos::null);
567 ro_ged->setOwnedVector(paramVec);
575 const RCP<panzer::ThyraObjContainer<Scalar> > thGlobalContainer =
578 TEUCHOS_ASSERT(!Teuchos::is_null(thGlobalContainer));
594 xContainer_->setOwnedVector(x);
599 xdotContainer_->setOwnedVector(x_dot);
609 for (
int i(0); i < num_param_vecs; ++i)
615 if (not parameters_[i]->is_distributed)
617 auto dxdp = rcp_const_cast<VectorBase<Scalar>>
618 (inArgs.get_p(vIndex + num_param_vecs));
619 if (not dxdp.is_null())
623 auto dxdpBlock = rcp_dynamic_cast<ProductVectorBase<Scalar>>(dxdp);
624 int numParams(parameters_[i]->scalar_value.size());
625 for (
int j(0); j < numParams; ++j)
627 RCP<ROVGED> dxdpContainer = lof_->buildReadOnlyDomainContainer();
628 dxdpContainer->
setOwnedVector(dxdpBlock->getNonconstVectorBlock(j));
629 string name(
"X TANGENT GATHER CONTAINER: " +
630 (*parameters_[i]->names)[j]);
634 if (build_transient_support_)
638 auto dxdotdp = rcp_const_cast<VectorBase<Scalar>>
639 (inArgs.get_p(vIndex + num_param_vecs + tangent_space_.size()));
640 if (not dxdotdp.is_null())
643 rcp_dynamic_cast<ProductVectorBase<Scalar>>(dxdotdp);
644 int numParams(parameters_[i]->scalar_value.size());
645 for (
int j(0); j < numParams; ++j)
647 RCP<ROVGED> dxdotdpContainer = lof_->buildReadOnlyDomainContainer();
648 dxdotdpContainer->setOwnedVector(
649 dxdotdpBlock->getNonconstVectorBlock(j));
650 string name(
"DXDT TANGENT GATHER CONTAINER: " +
651 (*parameters_[i]->names)[j]);
664 template <
typename Scalar>
665 Thyra::ModelEvaluatorBase::OutArgs<Scalar>
668 typedef Thyra::ModelEvaluatorBase MEB;
670 if(require_out_args_refresh_) {
671 MEB::OutArgsSetup<Scalar> outArgs;
672 outArgs.setModelEvalDescription(this->description());
673 outArgs.set_Np_Ng(num_me_parameters_, responses_.size());
674 outArgs.setSupports(MEB::OUT_ARG_f);
675 outArgs.setSupports(MEB::OUT_ARG_W_op);
678 for(std::size_t i=0;i<responses_.size();i++) {
682 Teuchos::RCP<panzer::ResponseBase> respJacBase
683 = responseLibrary_->getResponse<RespEvalT>(responses_[i]->name);
684 if(respJacBase!=Teuchos::null) {
686 Teuchos::RCP<panzer::ResponseMESupportBase<RespEvalT> > resp
690 if(resp->supportsDerivative()) {
691 outArgs.setSupports(MEB::OUT_ARG_DgDx,i,MEB::DerivativeSupport(MEB::DERIV_MV_GRADIENT_FORM));
694 for(std::size_t p=0;p<parameters_.size();p++) {
695 if(parameters_[p]->is_distributed && parameters_[p]->global_indexer!=Teuchos::null)
696 outArgs.setSupports(MEB::OUT_ARG_DgDp,i,p,MEB::DerivativeSupport(MEB::DERIV_MV_GRADIENT_FORM));
697 if(!parameters_[p]->is_distributed)
698 outArgs.setSupports(MEB::OUT_ARG_DgDp,i,p,MEB::DerivativeSupport(MEB::DERIV_MV_JACOBIAN_FORM));
705 for(std::size_t p=0;p<parameters_.size();p++) {
707 if(!parameters_[p]->is_distributed)
708 outArgs.setSupports(MEB::OUT_ARG_DfDp,p,MEB::DerivativeSupport(MEB::DERIV_MV_BY_COL));
709 else if(parameters_[p]->is_distributed && parameters_[p]->global_indexer!=Teuchos::null)
710 outArgs.setSupports(MEB::OUT_ARG_DfDp,p,MEB::DerivativeSupport(MEB::DERIV_LINEAR_OP));
713 prototypeOutArgs_ = outArgs;
717 require_out_args_refresh_ =
false;
719 return prototypeOutArgs_;
722 template <
typename Scalar>
723 Teuchos::RCP<Thyra::LinearOpBase<Scalar> >
727 Teuchos::RCP<const ThyraObjFactory<Scalar> > tof
733 template <
typename Scalar>
734 Teuchos::RCP<const Thyra::LinearOpWithSolveFactoryBase<Scalar> >
738 return solverFactory_;
741 template <
typename Scalar>
742 Teuchos::RCP<Thyra::LinearOpBase<Scalar> >
747 using Teuchos::rcp_dynamic_cast;
749 typedef Thyra::ModelEvaluatorBase MEB;
761 if(require_out_args_refresh_) {
762 this->createOutArgs();
765 TEUCHOS_ASSERT(0<=p && p<Teuchos::as<int>(parameters_.size()));
771 TEUCHOS_ASSERT(prototypeOutArgs_.supports(MEB::OUT_ARG_DfDp,p).supports(MEB::DERIV_LINEAR_OP));
775 RCP<Response_Residual<Traits::Jacobian> > response_jacobian
778 return response_jacobian->allocateJacobian();
781 TEUCHOS_ASSERT(prototypeOutArgs_.supports(MEB::OUT_ARG_DfDp,p).supports(MEB::DERIV_MV_BY_COL));
784 return Thyra::createMember(*get_f_space());
788 TEUCHOS_ASSERT(
false);
790 return Teuchos::null;
793 template <
typename Scalar>
797 Teuchos::Array<std::string> tmp_names;
798 tmp_names.push_back(name);
800 Teuchos::Array<Scalar> tmp_values;
801 tmp_values.push_back(initialValue);
803 return addParameter(tmp_names,tmp_values);
806 template <
typename Scalar>
809 const Teuchos::Array<Scalar> & initialValues)
813 using Teuchos::rcp_dynamic_cast;
814 using Teuchos::ptrFromRef;
816 TEUCHOS_ASSERT(names.size()==initialValues.size());
818 int parameter_index = parameters_.size();
821 RCP<ParameterObject> param = createScalarParameter(names,initialValues);
822 parameters_.push_back(param);
825 RCP< Thyra::VectorSpaceBase<double> > tan_space =
826 Thyra::multiVectorProductVectorSpace(x_space_, param->names->size());
827 tangent_space_.push_back(tan_space);
831 num_me_parameters_ += 2;
832 if (build_transient_support_)
833 ++num_me_parameters_;
835 require_in_args_refresh_ =
true;
836 require_out_args_refresh_ =
true;
837 this->resetDefaultBase();
839 return parameter_index;
842 template <
typename Scalar>
846 const Teuchos::RCP<GlobalEvaluationData> & ged,
848 const Teuchos::RCP<const GlobalIndexer> & ugi)
850 distrParamGlobalEvaluationData_.addDataObject(key,ged);
852 int parameter_index = parameters_.size();
853 parameters_.push_back(createDistributedParameter(key,vs,initial,ugi));
854 ++num_me_parameters_;
856 require_in_args_refresh_ =
true;
857 require_out_args_refresh_ =
true;
858 this->resetDefaultBase();
860 return parameter_index;
863 template <
typename Scalar>
866 const Teuchos::RCP<GlobalEvaluationData> & ged)
868 nonParamGlobalEvaluationData_.addDataObject(key,ged);
871 template <
typename Scalar>
874 const std::vector<WorksetDescriptor> & wkst_desc,
875 const Teuchos::RCP<ResponseMESupportBuilderBase> & builder)
878 builder->setDerivativeInformation(lof_);
880 int respIndex = addResponse(responseName,wkst_desc,*builder);
883 responses_[respIndex]->builder = builder;
889 template <
typename Scalar>
895 using Teuchos::ArrayRCP;
896 using Teuchos::Array;
897 using Teuchos::tuple;
898 using Teuchos::rcp_dynamic_cast;
901 if(Teuchos::is_null(ghostedContainer_)) {
902 ghostedContainer_ = lof_->buildGhostedLinearObjContainer();
908 ae_inargs.
container_ = lof_->buildLinearObjContainer();
910 ae_inargs.
alpha = 0.0;
911 ae_inargs.
beta = 1.0;
923 const RCP<panzer::ThyraObjContainer<Scalar> > thGlobalContainer =
926 TEUCHOS_ASSERT(!Teuchos::is_null(thGlobalContainer));
931 const RCP<panzer::ThyraObjContainer<Scalar> > thGhostedContainer =
933 Thyra::assign(thGhostedContainer->get_f_th().ptr(),0.0);
942 thGlobalContainer->set_x_th(x);
945 RCP<panzer::LinearObjContainer> counter
946 = ae_tm_.template getAsObject<panzer::Traits::Residual>()->evaluateOnlyDirichletBCs(ae_inargs);
949 RCP<panzer::LinearObjContainer>
result = lof_->buildLinearObjContainer();
953 thGlobalContainer->get_f_th());
959 lof_->applyDirichletBCs(*counter,*
result);
962 template <
typename Scalar>
965 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
969 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 972 setParameters(inArgs);
975 std::string responseName = responses_[respIndex]->name;
976 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Hessian> > resp
979 resp->setDerivative(D2gDx2);
986 setupAssemblyInArgs(inArgs,ae_inargs);
988 ae_inargs.
beta = 1.0;
990 auto deltaXContainer = lof_->buildReadOnlyDomainContainer();
991 deltaXContainer->setOwnedVector(delta_x);
1005 TEUCHOS_ASSERT(
false);
1009 template <
typename Scalar>
1013 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
1017 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1020 setParameters(inArgs);
1023 std::string responseName = responses_[respIndex]->name;
1024 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Hessian> > resp
1027 resp->setDerivative(D2gDxDp);
1034 setupAssemblyInArgs(inArgs,ae_inargs);
1036 ae_inargs.
beta = 1.0;
1039 auto deltaPContainer = parameters_[pIndex]->dfdp_rl->getLinearObjFactory()->buildReadOnlyDomainContainer();
1040 deltaPContainer->setOwnedVector(delta_p);
1055 TEUCHOS_ASSERT(
false);
1059 template <
typename Scalar>
1063 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
1067 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1070 setParameters(inArgs);
1075 std::string responseName = responses_[respIndex]->name;
1076 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Hessian> > resp
1079 resp->setDerivative(D2gDp2);
1086 setupAssemblyInArgs(inArgs,ae_inargs);
1093 auto deltaPContainer = parameters_[pIndex]->dfdp_rl->getLinearObjFactory()->buildReadOnlyDomainContainer();
1094 deltaPContainer->setOwnedVector(delta_p);
1109 TEUCHOS_ASSERT(
false);
1113 template <
typename Scalar>
1117 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
1121 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1124 setParameters(inArgs);
1129 std::string responseName = responses_[respIndex]->name;
1130 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Hessian> > resp
1133 resp->setDerivative(D2gDpDx);
1140 setupAssemblyInArgs(inArgs,ae_inargs);
1147 auto deltaXContainer = lof_->buildReadOnlyDomainContainer();
1148 deltaXContainer->setOwnedVector(delta_x);
1163 TEUCHOS_ASSERT(
false);
1167 template <
typename Scalar>
1173 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1176 using Teuchos::ArrayRCP;
1177 using Teuchos::Array;
1178 using Teuchos::tuple;
1179 using Teuchos::rcp_dynamic_cast;
1181 typedef Thyra::ModelEvaluatorBase MEB;
1186 bool is_transient =
false;
1187 if (inArgs.supports(MEB::IN_ARG_x_dot ))
1188 is_transient = !Teuchos::is_null(inArgs.get_x_dot());
1191 TEUCHOS_TEST_FOR_EXCEPTION(is_transient && !build_transient_support_, std::runtime_error,
1192 "ModelEvaluator was not built with transient support enabled!");
1197 const RCP<Thyra::LinearOpBase<Scalar> > W_out = D2fDx2;
1203 setupAssemblyInArgs(inArgs,ae_inargs);
1205 auto deltaXContainer = lof_->buildReadOnlyDomainContainer();
1206 deltaXContainer->setOwnedVector(delta_x);
1210 setParameters(inArgs);
1216 if(oneTimeDirichletBeta_on_) {
1220 oneTimeDirichletBeta_on_ =
false;
1224 const RCP<panzer::ThyraObjContainer<Scalar> > thGlobalContainer =
1226 const RCP<panzer::ThyraObjContainer<Scalar> > thGhostedContainer =
1230 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(D2fDx2)");
1233 RCP<Thyra::VectorBase<Scalar> > dummy_f = Thyra::createMember(f_space_);
1234 thGlobalContainer->set_f_th(dummy_f);
1235 thGlobalContainer->set_A_th(W_out);
1238 thGhostedContainer->initializeMatrix(0.0);
1240 ae_tm_.template getAsObject<panzer::Traits::Hessian>()->evaluate(ae_inargs);
1246 thGlobalContainer->set_A_th(Teuchos::null);
1252 thGlobalContainer->set_x_th(Teuchos::null);
1253 thGlobalContainer->set_dxdt_th(Teuchos::null);
1254 thGlobalContainer->set_f_th(Teuchos::null);
1255 thGlobalContainer->set_A_th(Teuchos::null);
1263 TEUCHOS_ASSERT(
false);
1267 template <
typename Scalar>
1270 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
1274 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1277 using Teuchos::ArrayRCP;
1278 using Teuchos::Array;
1279 using Teuchos::tuple;
1280 using Teuchos::rcp_dynamic_cast;
1282 typedef Thyra::ModelEvaluatorBase MEB;
1287 bool is_transient =
false;
1288 if (inArgs.supports(MEB::IN_ARG_x_dot ))
1289 is_transient = !Teuchos::is_null(inArgs.get_x_dot());
1292 TEUCHOS_TEST_FOR_EXCEPTION(is_transient && !build_transient_support_, std::runtime_error,
1293 "ModelEvaluator was not built with transient support enabled!");
1298 const RCP<Thyra::LinearOpBase<Scalar> > W_out = D2fDxDp;
1304 setupAssemblyInArgs(inArgs,ae_inargs);
1308 auto deltaPContainer = parameters_[pIndex]->dfdp_rl->getLinearObjFactory()->buildReadOnlyDomainContainer();
1309 deltaPContainer->setOwnedVector(delta_p);
1313 setParameters(inArgs);
1319 if(oneTimeDirichletBeta_on_) {
1323 oneTimeDirichletBeta_on_ =
false;
1327 const RCP<panzer::ThyraObjContainer<Scalar> > thGlobalContainer =
1329 const RCP<panzer::ThyraObjContainer<Scalar> > thGhostedContainer =
1333 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(D2fDxDp)");
1336 RCP<Thyra::VectorBase<Scalar> > dummy_f = Thyra::createMember(f_space_);
1337 thGlobalContainer->set_f_th(dummy_f);
1338 thGlobalContainer->set_A_th(W_out);
1341 thGhostedContainer->initializeMatrix(0.0);
1343 ae_tm_.template getAsObject<panzer::Traits::Hessian>()->evaluate(ae_inargs);
1349 thGlobalContainer->set_A_th(Teuchos::null);
1355 thGlobalContainer->set_x_th(Teuchos::null);
1356 thGlobalContainer->set_dxdt_th(Teuchos::null);
1357 thGlobalContainer->set_f_th(Teuchos::null);
1358 thGlobalContainer->set_A_th(Teuchos::null);
1367 TEUCHOS_ASSERT(
false);
1371 template <
typename Scalar>
1374 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
1378 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1380 using Teuchos::rcp_dynamic_cast;
1381 using Teuchos::null;
1384 TEUCHOS_ASSERT(parameters_[pIndex]->is_distributed);
1388 TEUCHOS_ASSERT(parameters_[pIndex]->dfdp_rl!=null);
1393 RCP<Response_Residual<Traits::Hessian> > response_hessian =
1395 response_hessian->setHessian(D2fDpDx);
1400 setupAssemblyInArgs(inArgs,ae_inargs);
1402 auto deltaXContainer = lof_->buildReadOnlyDomainContainer();
1403 deltaXContainer->setOwnedVector(delta_x);
1418 TEUCHOS_ASSERT(
false);
1422 template <
typename Scalar>
1425 const Thyra::ModelEvaluatorBase::InArgs<Scalar> & inArgs,
1429 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 1431 using Teuchos::rcp_dynamic_cast;
1432 using Teuchos::null;
1435 TEUCHOS_ASSERT(parameters_[pIndex]->is_distributed);
1439 TEUCHOS_ASSERT(parameters_[pIndex]->dfdp_rl!=null);
1444 RCP<Response_Residual<Traits::Hessian> > response_hessian =
1446 response_hessian->setHessian(D2fDp2);
1451 setupAssemblyInArgs(inArgs,ae_inargs);
1453 auto deltaPContainer = parameters_[pIndex]->dfdp_rl->getLinearObjFactory()->buildReadOnlyDomainContainer();
1454 deltaPContainer->setOwnedVector(delta_p);
1469 TEUCHOS_ASSERT(
false);
1473 template <
typename Scalar>
1476 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1478 evalModelImpl_basic(inArgs,outArgs);
1481 if(required_basic_g(outArgs))
1482 evalModelImpl_basic_g(inArgs,outArgs);
1485 if(required_basic_dgdx(outArgs))
1486 evalModelImpl_basic_dgdx(inArgs,outArgs);
1489 if(required_basic_dgdp_scalar(outArgs))
1490 evalModelImpl_basic_dgdp_scalar(inArgs,outArgs);
1493 if(required_basic_dgdp_distro(outArgs))
1494 evalModelImpl_basic_dgdp_distro(inArgs,outArgs);
1496 if(required_basic_dfdp_scalar(outArgs)) {
1498 evalModelImpl_basic_dfdp_scalar_fd(inArgs,outArgs);
1500 evalModelImpl_basic_dfdp_scalar(inArgs,outArgs);
1503 if(required_basic_dfdp_distro(outArgs))
1504 evalModelImpl_basic_dfdp_distro(inArgs,outArgs);
1507 template <
typename Scalar>
1510 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1513 using Teuchos::ArrayRCP;
1514 using Teuchos::Array;
1515 using Teuchos::tuple;
1516 using Teuchos::rcp_dynamic_cast;
1518 typedef Thyra::ModelEvaluatorBase MEB;
1523 bool is_transient =
false;
1524 if (inArgs.supports(MEB::IN_ARG_x_dot ))
1525 is_transient = !Teuchos::is_null(inArgs.get_x_dot());
1528 TEUCHOS_TEST_FOR_EXCEPTION(is_transient && !build_transient_support_, std::runtime_error,
1529 "ModelEvaluator was not built with transient support enabled!");
1534 const RCP<Thyra::VectorBase<Scalar> > f_out = outArgs.get_f();
1535 const RCP<Thyra::LinearOpBase<Scalar> > W_out = outArgs.get_W_op();
1538 if(Teuchos::is_null(f_out) && Teuchos::is_null(W_out) ) {
1546 setupAssemblyInArgs(inArgs,ae_inargs);
1549 setParameters(inArgs);
1555 if(oneTimeDirichletBeta_on_) {
1559 oneTimeDirichletBeta_on_ =
false;
1563 const RCP<panzer::ThyraObjContainer<Scalar> > thGlobalContainer =
1565 const RCP<panzer::ThyraObjContainer<Scalar> > thGhostedContainer =
1568 if (!Teuchos::is_null(f_out) && !Teuchos::is_null(W_out)) {
1569 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(f and J)");
1575 thGlobalContainer->set_f_th(f_out);
1576 thGlobalContainer->set_A_th(W_out);
1579 Thyra::assign(thGhostedContainer->get_f_th().ptr(),0.0);
1580 thGhostedContainer->initializeMatrix(0.0);
1582 ae_tm_.template getAsObject<panzer::Traits::Jacobian>()->evaluate(ae_inargs);
1584 else if(!Teuchos::is_null(f_out) && Teuchos::is_null(W_out)) {
1586 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(f)");
1591 thGlobalContainer->set_f_th(f_out);
1594 Thyra::assign(thGhostedContainer->get_f_th().ptr(),0.0);
1596 ae_tm_.template getAsObject<panzer::Traits::Residual>()->evaluate(ae_inargs);
1598 else if(Teuchos::is_null(f_out) && !Teuchos::is_null(W_out)) {
1600 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(J)");
1606 RCP<Thyra::VectorBase<Scalar> > dummy_f = Thyra::createMember(f_space_);
1607 thGlobalContainer->set_f_th(dummy_f);
1608 thGlobalContainer->set_A_th(W_out);
1611 thGhostedContainer->initializeMatrix(0.0);
1613 ae_tm_.template getAsObject<panzer::Traits::Jacobian>()->evaluate(ae_inargs);
1619 thGlobalContainer->set_A_th(Teuchos::null);
1625 thGlobalContainer->set_x_th(Teuchos::null);
1626 thGlobalContainer->set_dxdt_th(Teuchos::null);
1627 thGlobalContainer->set_f_th(Teuchos::null);
1628 thGlobalContainer->set_A_th(Teuchos::null);
1634 template <
typename Scalar>
1637 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1646 setupAssemblyInArgs(inArgs,ae_inargs);
1649 setParameters(inArgs);
1651 for(std::size_t i=0;i<responses_.size();i++) {
1652 Teuchos::RCP<Thyra::VectorBase<Scalar> > vec = outArgs.get_g(i);
1653 if(vec!=Teuchos::null) {
1654 std::string responseName = responses_[i]->name;
1655 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Residual> > resp
1658 resp->setVector(vec);
1670 template <
typename Scalar>
1674 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1676 typedef Thyra::ModelEvaluatorBase MEB;
1679 TEUCHOS_ASSERT(required_basic_dgdx(outArgs));
1682 setParameters(inArgs);
1684 for(std::size_t i=0;i<responses_.size();i++) {
1686 MEB::Derivative<Scalar> deriv = outArgs.get_DgDx(i);
1690 Teuchos::RCP<Thyra::MultiVectorBase<Scalar> > vec = deriv.getMultiVector();
1692 if(vec!=Teuchos::null) {
1694 std::string responseName = responses_[i]->name;
1695 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Jacobian> > resp
1698 resp->setDerivative(vec);
1706 setupAssemblyInArgs(inArgs,ae_inargs);
1716 template <
typename Scalar>
1720 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1724 using Teuchos::rcp_dynamic_cast;
1726 typedef Thyra::ModelEvaluatorBase MEB;
1729 TEUCHOS_ASSERT(required_basic_dgdp_scalar(outArgs));
1732 std::vector<std::string> activeParameterNames;
1733 std::vector<int> activeParameters;
1734 int totalParameterCount = 0;
1735 for(std::size_t j=0; j<parameters_.size(); j++) {
1738 if(parameters_[j]->is_distributed)
1741 bool is_active =
false;
1742 for(std::size_t i=0;i<responses_.size(); i++) {
1744 MEB::Derivative<Scalar> deriv = outArgs.get_DgDp(i,j);
1748 Teuchos::RCP<Thyra::MultiVectorBase<Scalar> > vec = deriv.getMultiVector();
1749 if(vec!=Teuchos::null) {
1751 std::string responseName = responses_[i]->name;
1752 RCP<panzer::ResponseMESupportBase<panzer::Traits::Tangent> > resp =
1755 resp->setVector(vec);
1761 for (std::size_t k=0; k<parameters_[j]->scalar_value.size(); k++) {
1762 std::string name =
"PARAMETER_SENSITIVIES: "+(*parameters_[j]->names)[k];
1763 activeParameterNames.push_back(name);
1764 totalParameterCount++;
1766 activeParameters.push_back(j);
1772 setupAssemblyInArgs(inArgs,ae_inargs);
1775 RCP<panzer::GlobalEvaluationData> ged_activeParameters =
1781 for (std::size_t ap=0; ap<activeParameters.size(); ++ap) {
1782 const int j = activeParameters[ap];
1783 for (
unsigned int k=0; k < parameters_[j]->scalar_value.size(); k++) {
1785 p.fastAccessDx(paramIndex) = 1.0;
1786 parameters_[j]->scalar_value[k].family->template setValue<panzer::Traits::Tangent>(p);
1792 TEUCHOS_ASSERT(paramIndex==totalParameterCount);
1795 if(totalParameterCount>0) {
1801 template <
typename Scalar>
1805 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1807 typedef Thyra::ModelEvaluatorBase MEB;
1810 TEUCHOS_ASSERT(required_basic_dgdp_distro(outArgs));
1816 for(std::size_t p=0;p<parameters_.size();p++) {
1820 if(!parameters_[p]->is_distributed)
1825 for(std::size_t r=0;r<responses_.size();r++) {
1827 MEB::Derivative<Scalar> deriv = outArgs.get_DgDp(r,p);
1831 Teuchos::RCP<Thyra::MultiVectorBase<Scalar> > vec = deriv.getMultiVector();
1833 if(vec!=Teuchos::null) {
1836 std::string responseName = responses_[r]->name;
1837 Teuchos::RCP<panzer::ResponseMESupportBase<panzer::Traits::Jacobian> > resp
1841 resp->setDerivative(vec);
1848 setupAssemblyInArgs(inArgs,ae_inargs);
1860 template <
typename Scalar>
1864 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1867 using Teuchos::rcp_dynamic_cast;
1869 typedef Thyra::ModelEvaluatorBase MEB;
1871 TEUCHOS_ASSERT(required_basic_dfdp_scalar(outArgs));
1877 setupAssemblyInArgs(inArgs,ae_inargs);
1883 std::vector<std::string> activeParameters;
1885 int totalParameterCount = 0;
1886 for(std::size_t i=0; i < parameters_.size(); i++) {
1888 if(parameters_[i]->is_distributed)
1892 MEB::Derivative<Scalar> deriv = outArgs.get_DfDp(i);
1897 Teuchos::RCP<Thyra::MultiVectorBase<Scalar> > mVec = deriv.getMultiVector();
1898 TEUCHOS_ASSERT(mVec->domain()->dim()==Teuchos::as<int>(parameters_[i]->scalar_value.size()));
1900 for (std::size_t j=0; j < parameters_[i]->scalar_value.size(); j++) {
1903 RCP<LOCPair_GlobalEvaluationData> loc_pair
1905 RCP<LinearObjContainer> globalContainer = loc_pair->getGlobalLOC();
1908 RCP<Thyra::VectorBase<Scalar> > vec = mVec->col(j);
1909 RCP<panzer::ThyraObjContainer<Scalar> > thGlobalContainer =
1914 std::string name =
"PARAMETER_SENSITIVIES: "+(*parameters_[i]->names)[j];
1918 activeParameters.push_back(name);
1919 totalParameterCount++;
1927 RCP<GlobalEvaluationData> ged_activeParameters
1936 for(std::size_t i=0; i < parameters_.size(); i++) {
1938 if(parameters_[i]->is_distributed)
1942 MEB::Derivative<Scalar> deriv = outArgs.get_DfDp(i);
1943 if(deriv.isEmpty()) {
1945 for (
unsigned int j=0; j < parameters_[i]->scalar_value.size(); j++) {
1947 parameters_[i]->scalar_value[j].baseValue);
1948 parameters_[i]->scalar_value[j].family->template setValue<panzer::Traits::Tangent>(p);
1954 for (
unsigned int j=0; j < parameters_[i]->scalar_value.size(); j++) {
1956 parameters_[i]->scalar_value[j].baseValue);
1957 p.fastAccessDx(paramIndex) = 1.0;
1958 parameters_[i]->scalar_value[j].family->template setValue<panzer::Traits::Tangent>(p);
1965 TEUCHOS_ASSERT(paramIndex==totalParameterCount);
1970 if(totalParameterCount>0) {
1971 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(df/dp)");
1976 template <
typename Scalar>
1980 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 1982 PANZER_FUNC_TIME_MONITOR(
"panzer::ModelEvaluator::evalModel(df/dp)");
1985 using Teuchos::rcp_dynamic_cast;
1987 typedef Thyra::ModelEvaluatorBase MEB;
1989 TEUCHOS_ASSERT(required_basic_dfdp_scalar(outArgs));
1994 MEB::OutArgs<Scalar> outArgs_base = this->createOutArgs();
1995 if (outArgs.get_f() == Teuchos::null)
1996 outArgs_base.set_f(Thyra::createMember(this->get_f_space()));
1998 outArgs_base.set_f(outArgs.get_f());
1999 outArgs_base.set_W_op(outArgs.get_W_op());
2000 this->evalModel(inArgs, outArgs_base);
2001 RCP<const Thyra::VectorBase<Scalar> > f = outArgs_base.get_f();
2002 RCP<const Thyra::VectorBase<Scalar> > x = inArgs.get_x();
2003 RCP<const Thyra::VectorBase<Scalar> > x_dot;
2004 if (inArgs.supports(MEB::IN_ARG_x_dot))
2005 x_dot = inArgs.get_x_dot();
2008 RCP<Thyra::VectorBase<Scalar> > fd = Thyra::createMember(this->get_f_space());
2009 MEB::OutArgs<Scalar> outArgs_fd = this->createOutArgs();
2010 outArgs_fd.set_f(fd);
2012 RCP<Thyra::VectorBase<Scalar> > xd = Thyra::createMember(this->get_x_space());
2013 RCP<Thyra::VectorBase<Scalar> > xd_dot;
2014 if (x_dot != Teuchos::null)
2015 xd_dot = Thyra::createMember(this->get_x_space());
2016 MEB::InArgs<Scalar> inArgs_fd = this->createInArgs();
2017 inArgs_fd.setArgs(inArgs);
2018 inArgs_fd.set_x(xd);
2019 if (x_dot != Teuchos::null)
2020 inArgs_fd.set_x_dot(xd_dot);
2022 const double h = fd_perturb_size_;
2023 for(std::size_t i=0; i < parameters_.size(); i++) {
2026 if(parameters_[i]->is_distributed)
2030 MEB::Derivative<Scalar> deriv = outArgs.get_DfDp(i);
2035 RCP<Thyra::MultiVectorBase<Scalar> > dfdp = deriv.getMultiVector();
2036 TEUCHOS_ASSERT(dfdp->domain()->dim()==Teuchos::as<int>(parameters_[i]->scalar_value.size()));
2039 RCP<const Thyra::VectorBase<Scalar> > p = inArgs.get_p(i);
2040 RCP<const Thyra::VectorBase<Scalar> > dx_v = inArgs.get_p(i+parameters_.size());
2041 RCP<const Thyra::MultiVectorBase<Scalar> > dx =
2042 rcp_dynamic_cast<
const Thyra::DefaultMultiVectorProductVector<Scalar> >(dx_v,
true)->
getMultiVector();
2043 RCP<const Thyra::VectorBase<Scalar> > dx_dot_v;
2044 RCP<const Thyra::MultiVectorBase<Scalar> > dx_dot;
2045 if (x_dot != Teuchos::null) {
2046 dx_dot_v =inArgs.get_p(i+parameters_.size()+tangent_space_.size());
2048 rcp_dynamic_cast<
const Thyra::DefaultMultiVectorProductVector<Scalar> >(dx_dot_v,
true)->
getMultiVector();
2052 RCP<Thyra::VectorBase<Scalar> > pd = Thyra::createMember(this->get_p_space(i));
2053 inArgs_fd.set_p(i,pd);
2055 for (std::size_t j=0; j < parameters_[i]->scalar_value.size(); j++) {
2058 Thyra::copy(*p, pd.ptr());
2059 Thyra::set_ele(j, Thyra::get_ele(*p,j)+h, pd.ptr());
2062 Thyra::V_VpStV(xd.ptr(), *x, h, *(dx)->col(j));
2063 if (x_dot != Teuchos::null)
2064 Thyra::V_VpStV(xd_dot.ptr(), *x_dot, h, *(dx_dot)->col(j));
2067 Thyra::assign(fd.ptr(), 0.0);
2068 this->evalModel(inArgs_fd, outArgs_fd);
2071 Thyra::V_StVpStV(dfdp->col(j).ptr(), 1.0/h, *fd, -1.0/h, *f);
2074 parameters_[i]->scalar_value[j].family->setRealValueForAllTypes(Thyra::get_ele(*p,j));
2080 template <
typename Scalar>
2084 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const 2087 using Teuchos::rcp_dynamic_cast;
2088 using Teuchos::null;
2090 typedef Thyra::ModelEvaluatorBase MEB;
2092 TEUCHOS_ASSERT(required_basic_dfdp_distro(outArgs));
2098 for(std::size_t p=0;p<parameters_.size();p++) {
2102 if(!parameters_[p]->is_distributed)
2107 if(parameters_[p]->dfdp_rl==null)
2111 MEB::Derivative<Scalar> deriv = outArgs.get_DfDp(p);
2118 RCP<Response_Residual<Traits::Jacobian> > response_jacobian =
2120 response_jacobian->setJacobian(deriv.getLinearOp());
2125 setupAssemblyInArgs(inArgs,ae_inargs);
2136 template <
typename Scalar>
2141 bool activeGArgs =
false;
2142 for(
int i=0;i<outArgs.Ng();i++)
2143 activeGArgs |= (outArgs.get_g(i)!=Teuchos::null);
2145 return activeGArgs | required_basic_dgdx(outArgs);
2148 template <
typename Scalar>
2152 typedef Thyra::ModelEvaluatorBase MEB;
2155 bool activeGArgs =
false;
2156 for(
int i=0;i<outArgs.Ng();i++) {
2158 if(outArgs.supports(MEB::OUT_ARG_DgDx,i).none())
2162 activeGArgs |= (!outArgs.get_DgDx(i).isEmpty());
2168 template <
typename Scalar>
2172 typedef Thyra::ModelEvaluatorBase MEB;
2175 bool activeGArgs =
false;
2176 for(
int i=0;i<outArgs.Ng();i++) {
2177 for(
int p=0;p<Teuchos::as<int>(parameters_.size());p++) {
2180 if(parameters_[p]->is_distributed)
2184 if(outArgs.supports(MEB::OUT_ARG_DgDp,i,p).none())
2187 activeGArgs |= (!outArgs.get_DgDp(i,p).isEmpty());
2194 template <
typename Scalar>
2198 typedef Thyra::ModelEvaluatorBase MEB;
2201 bool activeGArgs =
false;
2202 for(
int i=0;i<outArgs.Ng();i++) {
2203 for(
int p=0;p<Teuchos::as<int>(parameters_.size());p++) {
2206 if(!parameters_[p]->is_distributed)
2210 if(outArgs.supports(MEB::OUT_ARG_DgDp,i,p).none())
2213 activeGArgs |= (!outArgs.get_DgDp(i,p).isEmpty());
2220 template <
typename Scalar>
2224 typedef Thyra::ModelEvaluatorBase MEB;
2227 bool activeFPArgs =
false;
2228 for(
int i=0;i<Teuchos::as<int>(parameters_.size());i++) {
2231 if(parameters_[i]->is_distributed)
2235 if(outArgs.supports(MEB::OUT_ARG_DfDp,i).none())
2239 activeFPArgs |= (!outArgs.get_DfDp(i).isEmpty());
2242 return activeFPArgs;
2245 template <
typename Scalar>
2249 typedef Thyra::ModelEvaluatorBase MEB;
2252 bool activeFPArgs =
false;
2253 for(
int i=0;i<Teuchos::as<int>(parameters_.size());i++) {
2256 if(!parameters_[i]->is_distributed)
2260 if(outArgs.supports(MEB::OUT_ARG_DfDp,i).none())
2264 activeFPArgs |= (!outArgs.get_DfDp(i).isEmpty());
2267 return activeFPArgs;
2270 template <
typename Scalar>
2273 const Teuchos::RCP<panzer::WorksetContainer> & wc,
2274 const std::vector<Teuchos::RCP<panzer::PhysicsBlock> >& physicsBlocks,
2275 const std::vector<panzer::BC> & bcs,
2279 const Teuchos::ParameterList& closure_models,
2280 const Teuchos::ParameterList& user_data,
2281 const bool write_graphviz_file,
2282 const std::string& graphviz_file_prefix)
2286 using Teuchos::null;
2292 for(std::size_t p=0;p<parameters_.size();p++) {
2295 if(!parameters_[p]->is_distributed)
2300 if(parameters_[p]->global_indexer==null)
2306 parameters_[p]->global_indexer);
2309 RCP<ResponseLibrary<Traits> > rLibrary
2312 rLibrary->buildResidualResponseEvaluators(physicsBlocks,eqset_factory,bcs,bc_factory,
2313 cm_factory,closure_models,user_data,
2314 write_graphviz_file,graphviz_file_prefix);
2317 parameters_[p]->dfdp_rl = rLibrary;
2321 template <
typename Scalar>
2324 const Teuchos::RCP<panzer::WorksetContainer> & wc,
2325 const std::vector<Teuchos::RCP<panzer::PhysicsBlock> >& physicsBlocks,
2326 const std::vector<panzer::BC>& ,
2330 const Teuchos::ParameterList& closure_models,
2331 const Teuchos::ParameterList& user_data,
2332 const bool write_graphviz_file,
2333 const std::string& graphviz_file_prefix)
2337 using Teuchos::null;
2343 for(std::size_t p=0;p<parameters_.size();p++) {
2346 if(!parameters_[p]->is_distributed)
2351 if(parameters_[p]->global_indexer==null)
2356 RCP<const LinearObjFactory<Traits> > param_lof = parameters_[p]->dfdp_rl->getLinearObjFactory();
2357 RCP<const GlobalIndexer > param_ugi = parameters_[p]->global_indexer;
2360 RCP<ResponseLibrary<Traits> > rLibrary
2365 for(std::size_t r=0;r<responses_.size();r++) {
2367 if(responses_[r]->builder==Teuchos::null)
2372 responses_[r]->builder->setDerivativeInformation(param_lof);
2375 rLibrary->addResponse(responses_[r]->name,
2376 responses_[r]->wkst_desc,
2377 *responses_[r]->builder);
2380 rLibrary->buildResponseEvaluators(physicsBlocks,eqset_factory,
2381 cm_factory,closure_models,user_data,
2382 write_graphviz_file,graphviz_file_prefix);
2385 parameters_[p]->dgdp_rl = rLibrary;
2389 template <
typename Scalar>
2393 oneTimeDirichletBeta_on_ =
true;
2394 oneTimeDirichletBeta_ = beta;
2397 template <
typename Scalar>
2398 Teuchos::RCP<typename panzer::ModelEvaluator<Scalar>::ParameterObject>
2401 const Teuchos::Array<Scalar> & in_values)
const 2405 using Teuchos::rcp_dynamic_cast;
2406 using Teuchos::ptrFromRef;
2408 TEUCHOS_ASSERT(in_names.size()==in_values.size());
2419 paramObj->names = rcp(
new Teuchos::Array<std::string>(in_names));
2420 paramObj->is_distributed =
false;
2423 for(
int i=0;i<in_names.size();i++)
2431 Thyra::locallyReplicatedDefaultSpmdVectorSpace<Scalar>(
2435 Teuchos::ArrayRCP<Scalar> data;
2436 RCP<Thyra::VectorBase<Scalar> > initial_value = Thyra::createMember(paramObj->space);
2437 RCP<Thyra::SpmdVectorBase<Scalar> > vec = rcp_dynamic_cast<Thyra::SpmdVectorBase<Scalar> >(initial_value);
2438 vec->getNonconstLocalData(ptrFromRef(data));
2439 for (
unsigned int i=0; i < paramObj->scalar_value.size(); i++)
2440 data[i] = in_values[i];
2442 paramObj->initial_value = initial_value;
2447 template <
typename Scalar>
2448 Teuchos::RCP<typename panzer::ModelEvaluator<Scalar>::ParameterObject>
2453 const Teuchos::RCP<const GlobalIndexer> & ugi)
const 2460 paramObj->is_distributed =
true;
2461 paramObj->names = rcp(
new Teuchos::Array<std::string>());
2462 paramObj->names->push_back(key);
2463 paramObj->space = vs;
2464 paramObj->initial_value = initial;
2466 paramObj->global_indexer = ugi;
2471 template <
typename Scalar>
2476 for(std::size_t i=0; i < parameters_.size(); i++) {
2479 if(parameters_[i]->is_distributed)
2483 Teuchos::RCP<const Thyra::VectorBase<Scalar> > p = inArgs.get_p(i);
2484 if (p != Teuchos::null) {
2485 for (
unsigned int j=0; j < parameters_[i]->scalar_value.size(); j++) {
2486 parameters_[i]->scalar_value[j].family->setRealValueForAllTypes(Thyra::get_ele(*p,j));
2493 template <
typename Scalar>
2498 for(std::size_t i=0; i < parameters_.size(); i++) {
2501 if(parameters_[i]->is_distributed)
2505 for (
unsigned int j=0; j < parameters_[i]->scalar_value.size(); j++) {
2506 parameters_[i]->scalar_value[j].family->setRealValueForAllTypes(Thyra::get_ele(*(parameters_[i]->initial_value),j));
2512 #endif // __Panzer_ModelEvaluator_impl_hpp__ Teuchos::RCP< const LinearObjFactory< panzer::Traits > > cloneWithNewDomain(const LinearObjFactory< panzer::Traits > &lof, const Teuchos::RCP< const GlobalIndexer > &dUgi)
Clone a linear object factory, but using a different domain.
Interface for constructing a BCStrategy_TemplateManager.
virtual void evalModelImpl_basic_dfdp_scalar(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Allocates and initializes an equation set template manager.
bool evaluate_transient_terms
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int i) const override
Teuchos::RCP< panzer::ResponseLibrary< panzer::Traits > > dfdp_rl
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const override
Teuchos::RCP< const GlobalIndexer > global_indexer
void setParameters(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs) const
int addDistributedParameter(const std::string &name, const Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > &vs, const Teuchos::RCP< GlobalEvaluationData > &ged, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &initial, const Teuchos::RCP< const GlobalIndexer > &ugi=Teuchos::null)
void addResponsesToInArgs(panzer::AssemblyEngineInArgs &input_args) const
int addFlexibleResponse(const std::string &responseName, const std::vector< WorksetDescriptor > &wkst_desc, const Teuchos::RCP< ResponseMESupportBuilderBase > &builder)
const std::string & get_g_name(int i) const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_f_space() const override
Teuchos::RCP< ResponseBase > getResponse(const std::string &responseName) const
panzer::AssemblyEngine_TemplateManager< panzer::Traits > ae_tm_
void buildDistroParamDgDp_RL(const Teuchos::RCP< panzer::WorksetContainer > &wc, const std::vector< Teuchos::RCP< panzer::PhysicsBlock > > &physicsBlocks, const std::vector< panzer::BC > &bcs, const panzer::EquationSetFactory &eqset_factory, const panzer::BCStrategyFactory &bc_factory, const panzer::ClosureModelFactory_TemplateManager< panzer::Traits > &cm_factory, const Teuchos::ParameterList &closure_models, const Teuchos::ParameterList &user_data, const bool write_graphviz_file=false, const std::string &graphviz_file_prefix="")
virtual void set_x_th(const Teuchos::RCP< Thyra::VectorBase< ScalarT > > &in)=0
void evalModel_D2gDpDx(int rIndex, int pIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_x, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &D2gDpDx) const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_DfDp_op(int i) const override
void evalModel_D2gDx2(int rIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_x, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &D2gDx2) const
void buildVolumeFieldManagers(const bool value)
Thyra::ModelEvaluatorBase::InArgs< Scalar > createInArgs() const override
void setupAssemblyInArgs(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, panzer::AssemblyEngineInArgs &ae_inargs) const
virtual void evalModelImpl_basic_dgdp_distro(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_g_space(int i) const override
void evalModel_D2fDp2(int pIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_x, const Teuchos::RCP< Thyra::LinearOpBase< Scalar > > &D2fDp2) const
void initializeNominalValues() const
Initialize the nominal values with good starting conditions.
virtual void evalModelImpl_basic_dgdp_scalar(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Sacado::ScalarParameterVector< panzer::EvaluationTraits > ParamVec
Teuchos::RCP< panzer::LinearObjContainer > ghostedContainer_
Teuchos::RCP< Epetra_MultiVector > getMultiVector(const std::string &modelEvalDescription, const ModelEvaluator::Derivative &deriv, const std::string &derivName, const ModelEvaluator::EDerivativeMultiVectorOrientation mvOrientation)
std::vector< double > gather_seeds
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_W_op() const override
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > x_space_
void evaluate(const panzer::AssemblyEngineInArgs &input_args)
void evalModel_D2gDp2(int rIndex, int pIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_x, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &D2gDp2) const
Teuchos::RCP< panzer::LinearObjContainer > container_
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > f_space_
bool required_basic_dfdp_scalar(const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Are derivatives of the residual with respect to the scalar parameters in the out args? DfDp.
Teuchos::ArrayView< const std::string > get_g_names(int i) const override
Thyra::ModelEvaluatorBase::OutArgs< Scalar > createOutArgsImpl() const override
void addGlobalEvaluationData(const std::string &key, const Teuchos::RCP< GlobalEvaluationData > &ged)
Teuchos::RCP< const Thyra::LinearOpWithSolveFactoryBase< Scalar > > get_W_factory() const override
PHX::MDField< ScalarT, panzer::Cell, panzer::IP > result
A field that will be used to build up the result of the integral we're performing.
bool required_basic_dfdp_distro(const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Are derivatives of the residual with respect to the distributed parameters in the out args...
Teuchos::RCP< const panzer::LinearObjFactory< panzer::Traits > > lof_
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_p_space(int i) const override
virtual void set_f_th(const Teuchos::RCP< Thyra::VectorBase< ScalarT > > &in)=0
void applyDirichletBCs(const Teuchos::RCP< Thyra::VectorBase< Scalar > > &x, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &f) const
virtual void evalModelImpl(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const override
bool required_basic_g(const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Does this set of out args require a simple response?
bool required_basic_dgdx(const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Are their required responses in the out args? DgDx.
virtual Teuchos::RCP< Thyra::LinearOpBase< ScalarT > > getThyraMatrix() const =0
Get a matrix operator.
Teuchos::RCP< ParameterObject > createScalarParameter(const Teuchos::Array< std::string > &names, const Teuchos::Array< Scalar > &in_values) const
void buildDistroParamDfDp_RL(const Teuchos::RCP< panzer::WorksetContainer > &wc, const std::vector< Teuchos::RCP< panzer::PhysicsBlock > > &physicsBlocks, const std::vector< panzer::BC > &bcs, const panzer::EquationSetFactory &eqset_factory, const panzer::BCStrategyFactory &bc_factory, const panzer::ClosureModelFactory_TemplateManager< panzer::Traits > &cm_factory, const Teuchos::ParameterList &closure_models, const Teuchos::ParameterList &user_data, const bool write_graphviz_file=false, const std::string &graphviz_file_prefix="")
void setupModel(const Teuchos::RCP< panzer::WorksetContainer > &wc, const std::vector< Teuchos::RCP< panzer::PhysicsBlock > > &physicsBlocks, const std::vector< panzer::BC > &bcs, const panzer::EquationSetFactory &eqset_factory, const panzer::BCStrategyFactory &bc_factory, const panzer::ClosureModelFactory_TemplateManager< panzer::Traits > &volume_cm_factory, const panzer::ClosureModelFactory_TemplateManager< panzer::Traits > &bc_cm_factory, const Teuchos::ParameterList &closure_models, const Teuchos::ParameterList &user_data, bool writeGraph=false, const std::string &graphPrefix="", const Teuchos::ParameterList &me_params=Teuchos::ParameterList())
void evalModel_D2fDpDx(int pIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_x, const Teuchos::RCP< Thyra::LinearOpBase< Scalar > > &D2fDpDx) const
void evalModel_D2gDxDp(int rIndex, int pIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_p, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &D2gDxDp) const
void evalModel_D2fDx2(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_x, const Teuchos::RCP< Thyra::LinearOpBase< Scalar > > &D2fDx2) const
Thyra::ModelEvaluatorBase::InArgs< Scalar > getNominalValues() const override
Teuchos::RCP< ParameterObject > createDistributedParameter(const std::string &key, const Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > &vs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &initial, const Teuchos::RCP< const GlobalIndexer > &ugi) const
virtual void evalModelImpl_basic_g(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Construct a simple response dicatated by this set of out args.
std::string first_sensitivities_name
void setOneTimeDirichletBeta(const Scalar &beta) const
bool apply_dirichlet_beta
Teuchos::RCP< panzer::ResponseLibrary< panzer::Traits > > responseLibrary_
std::string second_sensitivities_name
void addNonParameterGlobalEvaluationData(const std::string &name, const Teuchos::RCP< GlobalEvaluationData > &ged)
void resetParameters() const
virtual void evalModelImpl_basic(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Evaluate a simple model, meaning a residual and a jacobian, no fancy stochastic galerkin or multipoin...
bool required_basic_dgdp_scalar(const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Are their required responses in the out args? DgDp.
virtual void evalModelImpl_basic_dfdp_distro(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
void registerScalarParameter(const std::string name, panzer::ParamLib &pl, double realValue)
virtual void evalModelImpl_basic_dgdx(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
bool required_basic_dgdp_distro(const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Are their required responses in the out args? DgDp.
void buildBCFieldManagers(const bool value)
int addParameter(const std::string &name, const Scalar &initial)
void evalModel_D2fDxDp(int pIndex, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &delta_p, const Teuchos::RCP< Thyra::LinearOpBase< Scalar > > &D2fDxDp) const
virtual void evalModelImpl_basic_dfdp_scalar_fd(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
virtual void setOwnedVector(const Teuchos::RCP< const Thyra::VectorBase< double > > &ownedVector)=0
Set the owned vector.