43 #ifndef PANZER_DIRICHLET_RESIDUAL_EDGEBASIS_IMPL_HPP
44 #define PANZER_DIRICHLET_RESIDUAL_EDGEBASIS_IMPL_HPP
50 #include "Intrepid2_Kernels.hpp"
51 #include "Intrepid2_CellTools.hpp"
52 #include "Intrepid2_OrientationTools.hpp"
54 #include "Phalanx_Print.hpp"
57 #include "Kokkos_ViewFactory.hpp"
62 template<
typename EvalT,
typename Traits>
67 std::string residual_name = p.
get<std::string>(
"Residual Name");
72 std::string field_name = p.
get<std::string>(
"DOF Name");
73 std::string dof_name = field_name+
"_"+pointRule->getName();
74 std::string value_name = p.
get<std::string>(
"Value Name");
82 TEUCHOS_ASSERT(basis_layout->extent(0)==vector_layout_dof->extent(0));
83 TEUCHOS_ASSERT(basis_layout->extent(1)==vector_layout_dof->extent(1));
84 TEUCHOS_ASSERT(Teuchos::as<unsigned>(basis->dimension())==vector_layout_dof->extent(2));
85 TEUCHOS_ASSERT(vector_layout_vector->extent(0)==vector_layout_dof->extent(0));
86 TEUCHOS_ASSERT(vector_layout_vector->extent(1)==vector_layout_dof->extent(1));
87 TEUCHOS_ASSERT(vector_layout_vector->extent(2)==vector_layout_dof->extent(2));
97 pointValues.setupArrays(pointRule);
100 this->addNonConstDependentField(pointValues.jac);
102 this->addEvaluatedField(residual);
103 this->addDependentField(dof);
104 this->addDependentField(value);
106 std::string n =
"Dirichlet Residual Edge Basis Evaluator";
111 template<
typename EvalT,
typename Traits>
119 this->utils.setFieldData(pointValues.jac,fm);
121 const auto cellTopo = *basis->getCellTopology();
122 const int edgeDim = 1;
123 const int faceDim = 2;
124 if(cellTopo.getDimension() > edgeDim)
125 edgeParam = Intrepid2::RefSubcellParametrization<Kokkos::HostSpace>::get(edgeDim, cellTopo.getKey());
127 if(cellTopo.getDimension() > faceDim)
128 faceParam = Intrepid2::RefSubcellParametrization<Kokkos::HostSpace>::get(faceDim, cellTopo.getKey());
132 template<
typename EvalT,
typename Traits>
142 residual.deep_copy(
ScalarT(0.0));
147 const int subcellOrd = this->wda(workset).subcell_index;
149 const auto cellTopo = *basis->getCellTopology();
150 const auto worksetJacobians = pointValues.jac.get_view();
152 const int cellDim = cellTopo.getDimension();
153 const int edgeDim = 1;
154 const int faceDim = 2;
156 auto intrepid_basis = basis->getIntrepid2Basis();
163 auto residual_h = Kokkos::create_mirror_view(residual.get_static_view());
164 auto dof_h = Kokkos::create_mirror_view(dof.get_static_view());
165 auto value_h = Kokkos::create_mirror_view(value.get_static_view());
166 Kokkos::deep_copy(dof_h, dof.get_static_view());
167 Kokkos::deep_copy(value_h, value.get_static_view());
168 auto worksetJacobians_h = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(),worksetJacobians);
169 switch (subcellDim) {
171 if (intrepid_basis->getDofCount(1, subcellOrd)) {
172 auto ortEdgeTan = Kokkos::subview(work, 0, Kokkos::ALL());
173 auto phyEdgeTan = Kokkos::subview(work, 1, Kokkos::ALL());
175 const int ndofsEdge = intrepid_basis->getDofCount(1, subcellOrd);
176 const int numEdges = cellTopo.getEdgeCount();
177 int edgeOrts[4] = {};
179 for(index_t c=0;c<workset.
num_cells;c++) {
180 orientations->at(details.
cell_local_ids[c]).getEdgeOrientation(edgeOrts, numEdges);
182 Intrepid2::Impl::OrientationTools::getRefSubcellTangents(work,
184 cellTopo.getKey(edgeDim,subcellOrd),
186 edgeOrts[subcellOrd]);
188 for (
int i=0;i<ndofsEdge;++i) {
189 const int b = intrepid_basis->getDofOrdinal(1, subcellOrd, i);
190 auto J = Kokkos::subview(worksetJacobians_h, c, b, Kokkos::ALL(), Kokkos::ALL());
191 Intrepid2::Kernels::Serial::matvec_product(phyEdgeTan, J, ortEdgeTan);
193 for(
int d=0;d<cellDim;d++) {
194 residual_h(c,b) += (dof_h(c,b,d)-value_h(c,b,d))*phyEdgeTan(d);
202 const int numEdges = cellTopo.getEdgeCount();
203 const int numFaces = cellTopo.getFaceCount();
207 auto ortEdgeTan = Kokkos::subview(work, 0, Kokkos::ALL());
208 auto phyEdgeTan = Kokkos::subview(work, 1, Kokkos::ALL());
210 const int numEdgesOfFace= cellTopo.getEdgeCount(2, subcellOrd);
212 int edgeOrts[12] = {};
213 for(index_t c=0;c<workset.
num_cells;c++) {
214 for (
int i=0;i<numEdgesOfFace;++i) {
216 const int edgeOrd = Intrepid2::Orientation::getEdgeOrdinalOfFace(i, subcellOrd, cellTopo);
217 const int b = edgeOrd;
218 orientations->at(details.
cell_local_ids[c]).getEdgeOrientation(edgeOrts, numEdges);
220 Intrepid2::Impl::OrientationTools::getRefSubcellTangents(work,
222 cellTopo.getKey(edgeDim,edgeOrd),
227 auto J = Kokkos::subview(worksetJacobians_h, c, b, Kokkos::ALL(), Kokkos::ALL());
228 Intrepid2::Kernels::Serial::matvec_product(phyEdgeTan, J, ortEdgeTan);
230 for(
int d=0;d<dof.extent_int(2);d++) {
231 residual_h(c,b) += (dof_h(c,b,d)-value_h(c,b,d))*phyEdgeTan(d);
238 if (intrepid_basis->getDofCount(2, subcellOrd)) {
239 auto ortFaceTanU = Kokkos::subview(work, 0, Kokkos::ALL());
240 auto ortFaceTanV = Kokkos::subview(work, 1, Kokkos::ALL());
241 auto phyFaceTanU = Kokkos::subview(work, 2, Kokkos::ALL());
242 auto phyFaceTanV = Kokkos::subview(work, 3, Kokkos::ALL());
244 int faceOrts[6] = {};
245 for(index_t c=0;c<workset.
num_cells;c++) {
246 orientations->at(details.
cell_local_ids[c]).getFaceOrientation(faceOrts, numFaces);
248 Intrepid2::Impl::OrientationTools::getRefSubcellTangents(work,
250 cellTopo.getKey(faceDim,subcellOrd),
252 faceOrts[subcellOrd]);
254 for(
int b=0;b<dof.extent_int(1);b++) {
255 auto J = Kokkos::subview(worksetJacobians, c, b, Kokkos::ALL(), Kokkos::ALL());
256 Intrepid2::Kernels::Serial::matvec_product(phyFaceTanU, J, ortFaceTanU);
257 Intrepid2::Kernels::Serial::matvec_product(phyFaceTanV, J, ortFaceTanV);
259 for(
int d=0;d<dof.extent_int(2);d++) {
260 residual_h(c,b) += (dof_h(c,b,d)-value_h(c,b,d))*phyFaceTanU(d);
261 residual_h(c,b) += (dof_h(c,b,d)-value_h(c,b,d))*phyFaceTanV(d);
270 Kokkos::deep_copy(residual.get_static_view(), residual_h);
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.
Teuchos::RCP< const std::vector< Intrepid2::Orientation > > orientations_
int num_cells
DEPRECATED - use: numCells()
T & get(const std::string &name, T def_value)
DirichletResidual_EdgeBasis(const Teuchos::ParameterList &p)
std::vector< size_t > cell_local_ids
int subcell_dim
DEPRECATED - use: getSubcellDimension()
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
void evaluateFields(typename Traits::EvalData d)
typename EvalT::ScalarT ScalarT
#define TEUCHOS_ASSERT(assertion_test)