GPU bug fixes: device-side results and device-callable helpers libeigen/eigen!3314 Co-authored-by: Rasmus Munk Larsen <rmlarsen@gmail.com>
diff --git a/Eigen/src/Core/MathFunctions.h b/Eigen/src/Core/MathFunctions.h index ad3118d..3bf3d6e 100644 --- a/Eigen/src/Core/MathFunctions.h +++ b/Eigen/src/Core/MathFunctions.h
@@ -1745,11 +1745,13 @@ return ::exp(x); } +// exp(x + 0i) = exp(x) + 0i even where exp(x) overflows, so a zero imaginary part is kept rather than computed as +// inf * 0. template <> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE std::complex<float> exp(const std::complex<float>& x) { float com = ::expf(x.real()); float res_real = com * ::cosf(x.imag()); - float res_imag = com * ::sinf(x.imag()); + float res_imag = x.imag() == 0.f ? x.imag() : com * ::sinf(x.imag()); return std::complex<float>(res_real, res_imag); } @@ -1757,7 +1759,7 @@ EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE std::complex<double> exp(const std::complex<double>& x) { double com = ::exp(x.real()); double res_real = com * ::cos(x.imag()); - double res_imag = com * ::sin(x.imag()); + double res_imag = x.imag() == 0.0 ? x.imag() : com * ::sin(x.imag()); return std::complex<double>(res_real, res_imag); } #endif
diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h index 8135342..ae17e0d 100644 --- a/Eigen/src/Core/MatrixBase.h +++ b/Eigen/src/Core/MatrixBase.h
@@ -163,13 +163,13 @@ const MatrixBase<OtherDerived>& other) const; template <typename OtherDerived> - Derived& operator*=(const EigenBase<OtherDerived>& other); + EIGEN_DEVICE_FUNC Derived& operator*=(const EigenBase<OtherDerived>& other); template <typename OtherDerived> - void applyOnTheLeft(const EigenBase<OtherDerived>& other); + EIGEN_DEVICE_FUNC void applyOnTheLeft(const EigenBase<OtherDerived>& other); template <typename OtherDerived> - void applyOnTheRight(const EigenBase<OtherDerived>& other); + EIGEN_DEVICE_FUNC void applyOnTheRight(const EigenBase<OtherDerived>& other); template <typename DiagonalDerived> EIGEN_DEVICE_FUNC const Product<Derived, DiagonalDerived, LazyProduct> operator*( @@ -530,7 +530,7 @@ */ template <typename Derived> template <typename OtherDerived> -inline Derived& MatrixBase<Derived>::operator*=(const EigenBase<OtherDerived>& other) { +EIGEN_DEVICE_FUNC inline Derived& MatrixBase<Derived>::operator*=(const EigenBase<OtherDerived>& other) { other.derived().applyThisOnTheRight(derived()); return derived(); } @@ -543,7 +543,7 @@ */ template <typename Derived> template <typename OtherDerived> -inline void MatrixBase<Derived>::applyOnTheRight(const EigenBase<OtherDerived>& other) { +EIGEN_DEVICE_FUNC inline void MatrixBase<Derived>::applyOnTheRight(const EigenBase<OtherDerived>& other) { other.derived().applyThisOnTheRight(derived()); } @@ -554,7 +554,7 @@ */ template <typename Derived> template <typename OtherDerived> -inline void MatrixBase<Derived>::applyOnTheLeft(const EigenBase<OtherDerived>& other) { +EIGEN_DEVICE_FUNC inline void MatrixBase<Derived>::applyOnTheLeft(const EigenBase<OtherDerived>& other) { other.derived().applyThisOnTheLeft(derived()); }
diff --git a/Eigen/src/Core/PlainObjectBase.h b/Eigen/src/Core/PlainObjectBase.h index efa4f1e..d758bea 100644 --- a/Eigen/src/Core/PlainObjectBase.h +++ b/Eigen/src/Core/PlainObjectBase.h
@@ -594,84 +594,92 @@ * \see class Map */ ///@{ - static inline ConstMapType Map(const Scalar* data) { return ConstMapType(data); } - static inline MapType Map(Scalar* data) { return MapType(data); } - static inline ConstMapType Map(const Scalar* data, Index size) { return ConstMapType(data, size); } - static inline MapType Map(Scalar* data, Index size) { return MapType(data, size); } - static inline ConstMapType Map(const Scalar* data, Index rows, Index cols) { return ConstMapType(data, rows, cols); } - static inline MapType Map(Scalar* data, Index rows, Index cols) { return MapType(data, rows, cols); } + EIGEN_DEVICE_FUNC static inline ConstMapType Map(const Scalar* data) { return ConstMapType(data); } + EIGEN_DEVICE_FUNC static inline MapType Map(Scalar* data) { return MapType(data); } + EIGEN_DEVICE_FUNC static inline ConstMapType Map(const Scalar* data, Index size) { return ConstMapType(data, size); } + EIGEN_DEVICE_FUNC static inline MapType Map(Scalar* data, Index size) { return MapType(data, size); } + EIGEN_DEVICE_FUNC static inline ConstMapType Map(const Scalar* data, Index rows, Index cols) { + return ConstMapType(data, rows, cols); + } + EIGEN_DEVICE_FUNC static inline MapType Map(Scalar* data, Index rows, Index cols) { + return MapType(data, rows, cols); + } - static inline ConstAlignedMapType MapAligned(const Scalar* data) { return ConstAlignedMapType(data); } - static inline AlignedMapType MapAligned(Scalar* data) { return AlignedMapType(data); } - static inline ConstAlignedMapType MapAligned(const Scalar* data, Index size) { + EIGEN_DEVICE_FUNC static inline ConstAlignedMapType MapAligned(const Scalar* data) { + return ConstAlignedMapType(data); + } + EIGEN_DEVICE_FUNC static inline AlignedMapType MapAligned(Scalar* data) { return AlignedMapType(data); } + EIGEN_DEVICE_FUNC static inline ConstAlignedMapType MapAligned(const Scalar* data, Index size) { return ConstAlignedMapType(data, size); } - static inline AlignedMapType MapAligned(Scalar* data, Index size) { return AlignedMapType(data, size); } - static inline ConstAlignedMapType MapAligned(const Scalar* data, Index rows, Index cols) { + EIGEN_DEVICE_FUNC static inline AlignedMapType MapAligned(Scalar* data, Index size) { + return AlignedMapType(data, size); + } + EIGEN_DEVICE_FUNC static inline ConstAlignedMapType MapAligned(const Scalar* data, Index rows, Index cols) { return ConstAlignedMapType(data, rows, cols); } - static inline AlignedMapType MapAligned(Scalar* data, Index rows, Index cols) { + EIGEN_DEVICE_FUNC static inline AlignedMapType MapAligned(Scalar* data, Index rows, Index cols) { return AlignedMapType(data, rows, cols); } template <int Outer, int Inner> - static inline typename StridedConstMapType<Stride<Outer, Inner>>::type Map(const Scalar* data, - const Stride<Outer, Inner>& stride) { + EIGEN_DEVICE_FUNC static inline typename StridedConstMapType<Stride<Outer, Inner>>::type Map( + const Scalar* data, const Stride<Outer, Inner>& stride) { return typename StridedConstMapType<Stride<Outer, Inner>>::type(data, stride); } template <int Outer, int Inner> - static inline typename StridedMapType<Stride<Outer, Inner>>::type Map(Scalar* data, - const Stride<Outer, Inner>& stride) { + EIGEN_DEVICE_FUNC static inline typename StridedMapType<Stride<Outer, Inner>>::type Map( + Scalar* data, const Stride<Outer, Inner>& stride) { return typename StridedMapType<Stride<Outer, Inner>>::type(data, stride); } template <int Outer, int Inner> - static inline typename StridedConstMapType<Stride<Outer, Inner>>::type Map(const Scalar* data, Index size, - const Stride<Outer, Inner>& stride) { + EIGEN_DEVICE_FUNC static inline typename StridedConstMapType<Stride<Outer, Inner>>::type Map( + const Scalar* data, Index size, const Stride<Outer, Inner>& stride) { return typename StridedConstMapType<Stride<Outer, Inner>>::type(data, size, stride); } template <int Outer, int Inner> - static inline typename StridedMapType<Stride<Outer, Inner>>::type Map(Scalar* data, Index size, - const Stride<Outer, Inner>& stride) { + EIGEN_DEVICE_FUNC static inline typename StridedMapType<Stride<Outer, Inner>>::type Map( + Scalar* data, Index size, const Stride<Outer, Inner>& stride) { return typename StridedMapType<Stride<Outer, Inner>>::type(data, size, stride); } template <int Outer, int Inner> - static inline typename StridedConstMapType<Stride<Outer, Inner>>::type Map(const Scalar* data, Index rows, Index cols, - const Stride<Outer, Inner>& stride) { + EIGEN_DEVICE_FUNC static inline typename StridedConstMapType<Stride<Outer, Inner>>::type Map( + const Scalar* data, Index rows, Index cols, const Stride<Outer, Inner>& stride) { return typename StridedConstMapType<Stride<Outer, Inner>>::type(data, rows, cols, stride); } template <int Outer, int Inner> - static inline typename StridedMapType<Stride<Outer, Inner>>::type Map(Scalar* data, Index rows, Index cols, - const Stride<Outer, Inner>& stride) { + EIGEN_DEVICE_FUNC static inline typename StridedMapType<Stride<Outer, Inner>>::type Map( + Scalar* data, Index rows, Index cols, const Stride<Outer, Inner>& stride) { return typename StridedMapType<Stride<Outer, Inner>>::type(data, rows, cols, stride); } template <int Outer, int Inner> - static inline typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type MapAligned( + EIGEN_DEVICE_FUNC static inline typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type MapAligned( const Scalar* data, const Stride<Outer, Inner>& stride) { return typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type(data, stride); } template <int Outer, int Inner> - static inline typename StridedAlignedMapType<Stride<Outer, Inner>>::type MapAligned( + EIGEN_DEVICE_FUNC static inline typename StridedAlignedMapType<Stride<Outer, Inner>>::type MapAligned( Scalar* data, const Stride<Outer, Inner>& stride) { return typename StridedAlignedMapType<Stride<Outer, Inner>>::type(data, stride); } template <int Outer, int Inner> - static inline typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type MapAligned( + EIGEN_DEVICE_FUNC static inline typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type MapAligned( const Scalar* data, Index size, const Stride<Outer, Inner>& stride) { return typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type(data, size, stride); } template <int Outer, int Inner> - static inline typename StridedAlignedMapType<Stride<Outer, Inner>>::type MapAligned( + EIGEN_DEVICE_FUNC static inline typename StridedAlignedMapType<Stride<Outer, Inner>>::type MapAligned( Scalar* data, Index size, const Stride<Outer, Inner>& stride) { return typename StridedAlignedMapType<Stride<Outer, Inner>>::type(data, size, stride); } template <int Outer, int Inner> - static inline typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type MapAligned( + EIGEN_DEVICE_FUNC static inline typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type MapAligned( const Scalar* data, Index rows, Index cols, const Stride<Outer, Inner>& stride) { return typename StridedConstAlignedMapType<Stride<Outer, Inner>>::type(data, rows, cols, stride); } template <int Outer, int Inner> - static inline typename StridedAlignedMapType<Stride<Outer, Inner>>::type MapAligned( + EIGEN_DEVICE_FUNC static inline typename StridedAlignedMapType<Stride<Outer, Inner>>::type MapAligned( Scalar* data, Index rows, Index cols, const Stride<Outer, Inner>& stride) { return typename StridedAlignedMapType<Stride<Outer, Inner>>::type(data, rows, cols, stride); }
diff --git a/Eigen/src/Core/arch/Default/Half.h b/Eigen/src/Core/arch/Default/Half.h index 807ee39..8139d5e 100644 --- a/Eigen/src/Core/arch/Default/Half.h +++ b/Eigen/src/Core/arch/Default/Half.h
@@ -476,7 +476,8 @@ return a; } EIGEN_STRONG_INLINE __device__ bool operator==(const half& a, const half& b) { return __heq(a, b); } -EIGEN_STRONG_INLINE __device__ bool operator!=(const half& a, const half& b) { return __hne(a, b); } +// __hne is the ordered comparison, false for a NaN operand; != must be true there. +EIGEN_STRONG_INLINE __device__ bool operator!=(const half& a, const half& b) { return !__heq(a, b); } EIGEN_STRONG_INLINE __device__ bool operator<(const half& a, const half& b) { return __hlt(a, b); } EIGEN_STRONG_INLINE __device__ bool operator<=(const half& a, const half& b) { return __hle(a, b); } EIGEN_STRONG_INLINE __device__ bool operator>(const half& a, const half& b) { return __hgt(a, b); }
diff --git a/Eigen/src/Geometry/Scaling.h b/Eigen/src/Geometry/Scaling.h index be50fde..53927a9 100644 --- a/Eigen/src/Geometry/Scaling.h +++ b/Eigen/src/Geometry/Scaling.h
@@ -55,25 +55,25 @@ public: /** Default constructor without initialization. */ - UniformScaling() {} + EIGEN_DEVICE_FUNC UniformScaling() {} /** Constructs and initialize a uniform scaling transformation */ - explicit inline UniformScaling(const Scalar& s) : m_factor(s) {} + EIGEN_DEVICE_FUNC explicit inline UniformScaling(const Scalar& s) : m_factor(s) {} - inline const Scalar& factor() const { return m_factor; } - inline Scalar& factor() { return m_factor; } + EIGEN_DEVICE_FUNC inline const Scalar& factor() const { return m_factor; } + EIGEN_DEVICE_FUNC inline Scalar& factor() { return m_factor; } /** Concatenates two uniform scaling */ - inline UniformScaling operator*(const UniformScaling& other) const { + EIGEN_DEVICE_FUNC inline UniformScaling operator*(const UniformScaling& other) const { return UniformScaling(m_factor * other.factor()); } /** Concatenates a uniform scaling and a translation */ template <int Dim> - inline Transform<Scalar, Dim, Affine> operator*(const Translation<Scalar, Dim>& t) const; + EIGEN_DEVICE_FUNC inline Transform<Scalar, Dim, Affine> operator*(const Translation<Scalar, Dim>& t) const; /** Concatenates a uniform scaling and an affine transformation */ template <int Dim, int Mode, int Options> - inline typename internal::uniformscaling_times_affine_returntype<Scalar, Dim, Mode>::type operator*( + EIGEN_DEVICE_FUNC inline typename internal::uniformscaling_times_affine_returntype<Scalar, Dim, Mode>::type operator*( const Transform<Scalar, Dim, Mode, Options>& t) const { typename internal::uniformscaling_times_affine_returntype<Scalar, Dim, Mode>::type res = t; res.prescale(factor()); @@ -83,17 +83,18 @@ /** Concatenates a uniform scaling and a linear transformation matrix */ // TODO: return an expression instead of a dense matrix. template <typename Derived> - inline typename Eigen::internal::plain_matrix_type<Derived>::type operator*(const MatrixBase<Derived>& other) const { + EIGEN_DEVICE_FUNC inline typename Eigen::internal::plain_matrix_type<Derived>::type operator*( + const MatrixBase<Derived>& other) const { return other * m_factor; } template <typename Derived, int Dim> - inline Matrix<Scalar, Dim, Dim> operator*(const RotationBase<Derived, Dim>& r) const { + EIGEN_DEVICE_FUNC inline Matrix<Scalar, Dim, Dim> operator*(const RotationBase<Derived, Dim>& r) const { return r.toRotationMatrix() * m_factor; } /** \returns the inverse scaling */ - inline UniformScaling inverse() const { return UniformScaling(Scalar(1) / m_factor); } + EIGEN_DEVICE_FUNC inline UniformScaling inverse() const { return UniformScaling(Scalar(1) / m_factor); } /** \returns \c *this with scalar type casted to \a NewScalarType * @@ -101,13 +102,13 @@ * then this function smartly returns a const reference to \c *this. */ template <typename NewScalarType> - inline UniformScaling<NewScalarType> cast() const { + EIGEN_DEVICE_FUNC inline UniformScaling<NewScalarType> cast() const { return UniformScaling<NewScalarType>(NewScalarType(m_factor)); } /** Copy constructor with scalar type conversion */ template <typename OtherScalarType> - inline explicit UniformScaling(const UniformScaling<OtherScalarType>& other) { + EIGEN_DEVICE_FUNC inline explicit UniformScaling(const UniformScaling<OtherScalarType>& other) { m_factor = Scalar(other.factor()); } @@ -130,29 +131,29 @@ // NOTE this operator is defined in MatrixBase and not as a friend function // of UniformScaling to fix an internal crash of Intel's ICC template <typename Derived, typename Scalar> -EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(Derived, Scalar, internal::scalar_product_op) -operator*(const MatrixBase<Derived>& matrix, const UniformScaling<Scalar>& s) { +EIGEN_DEVICE_FUNC EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(Derived, Scalar, internal::scalar_product_op) operator*( + const MatrixBase<Derived>& matrix, const UniformScaling<Scalar>& s) { return matrix.derived() * s.factor(); } /** Constructs a uniform scaling from scale factor \a s */ -inline UniformScaling<float> Scaling(float s) { return UniformScaling<float>(s); } +EIGEN_DEVICE_FUNC inline UniformScaling<float> Scaling(float s) { return UniformScaling<float>(s); } /** Constructs a uniform scaling from scale factor \a s */ -inline UniformScaling<double> Scaling(double s) { return UniformScaling<double>(s); } +EIGEN_DEVICE_FUNC inline UniformScaling<double> Scaling(double s) { return UniformScaling<double>(s); } /** Constructs a uniform scaling from scale factor \a s */ template <typename RealScalar> -inline UniformScaling<std::complex<RealScalar> > Scaling(const std::complex<RealScalar>& s) { +EIGEN_DEVICE_FUNC inline UniformScaling<std::complex<RealScalar> > Scaling(const std::complex<RealScalar>& s) { return UniformScaling<std::complex<RealScalar> >(s); } /** Constructs a 2D axis aligned scaling */ template <typename Scalar> -inline DiagonalMatrix<Scalar, 2> Scaling(const Scalar& sx, const Scalar& sy) { +EIGEN_DEVICE_FUNC inline DiagonalMatrix<Scalar, 2> Scaling(const Scalar& sx, const Scalar& sy) { return DiagonalMatrix<Scalar, 2>(sx, sy); } /** Constructs a 3D axis aligned scaling */ template <typename Scalar> -inline DiagonalMatrix<Scalar, 3> Scaling(const Scalar& sx, const Scalar& sy, const Scalar& sz) { +EIGEN_DEVICE_FUNC inline DiagonalMatrix<Scalar, 3> Scaling(const Scalar& sx, const Scalar& sy, const Scalar& sz) { return DiagonalMatrix<Scalar, 3>(sx, sy, sz); } @@ -160,13 +161,13 @@ * This is an alias for coeffs.asDiagonal() */ template <typename Derived> -inline const DiagonalWrapper<const Derived> Scaling(const MatrixBase<Derived>& coeffs) { +EIGEN_DEVICE_FUNC inline const DiagonalWrapper<const Derived> Scaling(const MatrixBase<Derived>& coeffs) { return coeffs.asDiagonal(); } /** Constructs an axis aligned scaling expression from vector \a coeffs when passed as an rvalue reference */ template <typename Derived> -inline typename DiagonalWrapper<const Derived>::PlainObject Scaling(MatrixBase<Derived>&& coeffs) { +EIGEN_DEVICE_FUNC inline typename DiagonalWrapper<const Derived>::PlainObject Scaling(MatrixBase<Derived>&& coeffs) { return typename DiagonalWrapper<const Derived>::PlainObject(std::move(coeffs.derived())); } @@ -182,7 +183,8 @@ template <typename Scalar> template <int Dim> -inline Transform<Scalar, Dim, Affine> UniformScaling<Scalar>::operator*(const Translation<Scalar, Dim>& t) const { +EIGEN_DEVICE_FUNC inline Transform<Scalar, Dim, Affine> UniformScaling<Scalar>::operator*( + const Translation<Scalar, Dim>& t) const { Transform<Scalar, Dim, Affine> res; res.matrix().setZero(); res.linear().diagonal().fill(factor());
diff --git a/contrib/Eigen/src/GPU/GpuFFT.h b/contrib/Eigen/src/GPU/GpuFFT.h index 6aac926..204474f 100644 --- a/contrib/Eigen/src/GPU/GpuFFT.h +++ b/contrib/Eigen/src/GPU/GpuFFT.h
@@ -157,7 +157,7 @@ const ComplexVector input(X.derived()); const int n = static_cast<int>(nfft); const int n_complex = n / 2 + 1; - eigen_assert(input.size() == n_complex); + eigen_assert(n == 0 || input.size() == n_complex); if (n == 0) return RealVector(0); ensure_buffers(n_complex * sizeof(Complex), n * sizeof(Scalar));
diff --git a/contrib/Eigen/src/GPU/GpuLLT.h b/contrib/Eigen/src/GPU/GpuLLT.h index 76ada76..fee4d8e 100644 --- a/contrib/Eigen/src/GPU/GpuLLT.h +++ b/contrib/Eigen/src/GPU/GpuLLT.h
@@ -161,6 +161,7 @@ // the sync — use info() explicitly when failure must be detected. eigen_assert(solver_ctx_.info() == Success && "LLT::solve called on a failed or uninitialized factorization"); eigen_assert(B.rows() == n_); + if (n_ == 0) return PlainMatrix(0, B.cols()); const Ref<const PlainMatrix> rhs(B.derived()); const int64_t nrhs = static_cast<int64_t>(rhs.cols()); @@ -245,8 +246,10 @@ constexpr cudaDataType_t dtype = internal::cusolver_data_type<Scalar>::value; constexpr cublasFillMode_t uplo = internal::cusolver_fill_mode<UpLo_>::value; - EIGEN_CUSOLVER_CHECK(cusolverDnXpotrs(solver_ctx_.cusolverHandle(), solver_ctx_.params_.p, uplo, n_, nrhs, dtype, - d_factor_.get(), lda_, dtype, d_x.get(), ldb, solver_ctx_.scratch_info())); + // cuSOLVER rejects the empty problem's leading dimensions. + if (n_ > 0) + EIGEN_CUSOLVER_CHECK(cusolverDnXpotrs(solver_ctx_.cusolverHandle(), solver_ctx_.params_.p, uplo, n_, nrhs, dtype, + d_factor_.get(), lda_, dtype, d_x.get(), ldb, solver_ctx_.scratch_info())); DeviceMatrix<Scalar> result = DeviceMatrix<Scalar>::adopt(static_cast<Scalar*>(d_x.release()), n_, static_cast<Index>(nrhs));
diff --git a/contrib/Eigen/src/GPU/GpuLU.h b/contrib/Eigen/src/GPU/GpuLU.h index b97c43d..dc57945 100644 --- a/contrib/Eigen/src/GPU/GpuLU.h +++ b/contrib/Eigen/src/GPU/GpuLU.h
@@ -160,6 +160,7 @@ // the sync — use info() explicitly when failure must be detected. eigen_assert(solver_ctx_.info() == Success && "LU::solve called on a failed or uninitialized factorization"); eigen_assert(B.rows() == n_); + if (n_ == 0) return PlainMatrix(0, B.cols()); const Ref<const PlainMatrix> rhs(B.derived()); const int64_t nrhs = static_cast<int64_t>(rhs.cols()); @@ -246,9 +247,11 @@ constexpr cudaDataType_t dtype = internal::cusolver_data_type<Scalar>::value; const cublasOperation_t trans = internal::to_cublas_op(op); - EIGEN_CUSOLVER_CHECK(cusolverDnXgetrs(solver_ctx_.cusolverHandle(), solver_ctx_.params_.p, trans, n_, nrhs, dtype, - d_lu_.get(), lda_, static_cast<const int64_t*>(d_ipiv_.get()), dtype, - d_x.get(), ldb, solver_ctx_.scratch_info())); + // cuSOLVER rejects the empty problem's leading dimensions. + if (n_ > 0) + EIGEN_CUSOLVER_CHECK(cusolverDnXgetrs(solver_ctx_.cusolverHandle(), solver_ctx_.params_.p, trans, n_, nrhs, dtype, + d_lu_.get(), lda_, static_cast<const int64_t*>(d_ipiv_.get()), dtype, + d_x.get(), ldb, solver_ctx_.scratch_info())); DeviceMatrix<Scalar> result = DeviceMatrix<Scalar>::adopt(static_cast<Scalar*>(d_x.release()), n_, static_cast<Index>(nrhs));
diff --git a/contrib/Eigen/src/GPU/GpuQR.h b/contrib/Eigen/src/GPU/GpuQR.h index 02a2dc4..bbc3278 100644 --- a/contrib/Eigen/src/GPU/GpuQR.h +++ b/contrib/Eigen/src/GPU/GpuQR.h
@@ -150,6 +150,8 @@ // the sync — use info() explicitly when failure must be detected. eigen_assert(solver_ctx_.info() == Success && "QR::solve called on a failed or uninitialized factorization"); eigen_assert(B.rows() == m_); + // An empty A has the zero minimum-norm solution, as in xGELS; cuSOLVER would reject 0 x 0's ldb = 0. + if (m_ == 0 || n_ == 0) return PlainMatrix::Zero(n_, B.cols()); const Ref<const PlainMatrix> rhs(B.derived()); const Index nrhs = rhs.cols(); @@ -165,6 +167,11 @@ DeviceMatrix<Scalar> solve(const DeviceMatrix<Scalar>& d_B) const { eigen_assert(solver_ctx_.info() == Success && "QR::solve called on a failed or uninitialized factorization"); eigen_assert(d_B.rows() == m_); + if (m_ == 0 || n_ == 0) { + DeviceMatrix<Scalar> X(n_, d_B.cols()); + X.setZero(solver_ctx_.stream()); + return X; + } d_B.waitReady(solver_ctx_.stream()); if (!transposed_) { @@ -215,13 +222,13 @@ bool begin_compute(const DeviceMatrix<Scalar>& d_A) { m_ = d_A.rows(); n_ = d_A.cols(); + transposed_ = (m_ < n_); + lda_ = (std::max)(int64_t(1), static_cast<int64_t>(transposed_ ? n_ : m_)); if (!solver_ctx_.begin_compute(m_ != 0 && n_ != 0)) { d_qr_ = internal::DeviceBuffer(); d_tau_ = internal::DeviceBuffer(); return false; } - transposed_ = (m_ < n_); - lda_ = static_cast<int64_t>(transposed_ ? n_ : m_); d_A.waitReady(solver_ctx_.stream()); return true; }
diff --git a/contrib/Eigen/src/GPU/GpuSparseContext.h b/contrib/Eigen/src/GPU/GpuSparseContext.h index 3fbb1c9..734b318 100644 --- a/contrib/Eigen/src/GPU/GpuSparseContext.h +++ b/contrib/Eigen/src/GPU/GpuSparseContext.h
@@ -670,8 +670,18 @@ exec_spmv(x_size, y_size, d_x_.get(), d_y_.get(), alpha, beta, op); - EIGEN_CUDA_RUNTIME_CHECK( - cudaMemcpyAsync(y.data(), d_y_.get(), y_size * sizeof(Scalar), cudaMemcpyDeviceToHost, stream_)); + // The spacing of y's elements. innerStride() is not it for a vector shaped only at run time: a 1 x n block of a + // column-major matrix has innerStride() == 1 with its elements outerStride() apart. + const Index incy = y.cols() == 1 ? y.rowStride() : y.colStride(); + if (incy == 1) { + EIGEN_CUDA_RUNTIME_CHECK( + cudaMemcpyAsync(y.data(), d_y_.get(), y_size * sizeof(Scalar), cudaMemcpyDeviceToHost, stream_)); + } else { + // y is a strided block, e.g. a row of a column-major matrix. + EIGEN_CUDA_RUNTIME_CHECK(cudaMemcpy2DAsync(y.data(), static_cast<size_t>(incy) * sizeof(Scalar), d_y_.get(), + sizeof(Scalar), sizeof(Scalar), static_cast<size_t>(y_size), + cudaMemcpyDeviceToHost, stream_)); + } EIGEN_CUDA_RUNTIME_CHECK(cudaStreamSynchronize(stream_)); }
diff --git a/contrib/Eigen/src/GPU/GpuSparseSolverBase.h b/contrib/Eigen/src/GPU/GpuSparseSolverBase.h index f3ba8f3..0d82072 100644 --- a/contrib/Eigen/src/GPU/GpuSparseSolverBase.h +++ b/contrib/Eigen/src/GPU/GpuSparseSolverBase.h
@@ -360,7 +360,7 @@ eigen_assert(info() == Success && "GpuSparseSolver::solve requires a successful factorization"); eigen_assert(B.rows() == n_); - if (n_ == 0) return DenseMatrix(0, B.cols()); + if (n_ == 0 || B.cols() == 0) return DenseMatrix(n_, B.cols()); const Ref<const DenseMatrix> rhs(B.derived()); const int64_t nrhs = static_cast<int64_t>(rhs.cols());
diff --git a/contrib/Eigen/src/Tensor/TensorDeviceGpu.h b/contrib/Eigen/src/Tensor/TensorDeviceGpu.h index 5443ff8..231156b 100644 --- a/contrib/Eigen/src/Tensor/TensorDeviceGpu.h +++ b/contrib/Eigen/src/Tensor/TensorDeviceGpu.h
@@ -208,12 +208,14 @@ EIGEN_GPU_RUNTIME_CHECK(gpuSetDevice(device_)); void* result = nullptr; EIGEN_GPU_RUNTIME_CHECK(gpuMalloc(&result, num_bytes)); - gpu_assert(result != nullptr); + // A zero-byte request (an empty tensor's temporary) succeeds with a null pointer. + gpu_assert(result != nullptr || num_bytes == 0); return result; } void deallocate(void* buffer) const override { + // The null pointer of a zero-byte allocate(); callers such as TensorForcedEval free it unconditionally. + if (buffer == nullptr) return; EIGEN_GPU_RUNTIME_CHECK(gpuSetDevice(device_)); - gpu_assert(buffer != nullptr); EIGEN_GPU_RUNTIME_CHECK(gpuFree(buffer)); }
diff --git a/contrib/Eigen/src/Tensor/TensorReductionGpu.h b/contrib/Eigen/src/Tensor/TensorReductionGpu.h index ef17207..cc4bb38 100644 --- a/contrib/Eigen/src/Tensor/TensorReductionGpu.h +++ b/contrib/Eigen/src/Tensor/TensorReductionGpu.h
@@ -464,8 +464,10 @@ static void run(const Self& self, Op& reducer, const GpuDevice& device, OutputType* output) { gpu_assert(HasOptimizedImplementation && "Should only be called on doubles, floats or half floats"); const Index num_coeffs = array_prod(self.m_impl.dimensions()); - // Don't crash when we're called with an input tensor of size 0. + // An empty input reduces to the reducer's identity, as on the other devices. if (num_coeffs == 0) { + const OutputType identity = reducer.finalize(reducer.initialize()); + device.fill(output, output + 1, identity); return; }
diff --git a/contrib/Eigen/src/Tensor/TensorScan.h b/contrib/Eigen/src/Tensor/TensorScan.h index 2406f94..9126239 100644 --- a/contrib/Eigen/src/Tensor/TensorScan.h +++ b/contrib/Eigen/src/Tensor/TensorScan.h
@@ -330,6 +330,8 @@ struct ScanLauncher<Self, Reducer, GpuDevice, Vectorize> { void operator()(const Self& self, typename Self::CoeffReturnType* data) const { Index total_size = internal::array_prod(self.dimensions()); + // Nothing to scan; also avoids dividing by an empty scan axis and launching an empty grid. + if (total_size == 0) return; Index num_blocks = (total_size / self.size() + 63) / 64; Index block_size = 64;
diff --git a/contrib/test/GPU/cudss_lu.cpp b/contrib/test/GPU/cudss_lu.cpp index b108674..8b188a0 100644 --- a/contrib/test/GPU/cudss_lu.cpp +++ b/contrib/test/GPU/cudss_lu.cpp
@@ -225,6 +225,10 @@ VERIFY_IS_EQUAL(lu.info(), Success); VERIFY_IS_EQUAL(lu.rows(), 0); VERIFY_IS_EQUAL(lu.cols(), 0); + + // No right-hand-side columns: an n x 0 solution, without a cuDSS solve. + gpu::SparseLU<Scalar> lu8(make_general<Scalar>(8)); + VERIFY_IS_EQUAL(lu8.solve(Matrix<Scalar, Dynamic, Dynamic>(8, 0)).rows(), 8); } // ---- Per-scalar driver ------------------------------------------------------
diff --git a/contrib/test/GPU/cufft.cpp b/contrib/test/GPU/cufft.cpp index abdc9d3..c81f6d3 100644 --- a/contrib/test/GPU/cufft.cpp +++ b/contrib/test/GPU/cufft.cpp
@@ -211,6 +211,10 @@ VERIFY_IS_EQUAL(X.size(), 0); Vec y = fft.inv(X); VERIFY_IS_EQUAL(y.size(), 0); + + // Real round trip: fwd of an empty real signal is empty, and invReal accepts it back. + const Matrix<Scalar, Dynamic, 1> r(0); + VERIFY_IS_EQUAL(fft.invReal(fft.fwd(r), 0).size(), 0); } // ---- Explicit-Context ctor: stream and cuBLAS plumbing ----------------------
diff --git a/contrib/test/GPU/cusolver_llt.cpp b/contrib/test/GPU/cusolver_llt.cpp index ea87072..5c7a533 100644 --- a/contrib/test/GPU/cusolver_llt.cpp +++ b/contrib/test/GPU/cusolver_llt.cpp
@@ -115,6 +115,14 @@ } } +// An empty factorization solves a 0 x nrhs right-hand side. +void test_empty_solve() { + const MatrixXd B(0, 3); + gpu::LLT<double> solver{MatrixXd(0, 0)}; + VERIFY_IS_EQUAL(solver.solve(B).cols(), 3); + VERIFY_IS_EQUAL(solver.solve(gpu::DeviceMatrix<double>::fromHost(B)).cols(), 3); +} + // Test that GpuLLT correctly detects a non-SPD matrix. void test_not_spd() { Eigen::MatrixXd h_A = -Eigen::MatrixXd::Identity(8, 8); // negative definite @@ -370,5 +378,6 @@ CALL_SUBTEST_3(test_scalar<std::complex<float>>()); CALL_SUBTEST_4(test_scalar<std::complex<double>>()); CALL_SUBTEST_5(test_not_spd()); + CALL_SUBTEST_5(test_empty_solve()); CALL_SUBTEST_5(test_not_spd_device_solve_asserts()); }
diff --git a/contrib/test/GPU/cusolver_lu.cpp b/contrib/test/GPU/cusolver_lu.cpp index c80379a..e5b24af 100644 --- a/contrib/test/GPU/cusolver_lu.cpp +++ b/contrib/test/GPU/cusolver_lu.cpp
@@ -80,6 +80,14 @@ } } +// An empty factorization solves a 0 x nrhs right-hand side. +void test_empty_solve() { + const MatrixXd B(0, 3); + gpu::LU<double> solver{MatrixXd(0, 0)}; + VERIFY_IS_EQUAL(solver.solve(B).cols(), 3); + VERIFY_IS_EQUAL(solver.solve(gpu::DeviceMatrix<double>::fromHost(B)).cols(), 3); +} + // ---- Singular matrix detection ---------------------------------------------- void test_singular() { @@ -244,5 +252,6 @@ CALL_SUBTEST_3(test_scalar<std::complex<float>>()); CALL_SUBTEST_4(test_scalar<std::complex<double>>()); CALL_SUBTEST_5(test_singular()); + CALL_SUBTEST_5(test_empty_solve()); CALL_SUBTEST_5(test_singular_device_solve_asserts()); }
diff --git a/contrib/test/GPU/cusolver_qr.cpp b/contrib/test/GPU/cusolver_qr.cpp index dc9764a..82fb511 100644 --- a/contrib/test/GPU/cusolver_qr.cpp +++ b/contrib/test/GPU/cusolver_qr.cpp
@@ -296,6 +296,19 @@ VERIFY_IS_EQUAL(qr.info(), Success); VERIFY_IS_EQUAL(qr.rows(), 0); VERIFY_IS_EQUAL(qr.cols(), 0); + VERIFY_IS_EQUAL(qr.solve(MatrixXd(0, 2)).cols(), 2); + VERIFY_IS_EQUAL(qr.solve(gpu::DeviceMatrix<double>::fromHost(MatrixXd(0, 2))).cols(), 2); + + // An empty factor still solves: X is cols() x nrhs, zero for the minimum-norm case. + for (Index m : {5, 0}) { + const MatrixXd B = MatrixXd::Random(m, 2); + qr.compute(MatrixXd(m, 5 - m)); + const MatrixXd X = qr.solve(B); + VERIFY_IS_EQUAL(X.rows(), 5 - m); + VERIFY_IS_EQUAL(X.cols(), 2); + VERIFY_IS_EQUAL(X.norm(), 0.0); + VERIFY_IS_EQUAL(qr.solve(gpu::DeviceMatrix<double>::fromHost(B)).toHost().rows(), 5 - m); + } } EIGEN_DECLARE_TEST(gpu_cusolver_qr) {
diff --git a/contrib/test/GPU/cusparse_spmv.cpp b/contrib/test/GPU/cusparse_spmv.cpp index f2c11ce..434b716 100644 --- a/contrib/test/GPU/cusparse_spmv.cpp +++ b/contrib/test/GPU/cusparse_spmv.cpp
@@ -85,6 +85,27 @@ VERIFY((y_gpu - y_cpu).norm() / (y_cpu.norm() + RealScalar(1)) < tol); } +// ---- SpMV into a strided destination: a row of a column-major matrix --------- + +template <typename Scalar> +void test_spmv_strided_dest(Index n) { + using Mat = Matrix<Scalar, Dynamic, Dynamic>; + const SparseMatrix<Scalar, ColMajor, int> A = make_sparse<Scalar>(n, n); + const Matrix<Scalar, Dynamic, 1> x = Matrix<Scalar, Dynamic, 1>::Random(n); + Mat M = Mat::Zero(3, n), M_cpu = Mat::Zero(3, n); + auto row = M.row(1); + gpu::SparseContext<Scalar> ctx; + ctx.multiply(A, x, row); + M_cpu.row(1) = A * x; + VERIFY_IS_APPROX(M, M_cpu); + + // A 1 x n block shaped at run time: innerStride() is 1, its elements are outerStride() apart. + M.setZero(); + auto block_row = M.middleRows(1, 1); + ctx.multiply(A, x, block_row); + VERIFY_IS_APPROX(M, M_cpu); +} + // ---- Transpose: y = A^T * x ------------------------------------------------ template <typename Scalar> @@ -573,6 +594,7 @@ CALL_SUBTEST(test_spmv<Scalar>(128, 64)); // non-square CALL_SUBTEST(test_spmv<Scalar>(64, 128)); // wide CALL_SUBTEST(test_spmv_alpha_beta<Scalar>(64)); + CALL_SUBTEST(test_spmv_strided_dest<Scalar>(64)); CALL_SUBTEST(test_spmv_transpose<Scalar>(128, 64)); // cuSPARSE < 12 cannot represent A^H * x for complex scalars with the // CSR-of-A^T trick used by SparseContext; SparseContext asserts in that
diff --git a/contrib/test/tensor_device.cu b/contrib/test/tensor_device.cu index 66c81ec..e149b42 100644 --- a/contrib/test/tensor_device.cu +++ b/contrib/test/tensor_device.cu
@@ -437,6 +437,16 @@ test_device_memory<OffByOneScalar<int>>(context.device()); } +// Forcing an empty tensor allocates a zero-byte temporary, null on the GPU, and frees it. +void test_gpu_empty_eval() { + Eigen::GpuStreamDevice stream; + Eigen::GpuDevice gpu_device(&stream); + Eigen::TensorMap<Eigen::Tensor<float, 1>> in(static_cast<float*>(nullptr), 0), out(static_cast<float*>(nullptr), 0); + out.device(gpu_device) = in.eval() * 2.f; + gpu_device.deallocate(gpu_device.allocate(0)); + gpu_device.synchronize(); +} + // The attribute cache has to agree with the properties struct on every field both carry, and the accessors built // on it have to describe the device the stream belongs to. void test_device_attributes() { @@ -500,6 +510,7 @@ EIGEN_DECLARE_TEST(tensor_device) { CALL_SUBTEST_1(test_cpu()); CALL_SUBTEST_2(test_gpu()); + CALL_SUBTEST_2(test_gpu_empty_eval()); CALL_SUBTEST_3(test_device_attributes()); CALL_SUBTEST_3(test_stream_device_attributes()); }
diff --git a/contrib/test/tensor_reduction_gpu.cu b/contrib/test/tensor_reduction_gpu.cu index e4f9ddf..a42c5ef 100644 --- a/contrib/test/tensor_reduction_gpu.cu +++ b/contrib/test/tensor_reduction_gpu.cu
@@ -52,6 +52,24 @@ gpu_device.deallocate(gpu_out_ptr); } +// A full reduction of an empty tensor writes the reducer's identity. +template <typename Type> +static void test_empty_full_reduction() { + Eigen::GpuStreamDevice stream; + Eigen::GpuDevice gpu_device(&stream); + Type* gpu_out_ptr = static_cast<Type*>(gpu_device.allocate(sizeof(Type))); + const Type sentinel(7); + gpu_device.memcpyHostToDevice(gpu_out_ptr, &sentinel, sizeof(Type)); + TensorMap<Tensor<Type, 2> > in_gpu(static_cast<Type*>(nullptr), 0, 3); + TensorMap<Tensor<Type, 0> > out_gpu(gpu_out_ptr); + out_gpu.device(gpu_device) = in_gpu.sum(); + Type result; + gpu_device.memcpyDeviceToHost(&result, gpu_out_ptr, sizeof(Type)); + gpu_device.synchronize(); + VERIFY_IS_EQUAL(result, Type(0)); + gpu_device.deallocate(gpu_out_ptr); +} + template <typename Type, int DataLayout> static void test_first_dim_reductions(int dim_x = 33, int dim_y = 1, int dim_z = 128) { Tensor<Type, 3, DataLayout> in(dim_x, dim_y, dim_z); @@ -127,6 +145,8 @@ EIGEN_DECLARE_TEST(tensor_reduction_gpu) { CALL_SUBTEST_1((test_full_reductions<float, ColMajor>())); CALL_SUBTEST_1((test_full_reductions<double, ColMajor>())); + CALL_SUBTEST_1((test_empty_full_reduction<float>())); + CALL_SUBTEST_1((test_empty_full_reduction<double>())); CALL_SUBTEST_2((test_full_reductions<float, RowMajor>())); CALL_SUBTEST_2((test_full_reductions<double, RowMajor>()));
diff --git a/contrib/test/tensor_scan_gpu.cu b/contrib/test/tensor_scan_gpu.cu index f586c08..fb8030b 100644 --- a/contrib/test/tensor_scan_gpu.cu +++ b/contrib/test/tensor_scan_gpu.cu
@@ -70,5 +70,7 @@ EIGEN_DECLARE_TEST(tensor_scan_gpu) { CALL_SUBTEST_1(test_gpu_cumsum<ColMajor>(128, 128, 128)); + CALL_SUBTEST_1(test_gpu_cumsum<ColMajor>(4, 0, 4)); CALL_SUBTEST_2(test_gpu_cumsum<RowMajor>(128, 128, 128)); + CALL_SUBTEST_2(test_gpu_cumsum<RowMajor>(0, 4, 4)); }
diff --git a/test/gpu_basic.cu b/test/gpu_basic.cu index 83a1694..ae183d7 100644 --- a/test/gpu_basic.cu +++ b/test/gpu_basic.cu
@@ -71,6 +71,35 @@ } }; +// The static Map() helpers of the plain types are device-callable. +template <typename T> +struct static_map { + EIGEN_DEVICE_FUNC void operator()(int i, const typename T::Scalar* in, typename T::Scalar* out) const { + T::Map(out + i * T::SizeAtCompileTime) = + 2 * T::Map(in + i) + T::Map(in + i, T::RowsAtCompileTime, T::ColsAtCompileTime); + } +}; + +// A transform times a diagonal matrix scales the linear part through MatrixBase::operator*=. +struct transform_times_diagonal { + EIGEN_DEVICE_FUNC void operator()(int i, const float* in, float* out) const { + using namespace Eigen; + Affine3f T; + T.matrix() = Matrix4f(in + i); + Map<Matrix4f>(out + i * 16) = (T * DiagonalMatrix<float, 3>(in[i], in[i + 1], in[i + 2])).matrix(); + } +}; + +// The Scaling() helpers and UniformScaling compose with transforms on the device. +struct transform_scaling { + EIGEN_DEVICE_FUNC void operator()(int i, const float* in, float* out) const { + using namespace Eigen; + Affine3f T = Translation3f(Vector3f(in + i)) * Scaling(in[i + 3]); + T *= Scaling(in[i + 4]).inverse(); + Map<Matrix4f>(out + i * 16) = (T * Scaling(in[i + 5], in[i + 6], in[i + 7])).matrix(); + } +}; + template <int Order> struct scaled_structured_product { EIGEN_DEVICE_FUNC void operator()(int i, const float* in, float* out) const { @@ -725,6 +754,54 @@ #endif } +// Eigen::half != is true when either operand is NaN, on the device as on the host. +struct half_nan_not_equal_test { + EIGEN_DEVICE_FUNC void operator()(int i, const float* in, float* out) const { + EIGEN_UNUSED_VARIABLE(i); + const Eigen::half nan(in[0]), one(in[1]); + out[0] = float(nan != nan); + out[1] = float(nan != one); + out[2] = float(one != nan); + out[3] = float(one != one); + } +}; + +void test_half_nan_not_equal() { + Eigen::ArrayXf in(2), out_ref(4), out_gpu(4); + in << std::numeric_limits<float>::quiet_NaN(), 1.f; + run_on_cpu(half_nan_not_equal_test(), 1, in, out_ref); + run_on_gpu(half_nan_not_equal_test(), 1, in, out_gpu); +#if !defined(EIGEN_GPU_COMPILE_PHASE) + VERIFY_IS_CWISE_EQUAL(out_gpu, Eigen::Array4f(1.f, 1.f, 1.f, 0.f)); + VERIFY_IS_CWISE_EQUAL(out_ref, out_gpu); +#endif +} + +// exp of a real-axis complex keeps the zero imaginary part, and its sign, where the real part overflows. +struct complex_exp_overflow_test { + EIGEN_DEVICE_FUNC void operator()(int i, const float* in, float* out) const { + EIGEN_UNUSED_VARIABLE(i); + const std::complex<float> e = numext::exp(std::complex<float>(in[0], 0.f)); + const std::complex<double> e2 = numext::exp(std::complex<double>(in[0] * 10, -0.0)); + const std::complex<double> p2 = numext::exp2(std::complex<double>(in[0] * 20, 0.0)); + out[0] = e.imag(); + out[1] = float(e2.imag()); + out[2] = float(p2.imag()); + } +}; + +void test_complex_exp_overflow() { + Eigen::ArrayXf in(1), out_ref(3), out_gpu(3); + in << 100.f; + run_on_cpu(complex_exp_overflow_test(), 1, in, out_ref); + run_on_gpu(complex_exp_overflow_test(), 1, in, out_gpu); +#if !defined(EIGEN_GPU_COMPILE_PHASE) + VERIFY_IS_CWISE_EQUAL(out_ref, Eigen::Array3f::Zero()); + VERIFY_IS_CWISE_EQUAL(out_gpu, Eigen::Array3f::Zero()); + VERIFY((std::signbit)(out_ref(1)) && (std::signbit)(out_gpu(1))); +#endif +} + template <typename Type1, typename Type2> bool verifyIsApproxWithInfsNans(const Type1& a, const Type2& b, typename Type1::Scalar* = 0) // Enabled for Eigen's type only @@ -861,6 +938,9 @@ CALL_SUBTEST(run_and_compare_to_gpu(coeff_wise<Vector3f>(), nthreads, in, out)); CALL_SUBTEST(run_and_compare_to_gpu(coeff_wise<Array44f>(), nthreads, in, out)); + CALL_SUBTEST(run_and_compare_to_gpu(static_map<Matrix3f>(), nthreads, in, out)); + CALL_SUBTEST(run_and_compare_to_gpu(transform_times_diagonal(), nthreads, in, out)); + CALL_SUBTEST(run_and_compare_to_gpu(transform_scaling(), nthreads, in, out)); CALL_SUBTEST(run_and_compare_to_gpu(make_householder_small_tail(), nthreads, in, out)); CALL_SUBTEST(run_and_compare_to_gpu(make_householder_complex_zero_tail(), nthreads, cfin, cfout)); @@ -933,6 +1013,8 @@ CALL_SUBTEST(test_custom_less_scalar_minmax()); CALL_SUBTEST(test_float_nan_minmax()); CALL_SUBTEST(test_bit_count()); + CALL_SUBTEST(test_half_nan_not_equal()); + CALL_SUBTEST(test_complex_exp_overflow()); // `preverse` on every GPU packet type. 32 elements cover several packets of // each, including `Packet4h2`, the widest at 8.