8 #ifndef Intrepid2_Data_h 9 #define Intrepid2_Data_h 12 #include "Intrepid2_ScalarView.hpp" 48 int blockPlusDiagonalLastNonDiagonal = -1;
52 KOKKOS_INLINE_FUNCTION
55 const int myNominalExtent = myData.logicalExtent;
56 #ifdef HAVE_INTREPID2_DEBUG 57 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(myNominalExtent != otherData.logicalExtent, std::invalid_argument,
"both Data objects must match in their logical extent in the specified dimension");
62 const int & myVariationModulus = myData.variationModulus;
63 const int & otherVariationModulus = otherData.variationModulus;
65 int myDataExtent = myData.dataExtent;
66 int otherDataExtent = otherData.dataExtent;
74 switch (otherVariation)
84 switch (otherVariation)
90 if (myVariationModulus == otherVariationModulus)
93 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(myDataExtent != otherDataExtent, std::logic_error,
"Unexpected data extent/modulus combination");
120 switch (otherVariation)
139 combinedDimensionInfo.blockPlusDiagonalLastNonDiagonal = min(myData.blockPlusDiagonalLastNonDiagonal, otherData.blockPlusDiagonalLastNonDiagonal);
143 switch (otherVariation)
163 template<
class DataScalar,
typename DeviceType>
166 static ScalarView<DataScalar,DeviceType> zeroView()
168 static ScalarView<DataScalar,DeviceType> zeroView = ScalarView<DataScalar,DeviceType>(
"zero",1);
169 static bool havePushedFinalizeHook =
false;
170 if (!havePushedFinalizeHook)
172 Kokkos::push_finalize_hook( [=] {
173 zeroView = ScalarView<DataScalar,DeviceType>();
175 havePushedFinalizeHook =
true;
198 template<
class DataScalar,
typename DeviceType>
201 using value_type = DataScalar;
202 using execution_space =
typename DeviceType::execution_space;
204 ordinal_type dataRank_;
205 Kokkos::View<DataScalar*, DeviceType> data1_;
206 Kokkos::View<DataScalar**, DeviceType> data2_;
207 Kokkos::View<DataScalar***, DeviceType> data3_;
208 Kokkos::View<DataScalar****, DeviceType> data4_;
209 Kokkos::View<DataScalar*****, DeviceType> data5_;
210 Kokkos::View<DataScalar******, DeviceType> data6_;
211 Kokkos::View<DataScalar*******, DeviceType> data7_;
212 Kokkos::Array<int,7> extents_;
213 Kokkos::Array<DataVariationType,7> variationType_;
214 Kokkos::Array<int,7> variationModulus_;
215 int blockPlusDiagonalLastNonDiagonal_ = -1;
217 bool hasNontrivialModulusUNUSED_;
218 bool underlyingMatchesLogical_;
219 Kokkos::Array<ordinal_type,7> activeDims_;
224 using reference_type =
typename ScalarView<DataScalar,DeviceType>::reference_type;
225 using const_reference_type =
typename ScalarView<const DataScalar,DeviceType>::reference_type;
227 using return_type = const_reference_type;
229 ScalarView<DataScalar,DeviceType> zeroView_;
232 KOKKOS_INLINE_FUNCTION
235 return (lastNondiagonal + 1) * (lastNondiagonal + 1);
239 KOKKOS_INLINE_FUNCTION
242 return i * (lastNondiagonal + 1) + j;
246 KOKKOS_INLINE_FUNCTION
249 return i - (lastNondiagonal + 1) + numNondiagonalEntries;
253 KOKKOS_INLINE_FUNCTION
258 case 1:
return data1_.extent_int(dim);
259 case 2:
return data2_.extent_int(dim);
260 case 3:
return data3_.extent_int(dim);
261 case 4:
return data4_.extent_int(dim);
262 case 5:
return data5_.extent_int(dim);
263 case 6:
return data6_.extent_int(dim);
264 case 7:
return data7_.extent_int(dim);
274 for (
int d=rank_; d<7; d++)
276 INTREPID2_TEST_FOR_EXCEPTION(extents_[d] > 1, std::invalid_argument,
"Nominal extents may not be > 1 in dimensions beyond the rank of the container");
280 int blockPlusDiagonalCount = 0;
281 underlyingMatchesLogical_ =
true;
282 for (ordinal_type i=0; i<7; i++)
284 if (variationType_[i] ==
GENERAL)
286 if (extents_[i] != 0)
288 variationModulus_[i] = extents_[i];
292 variationModulus_[i] = 1;
294 activeDims_[numActiveDims_] = i;
297 else if (variationType_[i] ==
MODULAR)
299 underlyingMatchesLogical_ =
false;
303 const int logicalExtent = extents_[i];
304 const int modulus = dataExtent;
306 INTREPID2_TEST_FOR_EXCEPTION( dataExtent * (logicalExtent / dataExtent) != logicalExtent, std::invalid_argument,
"data extent must evenly divide logical extent");
308 variationModulus_[i] = modulus;
312 variationModulus_[i] = extents_[i];
314 activeDims_[numActiveDims_] = i;
319 underlyingMatchesLogical_ =
false;
320 blockPlusDiagonalCount++;
321 if (blockPlusDiagonalCount == 1)
324 #ifdef HAVE_INTREPID2_DEBUG 327 const int logicalExtent = extents_[i];
328 const int numDiagonalEntries = logicalExtent - (blockPlusDiagonalLastNonDiagonal_ + 1);
329 const int expectedDataExtent = numNondiagonalEntries + numDiagonalEntries;
330 INTREPID2_TEST_FOR_EXCEPTION(dataExtent != expectedDataExtent, std::invalid_argument, (
"BLOCK_PLUS_DIAGONAL data extent of " + std::to_string(dataExtent) +
" does not match expected based on blockPlusDiagonalLastNonDiagonal setting of " + std::to_string(blockPlusDiagonalLastNonDiagonal_)).c_str());
333 activeDims_[numActiveDims_] = i;
337 INTREPID2_TEST_FOR_EXCEPTION(variationType_[i+1] !=
BLOCK_PLUS_DIAGONAL, std::invalid_argument,
"BLOCK_PLUS_DIAGONAL ranks must be contiguous");
339 variationModulus_[i] = 1;
340 INTREPID2_TEST_FOR_EXCEPTION(blockPlusDiagonalCount > 1, std::invalid_argument,
"BLOCK_PLUS_DIAGONAL can only apply to two ranks");
346 underlyingMatchesLogical_ =
false;
348 variationModulus_[i] = 1;
352 if (rank_ != dataRank_)
354 underlyingMatchesLogical_ =
false;
357 for (
int d=numActiveDims_; d<7; d++)
363 for (
int d=0; d<7; d++)
365 INTREPID2_TEST_FOR_EXCEPTION(variationModulus_[d] == 0, std::logic_error,
"variationModulus should not ever be 0");
373 template<
class ViewType,
class ...IntArgs>
374 static KOKKOS_INLINE_FUNCTION reference_type
get(
const ViewType &view,
const IntArgs&... intArgs)
376 return view.getWritableEntry(intArgs...);
380 template<
class BinaryOperator,
class ThisUnderlyingViewType,
class AUnderlyingViewType,
class BUnderlyingViewType,
381 class ArgExtractorThis,
class ArgExtractorA,
class ArgExtractorB,
bool includeInnerLoop=
false>
385 ThisUnderlyingViewType this_underlying_;
386 AUnderlyingViewType A_underlying_;
387 BUnderlyingViewType B_underlying_;
388 BinaryOperator binaryOperator_;
391 InPlaceCombinationFunctor(ThisUnderlyingViewType this_underlying, AUnderlyingViewType A_underlying, BUnderlyingViewType B_underlying,
392 BinaryOperator binaryOperator)
394 this_underlying_(this_underlying),
395 A_underlying_(A_underlying),
396 B_underlying_(B_underlying),
397 binaryOperator_(binaryOperator)
399 INTREPID2_TEST_FOR_EXCEPTION(includeInnerLoop,std::invalid_argument,
"If includeInnerLoop is true, must specify the size of the inner loop");
402 InPlaceCombinationFunctor(ThisUnderlyingViewType this_underlying, AUnderlyingViewType A_underlying, BUnderlyingViewType B_underlying,
403 BinaryOperator binaryOperator,
int innerLoopSize)
405 this_underlying_(this_underlying),
406 A_underlying_(A_underlying),
407 B_underlying_(B_underlying),
408 binaryOperator_(binaryOperator),
409 innerLoopSize_(innerLoopSize)
411 INTREPID2_TEST_FOR_EXCEPTION(includeInnerLoop,std::invalid_argument,
"If includeInnerLoop is true, must specify the size of the inner loop");
414 template<
class ...IntArgs,
bool M=includeInnerLoop>
415 KOKKOS_INLINE_FUNCTION
416 enable_if_t<!M, void>
417 operator()(
const IntArgs&... args)
const 419 auto & result = ArgExtractorThis::get( this_underlying_, args... );
420 const auto & A_val = ArgExtractorA::get( A_underlying_, args... );
421 const auto & B_val = ArgExtractorB::get( B_underlying_, args... );
423 result = binaryOperator_(A_val,B_val);
426 template<
class ...IntArgs,
bool M=includeInnerLoop>
427 KOKKOS_INLINE_FUNCTION
429 operator()(
const IntArgs&... args)
const 431 using int_type = std::tuple_element_t<0, std::tuple<IntArgs...>>;
432 for (int_type iFinal=0; iFinal<static_cast<int_type>(innerLoopSize_); iFinal++)
434 auto & result = ArgExtractorThis::get( this_underlying_, args..., iFinal );
435 const auto & A_val = ArgExtractorA::get( A_underlying_, args..., iFinal );
436 const auto & B_val = ArgExtractorB::get( B_underlying_, args..., iFinal );
438 result = binaryOperator_(A_val,B_val);
444 template<
class BinaryOperator,
class PolicyType,
class ThisUnderlyingViewType,
class AUnderlyingViewType,
class BUnderlyingViewType,
445 class ArgExtractorThis,
class ArgExtractorA,
class ArgExtractorB>
447 AUnderlyingViewType &A_underlying, BUnderlyingViewType &B_underlying,
448 BinaryOperator &binaryOperator, ArgExtractorThis argThis, ArgExtractorA argA, ArgExtractorB argB)
450 using Functor = InPlaceCombinationFunctor<BinaryOperator, ThisUnderlyingViewType, AUnderlyingViewType, BUnderlyingViewType, ArgExtractorThis, ArgExtractorA, ArgExtractorB>;
451 Functor functor(this_underlying, A_underlying, B_underlying, binaryOperator);
452 Kokkos::parallel_for(
"compute in-place", policy, functor);
456 template<
class BinaryOperator,
int rank>
457 enable_if_t<rank != 7, void>
460 auto policy = dataExtentRangePolicy<rank>();
479 const FullArgExtractorWritableData fullArgsWritable;
492 const auto & variationTypes = data.getVariationTypes();
493 for (
int d=0; d<
rank; d++)
495 if (variationTypes[d] ==
GENERAL)
505 auto thisAE = constArg;
508 auto & this_underlying = this->getUnderlyingView<1>();
511 storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, BAE);
513 else if (this_full && A_full && B_full)
515 auto thisAE = fullArgs;
519 auto & this_underlying = this->getUnderlyingView<rank>();
523 storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, BAE);
531 auto thisAE = fullArgs;
532 auto & this_underlying = this->getUnderlyingView<rank>();
538 storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, BAE);
542 auto BAE = fullArgsConst;
549 if (B_1D && (get1DArgIndex(B) != -1) )
552 const int argIndex = get1DArgIndex(B);
554 auto & this_underlying = this->getUnderlyingView<1>();
557 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg0, AAE, arg0);
break;
558 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg1, AAE, arg1);
break;
559 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg2, AAE, arg2);
break;
560 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg3, AAE, arg3);
break;
561 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg4, AAE, arg4);
break;
562 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg5, AAE, arg5);
break;
563 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
569 auto thisAE = fullArgsWritable;
570 auto BAE = fullArgsConst;
581 auto thisAE = fullArgs;
582 auto & this_underlying = this->getUnderlyingView<rank>();
588 storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, BAE);
593 auto AAE = fullArgsConst;
600 if (A_1D && (get1DArgIndex(A) != -1) )
603 const int argIndex = get1DArgIndex(A);
605 auto & this_underlying = this->getUnderlyingView<1>();
608 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg0, arg0, BAE);
break;
609 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg1, arg1, BAE);
break;
610 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg2, arg2, BAE);
break;
611 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg3, arg3, BAE);
break;
612 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg4, arg4, BAE);
break;
613 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg5, arg5, BAE);
break;
614 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
620 auto thisAE = fullArgsWritable;
621 auto AAE = fullArgsConst;
628 if (this_1D && (get1DArgIndex(thisData) != -1))
635 const int argThis = get1DArgIndex(thisData);
636 const int argA = get1DArgIndex(A);
637 const int argB = get1DArgIndex(B);
641 auto & this_underlying = this->getUnderlyingView<1>();
642 if ((argA != -1) && (argB != -1))
644 #ifdef INTREPID2_HAVE_DEBUG 645 INTREPID2_TEST_FOR_EXCEPTION(argA != argThis, std::logic_error,
"Unexpected 1D arg combination.");
646 INTREPID2_TEST_FOR_EXCEPTION(argB != argThis, std::logic_error,
"Unexpected 1D arg combination.");
650 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg0, arg0, arg0);
break;
651 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg1, arg1, arg1);
break;
652 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg2, arg2, arg2);
break;
653 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg3, arg3, arg3);
break;
654 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg4, arg4, arg4);
break;
655 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, arg5, arg5, arg5);
break;
656 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
664 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, arg0, arg0, fullArgsConst);
break;
665 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, arg1, arg1, fullArgsConst);
break;
666 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, arg2, arg2, fullArgsConst);
break;
667 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, arg3, arg3, fullArgsConst);
break;
668 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, arg4, arg4, fullArgsConst);
break;
669 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, arg5, arg5, fullArgsConst);
break;
670 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
678 case 0:
storeInPlaceCombination(policy, this_underlying, A, B_underlying, binaryOperator, arg0, fullArgsConst, arg0);
break;
679 case 1:
storeInPlaceCombination(policy, this_underlying, A, B_underlying, binaryOperator, arg1, fullArgsConst, arg1);
break;
680 case 2:
storeInPlaceCombination(policy, this_underlying, A, B_underlying, binaryOperator, arg2, fullArgsConst, arg2);
break;
681 case 3:
storeInPlaceCombination(policy, this_underlying, A, B_underlying, binaryOperator, arg3, fullArgsConst, arg3);
break;
682 case 4:
storeInPlaceCombination(policy, this_underlying, A, B_underlying, binaryOperator, arg4, fullArgsConst, arg4);
break;
683 case 5:
storeInPlaceCombination(policy, this_underlying, A, B_underlying, binaryOperator, arg5, fullArgsConst, arg5);
break;
684 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
691 auto & this_underlying = this->getUnderlyingView<rank>();
692 auto thisAE = fullArgs;
699 if (B_1D && (get1DArgIndex(B) != -1))
701 const int argIndex = get1DArgIndex(B);
705 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, arg0);
break;
706 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, arg1);
break;
707 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, arg2);
break;
708 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, arg3);
break;
709 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, arg4);
break;
710 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, AAE, arg5);
break;
711 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
718 auto BAE = fullArgsConst;
724 if (A_1D && (get1DArgIndex(A) != -1))
726 const int argIndex = get1DArgIndex(A);
734 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, arg0, BAE);
break;
735 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, arg1, BAE);
break;
736 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, arg2, BAE);
break;
737 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, arg3, BAE);
break;
738 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, arg4, BAE);
break;
739 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B_underlying, binaryOperator, thisAE, arg5, BAE);
break;
740 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
745 auto BAE = fullArgsConst;
748 case 0:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, thisAE, arg0, BAE);
break;
749 case 1:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, thisAE, arg1, BAE);
break;
750 case 2:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, thisAE, arg2, BAE);
break;
751 case 3:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, thisAE, arg3, BAE);
break;
752 case 4:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, thisAE, arg4, BAE);
break;
753 case 5:
storeInPlaceCombination(policy, this_underlying, A_underlying, B, binaryOperator, thisAE, arg5, BAE);
break;
754 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid/unexpected arg index");
761 auto AAE = fullArgsConst;
762 auto BAE = fullArgsConst;
770 auto thisAE = fullArgsWritable;
771 auto AAE = fullArgsConst;
772 auto BAE = fullArgsConst;
779 template<
class BinaryOperator,
int rank>
780 enable_if_t<rank == 7, void>
783 auto policy = dataExtentRangePolicy<rank>();
786 using ThisAE = FullArgExtractorWritableData;
791 const bool includeInnerLoop =
true;
792 using Functor = InPlaceCombinationFunctor<BinaryOperator, DataType, DataType, DataType, ThisAE, AAE, BAE, includeInnerLoop>;
793 Functor functor(*
this, A, B, binaryOperator, dim6);
794 Kokkos::parallel_for(
"compute in-place", policy, functor);
798 template<
class UnaryOperator>
801 using ExecutionSpace =
typename DeviceType::execution_space;
807 const int dataRank = 1;
808 auto view = getUnderlyingView<dataRank>();
811 Kokkos::RangePolicy<ExecutionSpace> policy(ExecutionSpace(),0,dataExtent);
812 Kokkos::parallel_for(
"apply operator in-place", policy,
813 KOKKOS_LAMBDA (
const int &i0) {
814 view(i0) = unaryOperator(view(i0));
821 const int dataRank = 2;
822 auto policy = dataExtentRangePolicy<dataRank>();
823 auto view = getUnderlyingView<dataRank>();
825 Kokkos::parallel_for(
"apply operator in-place", policy,
826 KOKKOS_LAMBDA (
const int &i0,
const int &i1) {
827 view(i0,i1) = unaryOperator(view(i0,i1));
833 const int dataRank = 3;
834 auto policy = dataExtentRangePolicy<dataRank>();
835 auto view = getUnderlyingView<dataRank>();
837 Kokkos::parallel_for(
"apply operator in-place", policy,
838 KOKKOS_LAMBDA (
const int &i0,
const int &i1,
const int &i2) {
839 view(i0,i1,i2) = unaryOperator(view(i0,i1,i2));
845 const int dataRank = 4;
846 auto policy = dataExtentRangePolicy<dataRank>();
847 auto view = getUnderlyingView<dataRank>();
849 Kokkos::parallel_for(
"apply operator in-place", policy,
850 KOKKOS_LAMBDA (
const int &i0,
const int &i1,
const int &i2,
const int &i3) {
851 view(i0,i1,i2,i3) = unaryOperator(view(i0,i1,i2,i3));
857 const int dataRank = 5;
858 auto policy = dataExtentRangePolicy<dataRank>();
859 auto view = getUnderlyingView<dataRank>();
861 Kokkos::parallel_for(
"apply operator in-place", policy,
862 KOKKOS_LAMBDA (
const int &i0,
const int &i1,
const int &i2,
const int &i3,
const int &i4) {
863 view(i0,i1,i2,i3,i4) = unaryOperator(view(i0,i1,i2,i3,i4));
869 const int dataRank = 6;
870 auto policy = dataExtentRangePolicy<dataRank>();
871 auto view = getUnderlyingView<dataRank>();
873 Kokkos::parallel_for(
"apply operator in-place", policy,
874 KOKKOS_LAMBDA (
const int &i0,
const int &i1,
const int &i2,
const int &i3,
const int &i4,
const int &i5) {
875 view(i0,i1,i2,i3,i4,i5) = unaryOperator(view(i0,i1,i2,i3,i4,i5));
881 const int dataRank = 7;
882 auto policy6 = dataExtentRangePolicy<6>();
883 auto view = getUnderlyingView<dataRank>();
885 const int dim_i6 = view.extent_int(6);
887 Kokkos::parallel_for(
"apply operator in-place", policy6,
888 KOKKOS_LAMBDA (
const int &i0,
const int &i1,
const int &i2,
const int &i3,
const int &i4,
const int &i5) {
889 for (
int i6=0; i6<dim_i6; i6++)
891 view(i0,i1,i2,i3,i4,i5,i6) = unaryOperator(view(i0,i1,i2,i3,i4,i5,i6));
897 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Unsupported data rank");
902 template<
class ...IntArgs>
903 KOKKOS_INLINE_FUNCTION
906 #ifdef INTREPID2_HAVE_DEBUG 907 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(numArgs != rank_, std::invalid_argument,
"getWritableEntry() should have the same number of arguments as the logical rank.");
909 constexpr
int numArgs =
sizeof...(intArgs);
910 if (underlyingMatchesLogical_)
913 return getUnderlyingView<numArgs>()(intArgs...);
917 using int_type = std::tuple_element_t<0, std::tuple<IntArgs...>>;
919 const Kokkos::Array<int_type, numArgs> args {intArgs...};
920 Kokkos::Array<int_type, 7> refEntry;
921 for (
int d=0; d<numArgs; d++)
923 switch (variationType_[d])
925 case CONSTANT: refEntry[d] = 0;
break;
926 case GENERAL: refEntry[d] = args[d];
break;
927 case MODULAR: refEntry[d] = args[d] % variationModulus_[d];
break;
932 const int_type &i = args[d];
935 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"BLOCK_PLUS_DIAGONAL must be present for two dimensions; here, encountered only one");
939 const int_type &j = args[d+1];
941 if ((i > static_cast<int_type>(blockPlusDiagonalLastNonDiagonal_)) || (j > static_cast<int_type>(blockPlusDiagonalLastNonDiagonal_)))
966 for (
int d=numArgs; d<7; d++)
973 return data1_(refEntry[activeDims_[0]]);
975 else if (dataRank_ == 2)
977 return data2_(refEntry[activeDims_[0]],refEntry[activeDims_[1]]);
979 else if (dataRank_ == 3)
981 return data3_(refEntry[activeDims_[0]],refEntry[activeDims_[1]],refEntry[activeDims_[2]]);
983 else if (dataRank_ == 4)
985 return data4_(refEntry[activeDims_[0]],refEntry[activeDims_[1]],refEntry[activeDims_[2]],refEntry[activeDims_[3]]);
987 else if (dataRank_ == 5)
989 return data5_(refEntry[activeDims_[0]],refEntry[activeDims_[1]],refEntry[activeDims_[2]],refEntry[activeDims_[3]],
990 refEntry[activeDims_[4]]);
992 else if (dataRank_ == 6)
994 return data6_(refEntry[activeDims_[0]],refEntry[activeDims_[1]],refEntry[activeDims_[2]],refEntry[activeDims_[3]],
995 refEntry[activeDims_[4]],refEntry[activeDims_[5]]);
999 return data7_(refEntry[activeDims_[0]],refEntry[activeDims_[1]],refEntry[activeDims_[2]],refEntry[activeDims_[3]],
1000 refEntry[activeDims_[4]],refEntry[activeDims_[5]],refEntry[activeDims_[6]]);
1005 template<
class ToContainer,
class FromContainer>
1009 auto policy = Kokkos::MDRangePolicy<execution_space,Kokkos::Rank<6>>({0,0,0,0,0,0},{from.extent_int(0),from.extent_int(1),from.extent_int(2), from.extent_int(3), from.extent_int(4), from.extent_int(5)});
1011 Kokkos::parallel_for(
"copyContainer", policy,
1012 KOKKOS_LAMBDA (
const int &i0,
const int &i1,
const int &i2,
const int &i3,
const int &i4,
const int &i5) {
1013 for (
int i6=0; i6<from.extent_int(6); i6++)
1015 to.access(i0,i1,i2,i3,i4,i5,i6) = from.access(i0,i1,i2,i3,i4,i5,i6);
1026 case 1: data1_ = Kokkos::View<DataScalar*, DeviceType>(
"Intrepid2 Data", data.extent_int(0));
break;
1027 case 2: data2_ = Kokkos::View<DataScalar**, DeviceType>(
"Intrepid2 Data", data.extent_int(0), data.extent_int(1));
break;
1028 case 3: data3_ = Kokkos::View<DataScalar***, DeviceType>(
"Intrepid2 Data", data.extent_int(0), data.extent_int(1), data.extent_int(2));
break;
1029 case 4: data4_ = Kokkos::View<DataScalar****, DeviceType>(
"Intrepid2 Data", data.extent_int(0), data.extent_int(1), data.extent_int(2), data.extent_int(3));
break;
1030 case 5: data5_ = Kokkos::View<DataScalar*****, DeviceType>(
"Intrepid2 Data", data.extent_int(0), data.extent_int(1), data.extent_int(2), data.extent_int(3), data.extent_int(4));
break;
1031 case 6: data6_ = Kokkos::View<DataScalar******, DeviceType>(
"Intrepid2 Data", data.extent_int(0), data.extent_int(1), data.extent_int(2), data.extent_int(3), data.extent_int(4), data.extent_int(5));
break;
1032 case 7: data7_ = Kokkos::View<DataScalar*******, DeviceType>(
"Intrepid2 Data", data.extent_int(0), data.extent_int(1), data.extent_int(2), data.extent_int(3), data.extent_int(4), data.extent_int(5), data.extent_int(6));
break;
1033 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1045 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1050 Data(std::vector<DimensionInfo> dimInfoVector)
1053 dataRank_(0), extents_({0,0,0,0,0,0,0}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(-1), rank_(dimInfoVector.size())
1057 if (dimInfoVector.size() != 0)
1059 std::vector<int> dataExtents;
1061 bool blockPlusDiagonalEncountered =
true;
1062 for (
int d=0; d<rank_; d++)
1064 const DimensionInfo & dimInfo = dimInfoVector[d];
1065 extents_[d] = dimInfo.logicalExtent;
1066 variationType_[d] = dimInfo.variationType;
1067 const bool isBlockPlusDiagonal = (variationType_[d] == BLOCK_PLUS_DIAGONAL);
1068 const bool isSecondBlockPlusDiagonal = isBlockPlusDiagonal && blockPlusDiagonalEncountered;
1069 if (isBlockPlusDiagonal)
1071 blockPlusDiagonalEncountered =
true;
1072 blockPlusDiagonalLastNonDiagonal_ = dimInfo.blockPlusDiagonalLastNonDiagonal;
1074 if ((variationType_[d] != CONSTANT) && (!isSecondBlockPlusDiagonal))
1076 dataExtents.push_back(dimInfo.dataExtent);
1079 if (dataExtents.size() == 0)
1082 dataExtents.push_back(1);
1084 dataRank_ = dataExtents.size();
1087 case 1: data1_ = Kokkos::View<DataScalar*, DeviceType>(
"Intrepid2 Data", dataExtents[0]);
break;
1088 case 2: data2_ = Kokkos::View<DataScalar**, DeviceType>(
"Intrepid2 Data", dataExtents[0], dataExtents[1]);
break;
1089 case 3: data3_ = Kokkos::View<DataScalar***, DeviceType>(
"Intrepid2 Data", dataExtents[0], dataExtents[1], dataExtents[2]);
break;
1090 case 4: data4_ = Kokkos::View<DataScalar****, DeviceType>(
"Intrepid2 Data", dataExtents[0], dataExtents[1], dataExtents[2], dataExtents[3]);
break;
1091 case 5: data5_ = Kokkos::View<DataScalar*****, DeviceType>(
"Intrepid2 Data", dataExtents[0], dataExtents[1], dataExtents[2], dataExtents[3], dataExtents[4]);
break;
1092 case 6: data6_ = Kokkos::View<DataScalar******, DeviceType>(
"Intrepid2 Data", dataExtents[0], dataExtents[1], dataExtents[2], dataExtents[3], dataExtents[4], dataExtents[5]);
break;
1093 case 7: data7_ = Kokkos::View<DataScalar*******, DeviceType>(
"Intrepid2 Data", dataExtents[0], dataExtents[1], dataExtents[2], dataExtents[3], dataExtents[4], dataExtents[5], dataExtents[6]);
break;
1094 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1110 template<typename OtherDeviceType, class = typename std::enable_if< std::is_same<typename DeviceType::memory_space, typename OtherDeviceType::memory_space>::value>::type,
1111 class = typename std::enable_if<!std::is_same<DeviceType,OtherDeviceType>::value>::type>
1112 Data(const Data<DataScalar,OtherDeviceType> &data)
1114 dataRank_(data.getUnderlyingViewRank()), extents_(data.getExtents()), variationType_(data.getVariationTypes()), blockPlusDiagonalLastNonDiagonal_(data.blockPlusDiagonalLastNonDiagonal()), rank_(data.rank())
1122 case 1: data1_ = data.getUnderlyingView1();
break;
1123 case 2: data2_ = data.getUnderlyingView2();
break;
1124 case 3: data3_ = data.getUnderlyingView3();
break;
1125 case 4: data4_ = data.getUnderlyingView4();
break;
1126 case 5: data5_ = data.getUnderlyingView5();
break;
1127 case 6: data6_ = data.getUnderlyingView6();
break;
1128 case 7: data7_ = data.getUnderlyingView7();
break;
1129 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1136 template<typename OtherDeviceType, class = typename std::enable_if<!std::is_same<typename DeviceType::memory_space, typename OtherDeviceType::memory_space>::value>::type>
1147 case 1: data1_ = Kokkos::View<DataScalar*, DeviceType>(
"Intrepid2 Data", view.extent_int(0));
break;
1148 case 2: data2_ = Kokkos::View<DataScalar**, DeviceType>(
"Intrepid2 Data", view.extent_int(0), view.extent_int(1));
break;
1149 case 3: data3_ = Kokkos::View<DataScalar***, DeviceType>(
"Intrepid2 Data", view.extent_int(0), view.extent_int(1), view.extent_int(2));
break;
1150 case 4: data4_ = Kokkos::View<DataScalar****, DeviceType>(
"Intrepid2 Data", view.extent_int(0), view.extent_int(1), view.extent_int(2), view.extent_int(3));
break;
1151 case 5: data5_ = Kokkos::View<DataScalar*****, DeviceType>(
"Intrepid2 Data", view.extent_int(0), view.extent_int(1), view.extent_int(2), view.extent_int(3), view.extent_int(4));
break;
1152 case 6: data6_ = Kokkos::View<DataScalar******, DeviceType>(
"Intrepid2 Data", view.extent_int(0), view.extent_int(1), view.extent_int(2), view.extent_int(3), view.extent_int(4), view.extent_int(5));
break;
1153 case 7: data7_ = Kokkos::View<DataScalar*******, DeviceType>(
"Intrepid2 Data", view.extent_int(0), view.extent_int(1), view.extent_int(2), view.extent_int(3), view.extent_int(4), view.extent_int(5), view.extent_int(6));
break;
1154 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1160 using MemorySpace =
typename DeviceType::memory_space;
1163 case 1: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView1());
copyContainer(data1_, dataViewMirror);}
break;
1164 case 2: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView2());
copyContainer(data2_, dataViewMirror);}
break;
1165 case 3: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView3());
copyContainer(data3_, dataViewMirror);}
break;
1166 case 4: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView4());
copyContainer(data4_, dataViewMirror);}
break;
1167 case 5: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView5());
copyContainer(data5_, dataViewMirror);}
break;
1168 case 6: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView6());
copyContainer(data6_, dataViewMirror);}
break;
1169 case 7: {
auto dataViewMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), data.
getUnderlyingView7());
copyContainer(data7_, dataViewMirror);}
break;
1170 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1218 Data(ScalarView<DataScalar,DeviceType> data)
1222 Kokkos::Array<int,7> {data.extent_int(0),data.extent_int(1),data.extent_int(2),data.extent_int(3),data.extent_int(4),data.extent_int(5),data.extent_int(6)},
1227 template<
size_t rank,
class ...DynRankViewProperties>
1228 Data(
const Kokkos::DynRankView<DataScalar,DeviceType, DynRankViewProperties...> &data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1230 dataRank_(data.
rank()), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1234 for (
unsigned d=0; d<
rank; d++)
1236 extents_[d] = extents[d];
1237 variationType_[d] = variationType[d];
1242 template<
size_t rank,
class ...ViewProperties>
1243 Data(Kokkos::View<DataScalar*,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1245 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1248 for (
unsigned d=0; d<
rank; d++)
1250 extents_[d] = extents[d];
1251 variationType_[d] = variationType[d];
1256 template<
size_t rank,
class ...ViewProperties>
1257 Data(Kokkos::View<DataScalar**,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1259 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1262 for (
unsigned d=0; d<
rank; d++)
1264 extents_[d] = extents[d];
1265 variationType_[d] = variationType[d];
1270 template<
size_t rank,
class ...ViewProperties>
1271 Data(Kokkos::View<DataScalar***,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1273 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1276 for (
unsigned d=0; d<
rank; d++)
1278 extents_[d] = extents[d];
1279 variationType_[d] = variationType[d];
1284 template<
size_t rank,
class ...ViewProperties>
1285 Data(Kokkos::View<DataScalar****,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1287 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1290 for (
unsigned d=0; d<
rank; d++)
1292 extents_[d] = extents[d];
1293 variationType_[d] = variationType[d];
1298 template<
size_t rank,
class ...ViewProperties>
1299 Data(Kokkos::View<DataScalar*****,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1301 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1304 for (
unsigned d=0; d<
rank; d++)
1306 extents_[d] = extents[d];
1307 variationType_[d] = variationType[d];
1312 template<
size_t rank,
class ...ViewProperties>
1313 Data(Kokkos::View<DataScalar******,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1315 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1318 for (
unsigned d=0; d<
rank; d++)
1320 extents_[d] = extents[d];
1321 variationType_[d] = variationType[d];
1326 template<
size_t rank,
class ...ViewProperties>
1327 Data(Kokkos::View<DataScalar*******,DeviceType, ViewProperties...> data, Kokkos::Array<int,rank> extents, Kokkos::Array<DataVariationType,rank> variationType,
const int blockPlusDiagonalLastNonDiagonal = -1)
1329 dataRank_(data.
rank), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(
blockPlusDiagonalLastNonDiagonal), rank_(
rank)
1332 for (
unsigned d=0; d<
rank; d++)
1334 extents_[d] = extents[d];
1335 variationType_[d] = variationType[d];
1341 template<
size_t rank>
1342 Data(DataScalar constantValue, Kokkos::Array<int,rank> extents)
1344 dataRank_(1), extents_({1,1,1,1,1,1,1}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(-1), rank_(
rank)
1346 data1_ = Kokkos::View<DataScalar*,DeviceType>(
"Constant Data",1);
1347 Kokkos::deep_copy(data1_, constantValue);
1348 for (
unsigned d=0; d<
rank; d++)
1350 extents_[d] = extents[d];
1358 dataRank_(0), extents_({0,0,0,0,0,0,0}), variationType_({
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT,
CONSTANT}), blockPlusDiagonalLastNonDiagonal_(-1), rank_(0)
1364 KOKKOS_INLINE_FUNCTION
1367 return blockPlusDiagonalLastNonDiagonal_;
1371 KOKKOS_INLINE_FUNCTION
1378 KOKKOS_INLINE_FUNCTION
1384 dimInfo.variationType = variationType_[dim];
1386 dimInfo.variationModulus = variationModulus_[dim];
1390 dimInfo.blockPlusDiagonalLastNonDiagonal = blockPlusDiagonalLastNonDiagonal_;
1396 KOKKOS_INLINE_FUNCTION
1407 KOKKOS_INLINE_FUNCTION
1408 enable_if_t<rank==1, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1411 #ifdef HAVE_INTREPID2_DEBUG 1412 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1419 KOKKOS_INLINE_FUNCTION
1420 enable_if_t<rank==2, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1423 #ifdef HAVE_INTREPID2_DEBUG 1424 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1431 KOKKOS_INLINE_FUNCTION
1432 enable_if_t<rank==3, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1435 #ifdef HAVE_INTREPID2_DEBUG 1436 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1443 KOKKOS_INLINE_FUNCTION
1444 enable_if_t<rank==4, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1447 #ifdef HAVE_INTREPID2_DEBUG 1448 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1455 KOKKOS_INLINE_FUNCTION
1456 enable_if_t<rank==5, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1459 #ifdef HAVE_INTREPID2_DEBUG 1460 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1467 KOKKOS_INLINE_FUNCTION
1468 enable_if_t<rank==6, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1471 #ifdef HAVE_INTREPID2_DEBUG 1472 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1479 KOKKOS_INLINE_FUNCTION
1480 enable_if_t<rank==7, const Kokkos::View<typename RankExpander<DataScalar, rank>::value_type, DeviceType> &>
1483 #ifdef HAVE_INTREPID2_DEBUG 1484 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(dataRank_ !=
rank, std::invalid_argument,
"getUnderlyingView() called for rank that does not match dataRank_");
1490 KOKKOS_INLINE_FUNCTION
1493 return getUnderlyingView<1>();
1497 KOKKOS_INLINE_FUNCTION
1500 return getUnderlyingView<2>();
1504 KOKKOS_INLINE_FUNCTION
1507 return getUnderlyingView<3>();
1511 KOKKOS_INLINE_FUNCTION
1514 return getUnderlyingView<4>();
1518 KOKKOS_INLINE_FUNCTION
1521 return getUnderlyingView<5>();
1525 KOKKOS_INLINE_FUNCTION
1528 return getUnderlyingView<6>();
1532 KOKKOS_INLINE_FUNCTION
1535 return getUnderlyingView<7>();
1539 KOKKOS_INLINE_FUNCTION
1546 KOKKOS_INLINE_FUNCTION
1553 KOKKOS_INLINE_FUNCTION
1560 KOKKOS_INLINE_FUNCTION
1567 KOKKOS_INLINE_FUNCTION
1574 KOKKOS_INLINE_FUNCTION
1581 KOKKOS_INLINE_FUNCTION
1592 case 1:
return data1_;
1593 case 2:
return data2_;
1594 case 3:
return data3_;
1595 case 4:
return data4_;
1596 case 5:
return data5_;
1597 case 6:
return data6_;
1598 case 7:
return data7_;
1599 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1604 KOKKOS_INLINE_FUNCTION
1611 KOKKOS_INLINE_FUNCTION
1614 ordinal_type size = 1;
1615 for (ordinal_type r=0; r<dataRank_; r++)
1629 case 3:
return getMatchingViewWithLabel(data3_,
"Intrepid2 Data", data3_.extent_int(0), data3_.extent_int(1), data3_.extent_int(2));
1630 case 4:
return getMatchingViewWithLabel(data4_,
"Intrepid2 Data", data4_.extent_int(0), data4_.extent_int(1), data4_.extent_int(2), data4_.extent_int(3));
1631 case 5:
return getMatchingViewWithLabel(data5_,
"Intrepid2 Data", data5_.extent_int(0), data5_.extent_int(1), data5_.extent_int(2), data5_.extent_int(3), data5_.extent_int(4));
1632 case 6:
return getMatchingViewWithLabel(data6_,
"Intrepid2 Data", data6_.extent_int(0), data6_.extent_int(1), data6_.extent_int(2), data6_.extent_int(3), data6_.extent_int(4), data6_.extent_int(5));
1633 case 7:
return getMatchingViewWithLabel(data7_,
"Intrepid2 Data", data7_.extent_int(0), data7_.extent_int(1), data7_.extent_int(2), data7_.extent_int(3), data7_.extent_int(4), data7_.extent_int(5), data7_.extent_int(6));
1634 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1639 template<
class ... DimArgs>
1651 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1660 case 1: Kokkos::deep_copy(data1_, 0.0);
break;
1661 case 2: Kokkos::deep_copy(data2_, 0.0);
break;
1662 case 3: Kokkos::deep_copy(data3_, 0.0);
break;
1663 case 4: Kokkos::deep_copy(data4_, 0.0);
break;
1664 case 5: Kokkos::deep_copy(data5_, 0.0);
break;
1665 case 6: Kokkos::deep_copy(data6_, 0.0);
break;
1666 case 7: Kokkos::deep_copy(data7_, 0.0);
break;
1667 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1684 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid data rank");
1691 for (
unsigned i=0; i<activeDims_.size(); i++)
1693 if (activeDims_[i] == d)
1697 else if (activeDims_[i] > d)
1716 KOKKOS_INLINE_FUNCTION
1719 return variationModulus_[d];
1723 KOKKOS_INLINE_FUNCTION
1726 return variationType_;
1730 template<
class ...IntArgs>
1731 KOKKOS_INLINE_FUNCTION
1739 template <
bool... v>
1742 template <
class ...IntArgs>
1743 using valid_args = all_true<std::is_integral<IntArgs>{}...>;
1745 static_assert(valid_args<int,long,unsigned>::value,
"valid args works");
1748 template <
class ...IntArgs>
1749 KOKKOS_INLINE_FUNCTION
1750 #ifndef __INTEL_COMPILER 1754 enable_if_t<valid_args<IntArgs...>::value && (
sizeof...(IntArgs) <= 7),return_type>
1763 KOKKOS_INLINE_FUNCTION
1769 template <
typename iType>
1770 KOKKOS_INLINE_FUNCTION constexpr
1771 typename std::enable_if<std::is_integral<iType>::value,
size_t>::type
1772 extent(
const iType& r)
const {
1779 if (blockPlusDiagonalLastNonDiagonal_ >= 1)
return false;
1780 int numBlockPlusDiagonalTypes = 0;
1781 for (
unsigned r = 0; r<variationType_.size(); r++)
1783 const auto &entryType = variationType_[r];
1787 if (numBlockPlusDiagonalTypes == 2)
return true;
1788 else if (numBlockPlusDiagonalTypes == 0)
return false;
1789 else INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::invalid_argument,
"Unexpected number of ranks marked as BLOCK_PLUS_DIAGONAL (should be 0 or 2)");
1800 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A.
rank() != B.
rank(), std::invalid_argument,
"A and B must have the same logical shape");
1802 std::vector<DimensionInfo> dimInfo(
rank);
1803 for (
int d=0; d<
rank; d++)
1805 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A.
extent_int(d) != B.
extent_int(d), std::invalid_argument,
"A and B must have the same logical shape");
1821 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A_MatData.
rank() != B_MatData.
rank(), std::invalid_argument,
"AmatData and BmatData have incompatible ranks");
1823 const int D1_DIM = A_MatData.
rank() - 2;
1824 const int D2_DIM = A_MatData.
rank() - 1;
1826 const int A_rows = A_MatData.
extent_int(D1_DIM);
1827 const int A_cols = A_MatData.
extent_int(D2_DIM);
1828 const int B_rows = B_MatData.
extent_int(D1_DIM);
1829 const int B_cols = B_MatData.
extent_int(D2_DIM);
1831 const int leftRows = transposeA ? A_cols : A_rows;
1832 const int leftCols = transposeA ? A_rows : A_cols;
1833 const int rightRows = transposeB ? B_cols : B_rows;
1834 const int rightCols = transposeB ? B_rows : B_cols;
1836 INTREPID2_TEST_FOR_EXCEPTION(leftCols != rightRows, std::invalid_argument,
"incompatible matrix dimensions");
1838 Kokkos::Array<int,7> resultExtents;
1839 Kokkos::Array<DataVariationType,7> resultVariationTypes;
1841 resultExtents[D1_DIM] = leftRows;
1842 resultExtents[D2_DIM] = rightCols;
1843 int resultBlockPlusDiagonalLastNonDiagonal = -1;
1852 const int resultRank = A_MatData.
rank();
1857 Kokkos::Array<int,7> resultActiveDims;
1858 Kokkos::Array<int,7> resultDataDims;
1859 int resultNumActiveDims = 0;
1861 for (
int i=0; i<resultRank-2; i++)
1863 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A_MatData.
extent_int(i) != B_MatData.
extent_int(i), std::invalid_argument,
"A and B extents must match in each non-matrix dimension");
1871 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A_VariationType ==
BLOCK_PLUS_DIAGONAL, std::invalid_argument,
"unsupported variationType");
1872 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(B_VariationType ==
BLOCK_PLUS_DIAGONAL, std::invalid_argument,
"unsupported variationType");
1876 if ((A_VariationType ==
GENERAL) || (B_VariationType ==
GENERAL))
1878 resultVariationType =
GENERAL;
1879 dataSize = resultExtents[i];
1886 else if ((B_VariationType ==
MODULAR) && (A_VariationType ==
CONSTANT))
1888 resultVariationType =
MODULAR;
1891 else if ((B_VariationType ==
CONSTANT) && (A_VariationType ==
MODULAR))
1893 resultVariationType =
MODULAR;
1901 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A_Modulus != B_Modulus, std::invalid_argument,
"If both matrices have variation type MODULAR, they must agree on the modulus");
1902 resultVariationType =
MODULAR;
1903 dataSize = A_Modulus;
1905 resultVariationTypes[i] = resultVariationType;
1907 if (resultVariationType !=
CONSTANT)
1909 resultActiveDims[resultNumActiveDims] = i;
1910 resultDataDims[resultNumActiveDims] = dataSize;
1911 resultNumActiveDims++;
1916 resultExtents[D1_DIM] = leftRows;
1917 resultExtents[D2_DIM] = rightCols;
1926 resultActiveDims[resultNumActiveDims] = resultRank - 2;
1928 const int numDiagonalEntries = leftRows - (resultBlockPlusDiagonalLastNonDiagonal + 1);
1929 const int numNondiagonalEntries = (resultBlockPlusDiagonalLastNonDiagonal + 1) * (resultBlockPlusDiagonalLastNonDiagonal + 1);
1931 resultDataDims[resultNumActiveDims] = numDiagonalEntries + numNondiagonalEntries;
1932 resultNumActiveDims++;
1937 resultVariationTypes[D1_DIM] =
GENERAL;
1938 resultVariationTypes[D2_DIM] =
GENERAL;
1940 resultActiveDims[resultNumActiveDims] = resultRank - 2;
1941 resultActiveDims[resultNumActiveDims+1] = resultRank - 1;
1943 resultDataDims[resultNumActiveDims] = leftRows;
1944 resultDataDims[resultNumActiveDims+1] = rightCols;
1945 resultNumActiveDims += 2;
1948 for (
int i=resultRank; i<7; i++)
1950 resultVariationTypes[i] =
CONSTANT;
1951 resultExtents[i] = 1;
1954 ScalarView<DataScalar,DeviceType> data;
1955 if (resultNumActiveDims == 1)
1960 else if (resultNumActiveDims == 2)
1965 else if (resultNumActiveDims == 3)
1968 data =
getMatchingViewWithLabel(viewToMatch,
"Data mat-mat result", resultDataDims[0], resultDataDims[1], resultDataDims[2]);
1970 else if (resultNumActiveDims == 4)
1973 data =
getMatchingViewWithLabel(viewToMatch,
"Data mat-mat result", resultDataDims[0], resultDataDims[1], resultDataDims[2],
1976 else if (resultNumActiveDims == 5)
1979 data =
getMatchingViewWithLabel(viewToMatch,
"Data mat-mat result", resultDataDims[0], resultDataDims[1], resultDataDims[2],
1980 resultDataDims[3], resultDataDims[4]);
1982 else if (resultNumActiveDims == 6)
1985 data =
getMatchingViewWithLabel(viewToMatch,
"Data mat-mat result", resultDataDims[0], resultDataDims[1], resultDataDims[2],
1986 resultDataDims[3], resultDataDims[4], resultDataDims[5]);
1991 data =
getMatchingViewWithLabel(viewToMatch,
"Data mat-mat result", resultDataDims[0], resultDataDims[1], resultDataDims[2],
1992 resultDataDims[3], resultDataDims[4], resultDataDims[5], resultDataDims[6]);
2003 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(matData.
rank() != vecData.
rank() + 1, std::invalid_argument,
"matData and vecData have incompatible ranks");
2008 const int resultRank = vecData.
rank();
2010 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(matCols != vecDim, std::invalid_argument,
"matData column count != vecData dimension");
2012 Kokkos::Array<int,7> resultExtents;
2013 Kokkos::Array<DataVariationType,7> resultVariationTypes;
2017 Kokkos::Array<int,7> resultActiveDims;
2018 Kokkos::Array<int,7> resultDataDims;
2019 int resultNumActiveDims = 0;
2021 for (
int i=0; i<resultRank-1; i++)
2024 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(vecData.
extent_int(i) != matData.
extent_int(i), std::invalid_argument,
"matData and vecData extents must match in each non-matrix/vector dimension");
2030 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(vecVariationType ==
BLOCK_PLUS_DIAGONAL, std::invalid_argument,
"unsupported variationType");
2031 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(matVariationType ==
BLOCK_PLUS_DIAGONAL, std::invalid_argument,
"unsupported variationType");
2035 if ((vecVariationType ==
GENERAL) || (matVariationType ==
GENERAL))
2037 resultVariationType =
GENERAL;
2038 dataSize = resultExtents[i];
2045 else if ((matVariationType ==
MODULAR) && (vecVariationType ==
CONSTANT))
2047 resultVariationType =
MODULAR;
2050 else if ((matVariationType ==
CONSTANT) && (vecVariationType ==
MODULAR))
2052 resultVariationType =
MODULAR;
2060 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(matModulus != vecModulus, std::invalid_argument,
"If matrix and vector both have variation type MODULAR, they must agree on the modulus");
2061 resultVariationType =
MODULAR;
2062 dataSize = matModulus;
2064 resultVariationTypes[i] = resultVariationType;
2066 if (resultVariationType !=
CONSTANT)
2068 resultActiveDims[resultNumActiveDims] = i;
2069 resultDataDims[resultNumActiveDims] = dataSize;
2070 resultNumActiveDims++;
2075 resultVariationTypes[resultNumActiveDims] =
GENERAL;
2076 resultActiveDims[resultNumActiveDims] = resultRank - 1;
2077 resultDataDims[resultNumActiveDims] = matRows;
2078 resultExtents[resultRank-1] = matRows;
2079 resultNumActiveDims++;
2081 for (
int i=resultRank; i<7; i++)
2083 resultVariationTypes[i] =
CONSTANT;
2084 resultExtents[i] = 1;
2087 ScalarView<DataScalar,DeviceType> data;
2088 if (resultNumActiveDims == 1)
2092 else if (resultNumActiveDims == 2)
2096 else if (resultNumActiveDims == 3)
2100 else if (resultNumActiveDims == 4)
2105 else if (resultNumActiveDims == 5)
2108 resultDataDims[3], resultDataDims[4]);
2110 else if (resultNumActiveDims == 6)
2113 resultDataDims[3], resultDataDims[4], resultDataDims[5]);
2118 resultDataDims[3], resultDataDims[4], resultDataDims[5], resultDataDims[6]);
2126 enable_if_t<(rank!=1) && (rank!=7), Kokkos::MDRangePolicy<typename DeviceType::execution_space,Kokkos::Rank<rank>> >
2129 using ExecutionSpace =
typename DeviceType::execution_space;
2130 Kokkos::Array<int,rank> startingOrdinals;
2131 Kokkos::Array<int,rank> extents;
2133 for (
int d=0; d<
rank; d++)
2135 startingOrdinals[d] = 0;
2138 auto policy = Kokkos::MDRangePolicy<ExecutionSpace,Kokkos::Rank<rank>>(startingOrdinals,extents);
2144 enable_if_t<rank==7, Kokkos::MDRangePolicy<typename DeviceType::execution_space,Kokkos::Rank<6>> >
2147 using ExecutionSpace =
typename DeviceType::execution_space;
2148 Kokkos::Array<int,6> startingOrdinals;
2149 Kokkos::Array<int,6> extents;
2151 for (
int d=0; d<6; d++)
2153 startingOrdinals[d] = 0;
2156 auto policy = Kokkos::MDRangePolicy<ExecutionSpace,Kokkos::Rank<6>>(startingOrdinals,extents);
2162 enable_if_t<rank==1, Kokkos::RangePolicy<typename DeviceType::execution_space> >
2165 using ExecutionSpace =
typename DeviceType::execution_space;
2166 Kokkos::RangePolicy<ExecutionSpace> policy(ExecutionSpace(),0,
getDataExtent(0));
2171 template<
class BinaryOperator>
2174 using ExecutionSpace =
typename DeviceType::execution_space;
2176 #ifdef INTREPID2_HAVE_DEBUG 2178 for (
int d=0; d<rank_; d++)
2180 INTREPID2_TEST_FOR_EXCEPTION(A.
extent_int(d) != this->
extent_int(d), std::invalid_argument,
"A, B, and this must agree on all logical extents");
2181 INTREPID2_TEST_FOR_EXCEPTION(B.
extent_int(d) != this->
extent_int(d), std::invalid_argument,
"A, B, and this must agree on all logical extents");
2192 Kokkos::RangePolicy<ExecutionSpace> policy(ExecutionSpace(),0,1);
2193 auto this_underlying = this->getUnderlyingView<1>();
2196 Kokkos::parallel_for(
"compute in-place", policy,
2197 KOKKOS_LAMBDA (
const int &i0) {
2198 auto & result = this_underlying(0);
2199 const auto & A_val = A_underlying(0);
2200 const auto & B_val = B_underlying(0);
2202 result = binaryOperator(A_val,B_val);
2209 case 1: storeInPlaceCombination<BinaryOperator, 1>(A, B, binaryOperator);
break;
2210 case 2: storeInPlaceCombination<BinaryOperator, 2>(A, B, binaryOperator);
break;
2211 case 3: storeInPlaceCombination<BinaryOperator, 3>(A, B, binaryOperator);
break;
2212 case 4: storeInPlaceCombination<BinaryOperator, 4>(A, B, binaryOperator);
break;
2213 case 5: storeInPlaceCombination<BinaryOperator, 5>(A, B, binaryOperator);
break;
2214 case 6: storeInPlaceCombination<BinaryOperator, 6>(A, B, binaryOperator);
break;
2215 case 7: storeInPlaceCombination<BinaryOperator, 7>(A, B, binaryOperator);
break;
2217 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::logic_error,
"unhandled rank in switch");
2225 auto sum = KOKKOS_LAMBDA(
const DataScalar &a,
const DataScalar &b) -> DataScalar
2235 auto product = KOKKOS_LAMBDA(
const DataScalar &a,
const DataScalar &b) -> DataScalar
2245 auto difference = KOKKOS_LAMBDA(
const DataScalar &a,
const DataScalar &b) -> DataScalar
2255 auto quotient = KOKKOS_LAMBDA(
const DataScalar &a,
const DataScalar &b) -> DataScalar
2274 using ExecutionSpace =
typename DeviceType::execution_space;
2280 Kokkos::parallel_for(
"compute mat-vec", policy,
2281 KOKKOS_LAMBDA (
const int &cellOrdinal,
const int &pointOrdinal,
const int &i) {
2284 for (
int j=0; j<matCols; j++)
2286 val_i += matData(cellOrdinal,pointOrdinal,i,j) * vecData(cellOrdinal,pointOrdinal,j);
2290 else if (rank_ == 2)
2293 auto policy = Kokkos::MDRangePolicy<ExecutionSpace,Kokkos::Rank<2>>({0,0},{
getDataExtent(0),matRows});
2294 Kokkos::parallel_for(
"compute mat-vec", policy,
2295 KOKKOS_LAMBDA (
const int &vectorOrdinal,
const int &i) {
2298 for (
int j=0; j<matCols; j++)
2300 val_i += matData(vectorOrdinal,i,j) * vecData(vectorOrdinal,j);
2304 else if (rank_ == 1)
2307 Kokkos::RangePolicy<ExecutionSpace> policy(0,matRows);
2308 Kokkos::parallel_for(
"compute mat-vec", policy,
2309 KOKKOS_LAMBDA (
const int &i) {
2312 for (
int j=0; j<matCols; j++)
2314 val_i += matData(i,j) * vecData(j);
2321 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::logic_error,
"rank not yet supported");
2337 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(A_MatData.
rank() != B_MatData.
rank(), std::invalid_argument,
"AmatData and BmatData have incompatible ranks");
2339 const int D1_DIM = A_MatData.
rank() - 2;
2340 const int D2_DIM = A_MatData.
rank() - 1;
2342 const int A_rows = A_MatData.
extent_int(D1_DIM);
2343 const int A_cols = A_MatData.
extent_int(D2_DIM);
2344 const int B_rows = B_MatData.
extent_int(D1_DIM);
2345 const int B_cols = B_MatData.
extent_int(D2_DIM);
2347 const int leftRows = transposeA ? A_cols : A_rows;
2348 const int leftCols = transposeA ? A_rows : A_cols;
2349 const int rightCols = transposeB ? B_rows : B_cols;
2351 #ifdef INTREPID2_HAVE_DEBUG 2352 const int rightRows = transposeB ? B_cols : B_rows;
2353 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(leftCols != rightRows, std::invalid_argument,
"inner dimensions do not match");
2359 using ExecutionSpace =
typename DeviceType::execution_space;
2361 const int diagonalStart = (variationType_[D1_DIM] ==
BLOCK_PLUS_DIAGONAL) ? blockPlusDiagonalLastNonDiagonal_ + 1 : leftRows;
2366 auto policy = Kokkos::RangePolicy<ExecutionSpace>(0,
getDataExtent(0));
2367 Kokkos::parallel_for(
"compute mat-mat", policy,
2368 KOKKOS_LAMBDA (
const int &matrixOrdinal) {
2369 for (
int i=0; i<diagonalStart; i++)
2371 for (
int j=0; j<rightCols; j++)
2375 for (
int k=0; k<leftCols; k++)
2377 const auto & left = transposeA ? A_MatData(matrixOrdinal,k,i) : A_MatData(matrixOrdinal,i,k);
2378 const auto & right = transposeB ? B_MatData(matrixOrdinal,j,k) : B_MatData(matrixOrdinal,k,j);
2379 val_ij += left * right;
2383 for (
int i=diagonalStart; i<leftRows; i++)
2386 const auto & left = A_MatData(matrixOrdinal,i,i);
2387 const auto & right = B_MatData(matrixOrdinal,i,i);
2388 val_ii = left * right;
2392 else if (rank_ == 4)
2396 if (underlyingMatchesLogical_)
2398 Kokkos::parallel_for(
"compute mat-mat", policy,
2399 KOKKOS_LAMBDA (
const int &cellOrdinal,
const int &pointOrdinal) {
2400 for (
int i=0; i<leftCols; i++)
2402 for (
int j=0; j<rightCols; j++)
2406 for (
int k=0; k<leftCols; k++)
2408 const auto & left = transposeA ? A_MatData(cellOrdinal,pointOrdinal,k,i) : A_MatData(cellOrdinal,pointOrdinal,i,k);
2409 const auto & right = transposeB ? B_MatData(cellOrdinal,pointOrdinal,j,k) : B_MatData(cellOrdinal,pointOrdinal,k,j);
2410 val_ij += left * right;
2418 Kokkos::parallel_for(
"compute mat-mat", policy,
2419 KOKKOS_LAMBDA (
const int &cellOrdinal,
const int &pointOrdinal) {
2420 for (
int i=0; i<diagonalStart; i++)
2422 for (
int j=0; j<rightCols; j++)
2426 for (
int k=0; k<leftCols; k++)
2428 const auto & left = transposeA ? A_MatData(cellOrdinal,pointOrdinal,k,i) : A_MatData(cellOrdinal,pointOrdinal,i,k);
2429 const auto & right = transposeB ? B_MatData(cellOrdinal,pointOrdinal,j,k) : B_MatData(cellOrdinal,pointOrdinal,k,j);
2430 val_ij += left * right;
2434 for (
int i=diagonalStart; i<leftRows; i++)
2437 const auto & left = A_MatData(cellOrdinal,pointOrdinal,i,i);
2438 const auto & right = B_MatData(cellOrdinal,pointOrdinal,i,i);
2439 val_ii = left * right;
2447 INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(
true, std::logic_error,
"rank not yet supported");
2452 KOKKOS_INLINE_FUNCTION constexpr
bool isValid()
const 2454 return extents_[0] > 0;
2458 KOKKOS_INLINE_FUNCTION
2470 void setExtent(
const ordinal_type &d,
const ordinal_type &newExtent)
2472 INTREPID2_TEST_FOR_EXCEPTION(variationType_[d] ==
BLOCK_PLUS_DIAGONAL, std::invalid_argument,
"setExtent is not supported for BLOCK_PLUS_DIAGONAL dimensions");
2474 if (variationType_[d] ==
MODULAR)
2476 bool dividesEvenly = ((newExtent / variationModulus_[d]) * variationModulus_[d] == newExtent);
2477 INTREPID2_TEST_FOR_EXCEPTION(!dividesEvenly, std::invalid_argument,
"when setExtent is called on dimenisions with MODULAR variation, the modulus must divide the new extent evenly");
2480 if ((newExtent != extents_[d]) && (variationType_[d] ==
GENERAL))
2483 std::vector<ordinal_type> newExtents(dataRank_,-1);
2484 for (
int r=0; r<dataRank_; r++)
2486 if (activeDims_[r] == d)
2489 newExtents[r] = newExtent;
2500 case 1: Kokkos::resize(data1_,newExtents[0]);
2502 case 2: Kokkos::resize(data2_,newExtents[0],newExtents[1]);
2504 case 3: Kokkos::resize(data3_,newExtents[0],newExtents[1],newExtents[2]);
2506 case 4: Kokkos::resize(data4_,newExtents[0],newExtents[1],newExtents[2],newExtents[3]);
2508 case 5: Kokkos::resize(data5_,newExtents[0],newExtents[1],newExtents[2],newExtents[3],newExtents[4]);
2510 case 6: Kokkos::resize(data6_,newExtents[0],newExtents[1],newExtents[2],newExtents[3],newExtents[4],newExtents[5]);
2512 case 7: Kokkos::resize(data7_,newExtents[0],newExtents[1],newExtents[2],newExtents[3],newExtents[4],newExtents[5],newExtents[6]);
2514 default: INTREPID2_TEST_FOR_EXCEPTION(
true, std::logic_error,
"Unexpected dataRank_ value");
2518 extents_[d] = newExtent;
2522 KOKKOS_INLINE_FUNCTION
2525 return underlyingMatchesLogical_;
enable_if_t< rank==7, Kokkos::MDRangePolicy< typename DeviceType::execution_space, Kokkos::Rank< 6 > > > dataExtentRangePolicy()
returns an MDRangePolicy over the first six underlying data extents (but with the logical shape)...
KOKKOS_INLINE_FUNCTION void setUnderlyingView4(Kokkos::View< DataScalar ****, DeviceType > &view) const
sets the View that stores the unique data. For rank-4 underlying containers.
KOKKOS_INLINE_FUNCTION bool isDiagonal() const
returns true for containers that have two dimensions marked as BLOCK_PLUS_DIAGONAL for which the non-...
void applyOperator(UnaryOperator unaryOperator)
applies the specified unary operator to each entry
KOKKOS_INLINE_FUNCTION DimensionInfo getDimensionInfo(const int &dim) const
Returns an object fully specifying the indicated dimension. This is used in determining appropriate s...
static KOKKOS_INLINE_FUNCTION int blockPlusDiagonalDiagonalEntryIndex(const int &lastNondiagonal, const int &numNondiagonalEntries, const int &i)
Returns flattened index of the specified (i,i) matrix entry, assuming that i > lastNondiagonal. Only applicable for BLOCK_PLUS_DIAGONAL DataVariationType.
ScalarView< DataScalar, DeviceType > getUnderlyingView() const
Returns a DynRankView constructed atop the same underlying data as the fixed-rank Kokkos::View used i...
KOKKOS_INLINE_FUNCTION enable_if_t< rank==4, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 4...
KOKKOS_INLINE_FUNCTION void setUnderlyingView6(Kokkos::View< DataScalar ******, DeviceType > &view) const
sets the View that stores the unique data. For rank-6 underlying containers.
static KOKKOS_INLINE_FUNCTION int blockPlusDiagonalBlockEntryIndex(const int &lastNondiagonal, const int &numNondiagonalEntries, const int &i, const int &j)
//! Returns flattened index of the specified (i,j) matrix entry, assuming that i,j ≤ lastNondiagonal...
static Data< DataScalar, DeviceType > allocateInPlaceCombinationResult(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B)
KOKKOS_INLINE_FUNCTION Kokkos::Array< int, 7 > getExtents() const
Returns an array containing the logical extents in each dimension.
static void copyContainer(ToContainer to, FromContainer from)
Generic data copying method to allow construction of Data object from DynRankViews for which deep_cop...
KOKKOS_INLINE_FUNCTION enable_if_t< rank==3, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 3...
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar ******, DeviceType > & getUnderlyingView6() const
returns the View that stores the unique data. For rank-6 underlying containers.
KOKKOS_INLINE_FUNCTION void setUnderlyingView7(Kokkos::View< DataScalar *******, DeviceType > &view) const
sets the View that stores the unique data. For rank-7 underlying containers.
KOKKOS_INLINE_FUNCTION ordinal_type getUnderlyingViewRank() const
returns the rank of the View that stores the unique data
enable_if_t< rank !=7, void > storeInPlaceCombination(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B, BinaryOperator binaryOperator)
storeInPlaceCombination with compile-time rank – implementation for rank < 7.
KOKKOS_INLINE_FUNCTION void setUnderlyingView2(Kokkos::View< DataScalar **, DeviceType > &view) const
sets the View that stores the unique data. For rank-2 underlying containers.
KOKKOS_INLINE_FUNCTION void setUnderlyingView3(Kokkos::View< DataScalar ***, DeviceType > &view) const
sets the View that stores the unique data. For rank-3 underlying containers.
KOKKOS_INLINE_FUNCTION DimensionInfo combinedDataDimensionInfo(const Data &otherData, const int &dim) const
Returns (DataVariationType, data extent) in the specified dimension for a Data container that combine...
void storeMatVec(const Data< DataScalar, DeviceType > &matData, const Data< DataScalar, DeviceType > &vecData)
Places the result of a matrix-vector multiply corresponding to the two provided containers into this ...
KOKKOS_INLINE_FUNCTION enable_if_t< rank==1, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 1...
KOKKOS_INLINE_FUNCTION enable_if_t< rank==7, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 7...
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 clear() const
Copies 0.0 to the underlying View.
KOKKOS_INLINE_FUNCTION constexpr bool isValid() const
returns true for containers that have data; false for those that don't (namely, those that have been ...
one of two dimensions in a matrix; bottom-right part of matrix is diagonal
KOKKOS_INLINE_FUNCTION enable_if_t< rank==5, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 5...
Wrapper around a Kokkos::View that allows data that is constant or repeating in various logical dimen...
varies according to modulus of the index
void storeInPlaceCombination(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B, BinaryOperator binaryOperator)
Places the result of an in-place combination (e.g., entrywise sum) into this data container...
Header function for Intrepid2::Util class and other utility functions.
KOKKOS_INLINE_FUNCTION bool underlyingMatchesLogical() const
Returns true if the underlying container has exactly the same rank and extents as the logical contain...
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar *, DeviceType > & getUnderlyingView1() const
returns the View that stores the unique data. For rank-1 underlying containers.
enable_if_t< rank==7, void > storeInPlaceCombination(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B, BinaryOperator binaryOperator)
storeInPlaceCombination with compile-time rank – implementation for rank of 7. (Not optimized; expec...
Data(std::vector< DimensionInfo > dimInfoVector)
Constructor in terms of DimensionInfo for each logical dimension; does not require a View to be speci...
Struct expressing all variation information about a Data object in a single dimension, including its logical extent and storage extent.
KOKKOS_INLINE_FUNCTION unsigned rank() const
Returns the logical rank of the Data container.
void storeInPlaceDifference(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B)
stores the in-place (entrywise) difference, A .- B, into this container.
KOKKOS_INLINE_FUNCTION enable_if_t< valid_args< IntArgs... >::value &&(sizeof...(IntArgs)<=7), return_type > operator()(const IntArgs &... intArgs) const
Returns a value corresponding to the specified logical data location.
Data(DataScalar constantValue, Kokkos::Array< int, rank > extents)
constructor for everywhere-constant data
KOKKOS_INLINE_FUNCTION enable_if_t< rank==2, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 2...
enable_if_t<(rank!=1) &&(rank!=7), Kokkos::MDRangePolicy< typename DeviceType::execution_space, Kokkos::Rank< rank > > > dataExtentRangePolicy()
returns an MDRangePolicy over the underlying data extents (but with the logical shape).
KOKKOS_INLINE_FUNCTION int extent_int(const int &r) const
Returns the logical extent in the specified dimension.
void copyDataFromDynRankViewMatchingUnderlying(const ScalarView< DataScalar, DeviceType > &dynRankView) const
Copies from the provided DynRankView into the underlying Kokkos::View container storing the unique da...
KOKKOS_INLINE_FUNCTION return_type getEntry(const IntArgs &... intArgs) const
Returns a (read-only) value corresponding to the specified logical data location. ...
static Data< DataScalar, DeviceType > allocateMatVecResult(const Data< DataScalar, DeviceType > &matData, const Data< DataScalar, DeviceType > &vecData)
void storeInPlaceQuotient(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B)
stores the in-place (entrywise) quotient, A ./ B, into this container.
void storeInPlaceSum(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B)
stores the in-place (entrywise) sum, A .+ B, into this container.
Data(const ScalarView< DataScalar, DeviceType > &data, int rank, Kokkos::Array< int, 7 > extents, Kokkos::Array< DataVariationType, 7 > variationType, const int blockPlusDiagonalLastNonDiagonal=-1)
DynRankView constructor. Will copy to a View of appropriate rank.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar ****, DeviceType > & getUnderlyingView4() const
returns the View that stores the unique data. For rank-4 underlying containers.
void storeMatMat(const bool transposeA, const Data< DataScalar, DeviceType > &A_MatData, const bool transposeB, const Data< DataScalar, DeviceType > &B_MatData)
void allocateAndCopyFromDynRankView(ScalarView< DataScalar, DeviceType > data)
allocate an underlying View that matches the provided DynRankView in dimensions, and copy...
KOKKOS_INLINE_FUNCTION DimensionInfo combinedDimensionInfo(const DimensionInfo &myData, const DimensionInfo &otherData)
Returns DimensionInfo for a Data container that combines (through multiplication, say...
static KOKKOS_INLINE_FUNCTION int blockPlusDiagonalNumNondiagonalEntries(const int &lastNondiagonal)
Returns the number of non-diagonal entries based on the last non-diagonal. Only applicable for BLOCK_...
KOKKOS_INLINE_FUNCTION enable_if_t< rank==6, const Kokkos::View< typename RankExpander< DataScalar, rank >::value_type, DeviceType > & > getUnderlyingView() const
Returns the underlying view. Throws an exception if the underlying view is not rank 6...
void setActiveDims()
class initialization method. Called by constructors.
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.
Data()
default constructor (empty data)
ScalarView< DataScalar, DeviceType > allocateDynRankViewMatchingUnderlying() const
Returns a DynRankView that matches the underlying Kokkos::View object in value_type, layout, and dimension.
void setExtent(const ordinal_type &d, const ordinal_type &newExtent)
sets the logical extent in the specified dimension. If needed, the underlying data container is resiz...
A singleton class for a DynRankView containing exactly one zero entry. (Technically, the entry is DataScalar(), the default value for the scalar type.) This allows View-wrapping classes to return a reference to zero, even when that zero is not explicitly stored in the wrapped views.
void storeInPlaceProduct(const Data< DataScalar, DeviceType > &A, const Data< DataScalar, DeviceType > &B)
stores the in-place (entrywise) product, A .* B, into this container.
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar *******, DeviceType > & getUnderlyingView7() const
returns the View that stores the unique data. For rank-7 underlying containers.
KOKKOS_INLINE_FUNCTION reference_type getWritableEntry(const IntArgs... intArgs) const
Returns an l-value reference to the specified logical entry in the underlying view. Note that for variation types other than GENERAL, multiple valid argument sets will refer to the same memory location. Intended for Intrepid2 developers and expert users only.
KOKKOS_INLINE_FUNCTION void setUnderlyingView1(Kokkos::View< DataScalar *, DeviceType > &view) const
sets the View that stores the unique data. For rank-1 underlying containers.
KOKKOS_INLINE_FUNCTION void setUnderlyingView5(Kokkos::View< DataScalar *****, DeviceType > &view) const
sets the View that stores the unique data. For rank-5 underlying containers.
void storeInPlaceCombination(PolicyType &policy, ThisUnderlyingViewType &this_underlying, AUnderlyingViewType &A_underlying, BUnderlyingViewType &B_underlying, BinaryOperator &binaryOperator, ArgExtractorThis argThis, ArgExtractorA argA, ArgExtractorB argB)
storeInPlaceCombination implementation for rank < 7, with compile-time underlying views and argument ...
KOKKOS_INLINE_FUNCTION int getVariationModulus(const int &d) const
Variation modulus accessor.
Data(ScalarView< DataScalar, DeviceType > data)
copy constructor modeled after the copy-like constructor above. Not as efficient as the implicit copy...
Data(const Kokkos::DynRankView< DataScalar, DeviceType, DynRankViewProperties... > &data, Kokkos::Array< int, rank > extents, Kokkos::Array< DataVariationType, rank > variationType, const int blockPlusDiagonalLastNonDiagonal=-1)
Constructor that accepts a DynRankView as an argument. The data belonging to the DynRankView will be ...
KOKKOS_INLINE_FUNCTION int getDataExtent(const ordinal_type &d) const
returns the true extent of the data corresponding to the logical dimension provided; if the data does...
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 > & getUnderlyingView5() const
returns the View that stores the unique data. For rank-5 underlying containers.
DataVariationType
Enumeration to indicate how data varies in a particular dimension of an Intrepid2::Data object...
ScalarView< DataScalar, DeviceType > allocateDynRankViewMatchingUnderlying(DimArgs... dims) const
Returns a DynRankView that matches the underlying Kokkos::View object value_type and layout...
KOKKOS_INLINE_FUNCTION const Kokkos::View< DataScalar ***, DeviceType > & getUnderlyingView3() const
returns the View that stores the unique data. For rank-3 underlying containers.
static Data< DataScalar, DeviceType > allocateMatMatResult(const bool transposeA, const Data< DataScalar, DeviceType > &A_MatData, const bool transposeB, const Data< DataScalar, DeviceType > &B_MatData)
KOKKOS_INLINE_FUNCTION ordinal_type getUnderlyingViewSize() const
returns the number of entries in the View that stores the unique data
Data(const Data< DataScalar, OtherDeviceType > &data)
copy-like constructor for differing execution spaces. This does a deep_copy of the underlying view...
KOKKOS_INLINE_FUNCTION int getUnderlyingViewExtent(const int &dim) const
Returns the extent of the underlying view in the specified dimension.
KOKKOS_INLINE_FUNCTION const Kokkos::Array< DataVariationType, 7 > & getVariationTypes() const
Returns an array with the variation types in each logical dimension.