43 #ifndef PANZER_NORMALS_IMPL_HPP 44 #define PANZER_NORMALS_IMPL_HPP 49 #include "Intrepid2_FunctionSpaceTools.hpp" 50 #include "Intrepid2_CellTools.hpp" 55 template<
typename EvalT,
typename Traits>
58 const Teuchos::ParameterList& p)
62 const std::string name = p.get<std::string>(
"Name");
63 side_id = p.get<
int>(
"Side ID");
64 Teuchos::RCP<panzer::IntegrationRule> quadRule
65 = p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR");
66 if(p.isParameter(
"Normalize"))
70 Teuchos::RCP<PHX::DataLayout> vector_dl = quadRule->dl_vector;
74 normals = PHX::MDField<ScalarT,Cell,Point,Dim>(name, vector_dl);
75 this->addEvaluatedField(
normals);
77 std::string n =
"Normals: " + name;
82 template<
typename EvalT,
typename Traits>
89 num_qp = normals.extent(1);
90 num_dim = normals.extent(2);
96 template<
typename EvalT,
typename Traits>
105 Intrepid2::CellTools<PHX::exec_space>::getPhysicalSideNormals(normals.get_view(),
106 this->wda(workset).int_rules[quad_index]->jac.get_view(),
107 side_id, *this->wda(workset).int_rules[quad_index]->int_rule->topology);
112 auto local_normals = normals;
113 auto local_num_qp = num_qp;
114 auto local_num_dim = num_dim;
115 Kokkos::parallel_for(
"normalize", workset.
num_cells, KOKKOS_LAMBDA (index_t c) {
116 for(std::size_t q=0;q<local_num_qp;q++) {
120 for(std::size_t d=0;d<local_num_dim;d++)
121 norm += local_normals(c,q,d)*local_normals(c,q,d);
125 for(std::size_t d=0;d<local_num_dim;d++)
126 local_normals(c,q,d) /= norm;
int num_cells
DEPRECATED - use: numCells()
Normals(const Teuchos::ParameterList &p)
std::vector< int >::size_type getIntegrationRuleIndex(int ir_degree, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
PHX::MDField< ScalarT, Cell, Point, Dim > normals
typename EvalT::ScalarT ScalarT
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
void evaluateFields(typename Traits::EvalData d)
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_