1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2014 Benoit Steiner <benoit.steiner.goog (at) gmail.com> 5 // Copyright (C) 2016 Mehdi Goli, Codeplay Software Ltd <eigen (at) codeplay.com> 6 // 7 // This Source Code Form is subject to the terms of the Mozilla 8 // Public License v. 2.0. If a copy of the MPL was not distributed 9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10 11 #ifndef EIGEN_CXX11_TENSOR_TENSOR_REDUCTION_H 12 #define EIGEN_CXX11_TENSOR_TENSOR_REDUCTION_H 13 14 namespace Eigen { 15 16 /** \class TensorReduction 17 * \ingroup CXX11_Tensor_Module 18 * 19 * \brief Tensor reduction class. 20 * 21 */ 22 23 namespace internal { 24 template<typename Op, typename Dims, typename XprType,template <class> class MakePointer_ > 25 struct traits<TensorReductionOp<Op, Dims, XprType, MakePointer_> > 26 : traits<XprType> 27 { 28 typedef traits<XprType> XprTraits; 29 typedef typename XprTraits::Scalar Scalar; 30 typedef typename XprTraits::StorageKind StorageKind; 31 typedef typename XprTraits::Index Index; 32 typedef typename XprType::Nested Nested; 33 static const int NumDimensions = XprTraits::NumDimensions - array_size<Dims>::value; 34 static const int Layout = XprTraits::Layout; 35 36 template <class T> struct MakePointer { 37 // Intermediate typedef to workaround MSVC issue. 38 typedef MakePointer_<T> MakePointerT; 39 typedef typename MakePointerT::Type Type; 40 }; 41 }; 42 43 template<typename Op, typename Dims, typename XprType, template <class> class MakePointer_> 44 struct eval<TensorReductionOp<Op, Dims, XprType, MakePointer_>, Eigen::Dense> 45 { 46 typedef const TensorReductionOp<Op, Dims, XprType, MakePointer_>& type; 47 }; 48 49 template<typename Op, typename Dims, typename XprType, template <class> class MakePointer_> 50 struct nested<TensorReductionOp<Op, Dims, XprType, MakePointer_>, 1, typename eval<TensorReductionOp<Op, Dims, XprType, MakePointer_> >::type> 51 { 52 typedef TensorReductionOp<Op, Dims, XprType, MakePointer_> type; 53 }; 54 55 56 template <typename OutputDims> struct DimInitializer { 57 template <typename InputDims, typename ReducedDims> EIGEN_DEVICE_FUNC 58 static void run(const InputDims& input_dims, 59 const array<bool, internal::array_size<InputDims>::value>& reduced, 60 OutputDims* output_dims, ReducedDims* reduced_dims) { 61 const int NumInputDims = internal::array_size<InputDims>::value; 62 int outputIndex = 0; 63 int reduceIndex = 0; 64 for (int i = 0; i < NumInputDims; ++i) { 65 if (reduced[i]) { 66 (*reduced_dims)[reduceIndex] = input_dims[i]; 67 ++reduceIndex; 68 } else { 69 (*output_dims)[outputIndex] = input_dims[i]; 70 ++outputIndex; 71 } 72 } 73 } 74 }; 75 76 template <> struct DimInitializer<Sizes<> > { 77 template <typename InputDims, typename Index, size_t Rank> EIGEN_DEVICE_FUNC 78 static void run(const InputDims& input_dims, const array<bool, Rank>&, 79 Sizes<>*, array<Index, Rank>* reduced_dims) { 80 const int NumInputDims = internal::array_size<InputDims>::value; 81 for (int i = 0; i < NumInputDims; ++i) { 82 (*reduced_dims)[i] = input_dims[i]; 83 } 84 } 85 }; 86 87 88 template <typename ReducedDims, int NumTensorDims, int Layout> 89 struct are_inner_most_dims { 90 static const bool value = false; 91 }; 92 template <typename ReducedDims, int NumTensorDims, int Layout> 93 struct preserve_inner_most_dims { 94 static const bool value = false; 95 }; 96 97 #if EIGEN_HAS_CONSTEXPR && EIGEN_HAS_VARIADIC_TEMPLATES 98 template <typename ReducedDims, int NumTensorDims> 99 struct are_inner_most_dims<ReducedDims, NumTensorDims, ColMajor>{ 100 static const bool tmp1 = indices_statically_known_to_increase<ReducedDims>(); 101 static const bool tmp2 = index_statically_eq<ReducedDims>(0, 0); 102 static const bool tmp3 = index_statically_eq<ReducedDims>(array_size<ReducedDims>::value-1, array_size<ReducedDims>::value-1); 103 static const bool value = tmp1 & tmp2 & tmp3; 104 }; 105 template <typename ReducedDims, int NumTensorDims> 106 struct are_inner_most_dims<ReducedDims, NumTensorDims, RowMajor>{ 107 static const bool tmp1 = indices_statically_known_to_increase<ReducedDims>(); 108 static const bool tmp2 = index_statically_eq<ReducedDims>(0, NumTensorDims - array_size<ReducedDims>::value); 109 static const bool tmp3 = index_statically_eq<ReducedDims>(array_size<ReducedDims>::value - 1, NumTensorDims - 1); 110 static const bool value = tmp1 & tmp2 & tmp3; 111 112 }; 113 template <typename ReducedDims, int NumTensorDims> 114 struct preserve_inner_most_dims<ReducedDims, NumTensorDims, ColMajor>{ 115 static const bool tmp1 = indices_statically_known_to_increase<ReducedDims>(); 116 static const bool tmp2 = index_statically_gt<ReducedDims>(0, 0); 117 static const bool value = tmp1 & tmp2; 118 119 }; 120 template <typename ReducedDims, int NumTensorDims> 121 struct preserve_inner_most_dims<ReducedDims, NumTensorDims, RowMajor>{ 122 static const bool tmp1 = indices_statically_known_to_increase<ReducedDims>(); 123 static const bool tmp2 = index_statically_lt<ReducedDims>(array_size<ReducedDims>::value - 1, NumTensorDims - 1); 124 static const bool value = tmp1 & tmp2; 125 }; 126 #endif 127 128 129 template <int DimIndex, typename Self, typename Op> 130 struct GenericDimReducer { 131 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self& self, typename Self::Index firstIndex, Op& reducer, typename Self::CoeffReturnType* accum) { 132 EIGEN_STATIC_ASSERT((DimIndex > 0), YOU_MADE_A_PROGRAMMING_MISTAKE); 133 for (int j = 0; j < self.m_reducedDims[DimIndex]; ++j) { 134 const typename Self::Index input = firstIndex + j * self.m_reducedStrides[DimIndex]; 135 GenericDimReducer<DimIndex-1, Self, Op>::reduce(self, input, reducer, accum); 136 } 137 } 138 }; 139 template <typename Self, typename Op> 140 struct GenericDimReducer<0, Self, Op> { 141 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self& self, typename Self::Index firstIndex, Op& reducer, typename Self::CoeffReturnType* accum) { 142 for (int j = 0; j < self.m_reducedDims[0]; ++j) { 143 const typename Self::Index input = firstIndex + j * self.m_reducedStrides[0]; 144 reducer.reduce(self.m_impl.coeff(input), accum); 145 } 146 } 147 }; 148 template <typename Self, typename Op> 149 struct GenericDimReducer<-1, Self, Op> { 150 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self& self, typename Self::Index index, Op& reducer, typename Self::CoeffReturnType* accum) { 151 reducer.reduce(self.m_impl.coeff(index), accum); 152 } 153 }; 154 155 template <typename Self, typename Op, bool Vectorizable = (Self::InputPacketAccess & Op::PacketAccess)> 156 struct InnerMostDimReducer { 157 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename Self::CoeffReturnType reduce(const Self& self, typename Self::Index firstIndex, typename Self::Index numValuesToReduce, Op& reducer) { 158 typename Self::CoeffReturnType accum = reducer.initialize(); 159 for (typename Self::Index j = 0; j < numValuesToReduce; ++j) { 160 reducer.reduce(self.m_impl.coeff(firstIndex + j), &accum); 161 } 162 return reducer.finalize(accum); 163 } 164 }; 165 166 template <typename Self, typename Op> 167 struct InnerMostDimReducer<Self, Op, true> { 168 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename Self::CoeffReturnType reduce(const Self& self, typename Self::Index firstIndex, typename Self::Index numValuesToReduce, Op& reducer) { 169 const int packetSize = internal::unpacket_traits<typename Self::PacketReturnType>::size; 170 const typename Self::Index VectorizedSize = (numValuesToReduce / packetSize) * packetSize; 171 typename Self::PacketReturnType p = reducer.template initializePacket<typename Self::PacketReturnType>(); 172 for (typename Self::Index j = 0; j < VectorizedSize; j += packetSize) { 173 reducer.reducePacket(self.m_impl.template packet<Unaligned>(firstIndex + j), &p); 174 } 175 typename Self::CoeffReturnType accum = reducer.initialize(); 176 for (typename Self::Index j = VectorizedSize; j < numValuesToReduce; ++j) { 177 reducer.reduce(self.m_impl.coeff(firstIndex + j), &accum); 178 } 179 return reducer.finalizeBoth(accum, p); 180 } 181 }; 182 183 template <int DimIndex, typename Self, typename Op, bool vectorizable = (Self::InputPacketAccess & Op::PacketAccess)> 184 struct InnerMostDimPreserver { 185 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self&, typename Self::Index, Op&, typename Self::PacketReturnType*) { 186 eigen_assert(false && "should never be called"); 187 } 188 }; 189 190 template <int DimIndex, typename Self, typename Op> 191 struct InnerMostDimPreserver<DimIndex, Self, Op, true> { 192 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self& self, typename Self::Index firstIndex, Op& reducer, typename Self::PacketReturnType* accum) { 193 EIGEN_STATIC_ASSERT((DimIndex > 0), YOU_MADE_A_PROGRAMMING_MISTAKE); 194 for (typename Self::Index j = 0; j < self.m_reducedDims[DimIndex]; ++j) { 195 const typename Self::Index input = firstIndex + j * self.m_reducedStrides[DimIndex]; 196 InnerMostDimPreserver<DimIndex-1, Self, Op>::reduce(self, input, reducer, accum); 197 } 198 } 199 }; 200 201 template <typename Self, typename Op> 202 struct InnerMostDimPreserver<0, Self, Op, true> { 203 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self& self, typename Self::Index firstIndex, Op& reducer, typename Self::PacketReturnType* accum) { 204 for (typename Self::Index j = 0; j < self.m_reducedDims[0]; ++j) { 205 const typename Self::Index input = firstIndex + j * self.m_reducedStrides[0]; 206 reducer.reducePacket(self.m_impl.template packet<Unaligned>(input), accum); 207 } 208 } 209 }; 210 template <typename Self, typename Op> 211 struct InnerMostDimPreserver<-1, Self, Op, true> { 212 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void reduce(const Self&, typename Self::Index, Op&, typename Self::PacketReturnType*) { 213 eigen_assert(false && "should never be called"); 214 } 215 }; 216 217 // Default full reducer 218 template <typename Self, typename Op, typename Device, bool Vectorizable = (Self::InputPacketAccess & Op::PacketAccess)> 219 struct FullReducer { 220 static const bool HasOptimizedImplementation = false; 221 222 static EIGEN_DEVICE_FUNC void run(const Self& self, Op& reducer, const Device&, typename Self::CoeffReturnType* output) { 223 const typename Self::Index num_coeffs = array_prod(self.m_impl.dimensions()); 224 *output = InnerMostDimReducer<Self, Op, Vectorizable>::reduce(self, 0, num_coeffs, reducer); 225 } 226 }; 227 228 229 #ifdef EIGEN_USE_THREADS 230 // Multithreaded full reducers 231 template <typename Self, typename Op, 232 bool Vectorizable = (Self::InputPacketAccess & Op::PacketAccess)> 233 struct FullReducerShard { 234 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(const Self& self, typename Self::Index firstIndex, 235 typename Self::Index numValuesToReduce, Op& reducer, 236 typename Self::CoeffReturnType* output) { 237 *output = InnerMostDimReducer<Self, Op, Vectorizable>::reduce( 238 self, firstIndex, numValuesToReduce, reducer); 239 } 240 }; 241 242 // Multithreaded full reducer 243 template <typename Self, typename Op, bool Vectorizable> 244 struct FullReducer<Self, Op, ThreadPoolDevice, Vectorizable> { 245 static const bool HasOptimizedImplementation = !Op::IsStateful; 246 static const int PacketSize = 247 unpacket_traits<typename Self::PacketReturnType>::size; 248 249 // launch one reducer per thread and accumulate the result. 250 static void run(const Self& self, Op& reducer, const ThreadPoolDevice& device, 251 typename Self::CoeffReturnType* output) { 252 typedef typename Self::Index Index; 253 const Index num_coeffs = array_prod(self.m_impl.dimensions()); 254 if (num_coeffs == 0) { 255 *output = reducer.finalize(reducer.initialize()); 256 return; 257 } 258 const TensorOpCost cost = 259 self.m_impl.costPerCoeff(Vectorizable) + 260 TensorOpCost(0, 0, internal::functor_traits<Op>::Cost, Vectorizable, 261 PacketSize); 262 const int num_threads = TensorCostModel<ThreadPoolDevice>::numThreads( 263 num_coeffs, cost, device.numThreads()); 264 if (num_threads == 1) { 265 *output = 266 InnerMostDimReducer<Self, Op, Vectorizable>::reduce(self, 0, num_coeffs, reducer); 267 return; 268 } 269 const Index blocksize = 270 std::floor<Index>(static_cast<float>(num_coeffs) / num_threads); 271 const Index numblocks = blocksize > 0 ? num_coeffs / blocksize : 0; 272 eigen_assert(num_coeffs >= numblocks * blocksize); 273 274 Barrier barrier(internal::convert_index<unsigned int>(numblocks)); 275 MaxSizeVector<typename Self::CoeffReturnType> shards(numblocks, reducer.initialize()); 276 for (Index i = 0; i < numblocks; ++i) { 277 device.enqueue_with_barrier(&barrier, &FullReducerShard<Self, Op, Vectorizable>::run, 278 self, i * blocksize, blocksize, reducer, 279 &shards[i]); 280 } 281 typename Self::CoeffReturnType finalShard; 282 if (numblocks * blocksize < num_coeffs) { 283 finalShard = InnerMostDimReducer<Self, Op, Vectorizable>::reduce( 284 self, numblocks * blocksize, num_coeffs - numblocks * blocksize, 285 reducer); 286 } else { 287 finalShard = reducer.initialize(); 288 } 289 barrier.Wait(); 290 291 for (Index i = 0; i < numblocks; ++i) { 292 reducer.reduce(shards[i], &finalShard); 293 } 294 *output = reducer.finalize(finalShard); 295 } 296 }; 297 298 #endif 299 300 301 // Default inner reducer 302 template <typename Self, typename Op, typename Device> 303 struct InnerReducer { 304 static const bool HasOptimizedImplementation = false; 305 306 EIGEN_DEVICE_FUNC static bool run(const Self&, Op&, const Device&, typename Self::CoeffReturnType*, typename Self::Index, typename Self::Index) { 307 eigen_assert(false && "Not implemented"); 308 return true; 309 } 310 }; 311 312 // Default outer reducer 313 template <typename Self, typename Op, typename Device> 314 struct OuterReducer { 315 static const bool HasOptimizedImplementation = false; 316 317 EIGEN_DEVICE_FUNC static bool run(const Self&, Op&, const Device&, typename Self::CoeffReturnType*, typename Self::Index, typename Self::Index) { 318 eigen_assert(false && "Not implemented"); 319 return true; 320 } 321 }; 322 323 324 #if defined(EIGEN_USE_GPU) && defined(__CUDACC__) 325 template <int B, int N, typename S, typename R, typename I> 326 __global__ void FullReductionKernel(R, const S, I, typename S::CoeffReturnType*, unsigned int*); 327 328 329 #ifdef EIGEN_HAS_CUDA_FP16 330 template <typename S, typename R, typename I> 331 __global__ void ReductionInitFullReduxKernelHalfFloat(R, const S, I, half2*); 332 template <int B, int N, typename S, typename R, typename I> 333 __global__ void FullReductionKernelHalfFloat(R, const S, I, half*, half2*); 334 template <int NPT, typename S, typename R, typename I> 335 __global__ void InnerReductionKernelHalfFloat(R, const S, I, I, half*); 336 337 #endif 338 339 template <int NPT, typename S, typename R, typename I> 340 __global__ void InnerReductionKernel(R, const S, I, I, typename S::CoeffReturnType*); 341 342 template <int NPT, typename S, typename R, typename I> 343 __global__ void OuterReductionKernel(R, const S, I, I, typename S::CoeffReturnType*); 344 #endif 345 346 } // end namespace internal 347 348 349 template <typename Op, typename Dims, typename XprType, template <class> class MakePointer_> 350 class TensorReductionOp : public TensorBase<TensorReductionOp<Op, Dims, XprType, MakePointer_>, ReadOnlyAccessors> { 351 public: 352 typedef typename Eigen::internal::traits<TensorReductionOp>::Scalar Scalar; 353 typedef typename Eigen::NumTraits<Scalar>::Real RealScalar; 354 typedef typename internal::remove_const<typename XprType::CoeffReturnType>::type CoeffReturnType; 355 typedef typename Eigen::internal::nested<TensorReductionOp>::type Nested; 356 typedef typename Eigen::internal::traits<TensorReductionOp>::StorageKind StorageKind; 357 typedef typename Eigen::internal::traits<TensorReductionOp>::Index Index; 358 359 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE 360 TensorReductionOp(const XprType& expr, const Dims& dims) : m_expr(expr), m_dims(dims) 361 { } 362 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE 363 TensorReductionOp(const XprType& expr, const Dims& dims, const Op& reducer) : m_expr(expr), m_dims(dims), m_reducer(reducer) 364 { } 365 366 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE 367 const XprType& expression() const { return m_expr; } 368 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE 369 const Dims& dims() const { return m_dims; } 370 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE 371 const Op& reducer() const { return m_reducer; } 372 373 protected: 374 typename XprType::Nested m_expr; 375 const Dims m_dims; 376 const Op m_reducer; 377 }; 378 379 380 // Eval as rvalue 381 template<typename Op, typename Dims, typename ArgType, template <class> class MakePointer_, typename Device> 382 struct TensorEvaluator<const TensorReductionOp<Op, Dims, ArgType, MakePointer_>, Device> 383 { 384 typedef TensorReductionOp<Op, Dims, ArgType, MakePointer_> XprType; 385 typedef typename XprType::Index Index; 386 typedef ArgType ChildType; 387 typedef typename TensorEvaluator<ArgType, Device>::Dimensions InputDimensions; 388 static const int NumInputDims = internal::array_size<InputDimensions>::value; 389 static const int NumReducedDims = internal::array_size<Dims>::value; 390 static const int NumOutputDims = NumInputDims - NumReducedDims; 391 typedef typename internal::conditional<NumOutputDims==0, Sizes<>, DSizes<Index, NumOutputDims> >::type Dimensions; 392 typedef typename XprType::Scalar Scalar; 393 typedef TensorEvaluator<const TensorReductionOp<Op, Dims, ArgType, MakePointer_>, Device> Self; 394 static const bool InputPacketAccess = TensorEvaluator<ArgType, Device>::PacketAccess; 395 typedef typename internal::remove_const<typename XprType::CoeffReturnType>::type CoeffReturnType; 396 typedef typename PacketType<CoeffReturnType, Device>::type PacketReturnType; 397 static const int PacketSize = internal::unpacket_traits<PacketReturnType>::size; 398 399 enum { 400 IsAligned = false, 401 PacketAccess = Self::InputPacketAccess && Op::PacketAccess, 402 Layout = TensorEvaluator<ArgType, Device>::Layout, 403 CoordAccess = false, // to be implemented 404 RawAccess = false 405 }; 406 407 static const bool ReducingInnerMostDims = internal::are_inner_most_dims<Dims, NumInputDims, Layout>::value; 408 static const bool PreservingInnerMostDims = internal::preserve_inner_most_dims<Dims, NumInputDims, Layout>::value; 409 static const bool RunningFullReduction = (NumOutputDims==0); 410 411 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE TensorEvaluator(const XprType& op, const Device& device) 412 : m_impl(op.expression(), device), m_reducer(op.reducer()), m_result(NULL), m_device(device), m_xpr_dims(op.dims()) 413 { 414 EIGEN_STATIC_ASSERT((NumInputDims >= NumReducedDims), YOU_MADE_A_PROGRAMMING_MISTAKE); 415 EIGEN_STATIC_ASSERT((!ReducingInnerMostDims | !PreservingInnerMostDims | (NumReducedDims == NumInputDims)), 416 YOU_MADE_A_PROGRAMMING_MISTAKE); 417 418 // Build the bitmap indicating if an input dimension is reduced or not. 419 for (int i = 0; i < NumInputDims; ++i) { 420 m_reduced[i] = false; 421 } 422 for (int i = 0; i < NumReducedDims; ++i) { 423 eigen_assert(op.dims()[i] >= 0); 424 eigen_assert(op.dims()[i] < NumInputDims); 425 m_reduced[op.dims()[i]] = true; 426 } 427 428 const typename TensorEvaluator<ArgType, Device>::Dimensions& input_dims = m_impl.dimensions(); 429 internal::DimInitializer<Dimensions>::run(input_dims, m_reduced, &m_dimensions, &m_reducedDims); 430 431 // Precompute output strides. 432 if (NumOutputDims > 0) { 433 if (static_cast<int>(Layout) == static_cast<int>(ColMajor)) { 434 m_outputStrides[0] = 1; 435 for (int i = 1; i < NumOutputDims; ++i) { 436 m_outputStrides[i] = m_outputStrides[i - 1] * m_dimensions[i - 1]; 437 } 438 } else { 439 m_outputStrides.back() = 1; 440 for (int i = NumOutputDims - 2; i >= 0; --i) { 441 m_outputStrides[i] = m_outputStrides[i + 1] * m_dimensions[i + 1]; 442 } 443 } 444 } 445 446 // Precompute input strides. 447 if (NumInputDims > 0) { 448 array<Index, NumInputDims> input_strides; 449 if (static_cast<int>(Layout) == static_cast<int>(ColMajor)) { 450 input_strides[0] = 1; 451 for (int i = 1; i < NumInputDims; ++i) { 452 input_strides[i] = input_strides[i-1] * input_dims[i-1]; 453 } 454 } else { 455 input_strides.back() = 1; 456 for (int i = NumInputDims - 2; i >= 0; --i) { 457 input_strides[i] = input_strides[i + 1] * input_dims[i + 1]; 458 } 459 } 460 461 int outputIndex = 0; 462 int reduceIndex = 0; 463 for (int i = 0; i < NumInputDims; ++i) { 464 if (m_reduced[i]) { 465 m_reducedStrides[reduceIndex] = input_strides[i]; 466 ++reduceIndex; 467 } else { 468 m_preservedStrides[outputIndex] = input_strides[i]; 469 ++outputIndex; 470 } 471 } 472 } 473 474 // Special case for full reductions 475 if (NumOutputDims == 0) { 476 m_preservedStrides[0] = internal::array_prod(input_dims); 477 } 478 } 479 480 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Dimensions& dimensions() const { return m_dimensions; } 481 482 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool evalSubExprsIfNeeded(typename MakePointer_<CoeffReturnType>::Type data) { 483 m_impl.evalSubExprsIfNeeded(NULL); 484 485 // Use the FullReducer if possible. 486 if ((RunningFullReduction && RunningOnSycl) ||(RunningFullReduction && 487 internal::FullReducer<Self, Op, Device>::HasOptimizedImplementation && 488 ((RunningOnGPU && (m_device.majorDeviceVersion() >= 3)) || 489 !RunningOnGPU))) { 490 bool need_assign = false; 491 if (!data) { 492 m_result = static_cast<CoeffReturnType*>(m_device.allocate(sizeof(CoeffReturnType))); 493 data = m_result; 494 need_assign = true; 495 } 496 Op reducer(m_reducer); 497 internal::FullReducer<Self, Op, Device>::run(*this, reducer, m_device, data); 498 return need_assign; 499 } 500 else if(RunningOnSycl){ 501 const Index num_values_to_reduce = internal::array_prod(m_reducedDims); 502 const Index num_coeffs_to_preserve = internal::array_prod(m_dimensions); 503 if (!data) { 504 data = static_cast<CoeffReturnType*>(m_device.allocate(sizeof(CoeffReturnType) * num_coeffs_to_preserve)); 505 m_result = data; 506 } 507 Op reducer(m_reducer); 508 internal::InnerReducer<Self, Op, Device>::run(*this, reducer, m_device, data, num_values_to_reduce, num_coeffs_to_preserve); 509 return (m_result != NULL); 510 } 511 512 // Attempt to use an optimized reduction. 513 else if (RunningOnGPU && (m_device.majorDeviceVersion() >= 3)) { 514 bool reducing_inner_dims = true; 515 for (int i = 0; i < NumReducedDims; ++i) { 516 if (static_cast<int>(Layout) == static_cast<int>(ColMajor)) { 517 reducing_inner_dims &= m_reduced[i]; 518 } else { 519 reducing_inner_dims &= m_reduced[NumInputDims - 1 - i]; 520 } 521 } 522 if (internal::InnerReducer<Self, Op, Device>::HasOptimizedImplementation && 523 (reducing_inner_dims || ReducingInnerMostDims)) { 524 const Index num_values_to_reduce = internal::array_prod(m_reducedDims); 525 const Index num_coeffs_to_preserve = internal::array_prod(m_dimensions); 526 if (!data) { 527 if (num_coeffs_to_preserve < 1024 && num_values_to_reduce > num_coeffs_to_preserve && num_values_to_reduce > 128) { 528 data = static_cast<CoeffReturnType*>(m_device.allocate(sizeof(CoeffReturnType) * num_coeffs_to_preserve)); 529 m_result = data; 530 } 531 else { 532 return true; 533 } 534 } 535 Op reducer(m_reducer); 536 if (internal::InnerReducer<Self, Op, Device>::run(*this, reducer, m_device, data, num_values_to_reduce, num_coeffs_to_preserve)) { 537 if (m_result) { 538 m_device.deallocate(m_result); 539 m_result = NULL; 540 } 541 return true; 542 } else { 543 return (m_result != NULL); 544 } 545 } 546 547 bool preserving_inner_dims = true; 548 for (int i = 0; i < NumReducedDims; ++i) { 549 if (static_cast<int>(Layout) == static_cast<int>(ColMajor)) { 550 preserving_inner_dims &= m_reduced[NumInputDims - 1 - i]; 551 } else { 552 preserving_inner_dims &= m_reduced[i]; 553 } 554 } 555 if (internal::OuterReducer<Self, Op, Device>::HasOptimizedImplementation && 556 preserving_inner_dims) { 557 const Index num_values_to_reduce = internal::array_prod(m_reducedDims); 558 const Index num_coeffs_to_preserve = internal::array_prod(m_dimensions); 559 if (!data) { 560 if (num_coeffs_to_preserve < 1024 && num_values_to_reduce > num_coeffs_to_preserve && num_values_to_reduce > 32) { 561 data = static_cast<CoeffReturnType*>(m_device.allocate(sizeof(CoeffReturnType) * num_coeffs_to_preserve)); 562 m_result = data; 563 } 564 else { 565 return true; 566 } 567 } 568 Op reducer(m_reducer); 569 if (internal::OuterReducer<Self, Op, Device>::run(*this, reducer, m_device, data, num_values_to_reduce, num_coeffs_to_preserve)) { 570 if (m_result) { 571 m_device.deallocate(m_result); 572 m_result = NULL; 573 } 574 return true; 575 } else { 576 return (m_result != NULL); 577 } 578 } 579 } 580 return true; 581 } 582 583 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void cleanup() { 584 m_impl.cleanup(); 585 if (m_result) { 586 m_device.deallocate(m_result); 587 m_result = NULL; 588 } 589 } 590 591 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE CoeffReturnType coeff(Index index) const 592 { 593 if ((RunningOnSycl || RunningFullReduction || RunningOnGPU) && m_result) { 594 return *(m_result + index); 595 } 596 Op reducer(m_reducer); 597 if (ReducingInnerMostDims || RunningFullReduction) { 598 const Index num_values_to_reduce = 599 (static_cast<int>(Layout) == static_cast<int>(ColMajor)) ? m_preservedStrides[0] : m_preservedStrides[NumPreservedStrides - 1]; 600 return internal::InnerMostDimReducer<Self, Op>::reduce(*this, firstInput(index), 601 num_values_to_reduce, reducer); 602 } else { 603 typename Self::CoeffReturnType accum = reducer.initialize(); 604 internal::GenericDimReducer<NumReducedDims-1, Self, Op>::reduce(*this, firstInput(index), reducer, &accum); 605 return reducer.finalize(accum); 606 } 607 } 608 609 // TODO(bsteiner): provide a more efficient implementation. 610 template<int LoadMode> 611 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE PacketReturnType packet(Index index) const 612 { 613 EIGEN_STATIC_ASSERT((PacketSize > 1), YOU_MADE_A_PROGRAMMING_MISTAKE) 614 eigen_assert(index + PacketSize - 1 < Index(internal::array_prod(dimensions()))); 615 616 if (RunningOnGPU && m_result) { 617 return internal::pload<PacketReturnType>(m_result + index); 618 } 619 620 EIGEN_ALIGN_MAX typename internal::remove_const<CoeffReturnType>::type values[PacketSize]; 621 if (ReducingInnerMostDims) { 622 const Index num_values_to_reduce = 623 (static_cast<int>(Layout) == static_cast<int>(ColMajor)) ? m_preservedStrides[0] : m_preservedStrides[NumPreservedStrides - 1]; 624 const Index firstIndex = firstInput(index); 625 for (Index i = 0; i < PacketSize; ++i) { 626 Op reducer(m_reducer); 627 values[i] = internal::InnerMostDimReducer<Self, Op>::reduce(*this, firstIndex + i * num_values_to_reduce, 628 num_values_to_reduce, reducer); 629 } 630 } else if (PreservingInnerMostDims) { 631 const Index firstIndex = firstInput(index); 632 const int innermost_dim = (static_cast<int>(Layout) == static_cast<int>(ColMajor)) ? 0 : NumOutputDims - 1; 633 // TBD: extend this the the n innermost dimensions that we preserve. 634 if (((firstIndex % m_dimensions[innermost_dim]) + PacketSize - 1) < m_dimensions[innermost_dim]) { 635 Op reducer(m_reducer); 636 typename Self::PacketReturnType accum = reducer.template initializePacket<typename Self::PacketReturnType>(); 637 internal::InnerMostDimPreserver<NumReducedDims-1, Self, Op>::reduce(*this, firstIndex, reducer, &accum); 638 return reducer.finalizePacket(accum); 639 } else { 640 for (int i = 0; i < PacketSize; ++i) { 641 values[i] = coeff(index + i); 642 } 643 } 644 } else { 645 for (int i = 0; i < PacketSize; ++i) { 646 values[i] = coeff(index + i); 647 } 648 } 649 PacketReturnType rslt = internal::pload<PacketReturnType>(values); 650 return rslt; 651 } 652 653 // Must be called after evalSubExprsIfNeeded(). 654 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE TensorOpCost costPerCoeff(bool vectorized) const { 655 if (RunningFullReduction && m_result) { 656 return TensorOpCost(sizeof(CoeffReturnType), 0, 0, vectorized, PacketSize); 657 } else { 658 const Index num_values_to_reduce = internal::array_prod(m_reducedDims); 659 const double compute_cost = num_values_to_reduce * internal::functor_traits<Op>::Cost; 660 return m_impl.costPerCoeff(vectorized) * num_values_to_reduce + 661 TensorOpCost(0, 0, compute_cost, vectorized, PacketSize); 662 } 663 } 664 665 EIGEN_DEVICE_FUNC typename MakePointer_<Scalar>::Type data() const { return m_result; } 666 /// required by sycl in order to extract the accessor 667 const TensorEvaluator<ArgType, Device>& impl() const { return m_impl; } 668 /// added for sycl in order to construct the buffer from the sycl device 669 const Device& device() const{return m_device;} 670 /// added for sycl in order to re-construct the reduction eval on the device for the sub-kernel 671 const Dims& xprDims() const {return m_xpr_dims;} 672 673 674 private: 675 template <int, typename, typename> friend struct internal::GenericDimReducer; 676 template <typename, typename, bool> friend struct internal::InnerMostDimReducer; 677 template <int, typename, typename, bool> friend struct internal::InnerMostDimPreserver; 678 template <typename S, typename O, typename D, bool V> friend struct internal::FullReducer; 679 #ifdef EIGEN_USE_THREADS 680 template <typename S, typename O, bool V> friend struct internal::FullReducerShard; 681 #endif 682 #if defined(EIGEN_USE_GPU) && defined(__CUDACC__) 683 template <int B, int N, typename S, typename R, typename I> friend void internal::FullReductionKernel(R, const S, I, typename S::CoeffReturnType*, unsigned int*); 684 #ifdef EIGEN_HAS_CUDA_FP16 685 template <typename S, typename R, typename I> friend void internal::ReductionInitFullReduxKernelHalfFloat(R, const S, I, half2*); 686 template <int B, int N, typename S, typename R, typename I> friend void internal::FullReductionKernelHalfFloat(R, const S, I, half*, half2*); 687 template <int NPT, typename S, typename R, typename I> friend void internal::InnerReductionKernelHalfFloat(R, const S, I, I, half*); 688 #endif 689 template <int NPT, typename S, typename R, typename I> friend void internal::InnerReductionKernel(R, const S, I, I, typename S::CoeffReturnType*); 690 691 template <int NPT, typename S, typename R, typename I> friend void internal::OuterReductionKernel(R, const S, I, I, typename S::CoeffReturnType*); 692 #endif 693 694 template <typename S, typename O, typename D> friend struct internal::InnerReducer; 695 696 // Returns the Index in the input tensor of the first value that needs to be 697 // used to compute the reduction at output index "index". 698 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Index firstInput(Index index) const { 699 if (ReducingInnerMostDims) { 700 if (static_cast<int>(Layout) == static_cast<int>(ColMajor)) { 701 return index * m_preservedStrides[0]; 702 } else { 703 return index * m_preservedStrides[NumPreservedStrides - 1]; 704 } 705 } 706 // TBD: optimize the case where we preserve the innermost dimensions. 707 Index startInput = 0; 708 if (static_cast<int>(Layout) == static_cast<int>(ColMajor)) { 709 for (int i = NumOutputDims - 1; i > 0; --i) { 710 // This is index_i in the output tensor. 711 const Index idx = index / m_outputStrides[i]; 712 startInput += idx * m_preservedStrides[i]; 713 index -= idx * m_outputStrides[i]; 714 } 715 if (PreservingInnerMostDims) { 716 eigen_assert(m_preservedStrides[0] == 1); 717 startInput += index; 718 } else { 719 startInput += index * m_preservedStrides[0]; 720 } 721 } else { 722 for (int i = 0; i < NumOutputDims - 1; ++i) { 723 // This is index_i in the output tensor. 724 const Index idx = index / m_outputStrides[i]; 725 startInput += idx * m_preservedStrides[i]; 726 index -= idx * m_outputStrides[i]; 727 } 728 if (PreservingInnerMostDims) { 729 eigen_assert(m_preservedStrides[NumPreservedStrides - 1] == 1); 730 startInput += index; 731 } else { 732 startInput += index * m_preservedStrides[NumPreservedStrides - 1]; 733 } 734 } 735 return startInput; 736 } 737 738 // Bitmap indicating if an input dimension is reduced or not. 739 array<bool, NumInputDims> m_reduced; 740 // Dimensions of the output of the operation. 741 Dimensions m_dimensions; 742 // Precomputed strides for the output tensor. 743 array<Index, NumOutputDims> m_outputStrides; 744 // Subset of strides of the input tensor for the non-reduced dimensions. 745 // Indexed by output dimensions. 746 static const int NumPreservedStrides = max_n_1<NumOutputDims>::size; 747 array<Index, NumPreservedStrides> m_preservedStrides; 748 749 // Subset of strides of the input tensor for the reduced dimensions. 750 // Indexed by reduced dimensions. 751 array<Index, NumReducedDims> m_reducedStrides; 752 // Size of the input dimensions that are reduced. 753 // Indexed by reduced dimensions. 754 array<Index, NumReducedDims> m_reducedDims; 755 756 // Evaluator for the input expression. 757 TensorEvaluator<ArgType, Device> m_impl; 758 759 // Operation to apply for computing the reduction. 760 Op m_reducer; 761 762 // For full reductions 763 #if defined(EIGEN_USE_GPU) && defined(__CUDACC__) 764 static const bool RunningOnGPU = internal::is_same<Device, Eigen::GpuDevice>::value; 765 static const bool RunningOnSycl = false; 766 #elif defined(EIGEN_USE_SYCL) 767 static const bool RunningOnSycl = internal::is_same<typename internal::remove_all<Device>::type, Eigen::SyclDevice>::value; 768 static const bool RunningOnGPU = false; 769 #else 770 static const bool RunningOnGPU = false; 771 static const bool RunningOnSycl = false; 772 #endif 773 typename MakePointer_<CoeffReturnType>::Type m_result; 774 775 const Device& m_device; 776 const Dims& m_xpr_dims; 777 }; 778 779 } // end namespace Eigen 780 781 #endif // EIGEN_CXX11_TENSOR_TENSOR_REDUCTION_H 782