50 #ifndef Intrepid2_TensorData_h
51 #define Intrepid2_TensorData_h
55 #include "Intrepid2_ScalarView.hpp"
63 template<
class Scalar,
typename DeviceType>
67 Kokkos::Array<ordinal_type, 7> extents_;
68 Kokkos::Array<Kokkos::Array<ordinal_type, Parameters::MaxTensorComponents>, 7> entryModulus_;
70 bool separateFirstComponent_ =
false;
71 ordinal_type numTensorComponents_ = 0;
78 ordinal_type maxComponentRank = -1;
79 for (ordinal_type r=0; r<numTensorComponents_; r++)
81 const ordinal_type componentRank = tensorComponents_[r].rank();
82 maxComponentRank = (maxComponentRank > componentRank) ? maxComponentRank : componentRank;
85 ordinal_type rd_start;
86 if ((maxComponentRank == 1) && separateFirstComponent_)
90 extents_[0] = tensorComponents_[0].extent_int(0);
91 entryModulus_[0][0] = extents_[0];
95 rank_ = maxComponentRank;
99 for (ordinal_type d=rd_start; d<7; d++)
102 for (ordinal_type r=rd_start; r<numTensorComponents_; r++)
104 extents_[d] *= tensorComponents_[r].extent_int(d-rd_start);
106 ordinal_type entryModulus = extents_[d];
107 for (ordinal_type r=rd_start; r<numTensorComponents_; r++)
109 entryModulus /= tensorComponents_[r].extent_int(d-rd_start);
110 entryModulus_[d][r] = entryModulus;
124 template<
size_t numTensorComponents>
132 tensorComponents_[r] = tensorComponents[r];
150 numTensorComponents_(tensorComponents.size())
152 for (ordinal_type r=0; r<numTensorComponents_; r++)
154 tensorComponents_[r] = tensorComponents[r];
194 TensorData(Kokkos::Array<
Data<Scalar,DeviceType>, 1>({tensorComponent}),
false)
204 extents_({0,0,0,0,0,0,0}),
217 numTensorComponents_(whichComps.size())
220 for (
const auto & componentOrdinal : whichComps)
229 template<typename OtherDeviceType, class = typename std::enable_if< std::is_same<typename DeviceType::memory_space, typename OtherDeviceType::memory_space>::value>::type,
230 class = typename std::enable_if<!std::is_same<DeviceType,OtherDeviceType>::value>::type>
231 TensorData(const TensorData<Scalar,OtherDeviceType> &tensorData)
233 if (tensorData.isValid())
236 for (ordinal_type r=0; r<numTensorComponents_; r++)
245 extents_ = Kokkos::Array<ordinal_type,7>{0,0,0,0,0,0,0};
253 template<typename OtherDeviceType, class = typename std::enable_if<!std::is_same<typename DeviceType::memory_space, typename OtherDeviceType::memory_space>::value>::type>
259 for (ordinal_type r=0; r<numTensorComponents_; r++)
268 extents_ = Kokkos::Array<ordinal_type,7>{0,0,0,0,0,0,0};
278 KOKKOS_INLINE_FUNCTION
281 return tensorComponents_[r];
289 template <
typename iType0>
290 KOKKOS_INLINE_FUNCTION
typename std::enable_if<std::is_integral<iType0>::value, Scalar>::type
292 #ifdef HAVE_INTREPID2_DEBUG
296 iType0 remainingEntryOrdinal = tensorEntryIndex;
297 for (ordinal_type r=0; r<numTensorComponents_; r++)
299 const ordinal_type componentEntryCount = tensorComponents_[r].extent_int(0);
300 const ordinal_type componentEntryOrdinal = remainingEntryOrdinal % componentEntryCount;
301 remainingEntryOrdinal /= componentEntryCount;
303 value *= tensorComponents_[r](componentEntryOrdinal);
314 template <
typename iType0, ordinal_type numTensorComponents>
315 KOKKOS_INLINE_FUNCTION
typename std::enable_if<std::is_integral<iType0>::value, Scalar>::type
316 operator()(
const Kokkos::Array<iType0,numTensorComponents>& entryComponents)
const {
317 #ifdef HAVE_INTREPID2_DEBUG
324 value *= tensorComponents_[r](entryComponents[r]);
339 template <
typename iType0,
typename iType1>
340 KOKKOS_INLINE_FUNCTION
typename std::enable_if<
341 (std::is_integral<iType0>::value && std::is_integral<iType1>::value),
343 operator()(
const iType0& tensorEntryIndex0,
const iType1& tensorEntryIndex1)
const {
344 #ifdef HAVE_INTREPID2_DEBUG
348 if (numTensorComponents_ == 1)
350 return tensorComponents_[0](tensorEntryIndex0,tensorEntryIndex1);
353 if (!separateFirstComponent_)
356 iType0 remainingEntryOrdinal0 = tensorEntryIndex0;
357 iType1 remainingEntryOrdinal1 = tensorEntryIndex1;
358 for (ordinal_type r=0; r<numTensorComponents_; r++)
360 auto & component = tensorComponents_[r];
361 const ordinal_type componentEntryCount0 = component.extent_int(0);
362 const ordinal_type componentEntryCount1 = component.extent_int(1);
363 const iType0 componentEntryOrdinal0 = remainingEntryOrdinal0 % componentEntryCount0;
364 const iType1 componentEntryOrdinal1 = remainingEntryOrdinal1 % componentEntryCount1;
365 remainingEntryOrdinal0 /= componentEntryCount0;
366 remainingEntryOrdinal1 /= componentEntryCount1;
368 const ordinal_type componentRank = component.rank();
370 if (componentRank == 2)
372 value *= component(componentEntryOrdinal0,componentEntryOrdinal1);
374 else if (componentRank == 1)
376 value *= component(componentEntryOrdinal0);
388 Scalar value = tensorComponents_[0](tensorEntryIndex0);
389 iType0 remainingEntryOrdinal = tensorEntryIndex1;
390 for (ordinal_type r=1; r<numTensorComponents_; r++)
392 const ordinal_type componentEntryCount = tensorComponents_[r].extent_int(0);
393 const ordinal_type componentEntryOrdinal = remainingEntryOrdinal % componentEntryCount;
394 remainingEntryOrdinal /= componentEntryCount;
396 value *= tensorComponents_[r](componentEntryOrdinal);
403 KOKKOS_INLINE_FUNCTION
404 ordinal_type
getTensorComponentIndex(
const ordinal_type &tensorComponent,
const ordinal_type &dim,
const ordinal_type &enumerationIndex)
const
406 ordinal_type remainingEntryOrdinal = enumerationIndex;
407 for (ordinal_type r=0; r<tensorComponent; r++)
409 const auto & component = tensorComponents_[r];
410 const ordinal_type & componentEntryCount = component.extent_int(dim);
412 remainingEntryOrdinal /= componentEntryCount;
414 return remainingEntryOrdinal % tensorComponents_[tensorComponent].extent_int(dim);
426 template <
typename iType0,
typename iType1,
typename iType2>
427 KOKKOS_INLINE_FUNCTION
typename std::enable_if<
428 (std::is_integral<iType0>::value && std::is_integral<iType1>::value && std::is_integral<iType2>::value),
430 operator()(
const iType0& tensorEntryIndex0,
const iType1& tensorEntryIndex1,
const iType2& tensorEntryIndex2)
const
432 #ifdef HAVE_INTREPID2_DEBUG
438 Kokkos::Array<ordinal_type,3> remainingEntryOrdinal {tensorEntryIndex0, tensorEntryIndex1, tensorEntryIndex2};
439 for (ordinal_type r=0; r<numTensorComponents_; r++)
441 auto & component = tensorComponents_[r];
442 const ordinal_type componentEntryCount0 = component.extent_int(0);
443 const ordinal_type componentEntryCount1 = component.extent_int(1);
444 const ordinal_type componentEntryCount2 = component.extent_int(2);
445 const ordinal_type componentEntryOrdinal0 = remainingEntryOrdinal[0] % componentEntryCount0;
446 const ordinal_type componentEntryOrdinal1 = remainingEntryOrdinal[1] % componentEntryCount1;
447 const ordinal_type componentEntryOrdinal2 = remainingEntryOrdinal[2] % componentEntryCount2;
448 remainingEntryOrdinal[0] /= componentEntryCount0;
449 remainingEntryOrdinal[1] /= componentEntryCount1;
450 remainingEntryOrdinal[2] /= componentEntryCount2;
452 const ordinal_type componentRank = component.rank();
454 if (componentRank == 3)
456 value *= component(componentEntryOrdinal0,componentEntryOrdinal1,componentEntryOrdinal2);
458 else if (componentRank == 2)
460 value *= component(componentEntryOrdinal0,componentEntryOrdinal1);
462 else if (componentRank == 1)
464 value *= component(componentEntryOrdinal0);
481 template <
typename iType0,
typename iType1, ordinal_type numTensorComponents>
482 KOKKOS_INLINE_FUNCTION
typename std::enable_if<
483 (std::is_integral<iType0>::value && std::is_integral<iType1>::value),
485 operator()(
const Kokkos::Array<iType0,numTensorComponents>& entryComponents0,
const Kokkos::Array<iType1,numTensorComponents>& entryComponents1)
const {
486 #ifdef HAVE_INTREPID2_DEBUG
492 auto & component = tensorComponents_[r];
493 const ordinal_type componentRank = component.rank();
494 if (componentRank == 2)
496 value *= component(entryComponents0[r],entryComponents1[r]);
498 else if (componentRank == 1)
500 value *= component(entryComponents0[r]);
515 template <
typename iType>
516 KOKKOS_INLINE_FUNCTION
517 typename std::enable_if<std::is_integral<iType>::value, ordinal_type>::type
527 template <
typename iType>
528 KOKKOS_INLINE_FUNCTION constexpr
529 typename std::enable_if<std::is_integral<iType>::value,
size_t>::type
535 KOKKOS_INLINE_FUNCTION constexpr
bool isValid()
const
537 return extents_[0] > 0;
541 KOKKOS_INLINE_FUNCTION
548 KOKKOS_INLINE_FUNCTION
551 return numTensorComponents_;
555 KOKKOS_INLINE_FUNCTION
558 return separateFirstComponent_;
564 INTREPID2_TEST_FOR_EXCEPTION(!separateFirstComponent_ && (numTensorComponents_ != 1), std::invalid_argument,
"setFirstComponentExtent() is only allowed when separateFirstComponent_ is true, or there is only one component");
565 tensorComponents_[0].setExtent(0,newExtent);
566 extents_[0] = newExtent;
TensorData(Kokkos::Array< Data< Scalar, DeviceType >, numTensorComponents > tensorComponents, bool separateFirstComponent=false)
Constructor with fixed-length Kokkos::Array argument.
#define INTREPID2_TEST_FOR_EXCEPTION_DEVICE_SAFE(test, x, msg)
Defines the Data class, a wrapper around a Kokkos::View that allows data that is constant or repeatin...
TensorData(const TensorData &first, const TensorData &second, bool separateFirstComponent=false)
Constructor to combine two other TensorData objects.
KOKKOS_INLINE_FUNCTION ordinal_type getTensorComponentIndex(const ordinal_type &tensorComponent, const ordinal_type &dim, const ordinal_type &enumerationIndex) const
return the index into the specified tensorial component in the dimension specified corresponding to t...
KOKKOS_INLINE_FUNCTION const Data< Scalar, DeviceType > & getTensorComponent(const ordinal_type &r) const
Returns the requested tensor component.
KOKKOS_INLINE_FUNCTION constexpr bool isValid() const
Returns true for containers that have data; false for those that don't (e.g., those that have been co...
void initialize()
Initialize members based on constructor parameters.
TensorData(const TensorData< Scalar, OtherDeviceType > &tensorData)
Copy-like constructor for differing execution spaces. This performs a deep copy of the underlying dat...
KOKKOS_INLINE_FUNCTION ordinal_type rank() const
Returns the rank of the container.
TensorData(TensorData otherTensorData, std::vector< int > whichComps)
Constructor that takes a subset of the tensorial components of another TensorData container...
TensorData(std::vector< Data< Scalar, DeviceType > > tensorComponents, bool separateFirstComponent=false)
Constructor with variable-length std::vector containing the components.
Contains definitions of custom data types in Intrepid2.
KOKKOS_INLINE_FUNCTION std::enable_if< std::is_integral< iType0 >::value, Scalar >::type operator()(const iType0 &tensorEntryIndex) const
Accessor for rank-1 objects.
TensorData()
Default constructor.
void setFirstComponentExtentInDimension0(const ordinal_type &newExtent)
Sets the extent of the first component. Only valid when either there is only one component, or when separateFirstComponent() returns true. The intended use case is when the 0 dimension in first component represents a cell index, and the container is resized to match a workset size that does not evenly divide the number of cells.
KOKKOS_INLINE_FUNCTION std::enable_if< std::is_integral< iType0 >::value, Scalar >::type operator()(const Kokkos::Array< iType0, numTensorComponents > &entryComponents) const
Accessor that accepts a fixed-length array with entries corresponding to component indices...
View-like interface to tensor data; tensor components are stored separately and multiplied together a...
KOKKOS_INLINE_FUNCTION std::enable_if< std::is_integral< iType >::value, ordinal_type >::type extent_int(const iType &d) const
Returns the logical extent in the requested dimension.
KOKKOS_INLINE_FUNCTION bool separateFirstComponent() const
Returns true if the first component is indexed separately; false if not.
KOKKOS_INLINE_FUNCTION ordinal_type numTensorComponents() const
Return the number of tensorial components.
TensorData(Data< Scalar, DeviceType > tensorComponent)
Simple constructor for the case of trivial tensor-product structure (single component) ...
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< std::is_integral< iType >::value, size_t >::type extent(const iType &d) const
Returns the logical extent in the requested dimension.
static constexpr ordinal_type MaxTensorComponents
Maximum number of tensor/Cartesian products that can be taken: this allows hypercube basis in 7D to b...