43 #ifndef __Panzer_GatherSolution_Epetra_Hessian_impl_hpp__ 44 #define __Panzer_GatherSolution_Epetra_Hessian_impl_hpp__ 47 #ifdef Panzer_BUILD_HESSIAN_SUPPORT 49 #include "Teuchos_Assert.hpp" 50 #include "Phalanx_DataLayout.hpp" 59 #include "Teuchos_FancyOStream.hpp" 61 #include "Epetra_Vector.h" 62 #include "Epetra_Map.h" 70 template<
typename TRAITS,
typename LO,
typename GO>
75 : globalIndexer_(indexer)
77 typedef std::vector< std::vector<std::string> > vvstring;
79 GatherSolution_Input
input;
80 input.setParameterList(p);
82 const std::vector<std::string> & names =
input.getDofNames();
85 indexerNames_ =
input.getIndexerNames();
86 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
87 globalDataKey_ =
input.getGlobalDataKey();
89 gatherSeedIndex_ =
input.getGatherSeedIndex();
90 sensitivitiesName_ =
input.getSensitivitiesName();
91 firstSensitivitiesAvailable_ =
input.firstSensitivitiesAvailable();
93 secondSensitivitiesAvailable_ =
input.secondSensitivitiesAvailable();
94 sensitivities2ndPrefix_ =
input.getSecondSensitivityDataKeyPrefix();
97 gatherFields_.resize(names.size());
98 for (std::size_t fd = 0; fd < names.size(); ++fd) {
100 gatherFields_[fd] = f;
101 this->addEvaluatedField(gatherFields_[fd]);
105 std::string firstName =
"<none>";
107 firstName = names[0];
110 if(!firstSensitivitiesAvailable_) {
111 std::string n =
"GatherSolution (Epetra, No First Sensitivities): "+firstName+
" (Hessian)";
115 std::string n =
"GatherSolution (Epetra): "+firstName+
" ("+PHX::typeAsString<EvalT>()+
") ";
121 template<
typename TRAITS,
typename LO,
typename GO>
129 fieldIds_.resize(gatherFields_.size());
131 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
133 const std::string& fieldName = indexerNames_[fd];
134 fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
137 if(fieldIds_[fd]==-1) {
139 "GatherSolution_Epetra<Hessian>: Could not find field \"" + fieldName +
"\" in the global indexer. ");
144 this->utils.setFieldData(gatherFields_[fd],fm);
147 indexerNames_.clear();
151 template<
typename TRAITS,
typename LO,
typename GO>
156 using Teuchos::rcp_dynamic_cast;
160 if(firstSensitivitiesAvailable_) {
161 if(d.first_sensitivities_name==sensitivitiesName_)
162 firstApplySensitivities_ =
true;
164 firstApplySensitivities_ =
false;
167 firstApplySensitivities_ =
false;
169 if(secondSensitivitiesAvailable_) {
170 if(d.second_sensitivities_name==sensitivitiesName_)
171 secondApplySensitivities_ =
true;
173 secondApplySensitivities_ =
false;
176 secondApplySensitivities_ =
false;
180 RCP<GlobalEvaluationData> ged;
183 std::string post = useTimeDerivativeSolutionVector_ ?
" - Xdot" :
" - X";
184 if(d.gedc.containsDataObject(globalDataKey_+post)) {
185 ged = d.gedc.getDataObject(globalDataKey_+post);
187 RCP<EpetraVector_ReadOnly_GlobalEvaluationData> ro_ged = rcp_dynamic_cast<EpetraVector_ReadOnly_GlobalEvaluationData>(ged,
true);
189 x_ = ro_ged->getGhostedVector_Epetra();
191 else if(d.gedc.containsDataObject(globalDataKey_)) {
193 ged = d.gedc.getDataObject(globalDataKey_);
196 RCP<EpetraVector_ReadOnly_GlobalEvaluationData> ro_ged = rcp_dynamic_cast<EpetraVector_ReadOnly_GlobalEvaluationData>(ged);
197 RCP<EpetraLinearObjContainer> epetraContainer = rcp_dynamic_cast<EpetraLinearObjContainer>(ged);
201 RCP<LOCPair_GlobalEvaluationData> loc_pair = rcp_dynamic_cast<LOCPair_GlobalEvaluationData>(ged);
203 if(loc_pair!=Teuchos::null) {
206 epetraContainer = rcp_dynamic_cast<EpetraLinearObjContainer>(loc);
210 if(ro_ged!=Teuchos::null) {
211 RCP<EpetraVector_ReadOnly_GlobalEvaluationData> ro_ged = rcp_dynamic_cast<EpetraVector_ReadOnly_GlobalEvaluationData>(ged,
true);
213 x_ = ro_ged->getGhostedVector_Epetra();
215 else if(epetraContainer!=Teuchos::null) {
216 if (useTimeDerivativeSolutionVector_)
217 x_ = epetraContainer->get_dxdt();
219 x_ = epetraContainer->get_x();
227 if(!secondApplySensitivities_) {
236 if(d.gedc.containsDataObject(sensitivities2ndPrefix_+globalDataKey_)) {
237 ged = d.gedc.getDataObject(sensitivities2ndPrefix_+globalDataKey_);
239 RCP<EpetraVector_ReadOnly_GlobalEvaluationData> ro_ged = rcp_dynamic_cast<EpetraVector_ReadOnly_GlobalEvaluationData>(ged,
true);
241 dx_ = ro_ged->getGhostedVector_Epetra();
246 "Cannot find sensitivity vector associated with \""+sensitivities2ndPrefix_+globalDataKey_+
"\" and \""+post+
"\"");
250 template<
typename TRAITS,
typename LO,
typename GO>
255 std::string blockId = this->wda(workset).block_id;
256 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
258 double seed_value = 0.0;
259 if(firstApplySensitivities_) {
260 if (useTimeDerivativeSolutionVector_ && gatherSeedIndex_<0) {
261 seed_value = workset.alpha;
263 else if (gatherSeedIndex_<0) {
264 seed_value = workset.beta;
266 else if(!useTimeDerivativeSolutionVector_) {
267 seed_value = workset.gather_seeds[gatherSeedIndex_];
279 if(!firstApplySensitivities_)
286 if (this->wda.getDetailsIndex() == 1) {
288 dos = globalIndexer_->getElementBlockGIDCount(workset.details(0).block_id);
292 for(std::size_t fieldIndex=0;
293 fieldIndex<gatherFields_.size();fieldIndex++) {
294 PHX::MDField<ScalarT,Cell,NODE> &
field = gatherFields_[fieldIndex];
295 int fieldNum = fieldIds_[fieldIndex];
296 const std::vector<int> & elmtOffset = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
299 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
300 std::size_t cellLocalId = localCellIds[worksetCellIndex];
302 const std::vector<int> & LIDs = globalIndexer_->getElementLIDs(cellLocalId);
304 if(!firstApplySensitivities_) {
327 if(secondApplySensitivities_) {
331 field(worksetCellIndex,
basis).val().fastAccessDx(0) = (*dx_)[lid];
Gathers solution values from the Newton solution vector into the nodal fields of the field manager...
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
PHX::MDField< const ScalarT > input
PHX::MDField< const ScalarT, Cell, IP > field
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
#define TEUCHOS_ASSERT(assertion_test)
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis> or <Cell,Basis>