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( cellCenter.rank() != 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>
101 getReferenceVertex( Kokkos::DynRankView<cellVertexValueType,cellVertexProperties...> cellVertex,
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( cellVertex.rank() != 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( subcellVertices.rank() != 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>
172 getReferenceNode( Kokkos::DynRankView<cellNodeValueType,cellNodeProperties...> cellNode,
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( cellNode.rank() != 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 subcellNodeValueType,
class ...subcellNodeProperties>
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( subcellNodes.rank() != 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 refEdgeTangentValueType,
class ...refEdgeTangentProperties>
244 getReferenceEdgeTangent( Kokkos::DynRankView<refEdgeTangentValueType,refEdgeTangentProperties...> refEdgeTangent,
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( refEdgeTangent.rank() != 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 refFaceTanValueType,
class ...refFaceTanProperties>
283 Kokkos::DynRankView<refFaceTanValueType,refFaceTanProperties...> 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( refFaceTanU.rank() != 1 || refFaceTanV.rank() != 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 refSideNormalValueType,
class ...refSideNormalProperties>
327 getReferenceSideNormal( Kokkos::DynRankView<refSideNormalValueType,refSideNormalProperties...> refSideNormal,
328 const ordinal_type sideOrd,
329 const shards::CellTopology parentCell ) {
330#ifdef HAVE_INTREPID2_DEBUG
331 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
332 parentCell.getDimension() != 3, std::invalid_argument,
333 ">>> ERROR (Intrepid2::CellTools::getReferenceSideNormal): two or three-dimensional parent cell required");
336 INTREPID2_TEST_FOR_EXCEPTION( sideOrd < 0 || sideOrd >=
static_cast<ordinal_type
>(parentCell.getSubcellCount(parentCell.getDimension() - 1)), std::invalid_argument,
337 ">>> ERROR (Intrepid2::CellTools::getReferenceSideNormal): side ordinal out of bounds");
339 constexpr bool is_accessible = Kokkos::Impl::MemorySpaceAccess<
341 typename decltype(refSideNormal)::memory_space>::accessible;
342 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceSideNormal(..): output view's memory space is not compatible with DeviceType");
344 const auto dim = parentCell.getDimension();
347 getReferenceEdgeTangent(refSideNormal, sideOrd, parentCell);
350 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,1),
351 KOKKOS_LAMBDA (
const int &) {
352 refSideNormalValueType tmp = refSideNormal(0);
353 refSideNormal(0) = refSideNormal(1);
354 refSideNormal(1) = -tmp;
358 getReferenceFaceNormal(refSideNormal, sideOrd, parentCell);
363 template<
typename DeviceType>
364 template<
typename refFaceNormalValueType,
class ...refFaceNormalProperties>
367 getReferenceFaceNormal( Kokkos::DynRankView<refFaceNormalValueType,refFaceNormalProperties...> refFaceNormal,
368 const ordinal_type faceOrd,
369 const shards::CellTopology parentCell ) {
370#ifdef HAVE_INTREPID2_DEBUG
371 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
372 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): three-dimensional parent cell required");
374 INTREPID2_TEST_FOR_EXCEPTION( faceOrd < 0 || faceOrd >=
static_cast<ordinal_type
>(parentCell.getSubcellCount(2)), std::invalid_argument,
375 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): face ordinal out of bounds");
377 INTREPID2_TEST_FOR_EXCEPTION( refFaceNormal.rank() != 1, std::invalid_argument,
378 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): rank = 1 required for output array");
380 INTREPID2_TEST_FOR_EXCEPTION( refFaceNormal.extent(0) != parentCell.getDimension(), std::invalid_argument,
381 ">>> ERROR (Intrepid2::CellTools::getReferenceFaceNormal): dim0 (spatial dim) must match that of parent cell");
383 constexpr bool is_accessible = Kokkos::Impl::MemorySpaceAccess<
385 typename decltype(refFaceNormal)::memory_space>::accessible;
386 static_assert(is_accessible,
"CellTools<DeviceType>::getReferenceFaceNormal(..): output view's memory space is not compatible with DeviceType");
389 const auto dim = parentCell.getDimension();
390 auto vcprop = Kokkos::common_view_alloc_prop(refFaceNormal);
391 using common_value_type =
typename decltype(vcprop)::value_type;
392 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanU ( Kokkos::view_alloc(
"CellTools::getReferenceFaceNormal::refFaceTanU", vcprop), dim );
393 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanV ( Kokkos::view_alloc(
"CellTools::getReferenceFaceNormal::refFaceTanV", vcprop), dim );
394 getReferenceFaceTangents(refFaceTanU, refFaceTanV, faceOrd, parentCell);
399 template<
typename DeviceType>
400 template<
typename edgeTangentValueType,
class ...edgeTangentProperties,
401 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
405 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
406 const ordinal_type worksetEdgeOrd,
407 const shards::CellTopology parentCell ) {
408#ifdef HAVE_INTREPID2_DEBUG
409 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3 &&
410 parentCell.getDimension() != 2, std::invalid_argument,
411 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): 2D or 3D parent cell required." );
414 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.rank() != 3, std::invalid_argument,
415 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents requires rank 3." );
416 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(2) != 2 &&
417 edgeTangents.extent(2) != 3, std::invalid_argument,
418 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension(2) must be 2 or 3." );
421 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.rank() != 4, std::invalid_argument,
422 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians requires rank 4." );
423 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 2 &&
424 worksetJacobians.extent(2) != 3, std::invalid_argument,
425 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) must be 2 or 3." );
426 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
427 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) and (3) must match each other." );
430 for (
auto i=0;i<3;++i) {
431 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
432 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension (i) does not match to worksetJacobians dimension(i)." );
435 constexpr bool are_accessible =
436 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
437 typename decltype(edgeTangents)::memory_space>::accessible &&
438 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
439 typename decltype(worksetJacobians)::memory_space>::accessible;
440 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalEdgeTangents(..): input/output views' memory spaces are not compatible with DeviceType");
444 const auto dim = parentCell.getDimension();
445 auto vcprop = Kokkos::common_view_alloc_prop(edgeTangents);
446 using common_value_type =
typename decltype(vcprop)::value_type;
447 Kokkos::DynRankView< common_value_type, DeviceType > refEdgeTan ( Kokkos::view_alloc(
"CellTools::getPhysicalEdgeTangents::refEdgeTan", vcprop), dim);
448 getReferenceEdgeTangent(refEdgeTan, worksetEdgeOrd, parentCell);
454 namespace FunctorCellTools {
456 template<
typename tangentViewType,
457 typename faceOrdinalViewType,
458 typename parametrizationViewType
461 tangentViewType refEdgeTan_;
462 const faceOrdinalViewType edgeOrdView_;
463 const parametrizationViewType edgeParametrization_;
465 KOKKOS_INLINE_FUNCTION
467 const faceOrdinalViewType edgeOrdView,
468 const parametrizationViewType edgeParametrization)
469 : refEdgeTan_(refEdgeTan), edgeOrdView_(edgeOrdView), edgeParametrization_(edgeParametrization){};
471 KOKKOS_INLINE_FUNCTION
472 void operator()(
const size_type ic)
const {
473 for (size_type d=0; d<refEdgeTan_.extent(1); d++) {
474 refEdgeTan_(ic,d) = edgeParametrization_(edgeOrdView_(ic), d, 1);
480 template<
typename DeviceType>
481 template<
typename edgeTangentValueType,
class ...edgeTangentProperties,
482 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
483 typename edgeOrdValueType,
class ...edgeOrdProperties>
487 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
488 const Kokkos::DynRankView<edgeOrdValueType,edgeOrdProperties...> worksetEdgeOrds,
489 const shards::CellTopology parentCell ) {
490#ifdef HAVE_INTREPID2_DEBUG
491 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3 &&
492 parentCell.getDimension() != 2, std::invalid_argument,
493 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): 2D or 3D parent cell required." );
496 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.rank() != 3, std::invalid_argument,
497 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents requires rank 3." );
498 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(2) != 2 &&
499 edgeTangents.extent(2) != 3, std::invalid_argument,
500 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension(2) must be 2 or 3." );
502 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(0) != worksetEdgeOrds.extent(0), std::invalid_argument,
503 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetEdgeOrds extent 0 should match that of edgeTangents." );
506 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.rank() != 4, std::invalid_argument,
507 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians requires rank 4." );
508 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 2 &&
509 worksetJacobians.extent(2) != 3, std::invalid_argument,
510 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) must be 2 or 3." );
511 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
512 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): worksetJacobians dimension(2) and (3) must match each other." );
515 for (
auto i=0;i<3;++i) {
516 INTREPID2_TEST_FOR_EXCEPTION( edgeTangents.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
517 ">>> ERROR (Intrepid2::CellTools::getPhysicalEdgeTangents): edgeTangents dimension (i) does not match to worksetJacobians dimension(i)." );
520 constexpr bool are_accessible =
521 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
522 typename decltype(edgeTangents)::memory_space>::accessible &&
523 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
524 typename decltype(worksetJacobians)::memory_space>::accessible &&
525 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
526 typename decltype(worksetEdgeOrds)::memory_space>::accessible;
527 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalEdgeTangents(..): input/output views' memory spaces are not compatible with DeviceType");
531 const ordinal_type dim = parentCell.getDimension();
532 auto vcprop = Kokkos::common_view_alloc_prop(edgeTangents);
533 using common_value_type =
typename decltype(vcprop)::value_type;
534 Kokkos::DynRankView< common_value_type, DeviceType > refEdgeTan ( Kokkos::view_alloc(
"CellTools::getPhysicalEdgeTangents::refEdgeTan", vcprop), edgeTangents.extent(0), dim);
539 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,refEdgeTan.extent(0)), FunctorType(refEdgeTan, worksetEdgeOrds, edgeMap) );
541 typename DeviceType::execution_space().fence();
545 template<
typename DeviceType>
546 template<
typename faceTanValueType,
class ...faceTanProperties,
547 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
551 Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanV,
552 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
553 const ordinal_type worksetFaceOrd,
554 const shards::CellTopology parentCell ) {
555#ifdef HAVE_INTREPID2_DEBUG
556 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
557 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): three-dimensional parent cell required");
560 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.rank() != 3 ||
561 faceTanV.rank() != 3, std::invalid_argument,
562 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V must have rank 3." );
564 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(2) != 3 ||
565 faceTanV.extent(2) != 3, std::invalid_argument,
566 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (2) must be 3." );
568 for (
auto i=0;i<3;++i) {
569 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != faceTanV.extent(i), std::invalid_argument,
570 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (i) must match each other." );
574 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.rank() != 4, std::invalid_argument,
575 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians must have rank 4." );
577 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
578 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) must be 3." );
580 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
581 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) and dimension(3) must match." );
584 for (
auto i=0;i<3;++i) {
585 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
586 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(i) and faceTan dimension (i) must match." );
589 constexpr bool are_accessible =
590 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
591 typename decltype(faceTanU)::memory_space>::accessible &&
592 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
593 typename decltype(worksetJacobians)::memory_space>::accessible;
594 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceTangents(..): input/output views' memory spaces are not compatible with DeviceType");
597 const auto dim = parentCell.getDimension();
599 auto vcprop = Kokkos::common_view_alloc_prop(faceTanU);
600 using common_value_type =
typename decltype(vcprop)::value_type;
601 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanU", vcprop), dim);
602 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanV", vcprop), dim);
604 getReferenceFaceTangents(refFaceTanU, refFaceTanV, worksetFaceOrd, parentCell);
610 namespace FunctorCellTools {
612 template<
typename tangentsViewType,
613 typename faceOrdinalViewType,
614 typename parametrizationViewType
617 tangentsViewType refFaceTanU_;
618 tangentsViewType refFaceTanV_;
619 const faceOrdinalViewType faceOrdView_;
620 const parametrizationViewType faceParametrization_;
622 KOKKOS_INLINE_FUNCTION
624 tangentsViewType refFaceTanV,
625 const faceOrdinalViewType faceOrdView,
626 const parametrizationViewType faceParametrization)
627 : refFaceTanU_(refFaceTanU), refFaceTanV_(refFaceTanV), faceOrdView_(faceOrdView), faceParametrization_(faceParametrization){};
629 KOKKOS_INLINE_FUNCTION
630 void operator()(
const size_type ic)
const {
631 for (size_type d=0; d<refFaceTanU_.extent(1); d++) {
632 refFaceTanU_(ic,d) = faceParametrization_(faceOrdView_(ic), d, 1);
633 refFaceTanV_(ic,d) = faceParametrization_(faceOrdView_(ic), d, 2);
640 template<
typename DeviceType>
641 template<
typename faceTanValueType,
class ...faceTanProperties,
642 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
643 typename faceOrdValueType,
class ...faceOrdProperties>
647 Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanV,
648 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
649 const Kokkos::DynRankView<faceOrdValueType,faceOrdProperties...> worksetFaceOrds,
650 const shards::CellTopology parentCell ) {
651#ifdef HAVE_INTREPID2_DEBUG
652 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
653 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): three-dimensional parent cell required");
656 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.rank() != 3 ||
657 faceTanV.rank() != 3, std::invalid_argument,
658 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V must have rank 3." );
660 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(2) != 3 ||
661 faceTanV.extent(2) != 3, std::invalid_argument,
662 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (2) must be 3." );
664 for (
auto i=0;i<3;++i) {
665 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != faceTanV.extent(i), std::invalid_argument,
666 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): faceTan U,V dimension (i) must match each other." );
669 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(0) != worksetFaceOrds.extent(0), std::invalid_argument,
670 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetFaceOrds extent 0 should match that of faceTanU." );
674 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.rank() != 4, std::invalid_argument,
675 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians must have rank 4." );
677 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
678 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) must be 3." );
680 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
681 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(2) and dimension(3) must match." );
684 for (
auto i=0;i<3;++i) {
685 INTREPID2_TEST_FOR_EXCEPTION( faceTanU.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
686 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceTangents): worksetJacobians dimension(i) and faceTan dimension (i) must match." );
689 constexpr bool are_accessible =
690 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
691 typename decltype(faceTanU)::memory_space>::accessible &&
692 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
693 typename decltype(worksetJacobians)::memory_space>::accessible &&
694 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
695 typename decltype(worksetFaceOrds)::memory_space>::accessible;
696 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceTangents(..): input/output views' memory spaces are not compatible with DeviceType");
699 const ordinal_type dim = parentCell.getDimension();
701 auto vcprop = Kokkos::common_view_alloc_prop(faceTanU);
702 using common_value_type =
typename decltype(vcprop)::value_type;
703 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanU", vcprop), faceTanU.extent(0), dim);
704 Kokkos::DynRankView< common_value_type, DeviceType > refFaceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceTangents::refFaceTanV", vcprop), faceTanV.extent(0), dim);
709 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,refFaceTanU.extent(0)), FunctorType(refFaceTanU, refFaceTanV, worksetFaceOrds, faceMap) );
711 typename DeviceType::execution_space().fence();
716 namespace FunctorCellTools {
718 template<
typename normalsViewType,
719 typename tangentsViewType>
721 normalsViewType edgeNormals_;
722 const tangentsViewType edgeTangents_;
724 KOKKOS_INLINE_FUNCTION
726 const tangentsViewType refEdgeTangents)
727 : edgeNormals_(edgeNormals), edgeTangents_(refEdgeTangents){};
729 KOKKOS_INLINE_FUNCTION
730 void operator()(
const size_type iter)
const {
732 unrollIndex( cell, pt, edgeTangents_.extent(0),
733 edgeTangents_.extent(1), iter );
734 edgeNormals_(cell,pt,0) = edgeTangents_(cell,pt,1);
735 edgeNormals_(cell,pt,1) = -edgeTangents_(cell,pt,0);
742 template<
typename DeviceType>
743 template<
typename sideNormalValueType,
class ...sideNormalProperties,
744 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
748 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
749 const ordinal_type worksetSideOrd,
750 const shards::CellTopology parentCell ) {
751#ifdef HAVE_INTREPID2_DEBUG
752 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
753 parentCell.getDimension() != 3, std::invalid_argument,
754 ">>> ERROR (Intrepid2::CellTools::getPhysicalSideNormals): two or three-dimensional parent cell required");
757 INTREPID2_TEST_FOR_EXCEPTION( worksetSideOrd < 0 ||
758 worksetSideOrd >=
static_cast<ordinal_type
>(parentCell.getSubcellCount(parentCell.getDimension() - 1)), std::invalid_argument,
759 ">>> ERROR (Intrepid2::CellTools::getPhysicalSideNormals): side ordinal out of bounds");
761 constexpr bool are_accessible =
762 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
763 typename decltype(sideNormals)::memory_space>::accessible &&
764 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
765 typename decltype(worksetJacobians)::memory_space>::accessible;
766 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalSideNormals(..): input/output views' memory spaces are not compatible with DeviceType");
768 const auto dim = parentCell.getDimension();
772 auto vcprop = Kokkos::common_view_alloc_prop(sideNormals);
773 using common_value_type =
typename decltype(vcprop)::value_type;
774 Kokkos::DynRankView< common_value_type, DeviceType > edgeTangents ( Kokkos::view_alloc(
"CellTools::getPhysicalSideNormals::edgeTan", vcprop),
775 sideNormals.extent(0),
776 sideNormals.extent(1),
777 sideNormals.extent(2));
778 getPhysicalEdgeTangents(edgeTangents, worksetJacobians, worksetSideOrd, parentCell);
782 const auto loopSize = edgeTangents.extent(0)*edgeTangents.extent(1);
783 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
784 Kokkos::parallel_for( policy, FunctorType(sideNormals, edgeTangents) );
787 getPhysicalFaceNormals(sideNormals, worksetJacobians, worksetSideOrd, parentCell);
792 template<
typename DeviceType>
793 template<
typename sideNormalValueType,
class ...sideNormalProperties,
794 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
795 typename edgeOrdValueType,
class ...edgeOrdProperties>
799 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
800 const Kokkos::DynRankView<edgeOrdValueType,edgeOrdProperties...> worksetSideOrds,
801 const shards::CellTopology parentCell ) {
802#ifdef HAVE_INTREPID2_DEBUG
803 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 2 &&
804 parentCell.getDimension() != 3, std::invalid_argument,
805 ">>> ERROR (Intrepid2::CellTools::getPhysicalSideNormals): two or three-dimensional parent cell required");
807 constexpr bool are_accessible =
808 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
809 typename decltype(sideNormals)::memory_space>::accessible &&
810 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
811 typename decltype(worksetJacobians)::memory_space>::accessible &&
812 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
813 typename decltype(worksetSideOrds)::memory_space>::accessible;
814 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalSideNormals(..): input/output views' memory spaces are not compatible with DeviceType");
816 const auto dim = parentCell.getDimension();
820 auto vcprop = Kokkos::common_view_alloc_prop(sideNormals);
821 using common_value_type =
typename decltype(vcprop)::value_type;
822 Kokkos::DynRankView< common_value_type, DeviceType > edgeTangents ( Kokkos::view_alloc(
"CellTools::getPhysicalSideNormals::edgeTan", vcprop),
823 sideNormals.extent(0),
824 sideNormals.extent(1),
825 sideNormals.extent(2));
826 getPhysicalEdgeTangents(edgeTangents, worksetJacobians, worksetSideOrds, parentCell);
830 const auto loopSize = edgeTangents.extent(0)*edgeTangents.extent(1);
831 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
832 Kokkos::parallel_for( policy, FunctorType(sideNormals, edgeTangents) );
835 getPhysicalFaceNormals(sideNormals, worksetJacobians, worksetSideOrds, parentCell);
840 template<
typename DeviceType>
841 template<
typename faceNormalValueType,
class ...faceNormalProperties,
842 typename worksetJacobianValueType,
class ...worksetJacobianProperties>
846 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
847 const ordinal_type worksetFaceOrd,
848 const shards::CellTopology parentCell ) {
849#ifdef HAVE_INTREPID2_DEBUG
850 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
851 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): three-dimensional parent cell required." );
854 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.rank() != 3, std::invalid_argument,
855 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals must have a rank 3." );
856 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(2) != 3, std::invalid_argument,
857 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (2) must be 3." );
860 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.rank() != 4, std::invalid_argument,
861 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians must have a rank 4." );
862 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
863 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must be 3." );
864 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
865 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must match to dimension (3)." );
868 for (
auto i=0;i<3;++i) {
869 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
870 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (i) must match to worksetJacobians dimension (i)." );
873 constexpr bool are_accessible =
874 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
875 typename decltype(faceNormals)::memory_space>::accessible &&
876 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
877 typename decltype(worksetJacobians)::memory_space>::accessible;
878 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceNormals(..): input/output views' memory spaces are not compatible with DeviceType");
883 const auto worksetSize = worksetJacobians.extent(0);
884 const auto facePtCount = worksetJacobians.extent(1);
885 const auto dim = parentCell.getDimension();
887 auto vcprop = Kokkos::common_view_alloc_prop(faceNormals);
888 using common_value_type =
typename decltype(vcprop)::value_type;
889 Kokkos::DynRankView< common_value_type, DeviceType > faceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanU", vcprop), worksetSize, facePtCount, dim);
890 Kokkos::DynRankView< common_value_type, DeviceType > faceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanV", vcprop), worksetSize, facePtCount, dim);
892 getPhysicalFaceTangents(faceTanU, faceTanV,
897 typename DeviceType::execution_space().fence();
901 template<
typename DeviceType>
902 template<
typename faceNormalValueType,
class ...faceNormalProperties,
903 typename worksetJacobianValueType,
class ...worksetJacobianProperties,
904 typename faceOrdValueType,
class ...faceOrdProperties>
908 const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
909 const Kokkos::DynRankView<faceOrdValueType,faceOrdProperties...> worksetFaceOrds,
910 const shards::CellTopology parentCell ) {
911#ifdef HAVE_INTREPID2_DEBUG
912 INTREPID2_TEST_FOR_EXCEPTION( parentCell.getDimension() != 3, std::invalid_argument,
913 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): three-dimensional parent cell required." );
916 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.rank() != 3, std::invalid_argument,
917 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals must have a rank 3." );
918 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(2) != 3, std::invalid_argument,
919 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (2) must be 3." );
920 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(0) != worksetFaceOrds.extent(0), std::invalid_argument,
921 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetFaceOrds extent 0 should match that of faceNormals." );
924 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.rank() != 4, std::invalid_argument,
925 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians must have a rank 4." );
926 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != 3, std::invalid_argument,
927 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must be 3." );
928 INTREPID2_TEST_FOR_EXCEPTION( worksetJacobians.extent(2) != worksetJacobians.extent(3), std::invalid_argument,
929 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): worksetJacobians dimension (2) must match to dimension (3)." );
932 for (
auto i=0;i<3;++i) {
933 INTREPID2_TEST_FOR_EXCEPTION( faceNormals.extent(i) != worksetJacobians.extent(i), std::invalid_argument,
934 ">>> ERROR (Intrepid2::CellTools::getPhysicalFaceNormals): faceNormals dimension (i) must match to worksetJacobians dimension (i)." );
937 constexpr bool are_accessible =
938 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
939 typename decltype(faceNormals)::memory_space>::accessible &&
940 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
941 typename decltype(worksetJacobians)::memory_space>::accessible &&
942 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
943 typename decltype(worksetFaceOrds)::memory_space>::accessible;
944 static_assert(are_accessible,
"CellTools<DeviceType>::getPhysicalFaceNormals(..): input/output views' memory spaces are not compatible with DeviceType");
949 const auto worksetSize = worksetJacobians.extent(0);
950 const auto facePtCount = worksetJacobians.extent(1);
951 const auto dim = parentCell.getDimension();
953 auto vcprop = Kokkos::common_view_alloc_prop(faceNormals);
954 using common_value_type =
typename decltype(vcprop)::value_type;
955 Kokkos::DynRankView< common_value_type, DeviceType > faceTanU ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanU", vcprop), worksetSize, facePtCount, dim);
956 Kokkos::DynRankView< common_value_type, DeviceType > faceTanV ( Kokkos::view_alloc(
"CellTools::getPhysicalFaceNormals::faceTanV", vcprop), worksetSize, facePtCount, dim);
958 getPhysicalFaceTangents(faceTanU, faceTanV,
963 typename DeviceType::execution_space().fence();
static ConstViewType get(const unsigned cellTopoKey)
Retrieves the Cartesian coordinates of a reference cell barycenter.
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...