43 #ifndef PANZER_EVALUATOR_DotProduct_IMPL_HPP 44 #define PANZER_EVALUATOR_DotProduct_IMPL_HPP 51 #include "Teuchos_RCP.hpp" 55 template <
typename EvalT,
typename TraitsT>
56 Teuchos::RCP<DotProduct<EvalT,TraitsT> >
59 const std::string & vecA,
60 const std::string & vecB,
62 const std::string & fieldMultiplier)
64 Teuchos::ParameterList pl;
65 pl.set(
"Result Name",resultName);
66 pl.set(
"Point Rule",Teuchos::rcpFromRef(pr));
67 pl.set(
"Vector A Name",vecA);
68 pl.set(
"Vector B Name",vecB);
70 pl.set(
"Field Multiplier",fieldMultiplier);
76 template<
typename EvalT,
typename Traits>
79 const Teuchos::ParameterList& p)
80 : multiplier_field_on(false)
82 std::string result_name = p.get<std::string>(
"Result Name");
83 std::string vec_a_name = p.get<std::string>(
"Vector A Name");
84 std::string vec_b_name = p.get<std::string>(
"Vector B Name");
86 std::string multiplier_name =
"";
87 if(p.isType<std::string>(
"Field Multiplier"))
88 multiplier_name = p.get<std::string>(
"Field Multiplier");
91 if(p.isType<
double>(
"Multiplier"))
94 const Teuchos::RCP<const panzer::PointRule> pr =
95 p.get< Teuchos::RCP<const panzer::PointRule> >(
"Point Rule");
98 vec_a = PHX::MDField<const ScalarT>(vec_a_name, pr->dl_vector);
99 vec_b = PHX::MDField<const ScalarT>(vec_b_name, pr->dl_vector);
101 if(multiplier_name!=
"") {
102 multiplier_field = PHX::MDField<const ScalarT>(multiplier_name,pr->dl_scalar);
108 this->addDependentField(
vec_a);
109 this->addDependentField(
vec_b);
111 std::string n =
"DotProduct: " + result_name +
" = " + vec_a_name +
" . " + vec_b_name;
116 template<
typename EvalT,
typename Traits>
123 num_pts = vec_a.extent(1);
124 num_dim = vec_a.extent(2);
126 TEUCHOS_ASSERT(vec_a.extent(1) == vec_b.extent(1));
127 TEUCHOS_ASSERT(vec_a.extent(2) == vec_b.extent(2));
131 template<
typename EvalT,
typename Traits>
138 auto vec_a_v = vec_a.get_static_view();
139 auto vec_b_v = vec_b.get_static_view();
140 auto vec_a_dot_vec_b_v = vec_a_dot_vec_b.get_static_view();
141 auto multiplier_field_v = multiplier_field.get_static_view();
143 int l_num_pts = num_pts, l_num_dim = num_dim;
144 auto l_multiplier_field_on = multiplier_field_on;
145 auto l_multiplier_value = multiplier_value;
147 Kokkos::parallel_for (workset.
num_cells, KOKKOS_LAMBDA (
const int cell) {
148 for (
int p = 0; p < l_num_pts; ++p) {
149 vec_a_dot_vec_b_v(cell,p) =
ScalarT(0.0);
150 for (
int dim = 0; dim < l_num_dim; ++dim)
151 vec_a_dot_vec_b_v(cell,p) += vec_a_v(cell,p,dim) * vec_b_v(cell,p,dim);
153 if(l_multiplier_field_on)
154 vec_a_dot_vec_b_v(cell,p) *= l_multiplier_value*multiplier_field_v(cell,p);
156 vec_a_dot_vec_b_v(cell,p) *= l_multiplier_value;
PHX::MDField< const ScalarT > vec_a
int num_cells
DEPRECATED - use: numCells()
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
typename EvalT::ScalarT ScalarT
Evaluates dot product at a set of points.
PHX::MDField< const ScalarT > vec_b
double multiplier
The scalar multiplier out in front of the integral ( ).
Teuchos::RCP< DotProduct< EvalT, TraitsT > > buildEvaluator_DotProduct(const std::string &resultName, const panzer::PointRule &pr, const std::string &vecA, const std::string &vecB, double multiplier=1, const std::string &fieldMultiplier="")
Build a dot product evaluator. Evaluates dot product at a set of points.
DotProduct(const Teuchos::ParameterList &p)
void evaluateFields(typename Traits::EvalData d)
PHX::MDField< const ScalarT > multiplier_field
PHX::MDField< ScalarT > vec_a_dot_vec_b