43 #ifndef PANZER_GATHER_TANGENT_BLOCKED_EPETRA_IMPL_HPP 44 #define PANZER_GATHER_TANGENT_BLOCKED_EPETRA_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 52 #include "Panzer_EpetraLinearObjFactory.hpp" 55 #include "Teuchos_FancyOStream.hpp" 57 #include "Thyra_SpmdVectorBase.hpp" 58 #include "Thyra_ProductVectorBase.hpp" 60 #include "Epetra_Map.h" 62 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
65 const Teuchos::ParameterList& p)
67 , useTimeDerivativeSolutionVector_(false)
68 , globalDataKey_(
"Tangent Gather Container")
70 const std::vector<std::string>& names =
71 *(p.get< Teuchos::RCP< std::vector<std::string> > >(
"DOF Names"));
73 indexerNames_ = p.get< Teuchos::RCP< std::vector<std::string> > >(
"Indexer Names");
75 Teuchos::RCP<panzer::PureBasis>
basis =
76 p.get< Teuchos::RCP<panzer::PureBasis> >(
"Basis");
79 for (std::size_t fd = 0; fd < names.size(); ++fd) {
81 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
85 if (p.isType<
bool>(
"Use Time Derivative Solution Vector"))
88 if (p.isType<std::string>(
"Global Data Key"))
91 this->setName(
"Gather Tangent");
95 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
100 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_->size());
102 indexerIds_.resize(gatherFields_.size());
103 subFieldIds_.resize(gatherFields_.size());
105 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
107 const std::string& fieldName = (*indexerNames_)[fd];
110 subFieldIds_[fd] = indexers_[indexerIds_[fd]]->getFieldNum(fieldName);
112 TEUCHOS_ASSERT(indexerIds_[fd]>=0);
115 this->utils.setFieldData(gatherFields_[fd],fm);
118 indexerNames_ = Teuchos::null;
122 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
131 if (d.gedc.containsDataObject(globalDataKey_)) {
132 Teuchos::RCP<const BLOC> blockedContainer = Teuchos::rcp_dynamic_cast<
const BLOC>(d.gedc.getDataObject(globalDataKey_),
true);
134 if (useTimeDerivativeSolutionVector_)
135 x_ = Teuchos::rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer->get_dxdt());
137 x_ = Teuchos::rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer->get_x());
142 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
147 using Teuchos::ArrayRCP;
148 using Teuchos::ptrFromRef;
149 using Teuchos::rcp_dynamic_cast;
152 using Thyra::SpmdVectorBase;
156 if (x_ == Teuchos::null)
160 std::string blockId = this->wda(workset).block_id;
161 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
164 Teuchos::ArrayRCP<const double> local_x;
165 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
167 int indexerId = indexerIds_[fieldIndex];
168 int subFieldNum = subFieldIds_[fieldIndex];
171 rcp_dynamic_cast<SpmdVectorBase<double> >(x_->getNonconstVectorBlock(indexerId))->getLocalData(ptrFromRef(local_x));
173 auto subRowIndexer = indexers_[indexerId];
174 const std::vector<int> & elmtOffset = subRowIndexer->getGIDFieldOffsets(blockId,subFieldNum);
177 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
178 LO cellLocalId = localCellIds[worksetCellIndex];
180 const std::vector<int> & LIDs = subRowIndexer->getElementLIDs(cellLocalId);
190 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) = local_x[lid];
int getFieldBlock(const std::string &fieldName, const std::vector< Teuchos::RCP< UniqueGlobalIndexer< LocalOrdinalT, GlobalOrdinalT > > > &ugis)
std::vector< PHX::MDField< ScalarT, Cell, NODE > > gatherFields_
bool useTimeDerivativeSolutionVector_
Teuchos::RCP< std::vector< std::string > > indexerNames_
PHX::MDField< ScalarT > vector
std::string globalDataKey_
void preEvaluate(typename TRAITS::PreEvalData d)
GatherTangent_BlockedEpetra()
void evaluateFields(typename TRAITS::EvalData d)
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)