Albany: a Trilinos-based PDE code: Albany_FieldManagerScalarResponseFunction.cpp Source File

Go to the documentation of this file.
00001 //*****************************************************************//
00002 //    Albany 2.0:  Copyright 2012 Sandia Corporation               //
00003 //    This Software is released under the BSD license detailed     //
00004 //    in the file "license.txt" in the top-level Albany directory  //
00005 //*****************************************************************//
00006 
00007 #include "Albany_FieldManagerScalarResponseFunction.hpp"
00008 #include <algorithm>
00009 
00010 Albany::FieldManagerScalarResponseFunction::
00011 FieldManagerScalarResponseFunction(
00012   const Teuchos::RCP<Albany::Application>& application_,
00013   const Teuchos::RCP<Albany::AbstractProblem>& problem_,
00014   const Teuchos::RCP<Albany::MeshSpecsStruct>&  meshSpecs_,
00015   const Teuchos::RCP<Albany::StateManager>& stateMgr_,
00016   Teuchos::ParameterList& responseParams) :
00017   ScalarResponseFunction(application_->getComm()),
00018   application(application_),
00019   problem(problem_),
00020   meshSpecs(meshSpecs_),
00021   stateMgr(stateMgr_)
00022 {
00023   setup(responseParams);
00024 }
00025 
00026 Albany::FieldManagerScalarResponseFunction::
00027 FieldManagerScalarResponseFunction(
00028   const Teuchos::RCP<Albany::Application>& application_,
00029   const Teuchos::RCP<Albany::AbstractProblem>& problem_,
00030   const Teuchos::RCP<Albany::MeshSpecsStruct>&  meshSpecs_,
00031   const Teuchos::RCP<Albany::StateManager>& stateMgr_) :
00032   ScalarResponseFunction(application_->getComm()),
00033   application(application_),
00034   problem(problem_),
00035   meshSpecs(meshSpecs_),
00036   stateMgr(stateMgr_)
00037 {
00038 }
00039 
00040 void
00041 Albany::FieldManagerScalarResponseFunction::
00042 setup(Teuchos::ParameterList& responseParams)
00043 {
00044   Teuchos::RCP<const Epetra_Comm> comm = application->getComm();
00045 
00046   // Create field manager
00047   rfm = Teuchos::rcp(new PHX::FieldManager<PHAL::AlbanyTraits>);
00048     
00049   // Create evaluators for field manager
00050   Teuchos::Array< Teuchos::RCP<const PHX::FieldTag> > tags = 
00051     problem->buildEvaluators(*rfm, *meshSpecs, *stateMgr, 
00052            BUILD_RESPONSE_FM,
00053            Teuchos::rcp(&responseParams,false));
00054   int rank = tags[0]->dataLayout().rank();
00055   num_responses = tags[0]->dataLayout().dimension(rank-1);
00056   if (num_responses == 0)
00057     num_responses = 1;
00058   
00059   // Do post-registration setup
00060   rfm->postRegistrationSetup("");
00061 
00062   // Visualize rfm graph -- get file name from name of response function
00063   // (with spaces replaced by _ and lower case)
00064   vis_response_graph = 
00065     responseParams.get("Phalanx Graph Visualization Detail", 0);
00066   vis_response_name = responseParams.get<std::string>("Name");
00067   std::replace(vis_response_name.begin(), vis_response_name.end(), ' ', '_');
00068   std::transform(vis_response_name.begin(), vis_response_name.end(), 
00069      vis_response_name.begin(), ::tolower);
00070 }
00071 
00072 Albany::FieldManagerScalarResponseFunction::
00073 ~FieldManagerScalarResponseFunction()
00074 {
00075 }
00076 
00077 unsigned int
00078 Albany::FieldManagerScalarResponseFunction::
00079 numResponses() const 
00080 {
00081   return num_responses;
00082 }
00083 
00084 void
00085 Albany::FieldManagerScalarResponseFunction::
00086 evaluateResponse(const double current_time,
00087      const Epetra_Vector* xdot,
00088      const Epetra_Vector* xdotdot,
00089      const Epetra_Vector& x,
00090      const Teuchos::Array<ParamVec>& p,
00091      Epetra_Vector& g)
00092 {
00093   visResponseGraph<PHAL::AlbanyTraits::Residual>("");
00094 
00095   // Set data in Workset struct
00096   PHAL::Workset workset;
00097   application->setupBasicWorksetInfo(workset, current_time, xdot, xdotdot, &x, p);
00098   workset.g = Teuchos::rcp(&g,false);
00099 
00100   // Perform fill via field manager
00101   int numWorksets = application->getNumWorksets();
00102   rfm->preEvaluate<PHAL::AlbanyTraits::Residual>(workset);
00103   for (int ws=0; ws < numWorksets; ws++) {
00104     application->loadWorksetBucketInfo<PHAL::AlbanyTraits::Residual>(
00105       workset, ws);
00106     rfm->evaluateFields<PHAL::AlbanyTraits::Residual>(workset);
00107   }
00108   rfm->postEvaluate<PHAL::AlbanyTraits::Residual>(workset);
00109 }
00110 
00111 void
00112 Albany::FieldManagerScalarResponseFunction::
00113 evaluateTangent(const double alpha, 
00114     const double beta,
00115     const double omega,
00116     const double current_time,
00117     bool sum_derivs,
00118     const Epetra_Vector* xdot,
00119     const Epetra_Vector* xdotdot,
00120     const Epetra_Vector& x,
00121     const Teuchos::Array<ParamVec>& p,
00122     ParamVec* deriv_p,
00123     const Epetra_MultiVector* Vxdot,
00124     const Epetra_MultiVector* Vxdotdot,
00125     const Epetra_MultiVector* Vx,
00126     const Epetra_MultiVector* Vp,
00127     Epetra_Vector* g,
00128     Epetra_MultiVector* gx,
00129     Epetra_MultiVector* gp)
00130 {
00131   visResponseGraph<PHAL::AlbanyTraits::Tangent>("_tangent");
00132 
00133   // Set data in Workset struct
00134   PHAL::Workset workset;
00135   application->setupTangentWorksetInfo(workset, sum_derivs, 
00136                current_time, xdot, xdotdot, &x, p, 
00137                deriv_p, Vxdot, Vxdotdot, Vx, Vp);
00138   workset.g = Teuchos::rcp(g, false);
00139   workset.dgdx = Teuchos::rcp(gx, false);
00140   workset.dgdp = Teuchos::rcp(gp, false);
00141   
00142   // Perform fill via field manager
00143   int numWorksets = application->getNumWorksets();
00144   rfm->preEvaluate<PHAL::AlbanyTraits::Tangent>(workset);
00145   for (int ws=0; ws < numWorksets; ws++) {
00146     application->loadWorksetBucketInfo<PHAL::AlbanyTraits::Tangent>(
00147       workset, ws);
00148     rfm->evaluateFields<PHAL::AlbanyTraits::Tangent>(workset);
00149   }
00150   rfm->postEvaluate<PHAL::AlbanyTraits::Tangent>(workset);
00151 }
00152 
00153 void
00154 Albany::FieldManagerScalarResponseFunction::
00155 evaluateGradient(const double current_time,
00156      const Epetra_Vector* xdot,
00157      const Epetra_Vector* xdotdot,
00158      const Epetra_Vector& x,
00159      const Teuchos::Array<ParamVec>& p,
00160      ParamVec* deriv_p,
00161      Epetra_Vector* g,
00162      Epetra_MultiVector* dg_dx,
00163      Epetra_MultiVector* dg_dxdot,
00164      Epetra_MultiVector* dg_dxdotdot,
00165      Epetra_MultiVector* dg_dp)
00166 {
00167   visResponseGraph<PHAL::AlbanyTraits::Jacobian>("_gradient");
00168 
00169   // Set data in Workset struct
00170   PHAL::Workset workset;
00171   application->setupBasicWorksetInfo(workset, current_time, xdot, xdotdot, &x, p);
00172   workset.g = Teuchos::rcp(g, false);
00173   
00174   // Perform fill via field manager (dg/dx)
00175   int numWorksets = application->getNumWorksets();
00176   if (dg_dx != NULL) {
00177     workset.m_coeff = 0.0;
00178     workset.j_coeff = 1.0;
00179     workset.n_coeff = 0.0;
00180     workset.dgdx = Teuchos::rcp(dg_dx, false);
00181     workset.overlapped_dgdx = 
00182       Teuchos::rcp(new Epetra_MultiVector(workset.x_importer->TargetMap(),
00183             dg_dx->NumVectors()));
00184     rfm->preEvaluate<PHAL::AlbanyTraits::Jacobian>(workset);
00185     for (int ws=0; ws < numWorksets; ws++) {
00186       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::Jacobian>(
00187   workset, ws);
00188       rfm->evaluateFields<PHAL::AlbanyTraits::Jacobian>(workset);
00189     }
00190     rfm->postEvaluate<PHAL::AlbanyTraits::Jacobian>(workset);
00191   }
00192 
00193   // Perform fill via field manager (dg/dxdot)
00194   if (dg_dxdot != NULL) {
00195     workset.m_coeff = 1.0;
00196     workset.j_coeff = 0.0;
00197     workset.n_coeff = 0.0;
00198     workset.dgdx = Teuchos::null;
00199     workset.dgdxdot = Teuchos::rcp(dg_dxdot, false);
00200     workset.overlapped_dgdxdot = 
00201       Teuchos::rcp(new Epetra_MultiVector(workset.x_importer->TargetMap(),
00202             dg_dxdot->NumVectors()));
00203     rfm->preEvaluate<PHAL::AlbanyTraits::Jacobian>(workset);
00204     for (int ws=0; ws < numWorksets; ws++) {
00205       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::Jacobian>(
00206   workset, ws);
00207       rfm->evaluateFields<PHAL::AlbanyTraits::Jacobian>(workset);
00208     }
00209     rfm->postEvaluate<PHAL::AlbanyTraits::Jacobian>(workset);
00210   }  
00211 
00212   // Perform fill via field manager (dg/dxdotdot)
00213   if (dg_dxdotdot != NULL) {
00214     workset.m_coeff = 0.0;
00215     workset.j_coeff = 0.0;
00216     workset.n_coeff = 1.0;
00217     workset.dgdx = Teuchos::null;
00218     workset.dgdxdotdot = Teuchos::rcp(dg_dxdotdot, false);
00219     workset.overlapped_dgdxdotdot = 
00220       Teuchos::rcp(new Epetra_MultiVector(workset.x_importer->TargetMap(),
00221             dg_dxdotdot->NumVectors()));
00222     rfm->preEvaluate<PHAL::AlbanyTraits::Jacobian>(workset);
00223     for (int ws=0; ws < numWorksets; ws++) {
00224       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::Jacobian>(
00225   workset, ws);
00226       rfm->evaluateFields<PHAL::AlbanyTraits::Jacobian>(workset);
00227     }
00228     rfm->postEvaluate<PHAL::AlbanyTraits::Jacobian>(workset);
00229   }  
00230 }
00231 #ifdef ALBANY_SG_MP
00232 void
00233 Albany::FieldManagerScalarResponseFunction::
00234 evaluateSGResponse(
00235   const double current_time,
00236   const Stokhos::EpetraVectorOrthogPoly* sg_xdot,
00237   const Stokhos::EpetraVectorOrthogPoly* sg_xdotdot,
00238   const Stokhos::EpetraVectorOrthogPoly& sg_x,
00239   const Teuchos::Array<ParamVec>& p,
00240   const Teuchos::Array<int>& sg_p_index,
00241   const Teuchos::Array< Teuchos::Array<SGType> >& sg_p_vals,
00242   Stokhos::EpetraVectorOrthogPoly& sg_g)
00243 {
00244   visResponseGraph<PHAL::AlbanyTraits::SGResidual>("_sg");
00245 
00246   // Set data in Workset struct
00247   PHAL::Workset workset;
00248   application->setupBasicWorksetInfo(workset, current_time, sg_xdot, sg_xdotdot, &sg_x, p,
00249              sg_p_index, sg_p_vals);
00250   workset.sg_g = Teuchos::rcp(&sg_g,false);
00251 
00252   // Perform fill via field manager
00253   int numWorksets = application->getNumWorksets();
00254   rfm->preEvaluate<PHAL::AlbanyTraits::SGResidual>(workset);
00255   for (int ws=0; ws < numWorksets; ws++) {
00256     application->loadWorksetBucketInfo<PHAL::AlbanyTraits::SGResidual>(
00257       workset, ws);
00258     rfm->evaluateFields<PHAL::AlbanyTraits::SGResidual>(workset);
00259   }
00260   rfm->postEvaluate<PHAL::AlbanyTraits::SGResidual>(workset);
00261 }
00262 
00263 void
00264 Albany::FieldManagerScalarResponseFunction::
00265 evaluateSGTangent(
00266   const double alpha, 
00267   const double beta, 
00268   const double omega, 
00269   const double current_time,
00270   bool sum_derivs,
00271   const Stokhos::EpetraVectorOrthogPoly* sg_xdot,
00272   const Stokhos::EpetraVectorOrthogPoly* sg_xdotdot,
00273   const Stokhos::EpetraVectorOrthogPoly& sg_x,
00274   const Teuchos::Array<ParamVec>& p,
00275   const Teuchos::Array<int>& sg_p_index,
00276   const Teuchos::Array< Teuchos::Array<SGType> >& sg_p_vals,
00277   ParamVec* deriv_p,
00278   const Epetra_MultiVector* Vx,
00279   const Epetra_MultiVector* Vxdot,
00280   const Epetra_MultiVector* Vxdotdot,
00281   const Epetra_MultiVector* Vp,
00282   Stokhos::EpetraVectorOrthogPoly* sg_g,
00283   Stokhos::EpetraMultiVectorOrthogPoly* sg_JV,
00284   Stokhos::EpetraMultiVectorOrthogPoly* sg_gp)
00285 {
00286   visResponseGraph<PHAL::AlbanyTraits::SGTangent>("_sg_tangent");
00287 
00288   // Set data in Workset struct
00289   PHAL::Workset workset;
00290   application->setupTangentWorksetInfo(workset, current_time, sum_derivs, 
00291                sg_xdot, sg_xdotdot, &sg_x, p, deriv_p, 
00292                sg_p_index, sg_p_vals,
00293                Vxdot, Vxdotdot, Vx, Vp);
00294   workset.sg_g = Teuchos::rcp(sg_g, false);
00295   workset.sg_dgdx = Teuchos::rcp(sg_JV, false);
00296   workset.sg_dgdp = Teuchos::rcp(sg_gp, false);
00297   
00298   // Perform fill via field manager
00299   int numWorksets = application->getNumWorksets();
00300   rfm->preEvaluate<PHAL::AlbanyTraits::SGTangent>(workset);
00301   for (int ws=0; ws < numWorksets; ws++) {
00302     application->loadWorksetBucketInfo<PHAL::AlbanyTraits::SGTangent>(
00303       workset, ws);
00304     rfm->evaluateFields<PHAL::AlbanyTraits::SGTangent>(workset);
00305   }
00306   rfm->postEvaluate<PHAL::AlbanyTraits::SGTangent>(workset);
00307 }
00308 
00309 void
00310 Albany::FieldManagerScalarResponseFunction::
00311 evaluateSGGradient(
00312   const double current_time,
00313   const Stokhos::EpetraVectorOrthogPoly* sg_xdot,
00314   const Stokhos::EpetraVectorOrthogPoly* sg_xdotdot,
00315   const Stokhos::EpetraVectorOrthogPoly& sg_x,
00316   const Teuchos::Array<ParamVec>& p,
00317   const Teuchos::Array<int>& sg_p_index,
00318   const Teuchos::Array< Teuchos::Array<SGType> >& sg_p_vals,
00319   ParamVec* deriv_p,
00320   Stokhos::EpetraVectorOrthogPoly* sg_g,
00321   Stokhos::EpetraMultiVectorOrthogPoly* sg_dg_dx,
00322   Stokhos::EpetraMultiVectorOrthogPoly* sg_dg_dxdot,
00323   Stokhos::EpetraMultiVectorOrthogPoly* sg_dg_dxdotdot,
00324   Stokhos::EpetraMultiVectorOrthogPoly* sg_dg_dp)
00325 {
00326   visResponseGraph<PHAL::AlbanyTraits::SGJacobian>("_sg_gradient");
00327 
00328   // Set data in Workset struct
00329   PHAL::Workset workset;
00330   application->setupBasicWorksetInfo(workset, current_time, sg_xdot, sg_xdotdot, &sg_x, p,
00331              sg_p_index, sg_p_vals);
00332   workset.sg_g = Teuchos::rcp(sg_g, false);
00333   
00334   // Perform fill via field manager (dg/dx)
00335   int numWorksets = application->getNumWorksets();
00336   if (sg_dg_dx != NULL) {
00337     workset.m_coeff = 0.0;
00338     workset.j_coeff = 1.0;
00339     workset.n_coeff = 0.0;
00340     workset.sg_dgdx = Teuchos::rcp(sg_dg_dx, false);
00341     workset.overlapped_sg_dgdx = 
00342       Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
00343          sg_dg_dx->basis(),
00344          sg_dg_dx->map(),
00345          Teuchos::rcp(&(workset.x_importer->TargetMap()),false),
00346          sg_dg_dx->productComm(),
00347          sg_dg_dx->numVectors()));
00348     rfm->preEvaluate<PHAL::AlbanyTraits::SGJacobian>(workset);
00349     for (int ws=0; ws < numWorksets; ws++) {
00350       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::SGJacobian>(
00351   workset, ws);
00352       rfm->evaluateFields<PHAL::AlbanyTraits::SGJacobian>(workset);
00353     }
00354     rfm->postEvaluate<PHAL::AlbanyTraits::SGJacobian>(workset);
00355   }
00356 
00357   // Perform fill via field manager (dg/dxdot)
00358   if (sg_dg_dxdot != NULL) {
00359     workset.m_coeff = 1.0;
00360     workset.j_coeff = 0.0;
00361     workset.n_coeff = 0.0;
00362     workset.sg_dgdx = Teuchos::null;
00363     workset.sg_dgdxdot = Teuchos::rcp(sg_dg_dxdot, false);
00364     workset.overlapped_sg_dgdx = Teuchos::null;
00365     workset.overlapped_sg_dgdxdot = 
00366       Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
00367          sg_dg_dxdot->basis(),
00368          sg_dg_dxdot->map(),
00369          Teuchos::rcp(&(workset.x_importer->TargetMap()),false),
00370          sg_dg_dxdot->productComm(),
00371          sg_dg_dxdot->numVectors()));
00372     rfm->preEvaluate<PHAL::AlbanyTraits::SGJacobian>(workset);
00373     for (int ws=0; ws < numWorksets; ws++) {
00374       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::SGJacobian>(
00375   workset, ws);
00376       rfm->evaluateFields<PHAL::AlbanyTraits::SGJacobian>(workset);
00377     }
00378     rfm->postEvaluate<PHAL::AlbanyTraits::SGJacobian>(workset);
00379   }  
00380 
00381   // Perform fill via field manager (dg/dxdotdot)
00382   if (sg_dg_dxdotdot != NULL) {
00383     workset.m_coeff = 0.0;
00384     workset.j_coeff = 0.0;
00385     workset.n_coeff = 1.0;
00386     workset.sg_dgdx = Teuchos::null;
00387     workset.sg_dgdxdotdot = Teuchos::rcp(sg_dg_dxdotdot, false);
00388     workset.overlapped_sg_dgdx = Teuchos::null;
00389     workset.overlapped_sg_dgdxdotdot = 
00390       Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
00391          sg_dg_dxdotdot->basis(),
00392          sg_dg_dxdotdot->map(),
00393          Teuchos::rcp(&(workset.x_importer->TargetMap()),false),
00394          sg_dg_dxdotdot->productComm(),
00395          sg_dg_dxdotdot->numVectors()));
00396     rfm->preEvaluate<PHAL::AlbanyTraits::SGJacobian>(workset);
00397     for (int ws=0; ws < numWorksets; ws++) {
00398       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::SGJacobian>(
00399   workset, ws);
00400       rfm->evaluateFields<PHAL::AlbanyTraits::SGJacobian>(workset);
00401     }
00402     rfm->postEvaluate<PHAL::AlbanyTraits::SGJacobian>(workset);
00403   }  
00404 }
00405 
00406 void
00407 Albany::FieldManagerScalarResponseFunction::
00408 evaluateMPResponse(
00409   const double current_time,
00410   const Stokhos::ProductEpetraVector* mp_xdot,
00411   const Stokhos::ProductEpetraVector* mp_xdotdot,
00412   const Stokhos::ProductEpetraVector& mp_x,
00413   const Teuchos::Array<ParamVec>& p,
00414   const Teuchos::Array<int>& mp_p_index,
00415   const Teuchos::Array< Teuchos::Array<MPType> >& mp_p_vals,
00416   Stokhos::ProductEpetraVector& mp_g)
00417 {
00418   visResponseGraph<PHAL::AlbanyTraits::MPResidual>("_mp");
00419 
00420   // Set data in Workset struct
00421   PHAL::Workset workset;
00422   application->setupBasicWorksetInfo(workset, current_time, mp_xdot, mp_xdotdot, &mp_x, p,
00423              mp_p_index, mp_p_vals);
00424   workset.mp_g = Teuchos::rcp(&mp_g,false);
00425 
00426   // Perform fill via field manager
00427   int numWorksets = application->getNumWorksets();
00428   rfm->preEvaluate<PHAL::AlbanyTraits::MPResidual>(workset);
00429   for (int ws=0; ws < numWorksets; ws++) {
00430     application->loadWorksetBucketInfo<PHAL::AlbanyTraits::MPResidual>(
00431       workset, ws);
00432     rfm->evaluateFields<PHAL::AlbanyTraits::MPResidual>(workset);
00433   }
00434   rfm->postEvaluate<PHAL::AlbanyTraits::MPResidual>(workset);
00435 }
00436 
00437 void
00438 Albany::FieldManagerScalarResponseFunction::
00439 evaluateMPTangent(
00440   const double alpha, 
00441   const double beta, 
00442   const double omega, 
00443   const double current_time,
00444   bool sum_derivs,
00445   const Stokhos::ProductEpetraVector* mp_xdot,
00446   const Stokhos::ProductEpetraVector* mp_xdotdot,
00447   const Stokhos::ProductEpetraVector& mp_x,
00448   const Teuchos::Array<ParamVec>& p,
00449   const Teuchos::Array<int>& mp_p_index,
00450   const Teuchos::Array< Teuchos::Array<MPType> >& mp_p_vals,
00451   ParamVec* deriv_p,
00452   const Epetra_MultiVector* Vx,
00453   const Epetra_MultiVector* Vxdot,
00454   const Epetra_MultiVector* Vxdotdot,
00455   const Epetra_MultiVector* Vp,
00456   Stokhos::ProductEpetraVector* mp_g,
00457   Stokhos::ProductEpetraMultiVector* mp_JV,
00458   Stokhos::ProductEpetraMultiVector* mp_gp)
00459 {
00460   visResponseGraph<PHAL::AlbanyTraits::MPTangent>("_mp_tangent");
00461 
00462   // Set data in Workset struct
00463   PHAL::Workset workset;
00464   application->setupTangentWorksetInfo(workset, current_time, sum_derivs, 
00465                mp_xdot, mp_xdotdot, &mp_x, p, deriv_p, 
00466                mp_p_index, mp_p_vals,
00467                Vxdot, Vxdotdot, Vx, Vp);
00468   workset.mp_g = Teuchos::rcp(mp_g, false);
00469   workset.mp_dgdx = Teuchos::rcp(mp_JV, false);
00470   workset.mp_dgdp = Teuchos::rcp(mp_gp, false);
00471   
00472   // Perform fill via field manager
00473   int numWorksets = application->getNumWorksets();
00474   rfm->preEvaluate<PHAL::AlbanyTraits::MPTangent>(workset);
00475   for (int ws=0; ws < numWorksets; ws++) {
00476     application->loadWorksetBucketInfo<PHAL::AlbanyTraits::MPTangent>(
00477       workset, ws);
00478     rfm->evaluateFields<PHAL::AlbanyTraits::MPTangent>(workset);
00479   }
00480   rfm->postEvaluate<PHAL::AlbanyTraits::MPTangent>(workset);
00481 }
00482 
00483 void
00484 Albany::FieldManagerScalarResponseFunction::
00485 evaluateMPGradient(
00486   const double current_time,
00487   const Stokhos::ProductEpetraVector* mp_xdot,
00488   const Stokhos::ProductEpetraVector* mp_xdotdot,
00489   const Stokhos::ProductEpetraVector& mp_x,
00490   const Teuchos::Array<ParamVec>& p,
00491   const Teuchos::Array<int>& mp_p_index,
00492   const Teuchos::Array< Teuchos::Array<MPType> >& mp_p_vals,
00493   ParamVec* deriv_p,
00494   Stokhos::ProductEpetraVector* mp_g,
00495   Stokhos::ProductEpetraMultiVector* mp_dg_dx,
00496   Stokhos::ProductEpetraMultiVector* mp_dg_dxdot,
00497   Stokhos::ProductEpetraMultiVector* mp_dg_dxdotdot,
00498   Stokhos::ProductEpetraMultiVector* mp_dg_dp)
00499 {
00500   visResponseGraph<PHAL::AlbanyTraits::MPJacobian>("_mp_gradient");
00501 
00502   // Set data in Workset struct
00503   PHAL::Workset workset;
00504   application->setupBasicWorksetInfo(workset, current_time, mp_xdot, mp_xdotdot, &mp_x, p,
00505              mp_p_index, mp_p_vals);
00506   workset.mp_g = Teuchos::rcp(mp_g, false);
00507   
00508   // Perform fill via field manager (dg/dx)
00509   int numWorksets = application->getNumWorksets();
00510   if (mp_dg_dx != NULL) {
00511     workset.m_coeff = 0.0;
00512     workset.j_coeff = 1.0;
00513     workset.n_coeff = 0.0;
00514     workset.mp_dgdx = Teuchos::rcp(mp_dg_dx, false);
00515     workset.overlapped_mp_dgdx = 
00516       Teuchos::rcp(new Stokhos::ProductEpetraMultiVector(
00517          mp_dg_dx->map(),
00518          Teuchos::rcp(&(workset.x_importer->TargetMap()),false),
00519          mp_dg_dx->productComm(),
00520          mp_dg_dx->numVectors()));
00521     rfm->preEvaluate<PHAL::AlbanyTraits::MPJacobian>(workset);
00522     for (int ws=0; ws < numWorksets; ws++) {
00523       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::MPJacobian>(
00524   workset, ws);
00525       rfm->evaluateFields<PHAL::AlbanyTraits::MPJacobian>(workset);
00526     }
00527     rfm->postEvaluate<PHAL::AlbanyTraits::MPJacobian>(workset);
00528   }
00529 
00530   // Perform fill via field manager (dg/dxdot)
00531   if (mp_dg_dxdot != NULL) {
00532     workset.m_coeff = 1.0;
00533     workset.j_coeff = 0.0;
00534     workset.n_coeff = 0.0;
00535     workset.mp_dgdx = Teuchos::null;
00536     workset.mp_dgdxdot = Teuchos::rcp(mp_dg_dxdot, false);
00537     workset.overlapped_mp_dgdx = Teuchos::null;
00538     workset.overlapped_mp_dgdxdot = 
00539       Teuchos::rcp(new Stokhos::ProductEpetraMultiVector(
00540          mp_dg_dxdot->map(),
00541          Teuchos::rcp(&(workset.x_importer->TargetMap()),false),
00542          mp_dg_dxdot->productComm(),
00543          mp_dg_dxdot->numVectors()));
00544     rfm->preEvaluate<PHAL::AlbanyTraits::MPJacobian>(workset);
00545     for (int ws=0; ws < numWorksets; ws++) {
00546       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::MPJacobian>(
00547   workset, ws);
00548       rfm->evaluateFields<PHAL::AlbanyTraits::MPJacobian>(workset);
00549     }
00550     rfm->postEvaluate<PHAL::AlbanyTraits::MPJacobian>(workset);
00551   }  
00552 
00553   // Perform fill via field manager (dg/dxdotdot)
00554   if (mp_dg_dxdotdot != NULL) {
00555     workset.m_coeff = 0.0;
00556     workset.j_coeff = 0.0;
00557     workset.n_coeff = 1.0;
00558     workset.mp_dgdx = Teuchos::null;
00559     workset.mp_dgdxdotdot = Teuchos::rcp(mp_dg_dxdotdot, false);
00560     workset.overlapped_mp_dgdx = Teuchos::null;
00561     workset.overlapped_mp_dgdxdotdot = 
00562       Teuchos::rcp(new Stokhos::ProductEpetraMultiVector(
00563          mp_dg_dxdotdot->map(),
00564          Teuchos::rcp(&(workset.x_importer->TargetMap()),false),
00565          mp_dg_dxdotdot->productComm(),
00566          mp_dg_dxdotdot->numVectors()));
00567     rfm->preEvaluate<PHAL::AlbanyTraits::MPJacobian>(workset);
00568     for (int ws=0; ws < numWorksets; ws++) {
00569       application->loadWorksetBucketInfo<PHAL::AlbanyTraits::MPJacobian>(
00570   workset, ws);
00571       rfm->evaluateFields<PHAL::AlbanyTraits::MPJacobian>(workset);
00572     }
00573     rfm->postEvaluate<PHAL::AlbanyTraits::MPJacobian>(workset);
00574   }  
00575 }
00576 #endif //ALBANY_SG_MP