43 #ifndef PANZER_GATHER_TANGENTS_IMPL_HPP 44 #define PANZER_GATHER_TANGENTS_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 49 #include "Intrepid2_Kernels.hpp" 50 #include "Intrepid2_OrientationTools.hpp" 53 #include "Kokkos_ViewFactory.hpp" 55 #include "Teuchos_FancyOStream.hpp" 57 template<
typename EvalT,
typename Traits>
60 const Teuchos::ParameterList& p)
62 dof_name = (p.get< std::string >(
"DOF Name"));
64 if(p.isType< Teuchos::RCP<PureBasis> >(
"Basis"))
65 basis = p.get< Teuchos::RCP<PureBasis> >(
"Basis");
67 basis = p.get< Teuchos::RCP<const PureBasis> >(
"Basis");
69 pointRule = p.get<Teuchos::RCP<const PointRule> >(
"Point Rule");
71 Teuchos::RCP<PHX::DataLayout> basis_layout = basis->functional;
72 Teuchos::RCP<PHX::DataLayout> vector_layout_vector = basis->functional_grad;
75 TEUCHOS_ASSERT(basis->isVectorBasis());
78 std::string orientationFieldName = basis->name() +
" Orientation";
84 constJac_ = pointValues.jac;
85 this->addDependentField(constJac_);
87 gatherFieldTangents = PHX::MDField<ScalarT,Cell,NODE,Dim>(dof_name+
"_Tangents",vector_layout_vector);
88 this->addEvaluatedField(gatherFieldTangents);
90 this->setName(
"Gather Tangents");
94 template<
typename EvalT,
typename Traits>
100 this->utils.setFieldData(pointValues.jac,fm);
101 const shards::CellTopology & parentCell = *basis->getCellTopology();
103 edgeParam = Intrepid2::RefSubcellParametrization<PHX::Device>::get(edgeDim, parentCell.getKey());
107 template<
typename EvalT,
typename Traits>
115 const shards::CellTopology & parentCell = *basis->getCellTopology();
116 int cellDim = parentCell.getDimension();
118 int numEdges = gatherFieldTangents.extent(1);
124 const auto worksetJacobians = pointValues.jac.get_view();
127 for(index_t c=0;c<workset.
num_cells;c++) {
129 int edgeOrts[12] = {};
130 orientations->at(
details.cell_local_ids[c]).getEdgeOrientation(edgeOrts, numEdges);
132 for(
int pt = 0; pt < numEdges; pt++) {
133 auto phyEdgeTan = Kokkos::subview(gatherFieldTangents.get_static_view(), c, pt, Kokkos::ALL());
134 auto ortEdgeTan = Kokkos::subview(workspace, 0, Kokkos::ALL());
137 Intrepid2::Impl::OrientationTools::getRefSubcellTangents(
140 parentCell.getKey(edgeDim,pt),
144 auto J = Kokkos::subview(worksetJacobians, c, pt, Kokkos::ALL(), Kokkos::ALL());
145 Intrepid2::Kernels::Serial::matvec_product(phyEdgeTan, J, ortEdgeTan);
Teuchos::RCP< const std::vector< Intrepid2::Orientation > > orientations_
int num_cells
DEPRECATED - use: numCells()
void setupArrays(const Teuchos::RCP< const panzer::PointRule > &pr)
Sizes/allocates memory for arrays.
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.
void evaluateFields(typename Traits::EvalData d)
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)