49 #ifndef __INTREPID2_CELLTOOLS_DEF_HPP__ 50 #define __INTREPID2_CELLTOOLS_DEF_HPP__ 53 #if defined (__clang__) && !defined (__INTEL_COMPILER) 54 #pragma clang system_header 68 namespace FunctorCellTools {
72 template<
typename jacobianViewType,
73 typename worksetCellType,
74 typename basisGradType>
76 jacobianViewType _jacobian;
77 const worksetCellType _worksetCells;
78 const basisGradType _basisGrads;
82 KOKKOS_INLINE_FUNCTION
84 worksetCellType worksetCells_,
85 basisGradType basisGrads_,
88 : _jacobian(jacobian_), _worksetCells(worksetCells_), _basisGrads(basisGrads_),
89 _startCell(startCell_), _endCell(endCell_) {}
91 KOKKOS_INLINE_FUNCTION
92 void operator()(
const ordinal_type cell,
93 const ordinal_type point)
const {
95 const ordinal_type dim = _jacobian.extent(2);
97 const ordinal_type gradRank = _basisGrads.rank();
100 const ordinal_type cardinality = _basisGrads.extent(0);
101 for (ordinal_type i=0;i<dim;++i)
102 for (ordinal_type j=0;j<dim;++j) {
103 _jacobian(cell, point, i, j) = 0;
104 for (ordinal_type bf=0;bf<cardinality;++bf)
105 _jacobian(cell, point, i, j) += _worksetCells(cell+_startCell, bf, i) * _basisGrads(bf, point, j);
110 const ordinal_type cardinality = _basisGrads.extent(1);
111 for (ordinal_type i=0;i<dim;++i)
112 for (ordinal_type j=0;j<dim;++j) {
113 _jacobian(cell, point, i, j) = 0;
114 for (ordinal_type bf=0;bf<cardinality;++bf)
115 _jacobian(cell, point, i, j) += _worksetCells(cell+_startCell, bf, i) * _basisGrads(cell, bf, point, j);
122 template<
typename DeviceType>
123 template<
class Po
intScalar>
128 const bool cellVaries = (variationTypes[0] !=
CONSTANT);
129 const bool pointVaries = (variationTypes[1] !=
CONSTANT);
136 if ( cellVaries && pointVaries )
142 else if (cellVaries || pointVaries)
156 template<
typename DeviceType>
157 template<
class Po
intScalar>
164 if ( jacDataRank == 4 )
167 auto invData =
getMatchingViewWithLabel(jacData,
"Jacobian inv data",jacData.extent(0),jacData.extent(1),jacData.extent(2),jacData.extent(3));
170 else if (jacDataRank == 3)
173 auto invData =
getMatchingViewWithLabel(jacData,
"Jacobian inv data",jacData.extent(0),jacData.extent(1),jacData.extent(2));
176 else if (jacDataRank == 2)
182 else if (jacDataRank == 1)
190 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"allocateJacobianInv requires jacobian to vary in *some* dimension…");
195 template<
typename DeviceType>
196 template<
class Po
intScalar>
201 const int CELL_DIM = 0;
202 const int POINT_DIM = 1;
203 const int D1_DIM = 2;
204 const bool cellVaries = (variationTypes[CELL_DIM] !=
CONSTANT);
205 const bool pointVaries = (variationTypes[POINT_DIM] !=
CONSTANT);
207 auto det2 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
const PointScalar &d) -> PointScalar
209 return a * d - b * c;
212 auto det3 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
213 const PointScalar &d,
const PointScalar &e,
const PointScalar &f,
214 const PointScalar &g,
const PointScalar &h,
const PointScalar &i) -> PointScalar
216 return a * det2(e,f,h,i) - b * det2(d,f,g,i) + c * det2(d,e,g,h);
219 auto det4 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
const PointScalar &d,
220 const PointScalar &e,
const PointScalar &f,
const PointScalar &g,
const PointScalar &h,
221 const PointScalar &i,
const PointScalar &j,
const PointScalar &k,
const PointScalar &l,
222 const PointScalar &m,
const PointScalar &n,
const PointScalar &o,
const PointScalar &p) -> PointScalar
224 return a * det3(f,g,h,j,k,l,n,o,p) - b * det3(e,g,h,i,k,l,m,o,p) + c * det3(e,f,h,i,j,l,m,n,p) - d * det3(e,f,g,i,j,k,m,n,o);
229 if (cellVaries && pointVaries)
234 Kokkos::parallel_for(
235 Kokkos::MDRangePolicy<ExecSpaceType,Kokkos::Rank<2>>({0,0},{data.extent_int(0),data.extent_int(1)}),
236 KOKKOS_LAMBDA (
int cellOrdinal,
int pointOrdinal) {
238 const int numDiagonals = data.extent_int(2) - blockWidth * blockWidth;
239 const int spaceDim = blockWidth + numDiagonals;
241 PointScalar det = 1.0;
249 det = data(cellOrdinal,pointOrdinal,0);
254 const auto & a = data(cellOrdinal,pointOrdinal,0);
255 const auto & b = data(cellOrdinal,pointOrdinal,1);
256 const auto & c = data(cellOrdinal,pointOrdinal,2);
257 const auto & d = data(cellOrdinal,pointOrdinal,3);
264 const auto & a = data(cellOrdinal,pointOrdinal,0);
265 const auto & b = data(cellOrdinal,pointOrdinal,1);
266 const auto & c = data(cellOrdinal,pointOrdinal,2);
267 const auto & d = data(cellOrdinal,pointOrdinal,3);
268 const auto & e = data(cellOrdinal,pointOrdinal,4);
269 const auto & f = data(cellOrdinal,pointOrdinal,5);
270 const auto & g = data(cellOrdinal,pointOrdinal,6);
271 const auto & h = data(cellOrdinal,pointOrdinal,7);
272 const auto & i = data(cellOrdinal,pointOrdinal,8);
273 det = det3(a,b,c,d,e,f,g,h,i);
279 const auto & a = data(cellOrdinal,pointOrdinal,0);
280 const auto & b = data(cellOrdinal,pointOrdinal,1);
281 const auto & c = data(cellOrdinal,pointOrdinal,2);
282 const auto & d = data(cellOrdinal,pointOrdinal,3);
283 const auto & e = data(cellOrdinal,pointOrdinal,4);
284 const auto & f = data(cellOrdinal,pointOrdinal,5);
285 const auto & g = data(cellOrdinal,pointOrdinal,6);
286 const auto & h = data(cellOrdinal,pointOrdinal,7);
287 const auto & i = data(cellOrdinal,pointOrdinal,8);
288 const auto & j = data(cellOrdinal,pointOrdinal,9);
289 const auto & k = data(cellOrdinal,pointOrdinal,10);
290 const auto & l = data(cellOrdinal,pointOrdinal,11);
291 const auto & m = data(cellOrdinal,pointOrdinal,12);
292 const auto & n = data(cellOrdinal,pointOrdinal,13);
293 const auto & o = data(cellOrdinal,pointOrdinal,14);
294 const auto & p = data(cellOrdinal,pointOrdinal,15);
295 det = det4(a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p);
300 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"blocks with block width > 4 not supported for BLOCK_PLUS_DIAGONAL");
303 const int offset = blockWidth * blockWidth;
305 for (
int d=blockWidth; d<spaceDim; d++)
307 const int index = d-blockWidth+offset;
308 det *= data(cellOrdinal,pointOrdinal,index);
310 detData(cellOrdinal,pointOrdinal) = det;
313 else if (cellVaries || pointVaries)
318 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,data.extent_int(0)),
319 KOKKOS_LAMBDA (
const int &cellPointOrdinal) {
321 const int numDiagonals = data.extent_int(2) - blockWidth * blockWidth;
322 const int spaceDim = blockWidth + numDiagonals;
324 PointScalar det = 1.0;
332 det = data(cellPointOrdinal,0);
337 const auto & a = data(cellPointOrdinal,0);
338 const auto & b = data(cellPointOrdinal,1);
339 const auto & c = data(cellPointOrdinal,2);
340 const auto & d = data(cellPointOrdinal,3);
347 const auto & a = data(cellPointOrdinal,0);
348 const auto & b = data(cellPointOrdinal,1);
349 const auto & c = data(cellPointOrdinal,2);
350 const auto & d = data(cellPointOrdinal,3);
351 const auto & e = data(cellPointOrdinal,4);
352 const auto & f = data(cellPointOrdinal,5);
353 const auto & g = data(cellPointOrdinal,6);
354 const auto & h = data(cellPointOrdinal,7);
355 const auto & i = data(cellPointOrdinal,8);
356 det = det3(a,b,c,d,e,f,g,h,i);
361 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"blocks with block width > 3 not supported for BLOCK_PLUS_DIAGONAL");
364 const int offset = blockWidth * blockWidth;
366 for (
int d=blockWidth; d<spaceDim; d++)
368 const int index = d-blockWidth+offset;
369 det *= data(cellPointOrdinal,index);
371 detData(cellPointOrdinal) = det;
378 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,1),
379 KOKKOS_LAMBDA (
const int &dummyArg) {
381 const int numDiagonals = jacobian.
extent_int(2) - blockWidth * blockWidth;
382 const int spaceDim = blockWidth + numDiagonals;
384 PointScalar det = 1.0;
397 const auto & a = data(0);
398 const auto & b = data(1);
399 const auto & c = data(2);
400 const auto & d = data(3);
407 const auto & a = data(0);
408 const auto & b = data(1);
409 const auto & c = data(2);
410 const auto & d = data(3);
411 const auto & e = data(4);
412 const auto & f = data(5);
413 const auto & g = data(6);
414 const auto & h = data(7);
415 const auto & i = data(8);
416 det = det3(a,b,c,d,e,f,g,h,i);
421 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"blocks with block width > 3 not supported for BLOCK_PLUS_DIAGONAL");
424 const int offset = blockWidth * blockWidth;
426 for (
int d=blockWidth; d<spaceDim; d++)
428 const int index = d-blockWidth+offset;
451 Kokkos::parallel_for(
"fill jacobian det", Kokkos::RangePolicy<ExecSpaceType>(0,1), KOKKOS_LAMBDA(
const int &i)
453 detData(0) = Intrepid2::Kernels::Serial::determinant(data);
458 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"jacobian's underlying view must have rank 2,3, or 4");
462 template<
typename DeviceType>
463 template<
class Po
intScalar>
468 const int CELL_DIM = 0;
469 const int POINT_DIM = 1;
470 const int D1_DIM = 2;
472 auto det2 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
const PointScalar &d) -> PointScalar
474 return a * d - b * c;
477 auto det3 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
478 const PointScalar &d,
const PointScalar &e,
const PointScalar &f,
479 const PointScalar &g,
const PointScalar &h,
const PointScalar &i) -> PointScalar
481 return a * det2(e,f,h,i) - b * det2(d,f,g,i) + c * det2(d,e,g,h);
486 const bool cellVaries = variationTypes[CELL_DIM] !=
CONSTANT;
487 const bool pointVaries = variationTypes[POINT_DIM] !=
CONSTANT;
488 if (cellVaries && pointVaries)
493 Kokkos::parallel_for(
494 Kokkos::MDRangePolicy<ExecSpaceType,Kokkos::Rank<2>>({0,0},{data.extent_int(0),data.extent_int(1)}),
495 KOKKOS_LAMBDA (
int cellOrdinal,
int pointOrdinal) {
497 const int numDiagonals = data.extent_int(2) - blockWidth * blockWidth;
498 const int spaceDim = blockWidth + numDiagonals;
506 invData(cellOrdinal,pointOrdinal,0) = 1.0 / data(cellOrdinal,pointOrdinal,0);
511 const auto & a = data(cellOrdinal,pointOrdinal,0);
512 const auto & b = data(cellOrdinal,pointOrdinal,1);
513 const auto & c = data(cellOrdinal,pointOrdinal,2);
514 const auto & d = data(cellOrdinal,pointOrdinal,3);
515 const PointScalar det = det2(a,b,c,d);
517 invData(cellOrdinal,pointOrdinal,0) = d/det;
518 invData(cellOrdinal,pointOrdinal,1) = - b/det;
519 invData(cellOrdinal,pointOrdinal,2) = - c/det;
520 invData(cellOrdinal,pointOrdinal,3) = a/det;
525 const auto & a = data(cellOrdinal,pointOrdinal,0);
526 const auto & b = data(cellOrdinal,pointOrdinal,1);
527 const auto & c = data(cellOrdinal,pointOrdinal,2);
528 const auto & d = data(cellOrdinal,pointOrdinal,3);
529 const auto & e = data(cellOrdinal,pointOrdinal,4);
530 const auto & f = data(cellOrdinal,pointOrdinal,5);
531 const auto & g = data(cellOrdinal,pointOrdinal,6);
532 const auto & h = data(cellOrdinal,pointOrdinal,7);
533 const auto & i = data(cellOrdinal,pointOrdinal,8);
534 const PointScalar det = det3(a,b,c,d,e,f,g,h,i);
537 const auto val0 = e*i - h*f;
538 const auto val1 = - d*i + g*f;
539 const auto val2 = d*h - g*e;
541 invData(cellOrdinal,pointOrdinal,0) = val0/det;
542 invData(cellOrdinal,pointOrdinal,1) = val1/det;
543 invData(cellOrdinal,pointOrdinal,2) = val2/det;
546 const auto val0 = h*c - b*i;
547 const auto val1 = a*i - g*c;
548 const auto val2 = - a*h + g*b;
550 invData(cellOrdinal,pointOrdinal,3) = val0/det;
551 invData(cellOrdinal,pointOrdinal,4) = val1/det;
552 invData(cellOrdinal,pointOrdinal,5) = val2/det;
555 const auto val0 = b*f - e*c;
556 const auto val1 = - a*f + d*c;
557 const auto val2 = a*e - d*b;
559 invData(cellOrdinal,pointOrdinal,6) = val0/det;
560 invData(cellOrdinal,pointOrdinal,7) = val1/det;
561 invData(cellOrdinal,pointOrdinal,8) = val2/det;
566 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"blocks with block width > 3 not supported for BLOCK_PLUS_DIAGONAL");
569 const int offset = blockWidth * blockWidth;
571 for (
int d=blockWidth; d<spaceDim; d++)
573 const int index = d-blockWidth+offset;
574 invData(cellOrdinal,pointOrdinal,index) = 1.0 / data(cellOrdinal,pointOrdinal,index);
578 else if (cellVaries || pointVaries)
583 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,data.extent_int(0)),
584 KOKKOS_LAMBDA (
const int &cellPointOrdinal) {
586 const int numDiagonals = data.extent_int(1) - blockWidth * blockWidth;
587 const int spaceDim = blockWidth + numDiagonals;
595 invData(cellPointOrdinal,0) = 1.0 / data(cellPointOrdinal,0);
600 const auto & a = data(cellPointOrdinal,0);
601 const auto & b = data(cellPointOrdinal,1);
602 const auto & c = data(cellPointOrdinal,2);
603 const auto & d = data(cellPointOrdinal,3);
604 const PointScalar det = det2(a,b,c,d);
606 invData(cellPointOrdinal,0) = d/det;
607 invData(cellPointOrdinal,1) = - b/det;
608 invData(cellPointOrdinal,2) = - c/det;
609 invData(cellPointOrdinal,3) = a/det;
614 const auto & a = data(cellPointOrdinal,0);
615 const auto & b = data(cellPointOrdinal,1);
616 const auto & c = data(cellPointOrdinal,2);
617 const auto & d = data(cellPointOrdinal,3);
618 const auto & e = data(cellPointOrdinal,4);
619 const auto & f = data(cellPointOrdinal,5);
620 const auto & g = data(cellPointOrdinal,6);
621 const auto & h = data(cellPointOrdinal,7);
622 const auto & i = data(cellPointOrdinal,8);
623 const PointScalar det = det3(a,b,c,d,e,f,g,h,i);
626 const auto val0 = e*i - h*f;
627 const auto val1 = - d*i + g*f;
628 const auto val2 = d*h - g*e;
630 invData(cellPointOrdinal,0) = val0/det;
631 invData(cellPointOrdinal,1) = val1/det;
632 invData(cellPointOrdinal,2) = val2/det;
635 const auto val0 = h*c - b*i;
636 const auto val1 = a*i - g*c;
637 const auto val2 = - a*h + g*b;
639 invData(cellPointOrdinal,3) = val0/det;
640 invData(cellPointOrdinal,4) = val1/det;
641 invData(cellPointOrdinal,5) = val2/det;
644 const auto val0 = b*f - e*c;
645 const auto val1 = - a*f + d*c;
646 const auto val2 = a*e - d*b;
648 invData(cellPointOrdinal,6) = val0/det;
649 invData(cellPointOrdinal,7) = val1/det;
650 invData(cellPointOrdinal,8) = val2/det;
655 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"blocks with block width > 3 not supported for BLOCK_PLUS_DIAGONAL in setJacobianInv()");
658 const int offset = blockWidth * blockWidth;
660 for (
int d=blockWidth; d<spaceDim; d++)
662 const int index = d-blockWidth+offset;
663 invData(cellPointOrdinal,index) = 1.0 / data(cellPointOrdinal,index);
672 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,1),
673 KOKKOS_LAMBDA (
const int &dummyArg) {
675 const int numDiagonals = data.extent_int(0) - blockWidth * blockWidth;
676 const int spaceDim = blockWidth + numDiagonals;
684 invData(0) = 1.0 / data(0);
689 const auto & a = data(0);
690 const auto & b = data(1);
691 const auto & c = data(2);
692 const auto & d = data(3);
693 const PointScalar det = det2(a,b,c,d);
696 invData(1) = - b/det;
697 invData(2) = - c/det;
703 const auto & a = data(0);
704 const auto & b = data(1);
705 const auto & c = data(2);
706 const auto & d = data(3);
707 const auto & e = data(4);
708 const auto & f = data(5);
709 const auto & g = data(6);
710 const auto & h = data(7);
711 const auto & i = data(8);
712 const PointScalar det = det3(a,b,c,d,e,f,g,h,i);
715 const auto val0 = e*i - h*f;
716 const auto val1 = - d*i + g*f;
717 const auto val2 = d*h - g*e;
719 invData(0) = val0/det;
720 invData(1) = val1/det;
721 invData(2) = val2/det;
724 const auto val0 = h*c - b*i;
725 const auto val1 = a*i - g*c;
726 const auto val2 = - a*h + g*b;
728 invData(3) = val0/det;
729 invData(4) = val1/det;
730 invData(5) = val2/det;
733 const auto val0 = b*f - e*c;
734 const auto val1 = - a*f + d*c;
735 const auto val2 = a*e - d*b;
737 invData(6) = val0/det;
738 invData(7) = val1/det;
739 invData(8) = val2/det;
744 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"blocks with block width > 3 not supported for BLOCK_PLUS_DIAGONAL in setJacobianInv()");
747 const int offset = blockWidth * blockWidth;
749 for (
int d=blockWidth; d<spaceDim; d++)
751 const int index = d-blockWidth+offset;
752 invData(index) = 1.0 / data(index);
776 Kokkos::parallel_for(
"fill jacobian inverse", Kokkos::RangePolicy<ExecSpaceType>(0,1), KOKKOS_LAMBDA(
const int &i)
778 Intrepid2::Kernels::Serial::inverse(invData,data);
783 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"jacobian's underlying view must have rank 2,3, or 4");
787 template<
typename DeviceType>
788 template<
typename jacobianValueType,
class ...jacobianProperties,
789 typename BasisGradientsType,
790 typename WorksetType>
793 setJacobian( Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian,
794 const WorksetType worksetCell,
795 const BasisGradientsType gradients,
const int startCell,
const int endCell)
797 constexpr
bool is_accessible =
798 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
799 typename decltype(jacobian)::memory_space>::accessible;
800 static_assert(is_accessible,
"CellTools<DeviceType>::setJacobian(..): output view's memory space is not compatible with DeviceType");
802 using JacobianViewType = Kokkos::DynRankView<jacobianValueType,jacobianProperties...>;
806 int endCellResolved = (endCell == -1) ? worksetCell.extent_int(0) : endCell;
808 using range_policy_type = Kokkos::Experimental::MDRangePolicy
809 < ExecSpaceType, Kokkos::Experimental::Rank<2>, Kokkos::IndexType<ordinal_type> >;
810 range_policy_type policy( { 0, 0 },
811 { jacobian.extent(0), jacobian.extent(1) } );
812 Kokkos::parallel_for( policy, FunctorType(jacobian, worksetCell, gradients, startCell, endCellResolved) );
815 template<
typename DeviceType>
816 template<
typename jacobianValueType,
class ...jacobianProperties,
817 typename pointValueType,
class ...pointProperties,
818 typename WorksetType,
819 typename HGradBasisType>
822 setJacobian( Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian,
823 const Kokkos::DynRankView<pointValueType,pointProperties...> points,
824 const WorksetType worksetCell,
825 const Teuchos::RCP<HGradBasisType> basis,
826 const int startCell,
const int endCell) {
827 constexpr
bool are_accessible =
828 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
829 typename decltype(jacobian)::memory_space>::accessible &&
830 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
831 typename decltype(points)::memory_space>::accessible;
832 static_assert(are_accessible,
"CellTools<DeviceType>::setJacobian(..): input/output views' memory spaces are not compatible with DeviceType");
834 #ifdef HAVE_INTREPID2_DEBUG 835 CellTools_setJacobianArgs(jacobian, points, worksetCell, basis->getBaseCellTopology(), startCell, endCell);
838 const auto cellTopo = basis->getBaseCellTopology();
839 const ordinal_type spaceDim = cellTopo.getDimension();
840 const ordinal_type numCells = jacobian.extent(0);
843 const ordinal_type pointRank = points.rank();
844 const ordinal_type numPoints = (pointRank == 2 ? points.extent(0) : points.extent(1));
845 const ordinal_type basisCardinality = basis->getCardinality();
851 auto vcprop = Kokkos::common_view_alloc_prop(points);
852 using GradViewType = Kokkos::DynRankView<typename decltype(vcprop)::value_type,DeviceType>;
859 grads = GradViewType(Kokkos::view_alloc(
"CellTools::setJacobian::grads", vcprop),basisCardinality, numPoints, spaceDim);
860 basis->getValues(grads,
867 grads = GradViewType(Kokkos::view_alloc(
"CellTools::setJacobian::grads", vcprop), numCells, basisCardinality, numPoints, spaceDim);
868 for (ordinal_type cell=0;cell<numCells;++cell)
869 basis->getValues(Kokkos::subview( grads, cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL() ),
870 Kokkos::subview( points, cell, Kokkos::ALL(), Kokkos::ALL() ),
876 setJacobian(jacobian, worksetCell, grads, startCell, endCell);
879 template<
typename DeviceType>
880 template<
typename jacobianInvValueType,
class ...jacobianInvProperties,
881 typename jacobianValueType,
class ...jacobianProperties>
884 setJacobianInv( Kokkos::DynRankView<jacobianInvValueType,jacobianInvProperties...> jacobianInv,
885 const Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian ) {
886 #ifdef HAVE_INTREPID2_DEBUG 892 template<
typename DeviceType>
893 template<
typename jacobianDetValueType,
class ...jacobianDetProperties,
894 typename jacobianValueType,
class ...jacobianProperties>
897 setJacobianDet( Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
898 const Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian ) {
899 #ifdef HAVE_INTREPID2_DEBUG
KOKKOS_INLINE_FUNCTION Kokkos::Array< int, 7 > getExtents() const
Returns an array containing the logical extents in each dimension.
KOKKOS_INLINE_FUNCTION ordinal_type getUnderlyingViewRank() const
returns the rank of the View that stores the unique data
KOKKOS_INLINE_FUNCTION const int & blockPlusDiagonalLastNonDiagonal() const
For a Data object containing data with variation type BLOCK_PLUS_DIAGONAL, returns the row and column...
one of two dimensions in a matrix; bottom-right part of matrix is diagonal
Wrapper around a Kokkos::View that allows data that is constant or repeating in various logical dimen...
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar *, DeviceType > & getUnderlyingView1() const
returns the View that stores the unique data. For rank-1 underlying containers.
KOKKOS_INLINE_FUNCTION int extent_int(const int &r) const
Returns the logical extent in the specified dimension.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar ****, DeviceType > & getUnderlyingView4() const
returns the View that stores the unique data. For rank-4 underlying containers.
Kokkos::DynRankView< typename ViewType::value_type, typename DeduceLayout< ViewType >::result_layout, typename ViewType::device_type > getMatchingViewWithLabel(const ViewType &view, const std::string &label, DimArgs... dims)
Creates and returns a view that matches the provided view in Kokkos Layout.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar **, DeviceType > & getUnderlyingView2() const
returns the View that stores the unique data. For rank-2 underlying containers.
Header file for small functions used in Intrepid2.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar ***, DeviceType > & getUnderlyingView3() const
returns the View that stores the unique data. For rank-3 underlying containers.
KOKKOS_INLINE_FUNCTION const Kokkos::Array< DataVariationType, 7 > & getVariationTypes() const
Returns an array with the variation types in each logical dimension.