49 #ifndef __INTREPID2_CELLTOOLS_DEF_NODE_INFO_HPP__
50 #define __INTREPID2_CELLTOOLS_DEF_NODE_INFO_HPP__
53 #if defined (__clang__) && !defined (__INTEL_COMPILER)
54 #pragma clang system_header
65 template<
typename DeviceType>
66 template<
typename cellCenterValueType,
class ...cellCenterProperties>
70 const shards::CellTopology cell ) {
71 #ifdef HAVE_INTREPID2_DEBUG
72 INTREPID2_TEST_FOR_EXCEPTION( !hasReferenceCell(cell), std::invalid_argument,
73 ">>> ERROR (Intrepid2::CellTools::getReferenceCellCenter): the specified cell topology does not have a reference cell." );
75 INTREPID2_TEST_FOR_EXCEPTION( rank(cellCenter) != 1, std::invalid_argument,
76 ">>> ERROR (Intrepid2::CellTools::getReferenceCellCenter): cellCenter must have rank 1." );
78 INTREPID2_TEST_FOR_EXCEPTION( cellCenter.extent(0) < cell.getDimension(), std::invalid_argument,
79 ">>> ERROR (Intrepid2::CellTools::getReferenceCellCenter): cellCenter must have dimension at least as large as cell.getDimension()." );
82 constexpr
bool is_accessible =
83 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
84 typename decltype(cellCenter)::memory_space>::accessible;
85 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceCellCenter(..): output view's memory space is not compatible with DeviceType");
87 const ordinal_type dim = cell.getDimension();
91 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,dim),
92 KOKKOS_LAMBDA (
const int &i) {cellCenter(i) = refCellCenter(i);}
97 template<
typename DeviceType>
98 template<
typename cellVertexValueType,
class ...cellVertexProperties>
102 const shards::CellTopology cell,
103 const ordinal_type vertexOrd ) {
104 #ifdef HAVE_INTREPID2_DEBUG
105 INTREPID2_TEST_FOR_EXCEPTION( !hasReferenceCell(cell), std::invalid_argument,
106 ">>> ERROR (Intrepid2::CellTools::getReferenceVertex): the specified cell topology does not have a reference cell." );
108 INTREPID2_TEST_FOR_EXCEPTION( (vertexOrd < 0) || vertexOrd > static_cast<ordinal_type>(cell.getVertexCount()), std::invalid_argument,
109 ">>> ERROR (Intrepid2::CellTools::getReferenceVertex): invalid node ordinal for the specified cell topology." );
111 INTREPID2_TEST_FOR_EXCEPTION( rank(cellVertex) != 1, std::invalid_argument,
112 ">>> ERROR (Intrepid2::CellTools::getReferenceNode): cellNodes must have rank 1." );
114 INTREPID2_TEST_FOR_EXCEPTION( cellVertex.extent(0) < cell.getDimension(), std::invalid_argument,
115 ">>> ERROR (Intrepid2::CellTools::getReferenceNode): cellNodes must have dimension at least as large as cell.getDimension()." );
118 constexpr
bool is_accessible =
119 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
120 typename decltype(cellVertex)::memory_space>::accessible;
121 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceVertex(..): output view's memory space is not compatible with DeviceType");
125 ordinal_type dim = cell.getDimension();
126 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,dim),
127 KOKKOS_LAMBDA (
const int &i) {cellVertex(i) = refNodes(vertexOrd,i);}
132 template<
typename DeviceType>
133 template<
typename subcellVertexValueType,
class ...subcellVertexProperties>
137 const ordinal_type subcellDim,
138 const ordinal_type subcellOrd,
139 const shards::CellTopology parentCell ) {
140 #ifdef HAVE_INTREPID2_DEBUG
141 INTREPID2_TEST_FOR_EXCEPTION( !hasReferenceCell(parentCell), std::invalid_argument,
142 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellVertices): the specified cell topology does not have a reference cell." );
144 INTREPID2_TEST_FOR_EXCEPTION( subcellDim > static_cast<ordinal_type>(parentCell.getDimension()), std::invalid_argument,
145 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellVertices): subcell dimension cannot exceed cell dimension." );
147 INTREPID2_TEST_FOR_EXCEPTION( subcellOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(subcellDim)), std::invalid_argument,
148 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellVertices): subcell ordinal cannot exceed subcell count." );
151 INTREPID2_TEST_FOR_EXCEPTION( rank(subcellVertices) != 2, std::invalid_argument,
152 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellVertices): subcellVertieces must have rank 2." );
155 INTREPID2_TEST_FOR_EXCEPTION( subcellVertices.extent(0) != parentCell.getNodeCount(subcellDim, subcellOrd), std::invalid_argument,
156 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellVertices): subcellVertieces dimension(0) must match to parent cell vertex count." );
158 INTREPID2_TEST_FOR_EXCEPTION( subcellVertices.extent(1) != parentCell.getDimension(), std::invalid_argument,
159 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellVertices): subcellVertieces dimension(1) must match to parent cell dimension." );
161 getReferenceSubcellNodes(subcellVertices,
168 template<
typename DeviceType>
169 template<
typename cellNodeValueType,
class ...cellNodeProperties>
173 const shards::CellTopology cell,
174 const ordinal_type nodeOrd ) {
175 #ifdef HAVE_INTREPID2_DEBUG
176 INTREPID2_TEST_FOR_EXCEPTION( nodeOrd >= static_cast<ordinal_type>(cell.getNodeCount()), std::invalid_argument,
177 ">>> ERROR (Intrepid2::CellTools::getReferenceNode): invalid node ordinal for the specified cell topology." );
179 INTREPID2_TEST_FOR_EXCEPTION( rank(cellNode) != 1, std::invalid_argument,
180 ">>> ERROR (Intrepid2::CellTools::getReferenceNode): cellNode must have rank 1." );
182 INTREPID2_TEST_FOR_EXCEPTION( cellNode.extent(0) < cell.getDimension(), std::invalid_argument,
183 ">>> ERROR (Intrepid2::CellTools::getReferenceNode): cellNode must have dimension at least as large as cell.getDimension()." );
186 constexpr
bool is_accessible =
187 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
188 typename decltype(cellNode)::memory_space>::accessible;
189 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceNode(..): output view's memory space is not compatible with DeviceType");
193 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,cell.getDimension()),
194 KOKKOS_LAMBDA (
const int &i) {cellNode(i) = refNodes(nodeOrd,i);}
198 template<
typename DeviceType>
199 template<
typename SubcellNodeViewType>
203 const ordinal_type subcellDim,
204 const ordinal_type subcellOrd,
205 const shards::CellTopology parentCell ) {
206 #ifdef HAVE_INTREPID2_DEBUG
207 INTREPID2_TEST_FOR_EXCEPTION( !hasReferenceCell(parentCell), std::invalid_argument,
208 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellNodes): the specified cell topology does not have a reference cell.");
210 INTREPID2_TEST_FOR_EXCEPTION( subcellDim > static_cast<ordinal_type>(parentCell.getDimension()), std::invalid_argument,
211 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellNodes): subcell dimension out of range.");
213 INTREPID2_TEST_FOR_EXCEPTION( subcellOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(subcellDim)), std::invalid_argument,
214 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellNodes): subcell ordinal out of range.");
217 INTREPID2_TEST_FOR_EXCEPTION( rank(subcellNodes) != 2, std::invalid_argument,
218 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellNodes): subcellNodes must have rank 2.");
220 INTREPID2_TEST_FOR_EXCEPTION( subcellNodes.extent(0) != parentCell.getNodeCount(subcellDim, subcellOrd), std::invalid_argument,
221 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellNodes): subcellNodes dimension (0) must match to node count in cell.");
223 INTREPID2_TEST_FOR_EXCEPTION( subcellNodes.extent(1) != parentCell.getDimension(), std::invalid_argument,
224 ">>> ERROR (Intrepid2::CellTools::getReferenceSubcellNodes): subcellNodes dimension (1) must match to cell dimension.");
228 const auto subcNodeCount = parentCell.getNodeCount(subcellDim, subcellOrd);
231 for (size_type subcNodeOrd=0;subcNodeOrd<subcNodeCount;++subcNodeOrd) {
233 const auto cellNodeOrd = parentCell.getNodeMap(subcellDim, subcellOrd, subcNodeOrd);
235 auto dst = Kokkos::subdynrankview(subcellNodes, subcNodeOrd, Kokkos::ALL());
236 getReferenceNode(dst, parentCell, cellNodeOrd);
240 template<
typename DeviceType>
241 template<
typename RefEdgeTangentViewType>
245 const ordinal_type edgeOrd,
246 const shards::CellTopology parentCell ) {
247 #ifdef HAVE_INTREPID2_DEBUG
248 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
249 parentCell.getDimension() != 3, std::invalid_argument,
250 ">>> ERROR (Intrepid2::CellTools::getReferenceEdgeTangent): two or three-dimensional parent cell required");
252 INTREPID2_TEST_FOR_EXCEPTION( rank(refEdgeTangent) != 1, std::invalid_argument,
253 ">>> ERROR (Intrepid2::CellTools::getReferenceEdgeTangent): rank = 1 required for output arrays");
255 INTREPID2_TEST_FOR_EXCEPTION( refEdgeTangent.extent(0) != parentCell.getDimension(), std::invalid_argument,
256 ">>> ERROR (Intrepid2::CellTools::getReferenceEdgeTangent): output array size is required to match space dimension");
258 INTREPID2_TEST_FOR_EXCEPTION( edgeOrd < 0 ||
259 edgeOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(1)), std::invalid_argument,
260 ">>> ERROR (Intrepid2::CellTools::getReferenceEdgeTangent): edge ordinal out of bounds");
263 constexpr
bool is_accessible =
264 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
265 typename decltype(refEdgeTangent)::memory_space>::accessible;
266 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceEdgeTangent(..): output view's memory space is not compatible with DeviceType");
272 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,parentCell.getDimension()),
273 KOKKOS_LAMBDA (
const int &i) {refEdgeTangent(i) = edgeMap(edgeOrd, i, 1);}
278 template<
typename DeviceType>
279 template<
typename RefFaceTanViewType>
283 RefFaceTanViewType refFaceTanV,
284 const ordinal_type faceOrd,
285 const shards::CellTopology parentCell ) {
286 #ifdef HAVE_INTREPID2_DEBUG
287 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
288 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceTangents): three-dimensional parent cell required");
290 INTREPID2_TEST_FOR_EXCEPTION( faceOrd < 0 || faceOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(2)), std::invalid_argument,
291 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceTangents): face ordinal out of bounds");
293 INTREPID2_TEST_FOR_EXCEPTION( rank(refFaceTanU) != 1 || rank(refFaceTanV) != 1, std::invalid_argument,
294 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceTangents): rank = 1 required for output arrays");
296 INTREPID2_TEST_FOR_EXCEPTION( refFaceTanU.extent(0) != parentCell.getDimension(), std::invalid_argument,
297 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceTangents): dim0 (spatial dim) must match that of parent cell");
299 INTREPID2_TEST_FOR_EXCEPTION( refFaceTanV.extent(0) != parentCell.getDimension(), std::invalid_argument,
300 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceTangents): dim0 (spatial dim) must match that of parent cell");
302 constexpr
bool is_accessible =
303 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
304 typename decltype(refFaceTanU)::memory_space>::accessible;
305 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceFaceTangents(..): output views' memory spaces are not compatible with DeviceType");
316 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,parentCell.getDimension()),
317 KOKKOS_LAMBDA (
const int &i) {
318 refFaceTanU(i) = faceMap(faceOrd, i, 1);
319 refFaceTanV(i) = faceMap(faceOrd, i, 2);
323 template<
typename DeviceType>
324 template<
typename RefS
ideNormalViewType>
328 const ordinal_type sideOrd,
329 const shards::CellTopology parentCell ) {
330 using refSideNormalValueType =
typename RefSideNormalViewType::non_const_value_type;
331 #ifdef HAVE_INTREPID2_DEBUG
332 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
333 parentCell.getDimension() != 3, std::invalid_argument,
334 ">>> ERROR (Intrepid2::CellTools::getReferenceSideNormal): two or three-dimensional parent cell required");
337 INTREPID2_TEST_FOR_EXCEPTION( sideOrd < 0 || sideOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(parentCell.getDimension() - 1)), std::invalid_argument,
338 ">>> ERROR (Intrepid2::CellTools::getReferenceSideNormal): side ordinal out of bounds");
340 constexpr
bool is_accessible = Kokkos::Impl::MemorySpaceAccess<
342 typename decltype(refSideNormal)::memory_space>::accessible;
343 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceSideNormal(..): output view's memory space is not compatible with DeviceType");
345 const auto dim = parentCell.getDimension();
348 getReferenceEdgeTangent(refSideNormal, sideOrd, parentCell);
351 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,1),
352 KOKKOS_LAMBDA (
const int &) {
353 refSideNormalValueType tmp = refSideNormal(0);
354 refSideNormal(0) = refSideNormal(1);
355 refSideNormal(1) = -tmp;
359 getReferenceFaceNormal(refSideNormal, sideOrd, parentCell);
364 template<
typename DeviceType>
365 template<
typename RefFaceNormalViewType>
369 const ordinal_type faceOrd,
370 const shards::CellTopology parentCell ) {
371 #ifdef HAVE_INTREPID2_DEBUG
372 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
373 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): three-dimensional parent cell required");
375 INTREPID2_TEST_FOR_EXCEPTION( faceOrd < 0 || faceOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(2)), std::invalid_argument,
376 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): face ordinal out of bounds");
378 INTREPID2_TEST_FOR_EXCEPTION( rank(refFaceNormal) != 1, std::invalid_argument,
379 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): rank = 1 required for output array");
381 INTREPID2_TEST_FOR_EXCEPTION( refFaceNormal.extent(0) != parentCell.getDimension(), std::invalid_argument,
382 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): dim0 (spatial dim) must match that of parent cell");
384 constexpr
bool is_accessible = Kokkos::Impl::MemorySpaceAccess<
386 typename decltype(refFaceNormal)::memory_space>::accessible;
387 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceFaceNormal(..): output view's memory space is not compatible with DeviceType");
390 const auto dim = parentCell.getDimension();
391 auto vcprop = Kokkos::common_view_alloc_prop(refFaceNormal);
392 using common_value_type =
typename decltype(vcprop)::value_type;
393 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanU ( Kokkos::view_alloc(
"CellTools::getReferenceFaceNormal::refFaceTanU", vcprop), dim );
394 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanV ( Kokkos::view_alloc(
"CellTools::getReferenceFaceNormal::refFaceTanV", vcprop), dim );
395 getReferenceFaceTangents(refFaceTanU, refFaceTanV, faceOrd, parentCell);
400 template<
typename DeviceType>
401 template<
typename edgeTangentValueType,
class ...edgeTangentProperties,
402 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
406 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
407 const ordinal_type worksetEdgeOrd,
408 const shards::CellTopology parentCell ) {
409 #ifdef HAVE_INTREPID2_DEBUG
410 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3 &&
411 parentCell.getDimension() != 2, std::invalid_argument,
412 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): 2D or 3D parent cell required." );
415 INTREPID2_TEST_FOR_EXCEPTION( rank(edgeTangents) != 3, std::invalid_argument,
416 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents requires rank 3." );
417 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(2) != 2 &&
418 edgeTangents.extent(2) != 3, std::invalid_argument,
419 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension(2) must be 2 or 3." );
422 INTREPID2_TEST_FOR_EXCEPTION( rank(worksetJacobians) != 4, std::invalid_argument,
423 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians requires rank 4." );
424 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 2 &&
425 worksetJacobians.extent(2) != 3, std::invalid_argument,
426 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) must be 2 or 3." );
427 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
428 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) and (3) must match each other." );
431 for (
auto i=0;i<3;++i) {
432 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
433 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension (i) does not match to worksetJacobians dimension(i)." );
436 constexpr
bool are_accessible =
437 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
438 typename decltype(edgeTangents)::memory_space>::accessible &&
439 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
440 typename decltype(worksetJacobians)::memory_space>::accessible;
441 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalEdgeTangents(..): input/output views' memory spaces are not compatible with DeviceType");
445 const auto dim = parentCell.getDimension();
446 auto vcprop = Kokkos::common_view_alloc_prop(edgeTangents);
447 using common_value_type =
typename decltype(vcprop)::value_type;
448 Kokkos::DynRankView< common_value_type, DeviceType > refEdgeTan ( Kokkos::view_alloc(
"CellTools::getPhysicalEdgeTangents::refEdgeTan", vcprop), dim);
449 getReferenceEdgeTangent(refEdgeTan, worksetEdgeOrd, parentCell);
455 namespace FunctorCellTools {
457 template<
typename tangentViewType,
458 typename faceOrdinalViewType,
459 typename parametrizationViewType
462 tangentViewType refEdgeTan_;
463 const faceOrdinalViewType edgeOrdView_;
464 const parametrizationViewType edgeParametrization_;
466 KOKKOS_INLINE_FUNCTION
468 const faceOrdinalViewType edgeOrdView,
469 const parametrizationViewType edgeParametrization)
470 : refEdgeTan_(refEdgeTan), edgeOrdView_(edgeOrdView), edgeParametrization_(edgeParametrization){};
472 KOKKOS_INLINE_FUNCTION
473 void operator()(
const size_type ic)
const {
474 for (size_type d=0; d<refEdgeTan_.extent(1); d++) {
475 refEdgeTan_(ic,d) = edgeParametrization_(edgeOrdView_(ic), d, 1);
481 template<
typename DeviceType>
482 template<
typename edgeTangentValueType,
class ...edgeTangentProperties,
483 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
484 typename edgeOrdValueType,
class ...edgeOrdProperties>
488 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
489 const Kokkos::DynRankView<edgeOrdValueType,edgeOrdProperties...> worksetEdgeOrds,
490 const shards::CellTopology parentCell ) {
491 #ifdef HAVE_INTREPID2_DEBUG
492 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3 &&
493 parentCell.getDimension() != 2, std::invalid_argument,
494 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): 2D or 3D parent cell required." );
497 INTREPID2_TEST_FOR_EXCEPTION( rank(edgeTangents) != 3, std::invalid_argument,
498 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents requires rank 3." );
499 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(2) != 2 &&
500 edgeTangents.extent(2) != 3, std::invalid_argument,
501 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension(2) must be 2 or 3." );
503 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(0) != worksetEdgeOrds.extent(0), std::invalid_argument,
504 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetEdgeOrds extent 0 should match that of edgeTangents." );
507 INTREPID2_TEST_FOR_EXCEPTION( rank(worksetJacobians) != 4, std::invalid_argument,
508 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians requires rank 4." );
509 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 2 &&
510 worksetJacobians.extent(2) != 3, std::invalid_argument,
511 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) must be 2 or 3." );
512 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
513 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) and (3) must match each other." );
516 for (
auto i=0;i<3;++i) {
517 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
518 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension (i) does not match to worksetJacobians dimension(i)." );
521 constexpr
bool are_accessible =
522 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
523 typename decltype(edgeTangents)::memory_space>::accessible &&
524 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
525 typename decltype(worksetJacobians)::memory_space>::accessible &&
526 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
527 typename decltype(worksetEdgeOrds)::memory_space>::accessible;
528 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalEdgeTangents(..): input/output views' memory spaces are not compatible with DeviceType");
532 const ordinal_type dim = parentCell.getDimension();
533 auto vcprop = Kokkos::common_view_alloc_prop(edgeTangents);
534 using common_value_type =
typename decltype(vcprop)::value_type;
535 Kokkos::DynRankView< common_value_type, DeviceType > refEdgeTan ( Kokkos::view_alloc(
"CellTools::getPhysicalEdgeTangents::refEdgeTan", vcprop), edgeTangents.extent(0), dim);
540 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,refEdgeTan.extent(0)), FunctorType(refEdgeTan, worksetEdgeOrds, edgeMap) );
542 typename DeviceType::execution_space().fence();
546 template<
typename DeviceType>
547 template<
typename faceTanValueType,
class ...faceTanProperties,
548 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
552 Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanV,
553 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
554 const ordinal_type worksetFaceOrd,
555 const shards::CellTopology parentCell ) {
556 #ifdef HAVE_INTREPID2_DEBUG
557 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
558 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): three-dimensional parent cell required");
561 INTREPID2_TEST_FOR_EXCEPTION( rank(faceTanU) != 3 ||
562 rank(faceTanV) != 3, std::invalid_argument,
563 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V must have rank 3." );
565 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(2) != 3 ||
566 faceTanV.extent(2) != 3, std::invalid_argument,
567 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (2) must be 3." );
569 for (
auto i=0;i<3;++i) {
570 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != faceTanV.extent(i), std::invalid_argument,
571 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (i) must match each other." );
575 INTREPID2_TEST_FOR_EXCEPTION( rank(worksetJacobians) != 4, std::invalid_argument,
576 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians must have rank 4." );
578 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
579 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) must be 3." );
581 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
582 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) and dimension(3) must match." );
585 for (
auto i=0;i<3;++i) {
586 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
587 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(i) and faceTan dimension (i) must match." );
590 constexpr
bool are_accessible =
591 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
592 typename decltype(faceTanU)::memory_space>::accessible &&
593 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
594 typename decltype(worksetJacobians)::memory_space>::accessible;
595 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceTangents(..): input/output views' memory spaces are not compatible with DeviceType");
598 const auto dim = parentCell.getDimension();
600 auto vcprop = Kokkos::common_view_alloc_prop(faceTanU);
601 using common_value_type =
typename decltype(vcprop)::value_type;
602 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanU", vcprop), dim);
603 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanV", vcprop), dim);
605 getReferenceFaceTangents(refFaceTanU, refFaceTanV, worksetFaceOrd, parentCell);
611 namespace FunctorCellTools {
613 template<
typename tangentsViewType,
614 typename faceOrdinalViewType,
615 typename parametrizationViewType
618 tangentsViewType refFaceTanU_;
619 tangentsViewType refFaceTanV_;
620 const faceOrdinalViewType faceOrdView_;
621 const parametrizationViewType faceParametrization_;
623 KOKKOS_INLINE_FUNCTION
625 tangentsViewType refFaceTanV,
626 const faceOrdinalViewType faceOrdView,
627 const parametrizationViewType faceParametrization)
628 : refFaceTanU_(refFaceTanU), refFaceTanV_(refFaceTanV), faceOrdView_(faceOrdView), faceParametrization_(faceParametrization){};
630 KOKKOS_INLINE_FUNCTION
631 void operator()(
const size_type ic)
const {
632 for (size_type d=0; d<refFaceTanU_.extent(1); d++) {
633 refFaceTanU_(ic,d) = faceParametrization_(faceOrdView_(ic), d, 1);
634 refFaceTanV_(ic,d) = faceParametrization_(faceOrdView_(ic), d, 2);
641 template<
typename DeviceType>
642 template<
typename faceTanValueType,
class ...faceTanProperties,
643 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
644 typename faceOrdValueType,
class ...faceOrdProperties>
648 Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanV,
649 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
650 const Kokkos::DynRankView<faceOrdValueType,faceOrdProperties...> worksetFaceOrds,
651 const shards::CellTopology parentCell ) {
652 #ifdef HAVE_INTREPID2_DEBUG
653 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
654 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): three-dimensional parent cell required");
657 INTREPID2_TEST_FOR_EXCEPTION( rank(faceTanU) != 3 ||
658 rank(faceTanV) != 3, std::invalid_argument,
659 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V must have rank 3." );
661 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(2) != 3 ||
662 faceTanV.extent(2) != 3, std::invalid_argument,
663 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (2) must be 3." );
665 for (
auto i=0;i<3;++i) {
666 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != faceTanV.extent(i), std::invalid_argument,
667 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (i) must match each other." );
670 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(0) != worksetFaceOrds.extent(0), std::invalid_argument,
671 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetFaceOrds extent 0 should match that of faceTanU." );
675 INTREPID2_TEST_FOR_EXCEPTION( rank(worksetJacobians) != 4, std::invalid_argument,
676 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians must have rank 4." );
678 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
679 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) must be 3." );
681 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
682 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) and dimension(3) must match." );
685 for (
auto i=0;i<3;++i) {
686 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
687 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(i) and faceTan dimension (i) must match." );
690 constexpr
bool are_accessible =
691 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
692 typename decltype(faceTanU)::memory_space>::accessible &&
693 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
694 typename decltype(worksetJacobians)::memory_space>::accessible &&
695 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
696 typename decltype(worksetFaceOrds)::memory_space>::accessible;
697 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceTangents(..): input/output views' memory spaces are not compatible with DeviceType");
700 const ordinal_type dim = parentCell.getDimension();
702 auto vcprop = Kokkos::common_view_alloc_prop(faceTanU);
703 using common_value_type =
typename decltype(vcprop)::value_type;
704 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanU", vcprop), faceTanU.extent(0), dim);
705 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanV", vcprop), faceTanV.extent(0), dim);
710 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,refFaceTanU.extent(0)), FunctorType(refFaceTanU, refFaceTanV, worksetFaceOrds, faceMap) );
712 typename DeviceType::execution_space().fence();
717 namespace FunctorCellTools {
719 template<
typename normalsViewType,
720 typename tangentsViewType>
722 normalsViewType edgeNormals_;
723 const tangentsViewType edgeTangents_;
725 KOKKOS_INLINE_FUNCTION
727 const tangentsViewType refEdgeTangents)
728 : edgeNormals_(edgeNormals), edgeTangents_(refEdgeTangents){};
730 KOKKOS_INLINE_FUNCTION
731 void operator()(
const size_type iter)
const {
733 unrollIndex( cell, pt, edgeTangents_.extent(0),
734 edgeTangents_.extent(1), iter );
735 edgeNormals_(cell,pt,0) = edgeTangents_(cell,pt,1);
736 edgeNormals_(cell,pt,1) = -edgeTangents_(cell,pt,0);
743 template<
typename DeviceType>
744 template<
typename sideNormalValueType,
class ...sideNormalProperties,
745 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
749 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
750 const ordinal_type worksetSideOrd,
751 const shards::CellTopology parentCell ) {
752 #ifdef HAVE_INTREPID2_DEBUG
753 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
754 parentCell.getDimension() != 3, std::invalid_argument,
755 ">>> ERROR (Intrepid2::CellTools::getPhysicalSideNormals): two or three-dimensional parent cell required");
758 INTREPID2_TEST_FOR_EXCEPTION( worksetSideOrd < 0 ||
759 worksetSideOrd >= static_cast<ordinal_type>(parentCell.getSubcellCount(parentCell.getDimension() - 1)), std::invalid_argument,
760 ">>> ERROR (Intrepid2::CellTools::getPhysicalSideNormals): side ordinal out of bounds");
762 constexpr
bool are_accessible =
763 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
764 typename decltype(sideNormals)::memory_space>::accessible &&
765 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
766 typename decltype(worksetJacobians)::memory_space>::accessible;
767 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalSideNormals(..): input/output views' memory spaces are not compatible with DeviceType");
769 const auto dim = parentCell.getDimension();
773 auto vcprop = Kokkos::common_view_alloc_prop(sideNormals);
774 using common_value_type =
typename decltype(vcprop)::value_type;
775 Kokkos::DynRankView< common_value_type, DeviceType > edgeTangents ( Kokkos::view_alloc(
"CellTools::getPhysicalSideNormals::edgeTan", vcprop),
776 sideNormals.extent(0),
777 sideNormals.extent(1),
778 sideNormals.extent(2));
779 getPhysicalEdgeTangents(edgeTangents, worksetJacobians, worksetSideOrd, parentCell);
783 const auto loopSize = edgeTangents.extent(0)*edgeTangents.extent(1);
784 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
785 Kokkos::parallel_for( policy, FunctorType(sideNormals, edgeTangents) );
788 getPhysicalFaceNormals(sideNormals, worksetJacobians, worksetSideOrd, parentCell);
793 template<
typename DeviceType>
794 template<
typename sideNormalValueType,
class ...sideNormalProperties,
795 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
796 typename edgeOrdValueType,
class ...edgeOrdProperties>
800 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
801 const Kokkos::DynRankView<edgeOrdValueType,edgeOrdProperties...> worksetSideOrds,
802 const shards::CellTopology parentCell ) {
803 #ifdef HAVE_INTREPID2_DEBUG
804 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
805 parentCell.getDimension() != 3, std::invalid_argument,
806 ">>> ERROR (Intrepid2::CellTools::getPhysicalSideNormals): two or three-dimensional parent cell required");
808 constexpr
bool are_accessible =
809 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
810 typename decltype(sideNormals)::memory_space>::accessible &&
811 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
812 typename decltype(worksetJacobians)::memory_space>::accessible &&
813 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
814 typename decltype(worksetSideOrds)::memory_space>::accessible;
815 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalSideNormals(..): input/output views' memory spaces are not compatible with DeviceType");
817 const auto dim = parentCell.getDimension();
821 auto vcprop = Kokkos::common_view_alloc_prop(sideNormals);
822 using common_value_type =
typename decltype(vcprop)::value_type;
823 Kokkos::DynRankView< common_value_type, DeviceType > edgeTangents ( Kokkos::view_alloc(
"CellTools::getPhysicalSideNormals::edgeTan", vcprop),
824 sideNormals.extent(0),
825 sideNormals.extent(1),
826 sideNormals.extent(2));
827 getPhysicalEdgeTangents(edgeTangents, worksetJacobians, worksetSideOrds, parentCell);
831 const auto loopSize = edgeTangents.extent(0)*edgeTangents.extent(1);
832 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
833 Kokkos::parallel_for( policy, FunctorType(sideNormals, edgeTangents) );
836 getPhysicalFaceNormals(sideNormals, worksetJacobians, worksetSideOrds, parentCell);
841 template<
typename DeviceType>
842 template<
typename faceNormalValueType,
class ...faceNormalProperties,
843 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
847 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
848 const ordinal_type worksetFaceOrd,
849 const shards::CellTopology parentCell ) {
850 #ifdef HAVE_INTREPID2_DEBUG
851 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
852 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): three-dimensional parent cell required." );
855 INTREPID2_TEST_FOR_EXCEPTION( rank(faceNormals) != 3, std::invalid_argument,
856 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals must have a rank 3." );
857 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(2) != 3, std::invalid_argument,
858 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (2) must be 3." );
861 INTREPID2_TEST_FOR_EXCEPTION( rank(worksetJacobians) != 4, std::invalid_argument,
862 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians must have a rank 4." );
863 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
864 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must be 3." );
865 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
866 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must match to dimension (3)." );
869 for (
auto i=0;i<3;++i) {
870 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
871 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (i) must match to worksetJacobians dimension (i)." );
874 constexpr
bool are_accessible =
875 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
876 typename decltype(faceNormals)::memory_space>::accessible &&
877 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
878 typename decltype(worksetJacobians)::memory_space>::accessible;
879 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceNormals(..): input/output views' memory spaces are not compatible with DeviceType");
884 const auto worksetSize = worksetJacobians.extent(0);
885 const auto facePtCount = worksetJacobians.extent(1);
886 const auto dim = parentCell.getDimension();
888 auto vcprop = Kokkos::common_view_alloc_prop(faceNormals);
889 using common_value_type =
typename decltype(vcprop)::value_type;
890 Kokkos::DynRankView< common_value_type, DeviceType > faceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanU", vcprop), worksetSize, facePtCount, dim);
891 Kokkos::DynRankView< common_value_type, DeviceType > faceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanV", vcprop), worksetSize, facePtCount, dim);
893 getPhysicalFaceTangents(faceTanU, faceTanV,
898 typename DeviceType::execution_space().fence();
902 template<
typename DeviceType>
903 template<
typename faceNormalValueType,
class ...faceNormalProperties,
904 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
905 typename faceOrdValueType,
class ...faceOrdProperties>
909 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
910 const Kokkos::DynRankView<faceOrdValueType,faceOrdProperties...> worksetFaceOrds,
911 const shards::CellTopology parentCell ) {
912 #ifdef HAVE_INTREPID2_DEBUG
913 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
914 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): three-dimensional parent cell required." );
917 INTREPID2_TEST_FOR_EXCEPTION( rank(faceNormals) != 3, std::invalid_argument,
918 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals must have a rank 3." );
919 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(2) != 3, std::invalid_argument,
920 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (2) must be 3." );
921 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(0) != worksetFaceOrds.extent(0), std::invalid_argument,
922 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetFaceOrds extent 0 should match that of faceNormals." );
925 INTREPID2_TEST_FOR_EXCEPTION( rank(worksetJacobians) != 4, std::invalid_argument,
926 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians must have a rank 4." );
927 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
928 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must be 3." );
929 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
930 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must match to dimension (3)." );
933 for (
auto i=0;i<3;++i) {
934 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
935 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (i) must match to worksetJacobians dimension (i)." );
938 constexpr
bool are_accessible =
939 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
940 typename decltype(faceNormals)::memory_space>::accessible &&
941 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
942 typename decltype(worksetJacobians)::memory_space>::accessible &&
943 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
944 typename decltype(worksetFaceOrds)::memory_space>::accessible;
945 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceNormals(..): input/output views' memory spaces are not compatible with DeviceType");
950 const auto worksetSize = worksetJacobians.extent(0);
951 const auto facePtCount = worksetJacobians.extent(1);
952 const auto dim = parentCell.getDimension();
954 auto vcprop = Kokkos::common_view_alloc_prop(faceNormals);
955 using common_value_type =
typename decltype(vcprop)::value_type;
956 Kokkos::DynRankView< common_value_type, DeviceType > faceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanU", vcprop), worksetSize, facePtCount, dim);
957 Kokkos::DynRankView< common_value_type, DeviceType > faceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanV", vcprop), worksetSize, facePtCount, dim);
959 getPhysicalFaceTangents(faceTanU, faceTanV,
964 typename DeviceType::execution_space().fence();
static ConstViewType get(const unsigned cellTopoKey)
Retrieves the Cartesian coordinates of reference cell nodes.
static ConstViewType get(const ordinal_type subcellDim, const unsigned parentCellKey)
Returns a Kokkos view with the coefficients of the parametrization maps for the edges or faces of a r...
static ConstViewType get(const unsigned cellTopoKey)
Retrieves the Cartesian coordinates of a reference cell barycenter.