43 #ifndef PANZER_GATHER_NORMALS_IMPL_HPP 44 #define PANZER_GATHER_NORMALS_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 51 #include "Kokkos_ViewFactory.hpp" 53 #include "Teuchos_FancyOStream.hpp" 55 template<
typename EvalT,
typename Traits>
58 const Teuchos::ParameterList& p)
60 dof_name = (p.get< std::string >(
"DOF Name"));
62 if(p.isType< Teuchos::RCP<PureBasis> >(
"Basis"))
63 basis = p.get< Teuchos::RCP<PureBasis> >(
"Basis");
65 basis = p.get< Teuchos::RCP<const PureBasis> >(
"Basis");
67 pointRule = p.get<Teuchos::RCP<const PointRule> >(
"Point Rule");
69 Teuchos::RCP<PHX::DataLayout> basis_layout =
basis->functional;
70 Teuchos::RCP<PHX::DataLayout> vector_layout_vector =
basis->functional_grad;
73 TEUCHOS_ASSERT(
basis->isVectorBasis());
76 std::string orientationFieldName =
basis->name() +
" Orientation";
77 dof_orientation = PHX::MDField<ScalarT,Cell,NODE>(orientationFieldName,basis_layout);
89 gatherFieldNormals = PHX::MDField<ScalarT,Cell,NODE,Dim>(dof_name+
"_Normals",vector_layout_vector);
90 this->addEvaluatedField(gatherFieldNormals);
92 this->setName(
"Gather Normals");
96 template<
typename EvalT,
typename Traits>
102 this->utils.setFieldData(gatherFieldNormals,fm);
108 gatherFieldNormals.dimension(0),
109 gatherFieldNormals.dimension(1),
110 gatherFieldNormals.dimension(2));
114 template<
typename EvalT,
typename Traits>
122 const shards::CellTopology & parentCell = *
basis->getCellTopology();
123 int cellDim = parentCell.getDimension();
124 int numFaces = gatherFieldNormals.dimension(1);
129 Kokkos::DynRankView<ScalarT,PHX::Device> refFaceTanU =
Kokkos::createDynRankView(gatherFieldNormals.get_static_view(),
"refFaceTanU",numFaces,cellDim);
130 Kokkos::DynRankView<ScalarT,PHX::Device> refFaceTanV =
Kokkos::createDynRankView(gatherFieldNormals.get_static_view(),
"refFaceTanV",numFaces,cellDim);
131 for(
int i=0;i<numFaces;i++) {
132 Kokkos::DynRankView<double,PHX::Device> refTanU = Kokkos::DynRankView<double,PHX::Device>(
"refTanU",cellDim);
133 Kokkos::DynRankView<double,PHX::Device> refTanV = Kokkos::DynRankView<double,PHX::Device>(
"refTanV",cellDim);
134 Intrepid2::CellTools<double>::getReferenceFaceTangents(refTanU, refTanV, i, parentCell);
135 for(
int d=0;d<cellDim;d++) {
136 refFaceTanU(i,d) = refTanU(d);
137 refFaceTanV(i,d) = refTanV(d);
146 for(index_t c=0;c<workset.
num_cells;c++) {
147 for(
int f = 0; f < numFaces; f++) {
149 std::vector<double> faceTanU(3);
150 std::vector<double> faceTanV(3);
151 for(
int i = 0; i < cellDim; i++) {
152 faceTanU[i] = Sacado::ScalarValue<ScalarT>::eval(
pointValues.
jac(c,f,i,0)*refFaceTanU(f,0))
153 + Sacado::ScalarValue<ScalarT>::eval(
pointValues.
jac(c,f,i,1)*refFaceTanU(f,1))
154 + Sacado::ScalarValue<ScalarT>::eval(
pointValues.
jac(c,f,i,2)*refFaceTanU(f,2));
155 faceTanV[i] = Sacado::ScalarValue<ScalarT>::eval(
pointValues.
jac(c,f,i,0)*refFaceTanV(f,0))
156 + Sacado::ScalarValue<ScalarT>::eval(
pointValues.
jac(c,f,i,1)*refFaceTanV(f,1))
157 + Sacado::ScalarValue<ScalarT>::eval(
pointValues.
jac(c,f,i,2)*refFaceTanV(f,2));
161 gatherFieldNormals(c,f,0) = (faceTanU[1]*faceTanV[2] - faceTanU[2]*faceTanV[1])*
dof_orientation(c,f);
162 gatherFieldNormals(c,f,1) = (faceTanU[2]*faceTanV[0] - faceTanU[0]*faceTanV[2])*
dof_orientation(c,f);
163 gatherFieldNormals(c,f,2) = (faceTanU[0]*faceTanV[1] - faceTanU[1]*faceTanV[0])*
dof_orientation(c,f);
PointValues2< ScalarT, PHX::MDField > pointValues
Interpolates basis DOF values to IP DOF values.
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.
PHX::MDField< const ScalarT, Cell, BASIS > dof_orientation
Teuchos::RCP< const panzer::PointRule > pointRule
void evaluateFields(typename Traits::EvalData d)
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
Array< Scalar, Cell, IP, Dim, Dim, void, void, void, void > jac
void setupArrays(const Teuchos::RCP< const panzer::PointRule > &pr, const ArrayFactory &af)
Sizes/allocates memory for arrays.
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)
Kokkos::DynRankView< ScalarT, PHX::Device > faceNormal