43 #ifndef PANZER_PROJECT_TO_EDGES_IMPL_HPP 44 #define PANZER_PROJECT_TO_EDGES_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 52 #include "Teuchos_FancyOStream.hpp" 54 template<
typename EvalT,
typename Traits>
57 const Teuchos::ParameterList& p)
59 dof_name = (p.get< std::string >(
"DOF Name"));
61 if(p.isType< Teuchos::RCP<PureBasis> >(
"Basis"))
62 basis = p.get< Teuchos::RCP<PureBasis> >(
"Basis");
64 basis = p.get< Teuchos::RCP<const PureBasis> >(
"Basis");
66 Teuchos::RCP<PHX::DataLayout> basis_layout =
basis->functional;
67 Teuchos::RCP<PHX::DataLayout> vector_layout =
basis->functional_grad;
70 TEUCHOS_ASSERT(
basis->isVectorBasis());
72 result = PHX::MDField<ScalarT,Cell,BASIS>(dof_name,basis_layout);
73 this->addEvaluatedField(result);
75 vector_values = PHX::MDField<ScalarT,Cell,BASIS,Dim>(dof_name+
"_Vector",vector_layout);
76 tangents = PHX::MDField<ScalarT,Cell,BASIS,Dim>(dof_name+
"_Tangents",vector_layout);
77 this->addDependentField(vector_values);
78 this->addDependentField(tangents);
80 this->setName(
"Project To Edges");
84 template<
typename EvalT,
typename Traits>
90 this->utils.setFieldData(result,fm);
91 this->utils.setFieldData(vector_values,fm);
92 this->utils.setFieldData(tangents,fm);
94 num_pts = vector_values.dimension(1);
95 num_dim = vector_values.dimension(2);
97 TEUCHOS_ASSERT(vector_values.dimension(1) == tangents.dimension(1));
98 TEUCHOS_ASSERT(vector_values.dimension(2) == tangents.dimension(2));
102 template<
typename EvalT,
typename Traits>
106 const shards::CellTopology & parentCell = *
basis->getCellTopology();
107 const int intDegree =
basis->order();
108 TEUCHOS_ASSERT(intDegree == 1);
109 Intrepid2::DefaultCubatureFactory<double,Kokkos::DynRankView<double,PHX::Device>,Kokkos::DynRankView<double,PHX::Device>> quadFactory;
110 Teuchos::RCP< Intrepid2::Cubature<double,Kokkos::DynRankView<double,PHX::Device>,Kokkos::DynRankView<double,PHX::Device>> > edgeQuad;
113 const unsigned num_edges = parentCell.getEdgeCount();
114 std::vector<double> refEdgeWt(num_edges, 0.0);
115 for (
unsigned e=0; e<num_edges; e++) {
116 edgeQuad = quadFactory.create(parentCell.getCellTopologyData(1,e), intDegree);
117 const int numQPoints = edgeQuad->getNumPoints();
118 Kokkos::DynRankView<double,PHX::Device> quadWts(
"quadWts",numQPoints);
119 Kokkos::DynRankView<double,PHX::Device> quadPts(
"quadPts",numQPoints,
num_dim);
120 edgeQuad->getCubature(quadPts,quadWts);
121 for (
int q=0; q<numQPoints; q++)
122 refEdgeWt[e] += quadWts(q);
126 for (index_t cell = 0; cell < workset.
num_cells; ++cell) {
127 for (
int p = 0; p <
num_pts; ++p) {
129 for (
int dim = 0; dim <
num_dim; ++dim)
130 result(cell,p) += vector_values(cell,p,dim) * tangents(cell,p,dim);
131 result(cell,p) *= refEdgeWt[p];
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)
void evaluateFields(typename Traits::EvalData d)
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.