50 #ifndef Intrepid2_TensorPoints_h
51 #define Intrepid2_TensorPoints_h
53 #include <Kokkos_Core.hpp>
59 template<
class Po
intScalar,
typename DeviceType>
63 ordinal_type numTensorComponents_;
64 ordinal_type totalPointCount_;
65 ordinal_type totalDimension_;
66 Kokkos::View<ordinal_type*, DeviceType> dimToComponent_;
67 Kokkos::View<ordinal_type*, DeviceType> dimToComponentDim_;
68 Kokkos::Array<ordinal_type, Parameters::MaxTensorComponents> pointModulus_;
69 Kokkos::Array<ordinal_type, Parameters::MaxTensorComponents> pointDivisor_;
72 using reference_type =
typename ScalarView<PointScalar,DeviceType>::reference_type;
74 void TEST_VALID_POINT_COMPONENTS()
76 #ifdef HAVE_INTREPID2_DEBUG
79 for (ordinal_type r=0; r<numTensorComponents_; r++)
81 const auto & pointComponent = pointTensorComponents_[r];
82 INTREPID2_TEST_FOR_EXCEPTION(2 != pointComponent.rank(), std::invalid_argument,
"Each component must have shape (P,D)");
83 INTREPID2_TEST_FOR_EXCEPTION(pointComponent.extent_int(0) <= 0, std::invalid_argument,
"Each component must have at least one point");
94 TEST_VALID_POINT_COMPONENTS();
98 for (ordinal_type r=0; r<numTensorComponents_; r++)
100 totalPointCount_ *= pointTensorComponents_[r].extent_int(0);
101 totalDimension_ += pointTensorComponents_[r].extent_int(1);
103 ordinal_type pointDivisor = 1;
104 for (ordinal_type r=0; r<numTensorComponents_; r++)
106 pointModulus_[r] = pointTensorComponents_[r].extent_int(0);
107 pointDivisor_[r] = pointDivisor;
108 pointDivisor *= pointTensorComponents_[r].extent_int(0);
110 dimToComponent_ = Kokkos::View<ordinal_type*, DeviceType>(
"dimToComponent_",totalDimension_);
111 dimToComponentDim_ = Kokkos::View<ordinal_type*, DeviceType>(
"dimToComponentDim_",totalDimension_);
113 ordinal_type dimsSoFar = 0;
115 auto dimToComponentHost = Kokkos::create_mirror_view(dimToComponent_);
116 auto dimToComponentDimHost = Kokkos::create_mirror_view(dimToComponentDim_);
117 for (ordinal_type r=0; r<numTensorComponents_; r++)
119 const int componentDim = pointTensorComponents_[r].extent_int(1);
120 for (
int i=0; i<componentDim; i++)
122 dimToComponentHost[d] = r;
123 dimToComponentDimHost[d] = d - dimsSoFar;
126 dimsSoFar += componentDim;
128 Kokkos::deep_copy(dimToComponent_,dimToComponentHost);
129 Kokkos::deep_copy(dimToComponentDim_,dimToComponentDimHost);
138 template<
size_t numTensorComponents>
146 pointTensorComponents_[r] = pointTensorComponents[r];
161 numTensorComponents_(whichDims.size()),
165 for (
const auto & componentOrdinal : whichDims)
167 pointTensorComponents_[r++] = otherPointsContainer.
getTensorComponent(componentOrdinal);
179 TensorPoints(std::vector< ScalarView<PointScalar,DeviceType>> pointTensorComponents)
181 numTensorComponents_(pointTensorComponents.size()),
184 for (ordinal_type r=0; r<numTensorComponents_; r++)
186 pointTensorComponents_[r] = pointTensorComponents[r];
200 numTensorComponents_(1),
203 pointTensorComponents_[0] = points;
212 template<
class OtherPo
intsContainer>
215 const int numPoints = fromPoints.extent_int(0);
216 const int numDims = fromPoints.extent_int(1);
217 using ExecutionSpace =
typename DeviceType::execution_space;
218 auto policy = Kokkos::MDRangePolicy<ExecutionSpace,Kokkos::Rank<2>>({0,0},{numPoints,numDims});
219 Kokkos::parallel_for(
"copy points", policy,
220 KOKKOS_LAMBDA (
const int &i0,
const int &i1) {
221 toPoints(i0,i1) = fromPoints(i0,i1);
226 template<typename OtherDeviceType, class = typename std::enable_if< std::is_same<typename DeviceType::memory_space, typename OtherDeviceType::memory_space>::value>::type,
227 class = typename std::enable_if<!std::is_same<DeviceType,OtherDeviceType>::value>::type>
228 TensorPoints(const TensorPoints<PointScalar,OtherDeviceType> &tensorPoints)
230 numTensorComponents_(tensorPoints.numTensorComponents()),
231 isValid_(tensorPoints.isValid())
235 for (ordinal_type r=0; r<numTensorComponents_; r++)
244 template<typename OtherDeviceType, class = typename std::enable_if<!std::is_same<typename DeviceType::memory_space, typename OtherDeviceType::memory_space>::value>::type>
248 isValid_(tensorPoints.
isValid())
252 for (ordinal_type r=0; r<numTensorComponents_; r++)
254 ScalarView<PointScalar,OtherDeviceType> otherPointComponent = tensorPoints.
getTensorComponent(r);
255 const int numPoints = otherPointComponent.extent_int(0);
256 const int numDims = otherPointComponent.extent_int(1);
257 pointTensorComponents_[r] = ScalarView<PointScalar,DeviceType>(
"Intrepid2 point component", numPoints, numDims);
259 using MemorySpace =
typename DeviceType::memory_space;
260 auto pointComponentMirror = Kokkos::create_mirror_view_and_copy(MemorySpace(), otherPointComponent);
279 return pointTensorComponents_[tensorComponentOrdinal].extent_int(0);
290 template <
typename iType0,
typename iType1>
291 KOKKOS_INLINE_FUNCTION
typename std::enable_if<
292 (std::is_integral<iType0>::value && std::is_integral<iType1>::value),
293 reference_type>::type
294 operator()(
const iType0& tensorPointIndex,
const iType1& dim)
const {
295 const ordinal_type component = dimToComponent_[dim];
296 const ordinal_type d = dimToComponentDim_[dim];
297 const ordinal_type componentPointOrdinal = (tensorPointIndex / pointDivisor_[component]) % pointModulus_[component];
298 return pointTensorComponents_[component](componentPointOrdinal,d);
309 template <
typename iType0,
typename iType1,
size_t numTensorComponents>
310 KOKKOS_INLINE_FUNCTION
typename std::enable_if<
311 (std::is_integral<iType0>::value && std::is_integral<iType1>::value),
312 reference_type>::type
313 operator()(
const Kokkos::Array<iType0,numTensorComponents>& pointOrdinalComponents,
const iType1& dim)
const {
314 const ordinal_type component = dimToComponent_[dim];
315 const ordinal_type d = dimToComponentDim_[dim];
316 const ordinal_type componentPointOrdinal = pointOrdinalComponents[component];
317 return pointTensorComponents_[component](componentPointOrdinal,d);
325 template <
typename iType>
326 KOKKOS_INLINE_FUNCTION
327 typename std::enable_if<std::is_integral<iType>::value,
int>::type
329 if (r == static_cast<iType>(0))
331 return totalPointCount_;
333 else if (r == static_cast<iType>(1))
335 return totalDimension_;
348 template <
typename iType>
349 KOKKOS_INLINE_FUNCTION constexpr
350 typename std::enable_if<std::is_integral<iType>::value,
size_t>::type
353 return (r == static_cast<iType>(0)) ? totalPointCount_
355 (r ==
static_cast<iType
>(1)) ? totalDimension_ : 1;
363 ScalarView<PointScalar,DeviceType> expandedRawPoints(
"expanded raw points from TensorPoints", numPoints, spaceDim);
365 using ExecutionSpace =
typename DeviceType::execution_space;
366 Kokkos::parallel_for(
367 Kokkos::MDRangePolicy<ExecutionSpace,Kokkos::Rank<2>>({0,0},{numPoints,spaceDim}),
368 KOKKOS_LAMBDA (
const int &pointOrdinal,
const int &d) {
369 expandedRawPoints(pointOrdinal,d) = tensorPoints(pointOrdinal,d);
371 return expandedRawPoints;
379 KOKKOS_INLINE_FUNCTION
382 return pointTensorComponents_[r];
386 KOKKOS_INLINE_FUNCTION
393 KOKKOS_INLINE_FUNCTION
396 return numTensorComponents_;
400 KOKKOS_INLINE_FUNCTION
401 constexpr ordinal_type
rank()
const
TensorPoints(ScalarView< PointScalar, DeviceType > points)
Constructor for point set with trivial tensor structure.
KOKKOS_INLINE_FUNCTION bool isValid() const
Returns true for containers that have data; false for those that don't (e.g., those that have been co...
View-like interface to tensor points; point components are stored separately; the appropriate coordin...
void initialize()
Initialize members based on constructor parameters.
TensorPoints(const TensorPoints< PointScalar, OtherDeviceType > &tensorPoints)
copy-like constructor for differing memory spaces. This does a deep_copy of the underlying view...
void copyPointsContainer(ScalarView< PointScalar, DeviceType > toPoints, OtherPointsContainer fromPoints)
Copy from one points container, which may be an arbitrary functor, to a DynRankView.
ScalarView< PointScalar, DeviceType > allocateAndFillExpandedRawPointView() const
This method is for compatibility with existing methods that take raw point views. Note that in genera...
KOKKOS_INLINE_FUNCTION constexpr ordinal_type rank() const
Return the rank of the container, which is 2.
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< std::is_integral< iType >::value, size_t >::type extent(const iType &r) const
Returns the logical extent in the requested dimension.
TensorPoints(TensorPoints otherPointsContainer, std::vector< int > whichDims)
Constructor that takes a subset of the tensorial components of another points container.
TensorPoints(Kokkos::Array< ScalarView< PointScalar, DeviceType >, numTensorComponents > pointTensorComponents)
Constructor with fixed-length Kokkos::Array argument.
KOKKOS_INLINE_FUNCTION ordinal_type numTensorComponents() const
Returns the number of tensorial components.
TensorPoints()
Default constructor. TensorPoints::isValid() will return false.
ordinal_type componentPointCount(const ordinal_type &tensorComponentOrdinal) const
Returns the number of points in the indicated component.
KOKKOS_INLINE_FUNCTION ScalarView< PointScalar, DeviceType > getTensorComponent(const ordinal_type &r) const
Returns the requested tensor component.
KOKKOS_INLINE_FUNCTION std::enable_if< std::is_integral< iType >::value, int >::type extent_int(const iType &r) const
Returns the logical extent in the requested dimension.
TensorPoints(std::vector< ScalarView< PointScalar, DeviceType >> pointTensorComponents)
Constructor with variable-length std::vector argument.
static constexpr ordinal_type MaxTensorComponents
Maximum number of tensor/Cartesian products that can be taken: this allows hypercube basis in 7D to b...