43 #ifndef PANZER_GATHER_NORMALS_IMPL_HPP
44 #define PANZER_GATHER_NORMALS_IMPL_HPP
46 #include "Teuchos_Assert.hpp"
47 #include "Phalanx_DataLayout.hpp"
50 #include "Kokkos_ViewFactory.hpp"
52 #include "Intrepid2_Kernels.hpp"
53 #include "Intrepid2_OrientationTools.hpp"
55 #include "Teuchos_FancyOStream.hpp"
57 template<
typename EvalT,
typename Traits>
62 dof_name_ = (p.
get< std::string >(
"DOF Name"));
78 std::string orientationFieldName = basis_->name() +
" Orientation";
84 constJac_ = pointValues_.jac;
85 this->addDependentField(constJac_);
88 this->addEvaluatedField(gatherFieldNormals_);
90 this->setName(
"Gather Normals");
94 template<
typename EvalT,
typename Traits>
101 auto orientations_host = Kokkos::create_mirror_view(orientations_);
102 for (
size_t i=0; i < orientations->size(); ++i)
103 orientations_host(i) = (*orientations)[i];
104 Kokkos::deep_copy(orientations_,orientations_host);
106 this->utils.setFieldData(pointValues_.jac,fm);
108 const shards::CellTopology & parentCell = *basis_->getCellTopology();
109 int sideDim = parentCell.getDimension()-1;
110 sideParam_ = Intrepid2::RefSubcellParametrization<PHX::Device>::get(sideDim, parentCell.getKey());
112 int numFaces = gatherFieldNormals_.extent(1);
114 auto keys_host = Kokkos::create_mirror_view(keys_);
115 for (
int i=0; i < numFaces; ++i)
116 keys_host(i) = parentCell.getKey(sideDim,i);
117 Kokkos::deep_copy(keys_,keys_host);
120 int cellDim = parentCell.getDimension();
121 refEdges_ = Kokkos::createDynRankViewWithType<Kokkos::DynRankView<ScalarT,PHX::Device>>(gatherFieldNormals_.get_static_view(),
"ref_edges", (*d.
worksets_)[0].num_cells, sideDim, cellDim);
122 phyEdges_ = Kokkos::createDynRankViewWithType<Kokkos::DynRankView<ScalarT,PHX::Device>>(gatherFieldNormals_.get_static_view(),
"phy_edges", (*d.
worksets_)[0].num_cells, sideDim, cellDim);
126 template<
typename EvalT,
typename Traits>
134 int numFaces = gatherFieldNormals_.extent(1);
135 const auto worksetJacobians = pointValues_.jac.get_view();
137 auto gatherFieldNormals = gatherFieldNormals_;
138 auto sideParam = sideParam_;
140 auto orientations = orientations_;
141 auto refEdges = refEdges_;
142 auto phyEdges = phyEdges_;
145 Kokkos::parallel_for(
"panzer::GatherNormals",workset.
num_cells,KOKKOS_LAMBDA(
const int c){
146 int faceOrts[6] = {};
147 orientations(cell_local_ids(c)).getFaceOrientation(faceOrts, numFaces);
149 for(
int pt = 0; pt < numFaces; pt++) {
150 auto ortEdgeTan_U = Kokkos::subview(refEdges, c, 0, Kokkos::ALL());
151 auto ortEdgeTan_V = Kokkos::subview(refEdges, c, 1, Kokkos::ALL());
153 auto tmpRefEdges = Kokkos::subview(refEdges, c, Kokkos::ALL(), Kokkos::ALL());
154 Intrepid2::Impl::OrientationTools::getRefSubcellTangents(tmpRefEdges,
160 auto phyEdgeTan_U = Kokkos::subview(phyEdges, c, 0, Kokkos::ALL());
161 auto phyEdgeTan_V = Kokkos::subview(phyEdges, c, 1, Kokkos::ALL());
162 auto J = Kokkos::subview(worksetJacobians, c, pt, Kokkos::ALL(), Kokkos::ALL());
164 Intrepid2::Kernels::Serial::matvec_product(phyEdgeTan_U, J, ortEdgeTan_U);
165 Intrepid2::Kernels::Serial::matvec_product(phyEdgeTan_V, J, ortEdgeTan_V);
168 gatherFieldNormals(c,pt,0) = (phyEdgeTan_U(1)*phyEdgeTan_V(2) - phyEdgeTan_U(2)*phyEdgeTan_V(1));
169 gatherFieldNormals(c,pt,1) = (phyEdgeTan_U(2)*phyEdgeTan_V(0) - phyEdgeTan_U(0)*phyEdgeTan_V(2));
170 gatherFieldNormals(c,pt,2) = (phyEdgeTan_U(0)*phyEdgeTan_V(1) - phyEdgeTan_U(1)*phyEdgeTan_V(0));
Teuchos::RCP< const std::vector< Intrepid2::Orientation > > orientations_
int num_cells
DEPRECATED - use: numCells()
T & get(const std::string &name, T def_value)
void setupArrays(const Teuchos::RCP< const panzer::PointRule > &pr)
Sizes/allocates memory for arrays.
void evaluateFields(typename Traits::EvalData d)
Kokkos::View< const int *, PHX::Device > getLocalCellIDs() const
Get the local cell IDs for the workset.
bool isType(const std::string &name) const
#define TEUCHOS_ASSERT(assertion_test)
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_