43 #ifndef __Panzer_SubcellSum_impl_hpp__ 44 #define __Panzer_SubcellSum_impl_hpp__ 50 #include "Phalanx_DataLayout_MDALayout.hpp" 55 template<
typename EvalT,
typename Traits>
58 const Teuchos::ParameterList& p)
59 : evaluateOnClosure_(false)
62 p.validateParameters(*valid_params);
64 const std::string inName = p.get<std::string>(
"Field Name");
65 const std::string outName = p.get<std::string>(
"Sum Name");
66 Teuchos::RCP<const PureBasis> basis = p.get< Teuchos::RCP<const PureBasis> >(
"Basis");
68 if(p.isType<
bool>(
"Evaluate On Closure"))
71 inField = PHX::MDField<const ScalarT,Cell,BASIS>( inName, basis->functional);
72 outField = PHX::MDField<ScalarT,Cell>( outName, basis->cell_data);
74 this->addDependentField(
inField);
80 std::string n =
"SubcellSum: " +
outField.fieldTag().name();
85 template<
typename EvalT,
typename Traits>
91 std::vector<int> indices;
96 if(evaluateOnClosure_)
97 fieldPattern_->getSubcellClosureIndices(workset.
subcell_dim,this->wda(workset).subcell_index,indices);
99 indices = fieldPattern_->getSubcellIndices(workset.
subcell_dim,this->wda(workset).subcell_index);
101 auto outField_h = Kokkos::create_mirror_view(outField.get_static_view());
102 auto inField_h = Kokkos::create_mirror_view(inField.get_static_view());
103 Kokkos::deep_copy(inField_h, inField.get_static_view());
104 for(index_t c=0;c<workset.
num_cells;c++) {
108 for(std::size_t i=0;i<indices.size();i++)
109 outField_h(c) += inField_h(c,indices[i]);
114 Kokkos::deep_copy(outField.get_static_view(), outField_h);
118 template<
typename EvalT,
typename TRAITS>
119 Teuchos::RCP<Teuchos::ParameterList>
122 Teuchos::RCP<Teuchos::ParameterList> p = Teuchos::rcp(
new Teuchos::ParameterList);
123 p->set<std::string>(
"Sum Name",
"?");
124 p->set<std::string>(
"Field Name",
"?");
125 p->set<
double>(
"Multiplier",1.0);
126 p->set<
bool>(
"Evaluate On Closure",
false);
128 Teuchos::RCP<const panzer::PureBasis> basis;
129 p->set(
"Basis", basis);
SubcellSum(const Teuchos::ParameterList &p)
int num_cells
DEPRECATED - use: numCells()
int subcell_dim
DEPRECATED - use: getSubcellDimension()
Teuchos::RCP< const panzer::FieldPattern > fieldPattern_
double multiplier
The scalar multiplier out in front of the integral ( ).
Teuchos::RCP< Teuchos::ParameterList > getValidParameters() const
PHX::MDField< const ScalarT, Cell, BASIS > inField
PHX::MDField< ScalarT, Cell > outField
void evaluateFields(typename Traits::EvalData d)