From 6fa6cdd2b988da98cbdd2b1a5fd2fd3b9d56a4b1 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Wed, 4 Jun 2014 09:21:48 -0700 Subject: [PATCH] Added support for tensor contractions Updated expression evaluation mechanism to also compute the size of the tensor result Misc fixes and improvements. --- unsupported/Eigen/CXX11/Tensor | 1 + .../CXX11/src/Core/util/EmulateCXX11Meta.h | 2 + unsupported/Eigen/CXX11/src/Tensor/Tensor.h | 38 ++- .../Eigen/CXX11/src/Tensor/TensorAssign.h | 5 +- .../Eigen/CXX11/src/Tensor/TensorBase.h | 50 ++-- .../CXX11/src/Tensor/TensorContraction.h | 229 ++++++++++++++++++ .../Eigen/CXX11/src/Tensor/TensorDevice.h | 2 +- .../Eigen/CXX11/src/Tensor/TensorDeviceType.h | 15 +- .../Eigen/CXX11/src/Tensor/TensorDimensions.h | 29 ++- .../Eigen/CXX11/src/Tensor/TensorEvaluator.h | 44 +++- .../Eigen/CXX11/src/Tensor/TensorExpr.h | 36 +-- .../Eigen/CXX11/src/Tensor/TensorFixedSize.h | 2 +- .../src/Tensor/TensorForwardDeclarations.h | 2 + .../Eigen/CXX11/src/Tensor/TensorStorage.h | 11 +- 14 files changed, 370 insertions(+), 96 deletions(-) create mode 100644 unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h diff --git a/unsupported/Eigen/CXX11/Tensor b/unsupported/Eigen/CXX11/Tensor index 323d9edff..d4e8d3a15 100644 --- a/unsupported/Eigen/CXX11/Tensor +++ b/unsupported/Eigen/CXX11/Tensor @@ -39,6 +39,7 @@ #include "unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h" #include "unsupported/Eigen/CXX11/src/Tensor/TensorExpr.h" +#include "unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h" #include "unsupported/Eigen/CXX11/src/Tensor/TensorAssign.h" #include "unsupported/Eigen/CXX11/src/Tensor/TensorDevice.h" diff --git a/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h b/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h index ab869177c..636063f9e 100644 --- a/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h +++ b/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h @@ -23,6 +23,8 @@ template class array { EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const T& operator[] (size_t index) const { return values[index]; } + static const std::size_t size = n; + T values[n]; EIGEN_DEVICE_FUNC diff --git a/unsupported/Eigen/CXX11/src/Tensor/Tensor.h b/unsupported/Eigen/CXX11/src/Tensor/Tensor.h index d8ff3f584..e034f8c03 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/Tensor.h +++ b/unsupported/Eigen/CXX11/src/Tensor/Tensor.h @@ -81,7 +81,7 @@ class Tensor : public TensorBase > typedef typename Base::PacketReturnType PacketReturnType; enum { - IsAligned = bool(EIGEN_ALIGN), + IsAligned = bool(EIGEN_ALIGN) & !(Options_&DontAlign), PacketAccess = true, }; @@ -94,11 +94,11 @@ class Tensor : public TensorBase > TensorStorage m_storage; public: - EIGEN_STRONG_INLINE Index dimension(std::size_t n) const { return m_storage.dimensions()[n]; } - EIGEN_STRONG_INLINE const DSizes& dimensions() const { return m_storage.dimensions(); } - EIGEN_STRONG_INLINE Index size() const { return m_storage.size(); } - EIGEN_STRONG_INLINE Scalar *data() { return m_storage.data(); } - EIGEN_STRONG_INLINE const Scalar *data() const { return m_storage.data(); } + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Index dimension(std::size_t n) const { return m_storage.dimensions()[n]; } + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const DSizes& dimensions() const { return m_storage.dimensions(); } + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Index size() const { return m_storage.size(); } + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar *data() { return m_storage.data(); } + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar *data() const { return m_storage.data(); } // This makes EIGEN_INITIALIZE_COEFFS_IF_THAT_OPTION_IS_ENABLED // work, because that uses base().coeffRef() - and we don't yet @@ -116,13 +116,13 @@ class Tensor : public TensorBase > } #endif - inline const Scalar& coeff(const array& indices) const + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar& coeff(const array& indices) const { eigen_internal_assert(checkIndexRange(indices)); return m_storage.data()[linearizedIndex(indices)]; } - inline const Scalar& coeff(Index index) const + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar& coeff(Index index) const { eigen_internal_assert(index >= 0 && index < size()); return m_storage.data()[index]; @@ -138,13 +138,13 @@ class Tensor : public TensorBase > } #endif - inline Scalar& coeffRef(const array& indices) + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar& coeffRef(const array& indices) { eigen_internal_assert(checkIndexRange(indices)); return m_storage.data()[linearizedIndex(indices)]; } - inline Scalar& coeffRef(Index index) + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar& coeffRef(Index index) { eigen_internal_assert(index >= 0 && index < size()); return m_storage.data()[index]; @@ -160,19 +160,19 @@ class Tensor : public TensorBase > } #endif - inline const Scalar& operator()(const array& indices) const + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar& operator()(const array& indices) const { eigen_assert(checkIndexRange(indices)); return coeff(indices); } - inline const Scalar& operator()(Index index) const + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar& operator()(Index index) const { eigen_internal_assert(index >= 0 && index < size()); return coeff(index); } - inline const Scalar& operator[](Index index) const + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar& operator[](Index index) const { // The bracket operator is only for vectors, use the parenthesis operator instead. EIGEN_STATIC_ASSERT(NumIndices == 1, YOU_MADE_A_PROGRAMMING_MISTAKE); @@ -189,19 +189,19 @@ class Tensor : public TensorBase > } #endif - inline Scalar& operator()(const array& indices) + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar& operator()(const array& indices) { eigen_assert(checkIndexRange(indices)); return coeffRef(indices); } - inline Scalar& operator()(Index index) + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar& operator()(Index index) { eigen_assert(index >= 0 && index < size()); return coeffRef(index); } - inline Scalar& operator[](Index index) + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar& operator[](Index index) { // The bracket operator is only for vectors, use the parenthesis operator instead EIGEN_STATIC_ASSERT(NumIndices == 1, YOU_MADE_A_PROGRAMMING_MISTAKE) @@ -223,11 +223,10 @@ class Tensor : public TensorBase > #ifdef EIGEN_HAS_VARIADIC_TEMPLATES template inline Tensor(Index firstDimension, IndexTypes... otherDimensions) - : m_storage() + : m_storage(internal::array_prod(array{{firstDimension, otherDimensions...}}), array{{firstDimension, otherDimensions...}}) { // The number of dimensions used to construct a tensor must be equal to the rank of the tensor. EIGEN_STATIC_ASSERT(sizeof...(otherDimensions) + 1 == NumIndices, YOU_MADE_A_PROGRAMMING_MISTAKE) - resize(array{{firstDimension, otherDimensions...}}); } #endif @@ -237,7 +236,6 @@ class Tensor : public TensorBase > EIGEN_INITIALIZE_COEFFS_IF_THAT_OPTION_IS_ENABLED } - template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Tensor& operator=(const OtherDerived& other) @@ -306,7 +304,7 @@ class Tensor : public TensorBase > array_zip_and_reduce(indices, m_storage.dimensions()); } - inline Index linearizedIndex(const array& indices) const + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Index linearizedIndex(const array& indices) const { if (Options&RowMajor) { return m_storage.dimensions().IndexOfRowMajor(indices); diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorAssign.h b/unsupported/Eigen/CXX11/src/Tensor/TensorAssign.h index e69ff6188..da1eb62cb 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorAssign.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorAssign.h @@ -53,7 +53,6 @@ template struct TensorAssign { typedef typename Derived1::Index Index; - EIGEN_DEVICE_FUNC static inline void run(Derived1& dst, const Derived2& src) { TensorEvaluator evalDst(dst); @@ -63,7 +62,7 @@ struct TensorAssign static const int LhsStoreMode = TensorEvaluator::IsAligned ? Aligned : Unaligned; static const int RhsLoadMode = TensorEvaluator::IsAligned ? Aligned : Unaligned; static const int PacketSize = unpacket_traits::PacketReturnType>::size; - static const int VectorizedSize = (size / PacketSize) * PacketSize; + const int VectorizedSize = (size / PacketSize) * PacketSize; for (Index i = 0; i < VectorizedSize; i += PacketSize) { evalDst.template writePacket(i, evalSrc.template packet(i)); @@ -148,7 +147,7 @@ struct TensorAssignMultiThreaded // GPU: the evaluation of the expressions is offloaded to a GPU. -#ifdef EIGEN_USE_GPU +#if defined(EIGEN_USE_GPU) && defined(__CUDACC__) template __global__ void EigenMetaKernelNoCheck(LhsEvaluator evalDst, const RhsEvaluator evalSrc) { const int index = blockIdx.x * blockDim.x + threadIdx.x; diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h b/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h index 8a88ba806..c5c711313 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h @@ -30,13 +30,16 @@ class TensorBase typedef Scalar CoeffReturnType; typedef typename internal::packet_traits::type PacketReturnType; - Derived& setZero() { + EIGEN_DEVICE_FUNC + EIGEN_STRONG_INLINE Derived& setZero() { return setConstant(Scalar(0)); } - Derived& setConstant(const Scalar& val) { + EIGEN_DEVICE_FUNC + EIGEN_STRONG_INLINE Derived& setConstant(const Scalar& val) { return derived() = constant(val); } - Derived& setRandom() { + EIGEN_DEVICE_FUNC + EIGEN_STRONG_INLINE Derived& setRandom() { return derived() = random(); } @@ -45,13 +48,13 @@ class TensorBase EIGEN_STRONG_INLINE const TensorCwiseNullaryOp, const Derived> constant(const Scalar& value) const { return TensorCwiseNullaryOp, const Derived> - (internal::scalar_constant_op(value)); + (derived(), internal::scalar_constant_op(value)); } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseNullaryOp, const Derived> random() const { - return TensorCwiseNullaryOp, const Derived>(); + return TensorCwiseNullaryOp, const Derived>(derived()); } // Coefficient-wise unary operators @@ -124,77 +127,86 @@ class TensorBase // Coefficient-wise binary operators. template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator+(const OtherDerived& other) const { + operator+(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator-(const OtherDerived& other) const { + operator-(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator*(const OtherDerived& other) const { + operator*(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator/(const OtherDerived& other) const { + operator/(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - cwiseMax(const OtherDerived& other) const { + cwiseMax(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - cwiseMin(const OtherDerived& other) const { + cwiseMin(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } // Comparisons and tests. template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator<(const OtherDerived& other) const { + operator<(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator<=(const OtherDerived& other) const { + operator<=(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator>(const OtherDerived& other) const { + operator>(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator>=(const OtherDerived& other) const { + operator>=(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator==(const OtherDerived& other) const { + operator==(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseBinaryOp, const Derived, const OtherDerived> - operator!=(const OtherDerived& other) const { + operator!=(const OtherDerived& other) const { return TensorCwiseBinaryOp, const Derived, const OtherDerived>(derived(), other.derived()); } + // Contractions. + typedef std::pair DimensionPair; + + template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE + const TensorContractionOp + contract(const OtherDerived& other, const Dimensions& dims) const { + return TensorContractionOp(derived(), other.derived(), dims); + } + // Coefficient-wise ternary operators. - template + template inline const TensorSelectOp - select(const ThenDerived& thenTensor, const ElseDerived& elseTensor) const{ + select(const ThenDerived& thenTensor, const ElseDerived& elseTensor) const { return TensorSelectOp(derived(), thenTensor.derived(), elseTensor.derived()); } diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h b/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h new file mode 100644 index 000000000..d424df36e --- /dev/null +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h @@ -0,0 +1,229 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2014 Benoit Steiner +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#ifndef EIGEN_CXX11_TENSOR_TENSOR_CONTRACTION_H +#define EIGEN_CXX11_TENSOR_TENSOR_CONTRACTION_H + +namespace Eigen { + +/** \class TensorContraction + * \ingroup CXX11_Tensor_Module + * + * \brief Tensor contraction class. + * + * + */ +namespace internal { +template +struct traits > +{ + // Type promotion to handle the case where the types of the lhs and the rhs are different. + typedef typename internal::promote_storage_type::ret Scalar; + typedef typename internal::packet_traits::type Packet; + typedef typename promote_storage_type::StorageKind, + typename traits::StorageKind>::ret StorageKind; + typedef typename promote_index_type::Index, + typename traits::Index>::type Index; + typedef typename LhsXprType::Nested LhsNested; + typedef typename RhsXprType::Nested RhsNested; + typedef typename remove_reference::type _LhsNested; + typedef typename remove_reference::type _RhsNested; +}; + +template +struct eval, Eigen::Dense> +{ + typedef const TensorContractionOp& type; +}; + +template +struct nested, 1, typename eval >::type> +{ + typedef TensorContractionOp type; +}; + +} // end namespace internal + + + +template +class TensorContractionOp : public TensorBase > +{ + public: + typedef typename Eigen::internal::traits::Scalar Scalar; + typedef typename Eigen::internal::traits::Packet Packet; + typedef typename Eigen::NumTraits::Real RealScalar; + typedef typename internal::promote_storage_type::ret CoeffReturnType; + typedef typename internal::promote_storage_type::ret PacketReturnType; + typedef typename Eigen::internal::nested::type Nested; + typedef typename Eigen::internal::traits::StorageKind StorageKind; + typedef typename Eigen::internal::traits::Index Index; + + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE TensorContractionOp(const LhsXprType& lhs, const RhsXprType& rhs, const Indices& dims) + : m_lhs_xpr(lhs), m_rhs_xpr(rhs), m_indices(dims) {} + + EIGEN_DEVICE_FUNC + const Indices& indices() const { return m_indices; } + + /** \returns the nested expressions */ + EIGEN_DEVICE_FUNC + const typename internal::remove_all::type& + lhsExpression() const { return m_lhs_xpr; } + + EIGEN_DEVICE_FUNC + const typename internal::remove_all::type& + rhsExpression() const { return m_rhs_xpr; } + + protected: + typename LhsXprType::Nested m_lhs_xpr; + typename RhsXprType::Nested m_rhs_xpr; + const Indices m_indices; +}; + + +template struct max_n_1 { + static const size_t size = n; +}; +template <> struct max_n_1<0> { + static const size_t size = 1; +}; + + +template +struct TensorEvaluator > +{ + typedef TensorContractionOp XprType; + + static const int NumDims = max_n_1::Dimensions::count + TensorEvaluator::Dimensions::count - 2 * Indices::size>::size; + typedef typename XprType::Index Index; + typedef DSizes Dimensions; + + enum { + IsAligned = TensorEvaluator::IsAligned & TensorEvaluator::IsAligned, + PacketAccess = /*TensorEvaluator::PacketAccess & TensorEvaluator::PacketAccess */ + false, + }; + + TensorEvaluator(const XprType& op) + : m_leftImpl(op.lhsExpression()), m_rightImpl(op.rhsExpression()) + { + Index index = 0; + Index stride = 1; + m_shiftright = 1; + + int skipped = 0; + const typename TensorEvaluator::Dimensions& left_dims = m_leftImpl.dimensions(); + for (int i = 0; i < TensorEvaluator::Dimensions::count; ++i) { + bool skip = false; + for (int j = 0; j < Indices::size; ++j) { + if (op.indices()[j].first == i) { + skip = true; + m_leftOffsets[2*skipped] = stride; + m_leftOffsets[2*skipped+1] = stride * left_dims[i]; + m_stitchsize[skipped] = left_dims[i]; + break; + } + } + if (!skip) { + m_dimensions[index++] = left_dims[i]; + m_shiftright *= left_dims[i]; + } else { + ++skipped; + } + stride *= left_dims[i]; + } + + stride = 1; + skipped = 0; + const typename TensorEvaluator::Dimensions& right_dims = m_rightImpl.dimensions(); + for (int i = 0; i < TensorEvaluator::Dimensions::count; ++i) { + bool skip = false; + for (int j = 0; j < Indices::size; ++j) { + if (op.indices()[j].second == i) { + skip = true; + m_rightOffsets[2*skipped] = stride; + m_rightOffsets[2*skipped+1] = stride * right_dims[i]; + break; + } + } + if (!skip) { + m_dimensions[index++] = right_dims[i]; + } else { + ++skipped; + } + stride *= right_dims[i]; + } + + // Scalar case + if (TensorEvaluator::Dimensions::count + TensorEvaluator::Dimensions::count == 2 * Indices::size) { + m_dimensions[0] = 1; + } + } + + // typedef typename XprType::Index Index; + typedef typename XprType::CoeffReturnType CoeffReturnType; + typedef typename XprType::PacketReturnType PacketReturnType; + + const Dimensions& dimensions() const { return m_dimensions; } + + void evalTo(typename XprType::Scalar* buffer) const { + for (int i = 0; i < dimensions().TotalSize(); ++i) { + buffer[i] += coeff(i); + } + } + + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE CoeffReturnType coeff(Index index) const + { + const Index startLeft = index % m_shiftright; + const Index startRight = index / m_shiftright; + CoeffReturnType result = CoeffReturnType(0); + partialStitch(startLeft, startRight, 0, result); + return result; + } + + /* TODO: vectorization + template + EIGEN_DEVICE_FUNC PacketReturnType packet(Index index) const + { + assert(false); + }*/ + + private: + EIGEN_DEVICE_FUNC void partialStitch(Index startLeft, Index startRight, int StitchIndex, CoeffReturnType& accum) const { + Index firstLeft = (startLeft / m_leftOffsets[2*StitchIndex]) * m_leftOffsets[2*StitchIndex+1] + (startLeft % m_leftOffsets[2*StitchIndex]); + Index firstRight = (startRight / m_rightOffsets[2*StitchIndex]) * m_rightOffsets[2*StitchIndex+1] + (startRight % m_rightOffsets[2*StitchIndex]); + + for (int j = 0; j < m_stitchsize[StitchIndex]; ++j) { + const Index left = firstLeft+j*m_leftOffsets[2*StitchIndex]; + const Index right = firstRight+j*m_rightOffsets[2*StitchIndex]; + if (StitchIndex < Indices::size-1) { + partialStitch(left, right, StitchIndex+1, accum); + } else { + accum += m_leftImpl.coeff(left) * m_rightImpl.coeff(right); + } + } + } + + private: + array m_leftOffsets; + array m_rightOffsets; + array m_stitchsize; + Index m_shiftright; + Dimensions m_dimensions; + TensorEvaluator m_leftImpl; + TensorEvaluator m_rightImpl; +}; + + +} // end namespace Eigen + +#endif // EIGEN_CXX11_TENSOR_TENSOR_CONTRACTION_H diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorDevice.h b/unsupported/Eigen/CXX11/src/Tensor/TensorDevice.h index 71890e187..dbe60a165 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorDevice.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorDevice.h @@ -59,7 +59,7 @@ template class TensorDevice class TensorDevice { public: diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorDeviceType.h b/unsupported/Eigen/CXX11/src/Tensor/TensorDeviceType.h index ded6ca604..d7f5ab7c9 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorDeviceType.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorDeviceType.h @@ -37,17 +37,14 @@ struct ThreadPoolDevice { // GPU offloading #ifdef EIGEN_USE_GPU struct GpuDevice { - // todo: support for multiple gpu; - GpuDevice() { - cudaStreamCreate(&stream_); - } - ~GpuDevice() { - cudaStreamDestroy(stream_); - } - const cudaStream_t& stream() const { return stream_; } + // The cudastream is not owned: the caller is responsible for its initialization and eventual destruction. + GpuDevice(const cudaStream_t* stream) : stream_(stream) { eigen_assert(stream); } + + const cudaStream_t& stream() const { return *stream_; } private: - cudaStream_t stream_; + // TODO: multigpu. + const cudaStream_t* stream_; }; #endif diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorDimensions.h b/unsupported/Eigen/CXX11/src/Tensor/TensorDimensions.h index 43e9d6550..c92b8c679 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorDimensions.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorDimensions.h @@ -35,14 +35,14 @@ namespace Eigen { namespace internal { template struct dget { - static const std::size_t value = internal::get::value; - }; + static const std::size_t value = get::value; +}; template struct fixed_size_tensor_index_linearization_helper { - template + template EIGEN_DEVICE_FUNC static inline Index run(array const& indices, const Dimensions& dimensions) { @@ -55,7 +55,7 @@ struct fixed_size_tensor_index_linearization_helper template struct fixed_size_tensor_index_linearization_helper { - template + template EIGEN_DEVICE_FUNC static inline Index run(array const& indices, const Dimensions&) { @@ -93,11 +93,11 @@ struct Sizes : internal::numeric_list { return *this; } - template + template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t IndexOfColMajor(const array& indices) const { return internal::fixed_size_tensor_index_linearization_helper::run(indices, *static_cast(this)); } - template + template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t IndexOfRowMajor(const array& indices) const { return internal::fixed_size_tensor_index_linearization_helper::run(indices, *static_cast(this)); } @@ -139,11 +139,11 @@ template + template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t IndexOfColMajor(const array& indices) const { return internal::fixed_size_tensor_index_linearization_helper::run(indices, *this); } - template + template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t IndexOfRowMajor(const array& indices) const { return internal::fixed_size_tensor_index_linearization_helper::run(indices, *this); } @@ -180,13 +180,18 @@ struct tensor_index_linearization_helper template struct DSizes : array { typedef array Base; + static const std::size_t count = NumDims; EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t TotalSize() const { return internal::array_prod(*static_cast(this)); } - DSizes() { } - explicit DSizes(const array& a) : Base(a) { } + EIGEN_DEVICE_FUNC DSizes() { + for (int i = 0 ; i < NumDims; ++i) { + (*this)[i] = 0; + } + } + EIGEN_DEVICE_FUNC explicit DSizes(const array& a) : Base(a) { } DSizes& operator = (const array& other) { *static_cast(this) = other; @@ -194,10 +199,10 @@ struct DSizes : array { } // A constexpr would be so much better here - size_t IndexOfColMajor(const array& indices) const { + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t IndexOfColMajor(const array& indices) const { return internal::tensor_index_linearization_helper::run(indices, *static_cast(this)); } - size_t IndexOfRowMajor(const array& indices) const { + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE size_t IndexOfRowMajor(const array& indices) const { return internal::tensor_index_linearization_helper::run(indices, *static_cast(this)); } }; diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h b/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h index e0c0863b7..ab2513cea 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h @@ -21,7 +21,6 @@ namespace Eigen { * * TODO: add support for more types of expressions, in particular expressions * leading to lvalues (slicing, reshaping, etc...) - * TODO: add support for vectorization */ template @@ -32,16 +31,19 @@ struct TensorEvaluator typedef typename Derived::Packet Packet; typedef typename Derived::Scalar CoeffReturnType; typedef typename Derived::Packet PacketReturnType; + typedef typename Derived::Dimensions Dimensions; enum { IsAligned = Derived::IsAligned, PacketAccess = Derived::PacketAccess, }; - TensorEvaluator(Derived& m) - : m_data(const_cast(m.data())) + EIGEN_DEVICE_FUNC TensorEvaluator(Derived& m) + : m_data(const_cast(m.data())), m_dims(m.dimensions()) { } + EIGEN_DEVICE_FUNC const Dimensions& dimensions() const { return m_dims; } + EIGEN_DEVICE_FUNC CoeffReturnType coeff(Index index) const { return m_data[index]; } @@ -64,29 +66,34 @@ struct TensorEvaluator protected: Scalar* m_data; + Dimensions m_dims; }; // -------------------- CwiseNullaryOp -------------------- -template -struct TensorEvaluator > +template +struct TensorEvaluator > { - typedef TensorCwiseNullaryOp XprType; + typedef TensorCwiseNullaryOp XprType; enum { IsAligned = true, PacketAccess = internal::functor_traits::PacketAccess, }; + EIGEN_DEVICE_FUNC TensorEvaluator(const XprType& op) - : m_functor(op.functor()) + : m_functor(op.functor()), m_argImpl(op.nestedExpression()) { } typedef typename XprType::Index Index; typedef typename XprType::CoeffReturnType CoeffReturnType; typedef typename XprType::PacketReturnType PacketReturnType; + typedef typename TensorEvaluator::Dimensions Dimensions; + + EIGEN_DEVICE_FUNC const Dimensions& dimensions() const { return m_argImpl.dimensions(); } EIGEN_DEVICE_FUNC CoeffReturnType coeff(Index index) const { @@ -101,6 +108,7 @@ struct TensorEvaluator > private: const NullaryOp m_functor; + TensorEvaluator m_argImpl; }; @@ -117,7 +125,7 @@ struct TensorEvaluator > PacketAccess = TensorEvaluator::PacketAccess & internal::functor_traits::PacketAccess, }; - TensorEvaluator(const XprType& op) + EIGEN_DEVICE_FUNC TensorEvaluator(const XprType& op) : m_functor(op.functor()), m_argImpl(op.nestedExpression()) { } @@ -125,6 +133,9 @@ struct TensorEvaluator > typedef typename XprType::Index Index; typedef typename XprType::CoeffReturnType CoeffReturnType; typedef typename XprType::PacketReturnType PacketReturnType; + typedef typename TensorEvaluator::Dimensions Dimensions; + + EIGEN_DEVICE_FUNC const Dimensions& dimensions() const { return m_argImpl.dimensions(); } EIGEN_DEVICE_FUNC CoeffReturnType coeff(Index index) const { @@ -156,7 +167,7 @@ struct TensorEvaluator::PacketAccess, }; - TensorEvaluator(const XprType& op) + EIGEN_DEVICE_FUNC TensorEvaluator(const XprType& op) : m_functor(op.functor()), m_leftImpl(op.lhsExpression()), m_rightImpl(op.rhsExpression()) @@ -165,6 +176,13 @@ struct TensorEvaluator::Dimensions Dimensions; + + EIGEN_DEVICE_FUNC const Dimensions& dimensions() const + { + // TODO: use right impl instead if right impl dimensions are known at compile time. + return m_leftImpl.dimensions(); + } EIGEN_DEVICE_FUNC CoeffReturnType coeff(Index index) const { @@ -196,7 +214,7 @@ struct TensorEvaluator TensorEvaluator::PacketAccess*/, }; - TensorEvaluator(const XprType& op) + EIGEN_DEVICE_FUNC TensorEvaluator(const XprType& op) : m_condImpl(op.ifExpression()), m_thenImpl(op.thenExpression()), m_elseImpl(op.elseExpression()) @@ -205,7 +223,13 @@ struct TensorEvaluator typedef typename XprType::Index Index; typedef typename XprType::CoeffReturnType CoeffReturnType; typedef typename XprType::PacketReturnType PacketReturnType; + typedef typename TensorEvaluator::Dimensions Dimensions; + EIGEN_DEVICE_FUNC const Dimensions& dimensions() const + { + // TODO: use then or else impl instead if they happen to be known at compile time. + return m_condImpl.dimensions(); + } EIGEN_DEVICE_FUNC CoeffReturnType coeff(Index index) const { return m_condImpl.coeff(index) ? m_thenImpl.coeff(index) : m_elseImpl.coeff(index); diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorExpr.h b/unsupported/Eigen/CXX11/src/Tensor/TensorExpr.h index 94cfae05c..60908ee94 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorExpr.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorExpr.h @@ -28,13 +28,13 @@ namespace Eigen { * */ namespace internal { -template -struct traits > - : traits +template +struct traits > + : traits { - typedef typename PlainObjectType::Packet Packet; - typedef typename PlainObjectType::Scalar Scalar; - typedef typename PlainObjectType::Nested XprTypeNested; + typedef typename XprType::Packet Packet; + typedef typename XprType::Scalar Scalar; + typedef typename XprType::Nested XprTypeNested; typedef typename remove_reference::type _XprTypeNested; }; @@ -42,27 +42,31 @@ struct traits > -template -class TensorCwiseNullaryOp : public TensorBase > +template +class TensorCwiseNullaryOp : public TensorBase > { public: typedef typename Eigen::internal::traits::Scalar Scalar; typedef typename Eigen::internal::traits::Packet Packet; typedef typename Eigen::NumTraits::Real RealScalar; - typedef typename PlainObjectType::CoeffReturnType CoeffReturnType; - typedef typename PlainObjectType::PacketReturnType PacketReturnType; - typedef TensorCwiseNullaryOp Nested; + typedef typename XprType::CoeffReturnType CoeffReturnType; + typedef typename XprType::PacketReturnType PacketReturnType; + typedef TensorCwiseNullaryOp Nested; typedef typename Eigen::internal::traits::StorageKind StorageKind; typedef typename Eigen::internal::traits::Index Index; - EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE TensorCwiseNullaryOp(const NullaryOp& func = NullaryOp()) - : m_functor(func) {} + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE TensorCwiseNullaryOp(const XprType& xpr, const NullaryOp& func = NullaryOp()) + : m_xpr(xpr), m_functor(func) {} + + EIGEN_DEVICE_FUNC + const typename internal::remove_all::type& + nestedExpression() const { return m_xpr; } EIGEN_DEVICE_FUNC const NullaryOp& functor() const { return m_functor; } protected: - // todo: add tensor dimension to be able to do some sanity checks + typename XprType::Nested m_xpr; const NullaryOp m_functor; }; @@ -71,7 +75,7 @@ class TensorCwiseNullaryOp : public TensorBase struct traits > - : traits + : traits { typedef typename result_of< UnaryOp(typename XprType::Scalar) @@ -207,7 +211,7 @@ class TensorCwiseBinaryOp : public TensorBase struct traits > - : traits + : traits { typedef typename traits::Scalar Scalar; typedef typename internal::packet_traits::type Packet; diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorFixedSize.h b/unsupported/Eigen/CXX11/src/Tensor/TensorFixedSize.h index dcc7ccd65..789c04238 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorFixedSize.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorFixedSize.h @@ -52,7 +52,7 @@ class TensorFixedSize : public TensorBase dimensions() const { return m_storage.dimensions(); } + EIGEN_STRONG_INLINE const Dimensions& dimensions() const { return m_storage.dimensions(); } EIGEN_STRONG_INLINE Index size() const { return m_storage.size(); } EIGEN_STRONG_INLINE Scalar *data() { return m_storage.data(); } EIGEN_STRONG_INLINE const Scalar *data() const { return m_storage.data(); } diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h b/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h index 03ac8d516..239b5cb67 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h @@ -21,6 +21,8 @@ template class TensorCwiseNullaryO template class TensorCwiseUnaryOp; template class TensorCwiseBinaryOp; template class TensorSelectOp; +template class TensorReductionOp; +template class TensorContractionOp; template class TensorDevice; diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h b/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h index 64098343e..c9d6517eb 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h @@ -53,7 +53,7 @@ class TensorStorage EIGEN_STRONG_INLINE const T *data() const { return m_data; } EIGEN_DEVICE_FUNC - EIGEN_STRONG_INLINE const FixedDimensions dimensions() const { return m_dimensions; } + EIGEN_STRONG_INLINE const FixedDimensions& dimensions() const { return m_dimensions; } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DenseIndex size() const { return m_dimensions.TotalSize(); } @@ -111,7 +111,8 @@ class TensorStorage(m_data, internal::array_prod(m_dimensions)); } void swap(Self_& other) { std::swap(m_data,other.m_data); std::swap(m_dimensions,other.m_dimensions); } - const DSizes& dimensions() const {return m_dimensions;} + + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const DSizes& dimensions() const {return m_dimensions;} void conservativeResize(DenseIndex size, const array& nbDimensions) { @@ -132,10 +133,10 @@ class TensorStorage