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)
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);
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);
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);
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>
466 setJacobianDet(jacobianDetInv, jacobian);
472 template<
typename DeviceType>
473 template<
class Po
intScalar>
478 const int CELL_DIM = 0;
479 const int POINT_DIM = 1;
480 const int D1_DIM = 2;
482 auto det2 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
const PointScalar &d) -> PointScalar
484 return a * d - b * c;
487 auto det3 = KOKKOS_LAMBDA (
const PointScalar &a,
const PointScalar &b,
const PointScalar &c,
488 const PointScalar &d,
const PointScalar &e,
const PointScalar &f,
489 const PointScalar &g,
const PointScalar &h,
const PointScalar &i) -> PointScalar
491 return a * det2(e,f,h,i) - b * det2(d,f,g,i) + c * det2(d,e,g,h);
496 const bool cellVaries = variationTypes[CELL_DIM] !=
CONSTANT;
497 const bool pointVaries = variationTypes[POINT_DIM] !=
CONSTANT;
498 if (cellVaries && pointVaries)
503 Kokkos::parallel_for(
504 Kokkos::MDRangePolicy<ExecSpaceType,Kokkos::Rank<2>>({0,0},{data.extent_int(0),data.extent_int(1)}),
505 KOKKOS_LAMBDA (
int cellOrdinal,
int pointOrdinal) {
507 const int numDiagonals = data.extent_int(2) - blockWidth * blockWidth;
508 const int spaceDim = blockWidth + numDiagonals;
516 invData(cellOrdinal,pointOrdinal,0) = 1.0 / data(cellOrdinal,pointOrdinal,0);
521 const auto & a = data(cellOrdinal,pointOrdinal,0);
522 const auto & b = data(cellOrdinal,pointOrdinal,1);
523 const auto & c = data(cellOrdinal,pointOrdinal,2);
524 const auto & d = data(cellOrdinal,pointOrdinal,3);
525 const PointScalar det = det2(a,b,c,d);
527 invData(cellOrdinal,pointOrdinal,0) = d/det;
528 invData(cellOrdinal,pointOrdinal,1) = - b/det;
529 invData(cellOrdinal,pointOrdinal,2) = - c/det;
530 invData(cellOrdinal,pointOrdinal,3) = a/det;
535 const auto & a = data(cellOrdinal,pointOrdinal,0);
536 const auto & b = data(cellOrdinal,pointOrdinal,1);
537 const auto & c = data(cellOrdinal,pointOrdinal,2);
538 const auto & d = data(cellOrdinal,pointOrdinal,3);
539 const auto & e = data(cellOrdinal,pointOrdinal,4);
540 const auto & f = data(cellOrdinal,pointOrdinal,5);
541 const auto & g = data(cellOrdinal,pointOrdinal,6);
542 const auto & h = data(cellOrdinal,pointOrdinal,7);
543 const auto & i = data(cellOrdinal,pointOrdinal,8);
544 const PointScalar det = det3(a,b,c,d,e,f,g,h,i);
547 const auto val0 = e*i - h*f;
548 const auto val1 = - d*i + g*f;
549 const auto val2 = d*h - g*e;
551 invData(cellOrdinal,pointOrdinal,0) = val0/det;
552 invData(cellOrdinal,pointOrdinal,1) = val1/det;
553 invData(cellOrdinal,pointOrdinal,2) = val2/det;
556 const auto val0 = h*c - b*i;
557 const auto val1 = a*i - g*c;
558 const auto val2 = - a*h + g*b;
560 invData(cellOrdinal,pointOrdinal,3) = val0/det;
561 invData(cellOrdinal,pointOrdinal,4) = val1/det;
562 invData(cellOrdinal,pointOrdinal,5) = val2/det;
565 const auto val0 = b*f - e*c;
566 const auto val1 = - a*f + d*c;
567 const auto val2 = a*e - d*b;
569 invData(cellOrdinal,pointOrdinal,6) = val0/det;
570 invData(cellOrdinal,pointOrdinal,7) = val1/det;
571 invData(cellOrdinal,pointOrdinal,8) = val2/det;
579 const int offset = blockWidth * blockWidth;
581 for (
int d=blockWidth; d<spaceDim; d++)
583 const int index = d-blockWidth+offset;
584 invData(cellOrdinal,pointOrdinal,index) = 1.0 / data(cellOrdinal,pointOrdinal,index);
588 else if (cellVaries || pointVaries)
593 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,data.extent_int(0)),
594 KOKKOS_LAMBDA (
const int &cellPointOrdinal) {
596 const int numDiagonals = data.extent_int(1) - blockWidth * blockWidth;
597 const int spaceDim = blockWidth + numDiagonals;
605 invData(cellPointOrdinal,0) = 1.0 / data(cellPointOrdinal,0);
610 const auto & a = data(cellPointOrdinal,0);
611 const auto & b = data(cellPointOrdinal,1);
612 const auto & c = data(cellPointOrdinal,2);
613 const auto & d = data(cellPointOrdinal,3);
614 const PointScalar det = det2(a,b,c,d);
616 invData(cellPointOrdinal,0) = d/det;
617 invData(cellPointOrdinal,1) = - b/det;
618 invData(cellPointOrdinal,2) = - c/det;
619 invData(cellPointOrdinal,3) = a/det;
624 const auto & a = data(cellPointOrdinal,0);
625 const auto & b = data(cellPointOrdinal,1);
626 const auto & c = data(cellPointOrdinal,2);
627 const auto & d = data(cellPointOrdinal,3);
628 const auto & e = data(cellPointOrdinal,4);
629 const auto & f = data(cellPointOrdinal,5);
630 const auto & g = data(cellPointOrdinal,6);
631 const auto & h = data(cellPointOrdinal,7);
632 const auto & i = data(cellPointOrdinal,8);
633 const PointScalar det = det3(a,b,c,d,e,f,g,h,i);
636 const auto val0 = e*i - h*f;
637 const auto val1 = - d*i + g*f;
638 const auto val2 = d*h - g*e;
640 invData(cellPointOrdinal,0) = val0/det;
641 invData(cellPointOrdinal,1) = val1/det;
642 invData(cellPointOrdinal,2) = val2/det;
645 const auto val0 = h*c - b*i;
646 const auto val1 = a*i - g*c;
647 const auto val2 = - a*h + g*b;
649 invData(cellPointOrdinal,3) = val0/det;
650 invData(cellPointOrdinal,4) = val1/det;
651 invData(cellPointOrdinal,5) = val2/det;
654 const auto val0 = b*f - e*c;
655 const auto val1 = - a*f + d*c;
656 const auto val2 = a*e - d*b;
658 invData(cellPointOrdinal,6) = val0/det;
659 invData(cellPointOrdinal,7) = val1/det;
660 invData(cellPointOrdinal,8) = val2/det;
668 const int offset = blockWidth * blockWidth;
670 for (
int d=blockWidth; d<spaceDim; d++)
672 const int index = d-blockWidth+offset;
673 invData(cellPointOrdinal,index) = 1.0 / data(cellPointOrdinal,index);
682 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,1),
683 KOKKOS_LAMBDA (
const int &dummyArg) {
685 const int numDiagonals = data.extent_int(0) - blockWidth * blockWidth;
686 const int spaceDim = blockWidth + numDiagonals;
694 invData(0) = 1.0 / data(0);
699 const auto & a = data(0);
700 const auto & b = data(1);
701 const auto & c = data(2);
702 const auto & d = data(3);
703 const PointScalar det = det2(a,b,c,d);
706 invData(1) = - b/det;
707 invData(2) = - c/det;
713 const auto & a = data(0);
714 const auto & b = data(1);
715 const auto & c = data(2);
716 const auto & d = data(3);
717 const auto & e = data(4);
718 const auto & f = data(5);
719 const auto & g = data(6);
720 const auto & h = data(7);
721 const auto & i = data(8);
722 const PointScalar det = det3(a,b,c,d,e,f,g,h,i);
725 const auto val0 = e*i - h*f;
726 const auto val1 = - d*i + g*f;
727 const auto val2 = d*h - g*e;
729 invData(0) = val0/det;
730 invData(1) = val1/det;
731 invData(2) = val2/det;
734 const auto val0 = h*c - b*i;
735 const auto val1 = a*i - g*c;
736 const auto val2 = - a*h + g*b;
738 invData(3) = val0/det;
739 invData(4) = val1/det;
740 invData(5) = val2/det;
743 const auto val0 = b*f - e*c;
744 const auto val1 = - a*f + d*c;
745 const auto val2 = a*e - d*b;
747 invData(6) = val0/det;
748 invData(7) = val1/det;
749 invData(8) = val2/det;
757 const int offset = blockWidth * blockWidth;
759 for (
int d=blockWidth; d<spaceDim; d++)
761 const int index = d-blockWidth+offset;
762 invData(index) = 1.0 / data(index);
786 Kokkos::parallel_for(
"fill jacobian inverse", Kokkos::RangePolicy<ExecSpaceType>(0,1), KOKKOS_LAMBDA(
const int &i)
788 Intrepid2::Kernels::Serial::inverse(invData,data);
793 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"jacobian's underlying view must have rank 2,3, or 4");
797 template<
typename DeviceType>
798 template<
typename jacobianValueType,
class ...jacobianProperties,
799 typename BasisGradientsType,
800 typename WorksetType>
803 setJacobian( Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian,
804 const WorksetType worksetCell,
805 const BasisGradientsType gradients,
const int startCell,
const int endCell)
807 constexpr bool is_accessible =
808 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
809 typename decltype(jacobian)::memory_space>::accessible;
810 static_assert(is_accessible,
"CellTools<DeviceType>::setJacobian(..): output view's memory space is not compatible with DeviceType");
812 using JacobianViewType = Kokkos::DynRankView<jacobianValueType,jacobianProperties...>;
816 int endCellResolved = (endCell == -1) ? worksetCell.extent_int(0) : endCell;
818 using range_policy_type = Kokkos::MDRangePolicy
819 < ExecSpaceType, Kokkos::Rank<2>, Kokkos::IndexType<ordinal_type> >;
820 range_policy_type policy( { 0, 0 },
821 { jacobian.extent(0), jacobian.extent(1) } );
822 Kokkos::parallel_for( policy, FunctorType(jacobian, worksetCell, gradients, startCell, endCellResolved) );
825 template<
typename DeviceType>
826 template<
typename jacobianValueType,
class ...jacobianProperties,
827 typename pointValueType,
class ...pointProperties,
828 typename WorksetType,
829 typename HGradBasisType>
832 setJacobian( Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian,
833 const Kokkos::DynRankView<pointValueType,pointProperties...> points,
834 const WorksetType worksetCell,
835 const Teuchos::RCP<HGradBasisType> basis,
836 const int startCell,
const int endCell) {
837 constexpr bool are_accessible =
838 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
839 typename decltype(jacobian)::memory_space>::accessible &&
840 Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
841 typename decltype(points)::memory_space>::accessible;
842 static_assert(are_accessible,
"CellTools<DeviceType>::setJacobian(..): input/output views' memory spaces are not compatible with DeviceType");
844#ifdef HAVE_INTREPID2_DEBUG
845 CellTools_setJacobianArgs(jacobian, points, worksetCell, basis->getBaseCellTopology(), startCell, endCell);
848 const auto cellTopo = basis->getBaseCellTopology();
849 const ordinal_type spaceDim = cellTopo.getDimension();
850 const ordinal_type numCells = jacobian.extent(0);
853 const ordinal_type pointRank = points.rank();
854 const ordinal_type numPoints = (pointRank == 2 ? points.extent(0) : points.extent(1));
855 const ordinal_type basisCardinality = basis->getCardinality();
861 auto vcprop = Kokkos::common_view_alloc_prop(points);
862 using GradViewType = Kokkos::DynRankView<
typename decltype(vcprop)::value_type,DeviceType>;
869 grads = GradViewType(Kokkos::view_alloc(
"CellTools::setJacobian::grads", vcprop),basisCardinality, numPoints, spaceDim);
870 basis->getValues(grads,
877 grads = GradViewType(Kokkos::view_alloc(
"CellTools::setJacobian::grads", vcprop), numCells, basisCardinality, numPoints, spaceDim);
878 for (ordinal_type cell=0;cell<numCells;++cell)
879 basis->getValues(Kokkos::subview( grads, cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL() ),
880 Kokkos::subview( points, cell, Kokkos::ALL(), Kokkos::ALL() ),
886 setJacobian(jacobian, worksetCell, grads, startCell, endCell);
889 template<
typename DeviceType>
890 template<
typename jacobianInvValueType,
class ...jacobianInvProperties,
891 typename jacobianValueType,
class ...jacobianProperties>
894 setJacobianInv( Kokkos::DynRankView<jacobianInvValueType,jacobianInvProperties...> jacobianInv,
895 const Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian ) {
896#ifdef HAVE_INTREPID2_DEBUG
902 template<
typename DeviceType>
903 template<
typename jacobianDetValueType,
class ...jacobianDetProperties,
904 typename jacobianValueType,
class ...jacobianProperties>
907 setJacobianDet( Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
908 const Kokkos::DynRankView<jacobianValueType,jacobianProperties...> jacobian ) {
909#ifdef HAVE_INTREPID2_DEBUG
@ BLOCK_PLUS_DIAGONAL
one of two dimensions in a matrix; bottom-right part of matrix is diagonal
Header file for small functions used in Intrepid2.
#define INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(test, x, msg)
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.
Wrapper around a Kokkos::View that allows data that is constant or repeating in various logical dimen...
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 > & getUnderlyingView2() const
returns the View that stores the unique data. For rank-2 underlying containers.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar ****, DeviceType > & getUnderlyingView4() const
returns the View that stores the unique data. For rank-4 underlying containers.
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 int & blockPlusDiagonalLastNonDiagonal() const
For a Data object containing data with variation type BLOCK_PLUS_DIAGONAL, returns the row and column...
void storeInPlaceQuotient(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B)
stores the in-place (entrywise) quotient, A ./ B, into this container.
KOKKOS_INLINE_FUNCTION const Kokkos::Array< DataVariationType, 7 > & getVariationTypes() const
Returns an array with the variation types in each logical dimension.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar *, DeviceType > & getUnderlyingView1() const
returns the View that stores the unique data. For rank-1 underlying containers.
Data< DataScalar, DeviceType > allocateConstantData(const DataScalar &value)
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