43 #ifndef PANZER_GATHER_TANGENT_TPETRA_IMPL_HPP 44 #define PANZER_GATHER_TANGENT_TPETRA_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 55 #include "Teuchos_FancyOStream.hpp" 57 #include "Tpetra_Vector.hpp" 58 #include "Tpetra_Map.hpp" 60 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO,
typename NodeT>
63 const Teuchos::RCP<const panzer::GlobalIndexer> & indexer,
64 const Teuchos::ParameterList& p)
65 : globalIndexer_(indexer)
66 , useTimeDerivativeSolutionVector_(false)
67 , globalDataKey_(
"Tangent Gather Container")
69 const std::vector<std::string>& names =
70 *(p.get< Teuchos::RCP< std::vector<std::string> > >(
"DOF Names"));
72 indexerNames_ = p.get< Teuchos::RCP< std::vector<std::string> > >(
"Indexer Names");
75 Teuchos::RCP<const panzer::PureBasis> basis;
76 if(p.isType< Teuchos::RCP<panzer::PureBasis> >(
"Basis"))
77 basis = p.get< Teuchos::RCP<panzer::PureBasis> >(
"Basis");
79 basis = p.get< Teuchos::RCP<const panzer::PureBasis> >(
"Basis");
82 for (std::size_t fd = 0; fd < names.size(); ++fd) {
84 PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
92 if (p.isType<
bool>(
"Use Time Derivative Solution Vector"))
95 if (p.isType<std::string>(
"Global Data Key"))
98 this->setName(
"Gather Tangent");
102 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO,
typename NodeT>
107 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_->size());
109 fieldIds_.resize(gatherFields_.size());
111 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
112 const std::string& fieldName = (*indexerNames_)[fd];
113 fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
116 indexerNames_ = Teuchos::null;
120 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO,
typename NodeT>
125 using Teuchos::rcp_dynamic_cast;
130 if (d.gedc->containsDataObject(globalDataKey_)) {
131 RCP<GlobalEvaluationData> ged = d.gedc->getDataObject(globalDataKey_);
132 RCP<LOCPair_GlobalEvaluationData> loc_pair =
135 if(loc_pair!=Teuchos::null) {
136 Teuchos::RCP<LinearObjContainer> loc = loc_pair->
getGhostedLOC();
137 tpetraContainer_ = rcp_dynamic_cast<LOC>(loc,
true);
140 if(tpetraContainer_==Teuchos::null) {
141 tpetraContainer_ = rcp_dynamic_cast<LOC>(ged,
true);
147 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO,
typename NodeT>
153 if (tpetraContainer_ == Teuchos::null)
158 std::string blockId = this->wda(workset).block_id;
159 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
161 Teuchos::RCP<typename LOC::VectorType> x;
162 if (useTimeDerivativeSolutionVector_)
165 x = tpetraContainer_->get_x();
167 Teuchos::ArrayRCP<const double> x_array = x->get1dView();
175 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
176 std::size_t cellLocalId = localCellIds[worksetCellIndex];
178 auto LIDs = globalIndexer_->getElementLIDs(cellLocalId);
181 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
182 int fieldNum = fieldIds_[fieldIndex];
183 const std::vector<int> & elmtOffset = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
186 for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
187 int offset = elmtOffset[basis];
188 LO lid = LIDs[offset];
189 (gatherFields_[fieldIndex])(worksetCellIndex,basis) = x_array[lid];
Teuchos::RCP< std::vector< std::string > > indexerNames_
void preEvaluate(typename TRAITS::PreEvalData d)
void evaluateFields(typename TRAITS::EvalData d)
std::vector< PHX::MDField< ScalarT, Cell, NODE > > gatherFields_
bool useTimeDerivativeSolutionVector_
virtual Teuchos::RCP< CloneableEvaluator > clone(const Teuchos::ParameterList &pl) const
std::string globalDataKey_
Teuchos::RCP< LinearObjContainer > getGhostedLOC() const
const Teuchos::RCP< VectorType > get_dxdt() const
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)