43 #ifndef PANZER_STK_GATHER_FIELDS_IMPL_HPP 44 #define PANZER_STK_GATHER_FIELDS_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 51 #include "Teuchos_FancyOStream.hpp" 57 template<
typename EvalT,
typename Traits>
59 GatherFields(
const Teuchos::RCP<const STK_Interface> & mesh,
const Teuchos::ParameterList& p)
66 const std::vector<std::string>& names =
67 *(p.get< Teuchos::RCP< std::vector<std::string> > >(
"Field Names"));
69 Teuchos::RCP<panzer::BasisIRLayout>
basis =
70 p.get< Teuchos::RCP<panzer::BasisIRLayout> >(
"Basis");
73 gatherFields_.resize(names.size());
75 for (std::size_t fd = 0; fd < names.size(); ++fd) {
77 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
78 this->addEvaluatedField(gatherFields_[fd]);
81 this->setName(
"Gather STK Fields");
85 template<
typename EvalT,
typename Traits>
90 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
91 std::string fieldName = gatherFields_[fd].fieldTag().name();
97 mesh_->printMetaData(ss);
98 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
99 "panzer_stk::GatherFields: STK field " <<
"\"" << fieldName <<
"\" " 100 "not found.\n STK meta data follows: \n\n" << ss.str());
105 this->utils.setFieldData(gatherFields_[fd],fm);
110 template<
typename EvalT,
typename Traits>
114 const std::vector<stk::mesh::Entity> & localElements = *
mesh_->getElementsOrderedByLID();
117 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
120 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
121 std::size_t cellLocalId = localCellIds[worksetCellIndex];
122 stk::mesh::Entity
const* relations =
mesh_->getBulkData()->begin_nodes(localElements[cellLocalId]);
125 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
128 std::size_t basisCnt = gatherFields_[fieldIndex].dimension(1);
132 (gatherFields_[fieldIndex])(worksetCellIndex,0) = *stk::mesh::field_data(*
field, localElements[cellLocalId]);
137 stk::mesh::Entity node = relations[
basis];
138 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) = *stk::mesh::field_data(*
field, node);
Teuchos::RCP< STK_Interface > mesh_
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)
PHX::MDField< const ScalarT, Cell, IP > field
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
std::vector< VariableField * > stkFields_
panzer_stk::STK_Interface::SolutionFieldType VariableField
void evaluateFields(typename Traits::EvalData d)