From 81f9e968fd390b003bda9b1373fefd2c4426c453 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Mon, 6 Jul 2015 13:32:38 -0700 Subject: [PATCH 01/10] Only attempt to use the texture path on GPUs when it's supported by CUDA --- .../Eigen/CXX11/src/Tensor/TensorEvaluator.h | 30 ++++++++++++------- 1 file changed, 20 insertions(+), 10 deletions(-) diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h b/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h index a38af84d5..36718e26f 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h @@ -109,6 +109,24 @@ struct TensorEvaluator const Device& m_device; }; +namespace { +template EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE +T loadConstant(const T* address) { + return *address; +} +// Use the texture cache on CUDA devices whenever possible +#if defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 350 +template <> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE +float loadConstant(const float* address) { + return __ldg(address); +} +template <> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE +double loadConstant(const double* address) { + return __ldg(address); +} +#endif +} + // Default evaluator for rvalues template @@ -150,11 +168,7 @@ struct TensorEvaluator EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE CoeffReturnType coeff(Index index) const { eigen_assert(m_data); -#ifdef __CUDA_ARCH__ - return __ldg(m_data+index); -#else - return m_data[index]; -#endif + return loadConstant(m_data+index); } template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE @@ -167,11 +181,7 @@ struct TensorEvaluator eigen_assert(m_data); const Index index = (static_cast(Layout) == static_cast(ColMajor)) ? m_dims.IndexOfColMajor(coords) : m_dims.IndexOfRowMajor(coords); -#ifdef __CUDA_ARCH__ - return __ldg(m_data+index); -#else - return m_data[index]; -#endif + return loadConstant(m_data+index); } EIGEN_DEVICE_FUNC const Scalar* data() const { return m_data; } From ebdacfc5ea84b5e5f1ca2c33419ac0fdfd71073e Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Mon, 6 Jul 2015 15:03:11 -0700 Subject: [PATCH 02/10] Fixed a compilation warning generated by clang --- .../Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h b/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h index b3bc16bc4..f73ce83f4 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorForwardDeclarations.h @@ -51,9 +51,9 @@ template class TensorForcedEvalOp; template class TensorDevice; template struct TensorEvaluator; -class DefaultDevice; -class ThreadPoolDevice; -class GpuDevice; +struct DefaultDevice; +struct ThreadPoolDevice; +struct GpuDevice; namespace internal { From 0485a2468da77b7b4f8ed3e83605e5919bb7f55e Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Mon, 6 Jul 2015 17:01:51 -0700 Subject: [PATCH 03/10] use Eigen smart_copy instead of std::copy --- unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h b/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h index 2244e40c2..0ae638fb9 100644 --- a/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h +++ b/unsupported/Eigen/CXX11/src/Core/util/EmulateCXX11Meta.h @@ -107,7 +107,7 @@ template class array { #ifdef EIGEN_HAS_VARIADIC_TEMPLATES array(std::initializer_list l) { eigen_assert(l.size() == n); - std::copy(l.begin(), l.end(), values); + internal::smart_copy(l.begin(), l.end(), values); } #endif }; From 3f2101b03b1fb96ef521dce3ae966ac18e90266d Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Mon, 6 Jul 2015 17:02:29 -0700 Subject: [PATCH 04/10] Use numext::swap instead of std::swap --- unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h | 2 +- unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h b/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h index fd2f3abc4..59ae4a2d0 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorContraction.h @@ -662,7 +662,7 @@ struct TensorContractionEvaluatorBase // If the layout is RowMajor, we need to reverse the m_dimensions if (static_cast(Layout) == static_cast(RowMajor)) { for (int i = 0, j = NumDims - 1; i < j; i++, j--) { - std::swap(m_dimensions[i], m_dimensions[j]); + numext::swap(m_dimensions[i], m_dimensions[j]); } } } diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h b/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h index f567b8c03..33670e36e 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorStorage.h @@ -95,7 +95,7 @@ class TensorStorage, Options_> EIGEN_DEVICE_FUNC ~TensorStorage() { internal::conditional_aligned_delete_auto(m_data, internal::array_prod(m_dimensions)); } EIGEN_DEVICE_FUNC void swap(Self& other) - { std::swap(m_data,other.m_data); std::swap(m_dimensions,other.m_dimensions); } + { numext::swap(m_data,other.m_data); numext::swap(m_dimensions,other.m_dimensions); } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Dimensions& dimensions() const {return m_dimensions;} From fa17358c4b2355cfc0fab48b4e1f5422f7fba9a7 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Tue, 7 Jul 2015 17:27:12 +0200 Subject: [PATCH 05/10] Rotation2D: fix slerp to take the shortest path, and add convenient method to get the angle in [-pi,pi] or [0,pi] --- Eigen/src/Geometry/Rotation2D.h | 20 +++++++++++++++++- test/geo_transformations.cpp | 37 ++++++++++++++++++++++++++++++--- 2 files changed, 53 insertions(+), 4 deletions(-) diff --git a/Eigen/src/Geometry/Rotation2D.h b/Eigen/src/Geometry/Rotation2D.h index 4feb3d4d2..65aa83be5 100644 --- a/Eigen/src/Geometry/Rotation2D.h +++ b/Eigen/src/Geometry/Rotation2D.h @@ -69,6 +69,20 @@ public: /** \returns a read-write reference to the rotation angle */ inline Scalar& angle() { return m_angle; } + + /** \returns the rotation angle in [0,2pi] */ + inline Scalar smallestPositiveAngle() const { + Scalar tmp = fmod(m_angle,Scalar(2)*EIGEN_PI); + return tmpScalar(EIGEN_PI)) tmp -= Scalar(2)*Scalar(EIGEN_PI); + else if(tmp<-Scalar(EIGEN_PI)) tmp += Scalar(2)*Scalar(EIGEN_PI); + return tmp; + } /** \returns the inverse rotation */ inline Rotation2D inverse() const { return Rotation2D(-m_angle); } @@ -93,7 +107,10 @@ public: * parameter \a t. It is in fact equivalent to a linear interpolation. */ inline Rotation2D slerp(const Scalar& t, const Rotation2D& other) const - { return Rotation2D(m_angle * (1-t) + other.angle() * t); } + { + Scalar dist = Rotation2D(other.m_angle-m_angle).smallestAngle(); + return Rotation2D(m_angle + dist*t); + } /** \returns \c *this with scalar type casted to \a NewScalarType * @@ -119,6 +136,7 @@ public: * \sa MatrixBase::isApprox() */ bool isApprox(const Rotation2D& other, const typename NumTraits::Real& prec = NumTraits::dummy_precision()) const { return internal::isApprox(m_angle,other.m_angle, prec); } + }; /** \ingroup Geometry_Module diff --git a/test/geo_transformations.cpp b/test/geo_transformations.cpp index c4025f32f..e296267cf 100644 --- a/test/geo_transformations.cpp +++ b/test/geo_transformations.cpp @@ -408,7 +408,24 @@ template void transformations() VERIFY_IS_APPROX(r2d1f.template cast(),r2d1); Rotation2D r2d1d = r2d1.template cast(); VERIFY_IS_APPROX(r2d1d.template cast(),r2d1); + + for(int k=0; k<100; ++k) + { + Scalar angle = internal::random(-100,100); + Rotation2D rot2(angle); + VERIFY( rot2.smallestPositiveAngle() >= 0 ); + VERIFY( rot2.smallestPositiveAngle() < Scalar(2)*Scalar(EIGEN_PI) ); + VERIFY_IS_APPROX( std::cos(rot2.smallestPositiveAngle()), std::cos(rot2.angle()) ); + VERIFY_IS_APPROX( std::sin(rot2.smallestPositiveAngle()), std::sin(rot2.angle()) ); + + VERIFY( rot2.smallestAngle() >= -Scalar(EIGEN_PI) ); + VERIFY( rot2.smallestAngle() <= Scalar(EIGEN_PI) ); + VERIFY_IS_APPROX( std::cos(rot2.smallestAngle()), std::cos(rot2.angle()) ); + VERIFY_IS_APPROX( std::sin(rot2.smallestAngle()), std::sin(rot2.angle()) ); + } + s0 = internal::random(-100,100); + s1 = internal::random(-100,100); Rotation2D R0(s0), R1(s1); t20 = Translation2(v20) * (R0 * Eigen::Scaling(s0)); @@ -420,9 +437,23 @@ template void transformations() VERIFY_IS_APPROX(t20,t21); VERIFY_IS_APPROX(s0, (R0.slerp(0, R1)).angle()); - VERIFY_IS_APPROX(s1, (R0.slerp(1, R1)).angle()); - VERIFY_IS_APPROX(s0, (R0.slerp(0.5, R0)).angle()); - VERIFY_IS_APPROX(Scalar(0), (R0.slerp(0.5, R0.inverse())).angle()); + VERIFY_IS_APPROX(R1.smallestPositiveAngle(), (R0.slerp(1, R1)).smallestPositiveAngle()); + VERIFY_IS_APPROX(R0.smallestPositiveAngle(), (R0.slerp(0.5, R0)).smallestPositiveAngle()); + + if(std::cos(s0)>0) + VERIFY_IS_MUCH_SMALLER_THAN((R0.slerp(0.5, R0.inverse())).smallestAngle(), Scalar(1)); + else + VERIFY_IS_APPROX(Scalar(EIGEN_PI), (R0.slerp(0.5, R0.inverse())).smallestPositiveAngle()); + + // Check path length + Scalar l = 0; + for(int k=0; k<100; ++k) + { + Scalar a1 = R0.slerp(Scalar(k)/Scalar(100), R1).angle(); + Scalar a2 = R0.slerp(Scalar(k+1)/Scalar(100), R1).angle(); + l += std::abs(a2-a1); + } + VERIFY(l<=EIGEN_PI); // check basic features { From 7fa6fe8d8c7f94b37c2b184615e283587a9eb068 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Tue, 7 Jul 2015 17:47:24 +0200 Subject: [PATCH 06/10] typo --- Eigen/src/SparseLU/SparseLU_gemm_kernel.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Eigen/src/SparseLU/SparseLU_gemm_kernel.h b/Eigen/src/SparseLU/SparseLU_gemm_kernel.h index 7420b4d17..1d456ee0c 100644 --- a/Eigen/src/SparseLU/SparseLU_gemm_kernel.h +++ b/Eigen/src/SparseLU/SparseLU_gemm_kernel.h @@ -165,7 +165,7 @@ void sparselu_gemm(Index m, Index n, Index d, const Scalar* A, Index lda, const Bc1 += RK; } // peeled loop on k } // peeled loop on the columns j - // process the last column (we now perform a matrux-vector product) + // process the last column (we now perform a matrix-vector product) if((n-n_end)>0) { const Scalar* Bc0 = B+(n-1)*ldb; From a93af659383002063513099ed35efa9fe177bec8 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Tue, 7 Jul 2015 08:52:14 -0700 Subject: [PATCH 07/10] Improved and cleaned up the 2d patch extraction code --- .../Eigen/CXX11/src/Tensor/TensorBase.h | 39 ++- .../Eigen/CXX11/src/Tensor/TensorImagePatch.h | 244 +++++++++++++----- .../Eigen/CXX11/src/Tensor/TensorIntDiv.h | 2 +- unsupported/test/cxx11_tensor_image_patch.cpp | 32 +-- 4 files changed, 212 insertions(+), 105 deletions(-) diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h b/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h index 1c31b9d28..165bd8a89 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorBase.h @@ -10,6 +10,8 @@ #ifndef EIGEN_CXX11_TENSOR_TENSOR_BASE_H #define EIGEN_CXX11_TENSOR_TENSOR_BASE_H +// clang-format off + namespace Eigen { /** \class TensorBase @@ -379,39 +381,28 @@ class TensorBase return TensorPatchOp(derived(), patch_dims); } - template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE - const TensorImagePatchOp - extract_image_patches() const { - return TensorImagePatchOp(derived(), Rows, Cols, 1, 1, PADDING_SAME); - } - - template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE - const TensorImagePatchOp - extract_image_patches(const PaddingType padding_type) const { - return TensorImagePatchOp(derived(), Rows, Cols, 1, 1, padding_type); - } - - template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE - const TensorImagePatchOp - extract_image_patches(const Index stride, const PaddingType padding_type) const { - return TensorImagePatchOp(derived(), Rows, Cols, stride, stride, padding_type); - } - EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorImagePatchOp - extract_image_patches(const Index patch_rows, const Index patch_cols, - const Index row_stride = 1, const Index col_stride = 1) const { + extract_image_patches(const Index patch_rows = 1, const Index patch_cols = 1, + const Index row_stride = 1, const Index col_stride = 1, + const Index in_row_stride = 1, const Index in_col_stride = 1, + const PaddingType padding_type = PADDING_SAME, const Scalar padding_value = Scalar(0)) const { return TensorImagePatchOp(derived(), patch_rows, patch_cols, row_stride, col_stride, - PADDING_SAME); + in_row_stride, in_col_stride, 1, 1, padding_type, padding_value); } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorImagePatchOp extract_image_patches(const Index patch_rows, const Index patch_cols, const Index row_stride, const Index col_stride, - const PaddingType padding_type) const { + const Index in_row_stride, const Index in_col_stride, + const Index row_inflate_stride, const Index col_inflate_stride, + const Index padding_top, const Index padding_bottom, + const Index padding_left,const Index padding_right, + const Scalar padding_value) const { return TensorImagePatchOp(derived(), patch_rows, patch_cols, row_stride, col_stride, - padding_type); + in_row_stride, in_col_stride, row_inflate_stride, col_inflate_stride, + padding_top, padding_bottom, padding_left, padding_right, padding_value); } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE @@ -481,7 +472,7 @@ class TensorBase return TensorStridingOp(derived(), strides); } - // Added support for custom unary and binary operations + // Support for custom unary and binary operations template EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCustomUnaryOp customOp(const CustomUnaryFunc& op) const { diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorImagePatch.h b/unsupported/Eigen/CXX11/src/Tensor/TensorImagePatch.h index 478696c65..f2c56a9ac 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorImagePatch.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorImagePatch.h @@ -30,7 +30,7 @@ namespace internal { template struct traits > : public traits { - typedef typename XprType::Scalar Scalar; + typedef typename internal::remove_const::type Scalar; typedef traits XprTraits; typedef typename packet_traits::type Packet; typedef typename XprTraits::StorageKind StorageKind; @@ -70,10 +70,30 @@ class TensorImagePatchOp : public TensorBase::type& @@ -96,10 +136,19 @@ class TensorImagePatchOp : public TensorBase struct TensorEvaluator, Device> @@ -109,7 +158,10 @@ struct TensorEvaluator, Device> static const int NumInputDims = internal::array_size::Dimensions>::value; static const int NumDims = NumInputDims + 1; typedef DSizes Dimensions; - typedef typename XprType::Scalar Scalar; + typedef typename internal::remove_const::type Scalar; + typedef TensorEvaluator, + Device> Self; + typedef TensorEvaluator Impl; enum { IsAligned = false, @@ -123,13 +175,17 @@ struct TensorEvaluator, Device> { EIGEN_STATIC_ASSERT(NumDims >= 4, YOU_MADE_A_PROGRAMMING_MISTAKE); + m_paddingValue = op.padding_value(); + const typename TensorEvaluator::Dimensions& input_dims = m_impl.dimensions(); // Caches a few variables. if (static_cast(Layout) == static_cast(ColMajor)) { + m_inputDepth = input_dims[0]; m_inputRows = input_dims[1]; m_inputCols = input_dims[2]; } else { + m_inputDepth = input_dims[NumInputDims-1]; m_inputRows = input_dims[NumInputDims-2]; m_inputCols = input_dims[NumInputDims-3]; } @@ -137,27 +193,57 @@ struct TensorEvaluator, Device> m_row_strides = op.row_strides(); m_col_strides = op.col_strides(); - // We only support same strides for both dimensions and square patches. - eigen_assert(m_row_strides == m_col_strides); + // Input strides and effective input/patch size + m_in_row_strides = op.in_row_strides(); + m_in_col_strides = op.in_col_strides(); + m_row_inflate_strides = op.row_inflate_strides(); + m_col_inflate_strides = op.col_inflate_strides(); + // The "effective" input rows and input cols are the input rows and cols + // after inflating them with zeros. + // For examples, a 2x3 matrix with row_inflate_strides and + // col_inflate_strides of 2 comes from: + // A B C + // D E F + // + // to a matrix is 3 x 5: + // + // A . B . C + // . . . . . + // D . E . F - switch (op.padding_type()) { - case PADDING_VALID: - m_outputRows = numext::ceil((m_inputRows - op.patch_rows() + 1.f) / static_cast(m_row_strides)); - m_outputCols = numext::ceil((m_inputCols - op.patch_cols() + 1.f) / static_cast(m_col_strides)); - // Calculate the padding - m_rowPaddingTop = ((m_outputRows - 1) * m_row_strides + op.patch_rows() - m_inputRows) / 2; - m_colPaddingLeft = ((m_outputCols - 1) * m_col_strides + op.patch_cols() - m_inputCols) / 2; - break; - case PADDING_SAME: - m_outputRows = numext::ceil(m_inputRows / static_cast(m_row_strides)); - m_outputCols = numext::ceil(m_inputCols / static_cast(m_col_strides)); - // Calculate the padding - m_rowPaddingTop = ((m_outputRows - 1) * m_row_strides + op.patch_rows() - m_inputRows) / 2; - m_colPaddingLeft = ((m_outputCols - 1) * m_col_strides + op.patch_cols() - m_inputCols) / 2; - break; - default: - eigen_assert(false && "unexpected padding"); + m_input_rows_eff = (m_inputRows - 1) * m_row_inflate_strides + 1; + m_input_cols_eff = (m_inputCols - 1) * m_col_inflate_strides + 1; + m_patch_rows_eff = op.patch_rows() + (op.patch_rows() - 1) * (m_in_row_strides - 1); + m_patch_cols_eff = op.patch_cols() + (op.patch_cols() - 1) * (m_in_col_strides - 1); + + if (op.padding_explicit()) { + m_outputRows = numext::ceil((m_input_rows_eff + op.padding_top() + op.padding_bottom() - m_patch_rows_eff + 1.f) / static_cast(m_row_strides)); + m_outputCols = numext::ceil((m_input_cols_eff + op.padding_left() + op.padding_right() - m_patch_cols_eff + 1.f) / static_cast(m_col_strides)); + m_rowPaddingTop = op.padding_top(); + m_colPaddingLeft = op.padding_left(); + } else { + // Computing padding from the type + switch (op.padding_type()) { + case PADDING_VALID: + m_outputRows = numext::ceil((m_input_rows_eff - m_patch_rows_eff + 1.f) / static_cast(m_row_strides)); + m_outputCols = numext::ceil((m_input_cols_eff - m_patch_cols_eff + 1.f) / static_cast(m_col_strides)); + // Calculate the padding + m_rowPaddingTop = ((m_outputRows - 1) * m_row_strides + m_patch_rows_eff - m_input_rows_eff) / 2; + m_colPaddingLeft = ((m_outputCols - 1) * m_col_strides + m_patch_cols_eff - m_input_cols_eff) / 2; + break; + case PADDING_SAME: + m_outputRows = numext::ceil(m_input_rows_eff / static_cast(m_row_strides)); + m_outputCols = numext::ceil(m_input_cols_eff / static_cast(m_col_strides)); + // Calculate the padding + m_rowPaddingTop = ((m_outputRows - 1) * m_row_strides + m_patch_rows_eff - m_input_rows_eff) / 2; + m_colPaddingLeft = ((m_outputCols - 1) * m_col_strides + m_patch_cols_eff - m_input_cols_eff) / 2; + break; + default: + eigen_assert(false && "unexpected padding"); } + } + eigen_assert(m_outputRows > 0); + eigen_assert(m_outputCols > 0); // Dimensions for result of extraction. if (static_cast(Layout) == static_cast(ColMajor)) { @@ -202,26 +288,24 @@ struct TensorEvaluator, Device> } // Strides for navigating through the input tensor. - if (static_cast(Layout) == static_cast(ColMajor)) { - m_rowInputStride = input_dims[0]; - m_colInputStride = input_dims[0] * input_dims[1]; - m_patchInputStride = input_dims[0] * input_dims[1] * input_dims[2]; - } else { - m_rowInputStride = input_dims[NumInputDims-1]; - m_colInputStride = input_dims[NumInputDims-1] * input_dims[NumInputDims-2]; - m_patchInputStride = input_dims[NumInputDims-1] * input_dims[NumInputDims-2] * input_dims[NumInputDims-3]; - } + m_rowInputStride = m_inputDepth; + m_colInputStride = m_inputDepth * m_inputRows; + m_patchInputStride = m_inputDepth * m_inputRows * m_inputCols; // Fast representations of different variables. m_fastOtherStride = internal::TensorIntDivisor(m_otherStride); m_fastPatchStride = internal::TensorIntDivisor(m_patchStride); m_fastColStride = internal::TensorIntDivisor(m_colStride); + m_fastInputRowStride = internal::TensorIntDivisor(m_row_inflate_strides); + m_fastInputColStride = internal::TensorIntDivisor(m_col_inflate_strides); + m_fastInputColsEff = internal::TensorIntDivisor(m_input_cols_eff); + // Number of patches in the width dimension. m_fastOutputRows = internal::TensorIntDivisor(m_outputRows); if (static_cast(Layout) == static_cast(ColMajor)) { - m_fastDimZero = internal::TensorIntDivisor(m_dimensions[0]); + m_fastOutputDepth = internal::TensorIntDivisor(m_dimensions[0]); } else { - m_fastDimZero = internal::TensorIntDivisor(m_dimensions[NumDims-1]); + m_fastOutputDepth = internal::TensorIntDivisor(m_dimensions[NumDims-1]); } } @@ -244,33 +328,36 @@ struct TensorEvaluator, Device> // Patch index corresponding to the passed in index. const Index patchIndex = index / m_fastPatchStride; // Find the offset of the element wrt the location of the first element. - const Index patchOffset = (index - patchIndex * m_patchStride) / m_fastDimZero; + const Index patchOffset = (index - patchIndex * m_patchStride) / m_fastOutputDepth; // Other ways to index this element. const Index otherIndex = (NumDims == 4) ? 0 : index / m_fastOtherStride; const Index patch2DIndex = (NumDims == 4) ? patchIndex : (index - otherIndex * m_otherStride) / m_fastPatchStride; + // Calculate col index in the input original tensor. const Index colIndex = patch2DIndex / m_fastOutputRows; const Index colOffset = patchOffset / m_fastColStride; - - // Calculate col index in the input original tensor. - const Index inputCol = colIndex * m_col_strides + colOffset - m_colPaddingLeft; - if (inputCol < 0 || inputCol >= m_inputCols) { - return Scalar(0); + const Index inputCol = colIndex * m_col_strides + colOffset * m_in_col_strides - m_colPaddingLeft; + const Index origInputCol = (m_col_inflate_strides == 1) ? inputCol : ((inputCol >= 0) ? (inputCol / m_fastInputColStride) : 0); + if (inputCol < 0 || inputCol >= m_input_cols_eff || + ((m_col_inflate_strides != 1) && (inputCol != origInputCol * m_col_inflate_strides))) { + return Scalar(m_paddingValue); } - const Index rowIndex = patch2DIndex - colIndex * m_outputRows; - const Index rowOffset = patchOffset - colOffset * m_colStride; // Calculate row index in the original input tensor. - const Index inputRow = rowIndex * m_row_strides + rowOffset - m_rowPaddingTop; - if (inputRow < 0 || inputRow >= m_inputRows) { - return Scalar(0); + const Index rowIndex = patch2DIndex - colIndex * m_outputRows; + const Index rowOffset = patchOffset - colOffset * m_colStride; + const Index inputRow = rowIndex * m_row_strides + rowOffset * m_in_row_strides - m_rowPaddingTop; + const Index origInputRow = (m_row_inflate_strides == 1) ? inputRow : ((inputRow >= 0) ? (inputRow / m_fastInputRowStride) : 0); + if (inputRow < 0 || inputRow >= m_input_rows_eff || + ((m_row_inflate_strides != 1) && (inputRow != origInputRow * m_row_inflate_strides))) { + return Scalar(m_paddingValue); } const int depth_index = static_cast(Layout) == static_cast(ColMajor) ? 0 : NumDims - 1; - const Index depth = index - (index / m_fastDimZero) * m_dimensions[depth_index]; + const Index depth = index - (index / m_fastOutputDepth) * m_dimensions[depth_index]; - const Index inputIndex = depth + inputRow * m_rowInputStride + inputCol * m_colInputStride + otherIndex * m_patchInputStride; + const Index inputIndex = depth + origInputRow * m_rowInputStride + origInputCol * m_colInputStride + otherIndex * m_patchInputStride; return m_impl.coeff(inputIndex); } @@ -281,6 +368,10 @@ struct TensorEvaluator, Device> EIGEN_STATIC_ASSERT(packetSize > 1, YOU_MADE_A_PROGRAMMING_MISTAKE) eigen_assert(index+packetSize-1 < dimensions().TotalSize()); + if (m_in_row_strides != 1 || m_in_col_strides != 1 || m_row_inflate_strides != 1 || m_col_inflate_strides != 1) { + return packetWithPossibleZero(index); + } + const Index indices[2] = {index, index + packetSize - 1}; const Index patchIndex = indices[0] / m_fastPatchStride; if (patchIndex != indices[1] / m_fastPatchStride) { @@ -290,8 +381,8 @@ struct TensorEvaluator, Device> eigen_assert(otherIndex == indices[1] / m_fastOtherStride); // Find the offset of the element wrt the location of the first element. - const Index patchOffsets[2] = {(indices[0] - patchIndex * m_patchStride) / m_fastDimZero, - (indices[1] - patchIndex * m_patchStride) / m_fastDimZero}; + const Index patchOffsets[2] = {(indices[0] - patchIndex * m_patchStride) / m_fastOutputDepth, + (indices[1] - patchIndex * m_patchStride) / m_fastOutputDepth}; const Index patch2DIndex = (NumDims == 4) ? patchIndex : (indices[0] - otherIndex * m_otherStride) / m_fastPatchStride; eigen_assert(patch2DIndex == (indices[1] - otherIndex * m_otherStride) / m_fastPatchStride); @@ -303,8 +394,7 @@ struct TensorEvaluator, Device> const Index inputCols[2] = {colIndex * m_col_strides + colOffsets[0] - m_colPaddingLeft, colIndex * m_col_strides + colOffsets[1] - m_colPaddingLeft}; if (inputCols[1] < 0 || inputCols[0] >= m_inputCols) { - // all zeros - return internal::pset1(Scalar(0)); + return internal::pset1(Scalar(m_paddingValue)); } if (inputCols[0] == inputCols[1]) { @@ -316,14 +406,13 @@ struct TensorEvaluator, Device> m_rowPaddingTop, rowIndex * m_row_strides + rowOffsets[1] - m_rowPaddingTop}; if (inputRows[1] < 0 || inputRows[0] >= m_inputRows) { - // all zeros - return internal::pset1(Scalar(0)); + return internal::pset1(Scalar(m_paddingValue)); } if (inputRows[0] >= 0 && inputRows[1] < m_inputRows) { // no padding const int depth_index = static_cast(Layout) == static_cast(ColMajor) ? 0 : NumDims - 1; - const Index depth = index - (index / m_fastDimZero) * m_dimensions[depth_index]; + const Index depth = index - (index / m_fastOutputDepth) * m_dimensions[depth_index]; const Index inputIndex = depth + inputRows[0] * m_rowInputStride + inputCols[0] * m_colInputStride + otherIndex * m_patchInputStride; return m_impl.template packet(inputIndex); } @@ -342,6 +431,10 @@ struct TensorEvaluator, Device> Index outputCols() const { return m_outputCols; } Index userRowStride() const { return m_row_strides; } Index userColStride() const { return m_col_strides; } + Index userInRowStride() const { return m_in_row_strides; } + Index userInColStride() const { return m_in_col_strides; } + Index rowInflateStride() const { return m_row_inflate_strides; } + Index colInflateStride() const { return m_col_inflate_strides; } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE CoeffReturnType coeff(const array& coords) const { @@ -350,24 +443,30 @@ struct TensorEvaluator, Device> // 0: d, 1: patch_rows, 2: patch_cols, 3: number of patches, 4: number of batches // RowMajor // 0: number of batches, 1: number of patches, 2: patch_cols , 3: patch_rows, 4: d - const Index patchIndex = coords[static_cast(Layout) == static_cast(ColMajor) ? 3 : 1]; + const Index patch2DIndex = coords[static_cast(Layout) == static_cast(ColMajor) ? 3 : 1]; array inputCoords; + Index input_col_idx = patch2DIndex / m_fastInputColsEff; + Index inputCol = input_col_idx + coords[1] * m_in_row_strides - m_rowPaddingTop; + Index inputRow = patch2DIndex - input_col_idx * m_input_cols_eff + coords[2] * m_in_col_strides - m_colPaddingLeft; + const Index origInputCol = (m_col_inflate_strides == 1) ? inputCol : ((inputCol >= 0) ? (inputCol / m_fastInputColStride) : 0); + const Index origInputRow = (m_row_inflate_strides == 1) ? inputRow : ((inputRow >= 0) ? (inputRow / m_fastInputRowStride) : 0); if (static_cast(Layout) == static_cast(ColMajor)) { inputCoords[0] = coords[0]; // depth - inputCoords[1] = patchIndex / m_inputCols + coords[1] - m_rowPaddingTop; - inputCoords[2] = patchIndex - patchIndex / m_inputCols * m_inputCols + coords[2] - m_colPaddingLeft; + inputCoords[1] = origInputCol; + inputCoords[2] = origInputRow; inputCoords[3] = coords[4]; // batch } else { inputCoords[3] = coords[4]; // depth - inputCoords[2] = patchIndex / m_inputCols + coords[3] - m_rowPaddingTop; - inputCoords[1] = patchIndex - patchIndex / m_inputCols * m_inputCols + coords[2] - m_colPaddingLeft; + inputCoords[2] = origInputCol; + inputCoords[1] = origInputRow; inputCoords[0] = coords[0]; // batch } // If the computed coordinates are outside the original image perimeter, return 0. - if (inputCoords[1] < 0 || inputCoords[1] >= m_inputRows || - inputCoords[2] < 0 || inputCoords[2] >= m_inputCols) { - return Scalar(0); + if (inputCol < 0 || inputCol >= m_input_cols_eff || inputRow < 0 || inputRow >= m_input_rows_eff || + ((m_col_inflate_strides != 1) && (inputCol != origInputCol * m_col_inflate_strides)) || + ((m_row_inflate_strides != 1) && (inputRow != origInputRow * m_row_inflate_strides))) { + return Scalar(m_paddingValue); } if (TensorEvaluator::CoordAccess) { return m_impl.coeff(inputCoords); @@ -409,14 +508,29 @@ struct TensorEvaluator, Device> Index m_colStride; Index m_row_strides; Index m_col_strides; + + Index m_in_row_strides; + Index m_in_col_strides; + Index m_row_inflate_strides; + Index m_col_inflate_strides; + + Index m_input_rows_eff; + Index m_input_cols_eff; + Index m_patch_rows_eff; + Index m_patch_cols_eff; + internal::TensorIntDivisor m_fastOtherStride; internal::TensorIntDivisor m_fastPatchStride; internal::TensorIntDivisor m_fastColStride; + internal::TensorIntDivisor m_fastInputRowStride; + internal::TensorIntDivisor m_fastInputColStride; + internal::TensorIntDivisor m_fastInputColsEff; Index m_rowInputStride; Index m_colInputStride; Index m_patchInputStride; + Index m_inputDepth; Index m_inputRows; Index m_inputCols; @@ -427,7 +541,9 @@ struct TensorEvaluator, Device> Index m_colPaddingLeft; internal::TensorIntDivisor m_fastOutputRows; - internal::TensorIntDivisor m_fastDimZero; + internal::TensorIntDivisor m_fastOutputDepth; + + Scalar m_paddingValue; TensorEvaluator m_impl; }; diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorIntDiv.h b/unsupported/Eigen/CXX11/src/Tensor/TensorIntDiv.h index 1de6ce3b4..4c5e784c9 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorIntDiv.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorIntDiv.h @@ -75,7 +75,7 @@ struct TensorIntDivisor { eigen_assert(numerator >= 0); eigen_assert(static_cast(numerator) <= (1ull<> 32; + uint32_t t1 = (multiplier * numerator) >> N; uint32_t t = (static_cast(numerator) - t1) >> shift1; return (t1 + t) >> shift2; } diff --git a/unsupported/test/cxx11_tensor_image_patch.cpp b/unsupported/test/cxx11_tensor_image_patch.cpp index e03e97316..5d6a49181 100644 --- a/unsupported/test/cxx11_tensor_image_patch.cpp +++ b/unsupported/test/cxx11_tensor_image_patch.cpp @@ -25,7 +25,7 @@ static void test_simple_patch() // Single pixel patch: ColMajor Tensor single_pixel_patch; - single_pixel_patch = tensor.extract_image_patches<1, 1>(); + single_pixel_patch = tensor.extract_image_patches(1, 1); VERIFY_IS_EQUAL(single_pixel_patch.dimension(0), 2); VERIFY_IS_EQUAL(single_pixel_patch.dimension(1), 1); VERIFY_IS_EQUAL(single_pixel_patch.dimension(2), 1); @@ -34,7 +34,7 @@ static void test_simple_patch() // Single pixel patch: RowMajor Tensor single_pixel_patch_row_major; - single_pixel_patch_row_major = tensor_row_major.extract_image_patches<1, 1>(); + single_pixel_patch_row_major = tensor_row_major.extract_image_patches(1, 1); VERIFY_IS_EQUAL(single_pixel_patch_row_major.dimension(0), 7); VERIFY_IS_EQUAL(single_pixel_patch_row_major.dimension(1), 3*5); VERIFY_IS_EQUAL(single_pixel_patch_row_major.dimension(2), 1); @@ -64,7 +64,7 @@ static void test_simple_patch() // Entire image patch: ColMajor Tensor entire_image_patch; - entire_image_patch = tensor.extract_image_patches<3, 5>(); + entire_image_patch = tensor.extract_image_patches(3, 5); VERIFY_IS_EQUAL(entire_image_patch.dimension(0), 2); VERIFY_IS_EQUAL(entire_image_patch.dimension(1), 3); VERIFY_IS_EQUAL(entire_image_patch.dimension(2), 5); @@ -73,7 +73,7 @@ static void test_simple_patch() // Entire image patch: RowMajor Tensor entire_image_patch_row_major; - entire_image_patch_row_major = tensor_row_major.extract_image_patches<3, 5>(); + entire_image_patch_row_major = tensor_row_major.extract_image_patches(3, 5); VERIFY_IS_EQUAL(entire_image_patch_row_major.dimension(0), 7); VERIFY_IS_EQUAL(entire_image_patch_row_major.dimension(1), 3*5); VERIFY_IS_EQUAL(entire_image_patch_row_major.dimension(2), 5); @@ -118,7 +118,7 @@ static void test_simple_patch() // 2D patch: ColMajor Tensor twod_patch; - twod_patch = tensor.extract_image_patches<2, 2>(); + twod_patch = tensor.extract_image_patches(2, 2); VERIFY_IS_EQUAL(twod_patch.dimension(0), 2); VERIFY_IS_EQUAL(twod_patch.dimension(1), 2); VERIFY_IS_EQUAL(twod_patch.dimension(2), 2); @@ -127,7 +127,7 @@ static void test_simple_patch() // 2D patch: RowMajor Tensor twod_patch_row_major; - twod_patch_row_major = tensor_row_major.extract_image_patches<2, 2>(); + twod_patch_row_major = tensor_row_major.extract_image_patches(2, 2); VERIFY_IS_EQUAL(twod_patch_row_major.dimension(0), 7); VERIFY_IS_EQUAL(twod_patch_row_major.dimension(1), 3*5); VERIFY_IS_EQUAL(twod_patch_row_major.dimension(2), 2); @@ -194,7 +194,7 @@ static void test_patch_padding_valid() tensor.data()[i] = i + 1; } // ColMajor - Tensor result = tensor.extract_image_patches(ksize, ksize, stride, stride, PADDING_VALID); + Tensor result = tensor.extract_image_patches(ksize, ksize, stride, stride, 1, 1, PADDING_VALID); VERIFY_IS_EQUAL(result.dimension(0), input_depth); // depth VERIFY_IS_EQUAL(result.dimension(1), ksize); // kernel rows @@ -209,7 +209,7 @@ static void test_patch_padding_valid() VERIFY_IS_EQUAL(tensor.dimension(2), tensor_row_major.dimension(1)); VERIFY_IS_EQUAL(tensor.dimension(3), tensor_row_major.dimension(0)); - Tensor result_row_major = tensor_row_major.extract_image_patches(ksize, ksize, stride, stride, PADDING_VALID); + Tensor result_row_major = tensor_row_major.extract_image_patches(ksize, ksize, stride, stride, 1, 1, PADDING_VALID); VERIFY_IS_EQUAL(result.dimension(0), result_row_major.dimension(4)); VERIFY_IS_EQUAL(result.dimension(1), result_row_major.dimension(3)); VERIFY_IS_EQUAL(result.dimension(2), result_row_major.dimension(2)); @@ -267,7 +267,7 @@ static void test_patch_padding_valid_same_value() // ColMajor Tensor tensor(input_depth, input_rows, input_cols, input_batches); tensor = tensor.constant(11.0f); - Tensor result = tensor.extract_image_patches(ksize, ksize, stride, stride, PADDING_VALID); + Tensor result = tensor.extract_image_patches(ksize, ksize, stride, stride, 1, 1, PADDING_VALID); VERIFY_IS_EQUAL(result.dimension(0), input_depth); // depth VERIFY_IS_EQUAL(result.dimension(1), ksize); // kernel rows @@ -282,7 +282,7 @@ static void test_patch_padding_valid_same_value() VERIFY_IS_EQUAL(tensor.dimension(2), tensor_row_major.dimension(1)); VERIFY_IS_EQUAL(tensor.dimension(3), tensor_row_major.dimension(0)); - Tensor result_row_major = tensor_row_major.extract_image_patches(ksize, ksize, stride, stride, PADDING_VALID); + Tensor result_row_major = tensor_row_major.extract_image_patches(ksize, ksize, stride, stride, 1, 1, PADDING_VALID); VERIFY_IS_EQUAL(result.dimension(0), result_row_major.dimension(4)); VERIFY_IS_EQUAL(result.dimension(1), result_row_major.dimension(3)); VERIFY_IS_EQUAL(result.dimension(2), result_row_major.dimension(2)); @@ -416,7 +416,7 @@ static void test_patch_no_extra_dim() // Single pixel patch: ColMajor Tensor single_pixel_patch; - single_pixel_patch = tensor.extract_image_patches<1, 1>(); + single_pixel_patch = tensor.extract_image_patches(1, 1); VERIFY_IS_EQUAL(single_pixel_patch.dimension(0), 2); VERIFY_IS_EQUAL(single_pixel_patch.dimension(1), 1); VERIFY_IS_EQUAL(single_pixel_patch.dimension(2), 1); @@ -424,7 +424,7 @@ static void test_patch_no_extra_dim() // Single pixel patch: RowMajor Tensor single_pixel_patch_row_major; - single_pixel_patch_row_major = tensor_row_major.extract_image_patches<1, 1>(); + single_pixel_patch_row_major = tensor_row_major.extract_image_patches(1, 1); VERIFY_IS_EQUAL(single_pixel_patch_row_major.dimension(0), 3*5); VERIFY_IS_EQUAL(single_pixel_patch_row_major.dimension(1), 1); VERIFY_IS_EQUAL(single_pixel_patch_row_major.dimension(2), 1); @@ -451,7 +451,7 @@ static void test_patch_no_extra_dim() // Entire image patch: ColMajor Tensor entire_image_patch; - entire_image_patch = tensor.extract_image_patches<3, 5>(); + entire_image_patch = tensor.extract_image_patches(3, 5); VERIFY_IS_EQUAL(entire_image_patch.dimension(0), 2); VERIFY_IS_EQUAL(entire_image_patch.dimension(1), 3); VERIFY_IS_EQUAL(entire_image_patch.dimension(2), 5); @@ -459,7 +459,7 @@ static void test_patch_no_extra_dim() // Entire image patch: RowMajor Tensor entire_image_patch_row_major; - entire_image_patch_row_major = tensor_row_major.extract_image_patches<3, 5>(); + entire_image_patch_row_major = tensor_row_major.extract_image_patches(3, 5); VERIFY_IS_EQUAL(entire_image_patch_row_major.dimension(0), 3*5); VERIFY_IS_EQUAL(entire_image_patch_row_major.dimension(1), 5); VERIFY_IS_EQUAL(entire_image_patch_row_major.dimension(2), 3); @@ -499,7 +499,7 @@ static void test_patch_no_extra_dim() // 2D patch: ColMajor Tensor twod_patch; - twod_patch = tensor.extract_image_patches<2, 2>(); + twod_patch = tensor.extract_image_patches(2, 2); VERIFY_IS_EQUAL(twod_patch.dimension(0), 2); VERIFY_IS_EQUAL(twod_patch.dimension(1), 2); VERIFY_IS_EQUAL(twod_patch.dimension(2), 2); @@ -507,7 +507,7 @@ static void test_patch_no_extra_dim() // 2D patch: RowMajor Tensor twod_patch_row_major; - twod_patch_row_major = tensor_row_major.extract_image_patches<2, 2>(); + twod_patch_row_major = tensor_row_major.extract_image_patches(2, 2); VERIFY_IS_EQUAL(twod_patch_row_major.dimension(0), 3*5); VERIFY_IS_EQUAL(twod_patch_row_major.dimension(1), 2); VERIFY_IS_EQUAL(twod_patch_row_major.dimension(2), 2); From 7b7df7b6b816c1752b99319639da4a6799b18125 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Tue, 7 Jul 2015 12:57:35 -0700 Subject: [PATCH 08/10] Updated internal::is_arithmetic::value to be true for complex numbers --- Eigen/src/Core/util/Meta.h | 1 + 1 file changed, 1 insertion(+) diff --git a/Eigen/src/Core/util/Meta.h b/Eigen/src/Core/util/Meta.h index 974f11516..385d0f650 100644 --- a/Eigen/src/Core/util/Meta.h +++ b/Eigen/src/Core/util/Meta.h @@ -67,6 +67,7 @@ template<> struct is_arithmetic { enum { value = true }; }; template<> struct is_arithmetic { enum { value = true }; }; template<> struct is_arithmetic { enum { value = true }; }; template<> struct is_arithmetic { enum { value = true }; }; +template struct is_arithmetic > { enum { value = true }; }; template struct add_const { typedef const T type; }; template struct add_const { typedef T& type; }; From 6de6fa94830ff6d2be0e1ceed4587cad88b11762 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Tue, 7 Jul 2015 15:23:56 -0700 Subject: [PATCH 09/10] Use NumTraits::RequireInitialization instead of internal::is_arithmetic::value to check whether it's possible to bypass the type constructor in the tensor code. --- Eigen/src/Core/util/Meta.h | 1 - unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h | 2 +- unsupported/Eigen/CXX11/src/Tensor/TensorForcedEval.h | 2 +- unsupported/Eigen/CXX11/src/Tensor/TensorMorphing.h | 2 +- 4 files changed, 3 insertions(+), 4 deletions(-) diff --git a/Eigen/src/Core/util/Meta.h b/Eigen/src/Core/util/Meta.h index 385d0f650..974f11516 100644 --- a/Eigen/src/Core/util/Meta.h +++ b/Eigen/src/Core/util/Meta.h @@ -67,7 +67,6 @@ template<> struct is_arithmetic { enum { value = true }; }; template<> struct is_arithmetic { enum { value = true }; }; template<> struct is_arithmetic { enum { value = true }; }; template<> struct is_arithmetic { enum { value = true }; }; -template struct is_arithmetic > { enum { value = true }; }; template struct add_const { typedef const T type; }; template struct add_const { typedef T& type; }; diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h b/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h index 36718e26f..e3bef8676 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorEvaluator.h @@ -157,7 +157,7 @@ struct TensorEvaluator EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Dimensions& dimensions() const { return m_dims; } EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE bool evalSubExprsIfNeeded(CoeffReturnType* data) { - if (internal::is_arithmetic::type>::value && data) { + if (!NumTraits::type>::RequireInitialization && data) { m_device.memcpy((void*)data, m_data, m_dims.TotalSize() * sizeof(Scalar)); return false; } diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorForcedEval.h b/unsupported/Eigen/CXX11/src/Tensor/TensorForcedEval.h index d253b70f3..d0202559a 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorForcedEval.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorForcedEval.h @@ -109,7 +109,7 @@ struct TensorEvaluator, Device> const Index numValues = m_impl.dimensions().TotalSize(); m_buffer = (CoeffReturnType*)m_device.allocate(numValues * sizeof(CoeffReturnType)); // Should initialize the memory in case we're dealing with non POD types. - if (!internal::is_arithmetic::value) { + if (NumTraits::RequireInitialization) { for (Index i = 0; i < numValues; ++i) { new(m_buffer+i) CoeffReturnType(); } diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorMorphing.h b/unsupported/Eigen/CXX11/src/Tensor/TensorMorphing.h index fc98c3d1e..29bc6ca21 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorMorphing.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorMorphing.h @@ -366,7 +366,7 @@ struct TensorEvaluator, Devi EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE bool evalSubExprsIfNeeded(CoeffReturnType* data) { m_impl.evalSubExprsIfNeeded(NULL); - if (internal::is_arithmetic::type>::value && data && m_impl.data()) { + if (!NumTraits::type>::RequireInitialization && data && m_impl.data()) { Index contiguous_values = 1; if (static_cast(Layout) == static_cast(ColMajor)) { for (int i = 0; i < NumDims; ++i) { From e6297741c9d5e6106b6fa4876afac9571e038161 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Tue, 7 Jul 2015 17:40:49 -0700 Subject: [PATCH 10/10] Added support for generation of random complex numbers on CUDA devices --- .../Eigen/CXX11/src/Tensor/TensorFunctors.h | 104 ++++++++++++++++++ unsupported/test/CMakeLists.txt | 1 + unsupported/test/cxx11_tensor_random_cuda.cpp | 35 ++++++ 3 files changed, 140 insertions(+) create mode 100644 unsupported/test/cxx11_tensor_random_cuda.cpp diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorFunctors.h b/unsupported/Eigen/CXX11/src/Tensor/TensorFunctors.h index 33e8c01c2..14ffd5c93 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorFunctors.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorFunctors.h @@ -387,6 +387,58 @@ template <> class UniformRandomGenerator { mutable curandStatePhilox4_32_10_t m_state; }; +template <> class UniformRandomGenerator > { + public: + static const bool PacketAccess = false; + + __device__ UniformRandomGenerator(bool deterministic = true) : m_deterministic(deterministic) { + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + __device__ UniformRandomGenerator(const UniformRandomGenerator& other) { + m_deterministic = other.m_deterministic; + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = m_deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + template + __device__ std::complex operator()(Index, Index = 0) const { + float4 vals = curand_uniform4(&m_state); + return std::complex(vals.x, vals.y); + } + + private: + bool m_deterministic; + mutable curandStatePhilox4_32_10_t m_state; +}; + +template <> class UniformRandomGenerator > { + public: + static const bool PacketAccess = false; + + __device__ UniformRandomGenerator(bool deterministic = true) : m_deterministic(deterministic) { + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + __device__ UniformRandomGenerator(const UniformRandomGenerator& other) { + m_deterministic = other.m_deterministic; + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = m_deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + template + __device__ std::complex operator()(Index, Index = 0) const { + double2 vals = curand_uniform2_double(&m_state); + return std::complex(vals.x, vals.y); + } + + private: + bool m_deterministic; + mutable curandStatePhilox4_32_10_t m_state; +}; + #endif @@ -489,6 +541,58 @@ template <> class NormalRandomGenerator { mutable curandStatePhilox4_32_10_t m_state; }; +template <> class NormalRandomGenerator > { + public: + static const bool PacketAccess = false; + + __device__ NormalRandomGenerator(bool deterministic = true) : m_deterministic(deterministic) { + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + __device__ NormalRandomGenerator(const NormalRandomGenerator& other) { + m_deterministic = other.m_deterministic; + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = m_deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + template + __device__ std::complex operator()(Index, Index = 0) const { + float4 vals = curand_normal4(&m_state); + return std::complex(vals.x, vals.y); + } + + private: + bool m_deterministic; + mutable curandStatePhilox4_32_10_t m_state; +}; + +template <> class NormalRandomGenerator > { + public: + static const bool PacketAccess = false; + + __device__ NormalRandomGenerator(bool deterministic = true) : m_deterministic(deterministic) { + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + __device__ NormalRandomGenerator(const NormalRandomGenerator& other) { + m_deterministic = other.m_deterministic; + const int tid = blockIdx.x * blockDim.x + threadIdx.x; + const int seed = m_deterministic ? 0 : get_random_seed(); + curand_init(seed, tid, 0, &m_state); + } + template + __device__ std::complex operator()(Index, Index = 0) const { + double2 vals = curand_normal2_double(&m_state); + return std::complex(vals.x, vals.y); + } + + private: + bool m_deterministic; + mutable curandStatePhilox4_32_10_t m_state; +}; + #else template class NormalRandomGenerator { diff --git a/unsupported/test/CMakeLists.txt b/unsupported/test/CMakeLists.txt index 155bfcd76..845cda8ce 100644 --- a/unsupported/test/CMakeLists.txt +++ b/unsupported/test/CMakeLists.txt @@ -144,5 +144,6 @@ if(EIGEN_TEST_CXX11) # ei_add_test(cxx11_tensor_cuda "-std=c++0x") # ei_add_test(cxx11_tensor_contract_cuda "-std=c++0x") # ei_add_test(cxx11_tensor_reduction_cuda "-std=c++0x") +# ei_add_test(cxx11_tensor_random_cuda "-std=c++0x") endif() diff --git a/unsupported/test/cxx11_tensor_random_cuda.cpp b/unsupported/test/cxx11_tensor_random_cuda.cpp new file mode 100644 index 000000000..5d091de15 --- /dev/null +++ b/unsupported/test/cxx11_tensor_random_cuda.cpp @@ -0,0 +1,35 @@ +// 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/. + +#define EIGEN_TEST_NO_LONGDOUBLE +#define EIGEN_TEST_NO_COMPLEX +#define EIGEN_TEST_FUNC cxx11_tensor_random_cuda +#define EIGEN_DEFAULT_DENSE_INDEX_TYPE int +#define EIGEN_USE_GPU + +#include "main.h" +#include + +static void test_default() +{ + Tensor, 1> vec(6); + vec.setRandom(); + + // Fixme: we should check that the generated numbers follow a uniform + // distribution instead. + for (int i = 1; i < 6; ++i) { + VERIFY_IS_NOT_EQUAL(vec(i), vec(i-1)); + } +} + + +void test_cxx11_tensor_random_cuda() +{ + CALL_SUBTEST(test_default()); +}