43 #ifndef PANZER_GATHER_TANGENT_BLOCKED_TPETRA_IMPL_HPP 44 #define PANZER_GATHER_TANGENT_BLOCKED_TPETRA_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 52 #include "Panzer_TpetraLinearObjFactory.hpp" 55 #include "Teuchos_FancyOStream.hpp" 57 #include "Thyra_SpmdVectorBase.hpp" 58 #include "Thyra_ProductVectorBase.hpp" 60 #include "Tpetra_Map.hpp" 62 template <
typename EvalT,
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
66 const Teuchos::ParameterList& p)
67 : gidIndexer_(indexer)
68 , useTimeDerivativeSolutionVector_(false)
69 , globalDataKey_(
"Tangent Gather Container")
71 const std::vector<std::string>& names =
72 *(p.get< Teuchos::RCP< std::vector<std::string> > >(
"DOF Names"));
74 indexerNames_ = p.get< Teuchos::RCP< std::vector<std::string> > >(
"Indexer Names");
76 Teuchos::RCP<panzer::PureBasis>
basis =
77 p.get< Teuchos::RCP<panzer::PureBasis> >(
"Basis");
80 for (std::size_t fd = 0; fd < names.size(); ++fd) {
82 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
86 if (p.isType<
bool>(
"Use Time Derivative Solution Vector"))
89 if (p.isType<std::string>(
"Global Data Key"))
92 this->setName(
"Gather Tangent");
96 template <
typename EvalT,
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
101 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_->size());
103 fieldIds_.resize(gatherFields_.size());
105 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
107 const std::string& fieldName = (*indexerNames_)[fd];
108 fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
111 this->utils.setFieldData(gatherFields_[fd],fm);
114 indexerNames_ = Teuchos::null;
118 template <
typename EvalT,
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
123 if (d.gedc.containsDataObject(globalDataKey_)) {
124 blockedContainer_ = Teuchos::rcp_dynamic_cast<
const ContainerType>(d.gedc.getDataObject(globalDataKey_),
true);
129 template <
typename EvalT,
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
134 using Teuchos::ArrayRCP;
135 using Teuchos::ptrFromRef;
136 using Teuchos::rcp_dynamic_cast;
139 using Thyra::SpmdVectorBase;
144 if (blockedContainer_ == Teuchos::null)
147 Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
148 out.setShowProcRank(
true);
149 out.setOutputToRootOnly(-1);
151 std::vector<std::pair<int,GO> > GIDs;
152 std::vector<LO> LIDs;
155 std::string blockId = this->wda(workset).block_id;
156 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
158 Teuchos::RCP<ProductVectorBase<double> > x;
159 if (useTimeDerivativeSolutionVector_)
160 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
162 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
165 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
166 LO cellLocalId = localCellIds[worksetCellIndex];
168 gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
171 LIDs.resize(GIDs.size());
172 for(std::size_t i=0;i<GIDs.size();i++) {
174 RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
176 LIDs[i] = x_map->getLocalElement(GIDs[i].second);
180 Teuchos::ArrayRCP<const double> local_x;
181 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
182 int fieldNum = fieldIds_[fieldIndex];
183 int indexerId = gidIndexer_->getFieldBlock(fieldNum);
186 RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
187 block_x->getLocalData(ptrFromRef(local_x));
189 const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
196 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) = local_x[lid];
Teuchos::RCP< std::vector< std::string > > indexerNames_
GatherTangent_BlockedTpetra()
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)
bool useTimeDerivativeSolutionVector_
std::string globalDataKey_
void evaluateFields(typename TRAITS::EvalData d)
std::vector< PHX::MDField< ScalarT, Cell, NODE > > gatherFields_
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
void preEvaluate(typename TRAITS::PreEvalData d)