43 #ifndef PANZER_GATHER_SOLUTION_BLOCKED_EPETRA_IMPL_HPP 44 #define PANZER_GATHER_SOLUTION_BLOCKED_EPETRA_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 52 #include "Panzer_EpetraLinearObjFactory.hpp" 57 #include "Teuchos_FancyOStream.hpp" 59 #include "Thyra_SpmdVectorBase.hpp" 60 #include "Thyra_ProductVectorBase.hpp" 66 template<
typename TRAITS,
typename LO,
typename GO>
69 const Teuchos::ParameterList& p)
71 , has_tangent_fields_(false)
73 typedef std::vector< std::vector<std::string> > vvstring;
76 input.setParameterList(p);
78 const std::vector<std::string> & names =
input.getDofNames();
79 Teuchos::RCP<const panzer::PureBasis>
basis =
input.getBasis();
80 const vvstring & tangent_field_names =
input.getTangentNames();
82 indexerNames_ =
input.getIndexerNames();
83 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
84 globalDataKey_ =
input.getGlobalDataKey();
87 gatherFields_.resize(names.size());
88 for (std::size_t fd = 0; fd < names.size(); ++fd) {
90 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
91 this->addEvaluatedField(gatherFields_[fd]);
95 if (tangent_field_names.size()>0) {
96 TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
98 has_tangent_fields_ =
true;
99 tangentFields_.resize(gatherFields_.size());
100 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
101 tangentFields_[fd].resize(tangent_field_names[fd].size());
102 for (std::size_t i=0; i<tangent_field_names[fd].size(); ++i) {
103 tangentFields_[fd][i] =
104 PHX::MDField<const ScalarT,Cell,NODE>(tangent_field_names[fd][i],
basis->functional);
105 this->addDependentField(tangentFields_[fd][i]);
111 std::string firstName =
"<none>";
113 firstName = names[0];
115 std::string n =
"GatherSolution (BlockedEpetra): "+firstName+
" ()";
120 template<
typename TRAITS,
typename LO,
typename GO>
125 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
127 indexerIds_.resize(gatherFields_.size());
128 subFieldIds_.resize(gatherFields_.size());
130 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
132 const std::string& fieldName = indexerNames_[fd];
135 subFieldIds_[fd] = indexers_[indexerIds_[fd]]->getFieldNum(fieldName);
137 TEUCHOS_ASSERT(indexerIds_[fd]>=0);
140 this->utils.setFieldData(gatherFields_[fd],fm);
143 if (has_tangent_fields_) {
144 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd)
145 for (std::size_t i=0; i<tangentFields_[fd].size(); ++i)
146 this->utils.setFieldData(tangentFields_[fd][i],fm);
149 indexerNames_.clear();
153 template<
typename TRAITS,
typename LO,
typename GO>
160 using Teuchos::rcp_dynamic_cast;
162 RCP<GlobalEvaluationData> ged;
165 std::string post = useTimeDerivativeSolutionVector_ ?
" - Xdot" :
" - X";
166 if(d.gedc.containsDataObject(globalDataKey_+post)) {
167 ged = d.gedc.getDataObject(globalDataKey_+post);
174 TEUCHOS_TEST_FOR_EXCEPTION(x_==Teuchos::null,std::logic_error,
175 "Gather Residual: Can't find x vector in GEDC \"" << globalDataKey_ <<
"\" (" << post <<
"). " 176 "A cast failed for " << ged <<
". Type is " << Teuchos::typeName(*ged));
181 ged = d.gedc.getDataObject(globalDataKey_);
185 RCP<const BlockedEpetraLinearObjContainer> blockedContainer = rcp_dynamic_cast<
const BLOC>(ged);
187 if(ro_ged!=Teuchos::null) {
192 else if(blockedContainer!=Teuchos::null) {
193 if (useTimeDerivativeSolutionVector_)
200 TEUCHOS_TEST_FOR_EXCEPTION(x_==Teuchos::null,std::logic_error,
201 "Gather Residual: Can't find x vector in GEDC \"" << globalDataKey_ <<
"\" (" << post <<
"). " 202 "A cast failed for " << ged <<
". Type is " << Teuchos::typeName(*ged));
207 template<
typename TRAITS,
typename LO,
typename GO>
212 using Teuchos::ArrayRCP;
213 using Teuchos::ptrFromRef;
214 using Teuchos::rcp_dynamic_cast;
217 using Thyra::SpmdVectorBase;
221 std::string blockId = this->wda(workset).block_id;
222 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
225 Teuchos::ArrayRCP<const double> local_x;
226 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
228 PHX::MDField<ScalarT,Cell,NODE> &
field = gatherFields_[fieldIndex];
230 int indexerId = indexerIds_[fieldIndex];
231 int subFieldNum = subFieldIds_[fieldIndex];
234 Teuchos::ArrayRCP<const double> local_x;
235 rcp_dynamic_cast<SpmdVectorBase<double> >(x_->getNonconstVectorBlock(indexerId))->getLocalData(ptrFromRef(local_x));
237 auto subRowIndexer = indexers_[indexerId];
238 const std::vector<int> & elmtOffset = subRowIndexer->getGIDFieldOffsets(blockId,subFieldNum);
241 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
242 LO cellLocalId = localCellIds[worksetCellIndex];
244 const std::vector<int> & LIDs = subRowIndexer->getElementLIDs(cellLocalId);
264 template<
typename TRAITS,
typename LO,
typename GO>
267 const Teuchos::ParameterList& p)
268 : indexers_(indexers)
269 , has_tangent_fields_(false)
271 typedef std::vector< std::vector<std::string> > vvstring;
274 input.setParameterList(p);
276 const std::vector<std::string> & names =
input.getDofNames();
277 Teuchos::RCP<const panzer::PureBasis>
basis =
input.getBasis();
278 const vvstring & tangent_field_names =
input.getTangentNames();
280 indexerNames_ =
input.getIndexerNames();
281 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
282 globalDataKey_ =
input.getGlobalDataKey();
285 gatherFields_.resize(names.size());
286 for (std::size_t fd = 0; fd < names.size(); ++fd) {
288 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
289 this->addEvaluatedField(gatherFields_[fd]);
293 if (tangent_field_names.size()>0) {
294 TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
296 has_tangent_fields_ =
true;
297 tangentFields_.resize(gatherFields_.size());
298 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
299 tangentFields_[fd].resize(tangent_field_names[fd].size());
300 for (std::size_t i=0; i<tangent_field_names[fd].size(); ++i) {
301 tangentFields_[fd][i] =
302 PHX::MDField<const ScalarT,Cell,NODE>(tangent_field_names[fd][i],
basis->functional);
303 this->addDependentField(tangentFields_[fd][i]);
309 std::string firstName =
"<none>";
311 firstName = names[0];
313 std::string n =
"GatherSolution Tangent (BlockedEpetra): "+firstName+
" ()";
318 template<
typename TRAITS,
typename LO,
typename GO>
323 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
325 indexerIds_.resize(gatherFields_.size());
326 subFieldIds_.resize(gatherFields_.size());
328 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
330 const std::string& fieldName = indexerNames_[fd];
333 subFieldIds_[fd] = indexers_[indexerIds_[fd]]->getFieldNum(fieldName);
335 TEUCHOS_ASSERT(indexerIds_[fd]>=0);
338 this->utils.setFieldData(gatherFields_[fd],fm);
341 if (has_tangent_fields_) {
342 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd)
343 for (std::size_t i=0; i<tangentFields_[fd].size(); ++i)
344 this->utils.setFieldData(tangentFields_[fd][i],fm);
347 indexerNames_.clear();
351 template<
typename TRAITS,
typename LO,
typename GO>
358 using Teuchos::rcp_dynamic_cast;
360 RCP<GlobalEvaluationData> ged;
363 std::string post = useTimeDerivativeSolutionVector_ ?
" - Xdot" :
" - X";
364 if(d.gedc.containsDataObject(globalDataKey_+post)) {
365 ged = d.gedc.getDataObject(globalDataKey_+post);
374 ged = d.gedc.getDataObject(globalDataKey_);
378 RCP<const BlockedEpetraLinearObjContainer> blockedContainer = rcp_dynamic_cast<
const BLOC>(ged);
380 if(ro_ged!=Teuchos::null) {
385 else if(blockedContainer!=Teuchos::null) {
386 if (useTimeDerivativeSolutionVector_)
394 TEUCHOS_ASSERT(x_!=Teuchos::null);
398 template<
typename TRAITS,
typename LO,
typename GO>
403 using Teuchos::ArrayRCP;
404 using Teuchos::ptrFromRef;
405 using Teuchos::rcp_dynamic_cast;
408 using Thyra::SpmdVectorBase;
411 std::vector<std::pair<int,GO> > GIDs;
412 std::vector<int> LIDs;
415 std::string blockId = this->wda(workset).block_id;
416 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
419 Teuchos::ArrayRCP<const double> local_x;
420 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
422 PHX::MDField<ScalarT,Cell,NODE> &
field = gatherFields_[fieldIndex];
424 int indexerId = indexerIds_[fieldIndex];
425 int subFieldNum = subFieldIds_[fieldIndex];
428 Teuchos::ArrayRCP<const double> local_x;
429 rcp_dynamic_cast<SpmdVectorBase<double> >(x_->getNonconstVectorBlock(indexerId))->getLocalData(ptrFromRef(local_x));
431 auto subRowIndexer = indexers_[indexerId];
432 const std::vector<int> & elmtOffset = subRowIndexer->getGIDFieldOffsets(blockId,subFieldNum);
435 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
436 LO cellLocalId = localCellIds[worksetCellIndex];
438 const std::vector<int> & LIDs = subRowIndexer->getElementLIDs(cellLocalId);
448 if (!has_tangent_fields_)
451 field(worksetCellIndex,
basis).val() = local_x[lid];
452 for (std::size_t i=0; i<tangentFields_[fieldIndex].size(); ++i)
453 field(worksetCellIndex,
basis).fastAccessDx(i) =
454 tangentFields_[fieldIndex][i](worksetCellIndex,
basis).val();
465 template<
typename TRAITS,
typename LO,
typename GO>
468 const Teuchos::ParameterList& p)
469 : indexers_(indexers)
474 input.setParameterList(p);
476 const std::vector<std::string> & names =
input.getDofNames();
477 Teuchos::RCP<const panzer::PureBasis>
basis =
input.getBasis();
480 indexerNames_ =
input.getIndexerNames();
481 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
482 globalDataKey_ =
input.getGlobalDataKey();
484 gatherSeedIndex_ =
input.getGatherSeedIndex();
485 sensitivitiesName_ =
input.getSensitivitiesName();
486 disableSensitivities_ = !
input.firstSensitivitiesAvailable();
489 gatherFields_.resize(names.size());
490 for (std::size_t fd = 0; fd < names.size(); ++fd) {
491 PHX::MDField<ScalarT,Cell,NODE> f(names[fd],
basis->functional);
492 gatherFields_[fd] = f;
493 this->addEvaluatedField(gatherFields_[fd]);
497 std::string firstName =
"<none>";
499 firstName = names[0];
502 if(disableSensitivities_) {
503 std::string n =
"GatherSolution (BlockedEpetra, No Sensitivities): "+firstName+
" ("+PHX::typeAsString<EvalT>()+
")";
507 std::string n =
"GatherSolution (BlockedEpetra): "+firstName+
" ("+PHX::typeAsString<EvalT>()+
") ";
513 template<
typename TRAITS,
typename LO,
typename GO>
518 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
520 indexerIds_.resize(gatherFields_.size());
521 subFieldIds_.resize(gatherFields_.size());
523 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
525 const std::string& fieldName = indexerNames_[fd];
528 subFieldIds_[fd] = indexers_[indexerIds_[fd]]->getFieldNum(fieldName);
530 TEUCHOS_ASSERT(indexerIds_[fd]>=0);
533 this->utils.setFieldData(gatherFields_[fd],fm);
536 indexerNames_.clear();
539 template<
typename TRAITS,
typename LO,
typename GO>
546 using Teuchos::rcp_dynamic_cast;
550 if(!disableSensitivities_) {
551 if(d.first_sensitivities_name==sensitivitiesName_)
552 applySensitivities_ =
true;
554 applySensitivities_ =
false;
557 applySensitivities_ =
false;
561 RCP<GlobalEvaluationData> ged;
564 std::string post = useTimeDerivativeSolutionVector_ ?
" - Xdot" :
" - X";
565 if(d.gedc.containsDataObject(globalDataKey_+post)) {
566 ged = d.gedc.getDataObject(globalDataKey_+post);
573 TEUCHOS_TEST_FOR_EXCEPTION(x_==Teuchos::null,std::logic_error,
574 "Gather Jacobian: Can't find x vector in GEDC \"" << globalDataKey_ << post <<
"\"" 575 "A cast failed for " << ged <<
". Type is " << Teuchos::typeName(*ged) << std::endl);
578 ged = d.gedc.getDataObject(globalDataKey_);
582 RCP<const BlockedEpetraLinearObjContainer> blockedContainer = rcp_dynamic_cast<
const BLOC>(ged);
584 if(ro_ged!=Teuchos::null) {
589 else if(blockedContainer!=Teuchos::null) {
590 if (useTimeDerivativeSolutionVector_)
597 TEUCHOS_TEST_FOR_EXCEPTION(x_==Teuchos::null,std::logic_error,
598 "Gather Jacobian: Can't find x vector in GEDC \"" << globalDataKey_ <<
"\" (" << post <<
"). " 599 "A cast failed for " << ged <<
". Type is " << Teuchos::typeName(*ged));
604 template<
typename TRAITS,
typename LO,
typename GO>
609 using Teuchos::ArrayRCP;
610 using Teuchos::ptrFromRef;
611 using Teuchos::rcp_dynamic_cast;
614 using Thyra::SpmdVectorBase;
618 std::string blockId = this->wda(workset).block_id;
619 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
621 double seed_value = 0.0;
622 if(applySensitivities_) {
623 if (useTimeDerivativeSolutionVector_ && gatherSeedIndex_<0) {
624 seed_value = workset.alpha;
626 else if (gatherSeedIndex_<0) {
627 seed_value = workset.beta;
629 else if(!useTimeDerivativeSolutionVector_) {
630 seed_value = workset.gather_seeds[gatherSeedIndex_];
633 TEUCHOS_ASSERT(
false);
640 if(!applySensitivities_)
643 std::vector<int> blockOffsets;
651 int numDerivs = blockOffsets[blockOffsets.size()-1];
654 for(std::size_t fieldIndex=0;
655 fieldIndex<gatherFields_.size();fieldIndex++) {
657 PHX::MDField<ScalarT,Cell,NODE> &
field = gatherFields_[fieldIndex];
659 int indexerId = indexerIds_[fieldIndex];
660 int subFieldNum = subFieldIds_[fieldIndex];
663 Teuchos::ArrayRCP<const double> local_x;
664 rcp_dynamic_cast<SpmdVectorBase<double> >(x_->getNonconstVectorBlock(indexerId))->getLocalData(ptrFromRef(local_x));
666 auto subRowIndexer = indexers_[indexerId];
667 const std::vector<int> & elmtOffset = subRowIndexer->getGIDFieldOffsets(blockId,subFieldNum);
669 int startBlkOffset = blockOffsets[indexerId];
672 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
673 LO cellLocalId = localCellIds[worksetCellIndex];
675 const std::vector<int> & LIDs = subRowIndexer->getElementLIDs(cellLocalId);
684 field(worksetCellIndex,
basis).fastAccessDx(startBlkOffset+
offset) = seed_value;
int getFieldBlock(const std::string &fieldName, const std::vector< Teuchos::RCP< UniqueGlobalIndexer< LocalOrdinalT, GlobalOrdinalT > > > &ugis)
void computeBlockOffsets(const std::string &blockId, const std::vector< Teuchos::RCP< UniqueGlobalIndexer< LocalOrdinalT, GlobalOrdinalT > > > &ugis, std::vector< int > &blockOffsets)
PHX::MDField< const ScalarT > input
PHX::MDField< ScalarT > vector
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)
GatherSolution_BlockedEpetra(const Teuchos::ParameterList &p)
PHX::MDField< const ScalarT, Cell, IP > field
Gathers solution values from the Newton solution vector into the nodal fields of the field manager...
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
void evaluateFields(typename TRAITS::EvalData d)
panzer::Traits::Jacobian::ScalarT ScalarT