43 #ifndef PANZER_EVALUATOR_SCALAR_IMPL_HPP 44 #define PANZER_EVALUATOR_SCALAR_IMPL_HPP 46 #include "Intrepid2_FunctionSpaceTools.hpp" 49 #include "Kokkos_ViewFactory.hpp" 50 #include "Phalanx_DataLayout_MDALayout.hpp" 55 template<
typename EvalT,
typename Traits>
58 const Teuchos::ParameterList& p) : quad_index(0)
61 p.validateParameters(*valid_params);
63 Teuchos::RCP<panzer::IntegrationRule> ir = p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR");
66 Teuchos::RCP<PHX::DataLayout> dl_cell = Teuchos::rcp(
new PHX::MDALayout<Cell>(ir->dl_scalar->extent(0)));
67 integral = PHX::MDField<ScalarT>( p.get<std::string>(
"Integral Name"), dl_cell);
68 scalar = PHX::MDField<const ScalarT,Cell,IP>( p.get<std::string>(
"Integrand Name"), ir->dl_scalar);
71 this->addDependentField(
scalar);
74 if(p.isType<
double>(
"Multiplier"))
77 if (p.isType<Teuchos::RCP<
const std::vector<std::string> > >(
"Field Multipliers")) {
78 const std::vector<std::string>& field_multiplier_names =
79 *(p.get<Teuchos::RCP<const std::vector<std::string> > >(
"Field Multipliers"));
81 for (std::vector<std::string>::const_iterator name =
82 field_multiplier_names.begin();
83 name != field_multiplier_names.end(); ++name) {
84 PHX::MDField<const ScalarT,Cell,IP> tmp_field(*name, p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar);
89 for (
typename std::vector<PHX::MDField<const ScalarT,Cell,IP> >::iterator
field =
field_multipliers.begin();
91 this->addDependentField(*
field);
93 std::string n =
"Integrator_Scalar: " +
integral.fieldTag().name();
98 template<
typename EvalT,
typename Traits>
105 num_qp = scalar.extent(1);
111 template<
typename EvalT,
typename Traits>
122 for (index_t cell = 0; cell < workset.
num_cells; ++cell) {
123 for (std::size_t qp = 0; qp < num_qp; ++qp) {
125 for (
typename std::vector<PHX::MDField<const ScalarT,Cell,IP> >::iterator
field = field_multipliers.begin();
127 tmp(cell,qp) = tmp(cell,qp) * (*field)(cell,qp);
152 for(index_t cl = 0; cl < workset.
num_cells; cl++) {
161 template<
typename EvalT,
typename TRAITS>
162 Teuchos::RCP<Teuchos::ParameterList>
165 Teuchos::RCP<Teuchos::ParameterList> p = Teuchos::rcp(
new Teuchos::ParameterList);
166 p->set<std::string>(
"Integral Name",
"?");
167 p->set<std::string>(
"Integrand Name",
"?");
169 Teuchos::RCP<panzer::IntegrationRule> ir;
171 p->set<
double>(
"Multiplier", 1.0);
173 Teuchos::RCP<const std::vector<std::string> > fms;
174 p->set(
"Field Multipliers", fms);
Teuchos::RCP< Teuchos::ParameterList > getValidParameters() const
std::vector< std::string >::size_type getIntegrationRuleIndex(int ir_degree, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
double multiplier
The scalar multiplier out in front of the integral ( ).
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack... dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
Integrator_Scalar(const Teuchos::ParameterList &p)
PHX::MDField< const ScalarT, Cell, IP > scalar
Array_CellIP weighted_measure
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we'll contribute, or in which we'll store, the result of computing this integral...
PHX::MDField< ScalarT > integral
std::vector< PHX::MDField< const ScalarT, Cell, IP > > field_multipliers
void evaluateFields(typename Traits::EvalData d)
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_