Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-01-19 09:51:38 

0001 // This file is part of Eigen, a lightweight C++ template library
0002 // for linear algebra.
0003 //
0004 // Copyright (C) 2014 Benoit Steiner <benoit.steiner.goog@gmail.com>
0005 // Copyright (C) 2021 C. Antonio Sanchez <cantonios@google.com>
0006 //
0007 // This Source Code Form is subject to the terms of the Mozilla
0008 // Public License v. 2.0. If a copy of the MPL was not distributed
0009 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
0010 
0011 #ifndef EIGEN_COMPLEX_CUDA_H
0012 #define EIGEN_COMPLEX_CUDA_H
0013 
0014 // clang-format off
0015 // Many std::complex methods such as operator+, operator-, operator* and
0016 // operator/ are not constexpr. Due to this, GCC and older versions of clang do
0017 // not treat them as device functions and thus Eigen functors making use of
0018 // these operators fail to compile. Here, we manually specialize these
0019 // operators and functors for complex types when building for CUDA to enable
0020 // their use on-device.
0021 
0022 #if defined(EIGEN_CUDACC) && defined(EIGEN_GPU_COMPILE_PHASE)
0023     
0024 // ICC already specializes std::complex<float> and std::complex<double>
0025 // operators, preventing us from making them device functions here.
0026 // This will lead to silent runtime errors if the operators are used on device.
0027 //
0028 // To allow std::complex operator use on device, define _OVERRIDE_COMPLEX_SPECIALIZATION_
0029 // prior to first inclusion of <complex>.  This prevents ICC from adding
0030 // its own specializations, so our custom ones below can be used instead.
0031 #if !(defined(EIGEN_COMP_ICC) && defined(_USE_COMPLEX_SPECIALIZATION_))
0032 
0033 // Import Eigen's internal operator specializations.
0034 #define EIGEN_USING_STD_COMPLEX_OPERATORS           \
0035   using Eigen::complex_operator_detail::operator+;  \
0036   using Eigen::complex_operator_detail::operator-;  \
0037   using Eigen::complex_operator_detail::operator*;  \
0038   using Eigen::complex_operator_detail::operator/;  \
0039   using Eigen::complex_operator_detail::operator+=; \
0040   using Eigen::complex_operator_detail::operator-=; \
0041   using Eigen::complex_operator_detail::operator*=; \
0042   using Eigen::complex_operator_detail::operator/=; \
0043   using Eigen::complex_operator_detail::operator==; \
0044   using Eigen::complex_operator_detail::operator!=;
0045 
0046 namespace Eigen {
0047 
0048 // Specialized std::complex overloads.
0049 namespace complex_operator_detail {
0050 
0051 template<typename T>
0052 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
0053 std::complex<T> complex_multiply(const std::complex<T>& a, const std::complex<T>& b) {
0054   const T a_real = numext::real(a);
0055   const T a_imag = numext::imag(a);
0056   const T b_real = numext::real(b);
0057   const T b_imag = numext::imag(b);
0058   return std::complex<T>(
0059       a_real * b_real - a_imag * b_imag,
0060       a_imag * b_real + a_real * b_imag);
0061 }
0062 
0063 template<typename T>
0064 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
0065 std::complex<T> complex_divide_fast(const std::complex<T>& a, const std::complex<T>& b) {
0066   const T a_real = numext::real(a);
0067   const T a_imag = numext::imag(a);
0068   const T b_real = numext::real(b);
0069   const T b_imag = numext::imag(b);
0070   const T norm = (b_real * b_real + b_imag * b_imag);
0071   return std::complex<T>((a_real * b_real + a_imag * b_imag) / norm,
0072                           (a_imag * b_real - a_real * b_imag) / norm);
0073 }
0074 
0075 template<typename T>
0076 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
0077 std::complex<T> complex_divide_stable(const std::complex<T>& a, const std::complex<T>& b) {
0078   const T a_real = numext::real(a);
0079   const T a_imag = numext::imag(a);
0080   const T b_real = numext::real(b);
0081   const T b_imag = numext::imag(b);
0082   // Smith's complex division (https://arxiv.org/pdf/1210.4539.pdf),
0083   // guards against over/under-flow.
0084   const bool scale_imag = numext::abs(b_imag) <= numext::abs(b_real);
0085   const T rscale = scale_imag ? T(1) : b_real / b_imag;
0086   const T iscale = scale_imag ? b_imag / b_real : T(1);
0087   const T denominator = b_real * rscale + b_imag * iscale;
0088   return std::complex<T>((a_real * rscale + a_imag * iscale) / denominator, 
0089                          (a_imag * rscale - a_real * iscale) / denominator);
0090 }
0091 
0092 template<typename T>
0093 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
0094 std::complex<T> complex_divide(const std::complex<T>& a, const std::complex<T>& b) {
0095 #if EIGEN_FAST_MATH
0096   return complex_divide_fast(a, b);
0097 #else
0098   return complex_divide_stable(a, b);
0099 #endif
0100 }
0101 
0102 // NOTE: We cannot specialize compound assignment operators with Scalar T,
0103 //         (i.e.  operator@=(const T&), for @=+,-,*,/)
0104 //       since they are already specialized for float/double/long double within
0105 //       the standard <complex> header. We also do not specialize the stream
0106 //       operators.
0107 #define EIGEN_CREATE_STD_COMPLEX_OPERATOR_SPECIALIZATIONS(T)                                    \
0108                                                                                                 \
0109 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0110 std::complex<T> operator+(const std::complex<T>& a) { return a; }                               \
0111                                                                                                 \
0112 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0113 std::complex<T> operator-(const std::complex<T>& a) {                                           \
0114   return std::complex<T>(-numext::real(a), -numext::imag(a));                                   \
0115 }                                                                                               \
0116                                                                                                 \
0117 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0118 std::complex<T> operator+(const std::complex<T>& a, const std::complex<T>& b) {                 \
0119   return std::complex<T>(numext::real(a) + numext::real(b), numext::imag(a) + numext::imag(b)); \
0120 }                                                                                               \
0121                                                                                                 \
0122 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0123 std::complex<T> operator+(const std::complex<T>& a, const T& b) {                               \
0124   return std::complex<T>(numext::real(a) + b, numext::imag(a));                                 \
0125 }                                                                                               \
0126                                                                                                 \
0127 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0128 std::complex<T> operator+(const T& a, const std::complex<T>& b) {                               \
0129   return std::complex<T>(a + numext::real(b), numext::imag(b));                                 \
0130 }                                                                                               \
0131                                                                                                 \
0132 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0133 std::complex<T> operator-(const std::complex<T>& a, const std::complex<T>& b) {                 \
0134   return std::complex<T>(numext::real(a) - numext::real(b), numext::imag(a) - numext::imag(b)); \
0135 }                                                                                               \
0136                                                                                                 \
0137 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0138 std::complex<T> operator-(const std::complex<T>& a, const T& b) {                               \
0139   return std::complex<T>(numext::real(a) - b, numext::imag(a));                                 \
0140 }                                                                                               \
0141                                                                                                 \
0142 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0143 std::complex<T> operator-(const T& a, const std::complex<T>& b) {                               \
0144   return std::complex<T>(a - numext::real(b), -numext::imag(b));                                \
0145 }                                                                                               \
0146                                                                                                 \
0147 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0148 std::complex<T> operator*(const std::complex<T>& a, const std::complex<T>& b) {                 \
0149   return complex_multiply(a, b);                                                                \
0150 }                                                                                               \
0151                                                                                                 \
0152 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0153 std::complex<T> operator*(const std::complex<T>& a, const T& b) {                               \
0154   return std::complex<T>(numext::real(a) * b, numext::imag(a) * b);                             \
0155 }                                                                                               \
0156                                                                                                 \
0157 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0158 std::complex<T> operator*(const T& a, const std::complex<T>& b) {                               \
0159   return std::complex<T>(a * numext::real(b), a * numext::imag(b));                             \
0160 }                                                                                               \
0161                                                                                                 \
0162 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0163 std::complex<T> operator/(const std::complex<T>& a, const std::complex<T>& b) {                 \
0164   return complex_divide(a, b);                                                                  \
0165 }                                                                                               \
0166                                                                                                 \
0167 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0168 std::complex<T> operator/(const std::complex<T>& a, const T& b) {                               \
0169   return std::complex<T>(numext::real(a) / b, numext::imag(a) / b);                             \
0170 }                                                                                               \
0171                                                                                                 \
0172 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0173 std::complex<T> operator/(const T& a, const std::complex<T>& b) {                               \
0174   return complex_divide(std::complex<T>(a, 0), b);                                              \
0175 }                                                                                               \
0176                                                                                                 \
0177 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0178 std::complex<T>& operator+=(std::complex<T>& a, const std::complex<T>& b) {                     \
0179   numext::real_ref(a) += numext::real(b);                                                       \
0180   numext::imag_ref(a) += numext::imag(b);                                                       \
0181   return a;                                                                                     \
0182 }                                                                                               \
0183                                                                                                 \
0184 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0185 std::complex<T>& operator-=(std::complex<T>& a, const std::complex<T>& b) {                     \
0186   numext::real_ref(a) -= numext::real(b);                                                       \
0187   numext::imag_ref(a) -= numext::imag(b);                                                       \
0188   return a;                                                                                     \
0189 }                                                                                               \
0190                                                                                                 \
0191 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0192 std::complex<T>& operator*=(std::complex<T>& a, const std::complex<T>& b) {                     \
0193   a = complex_multiply(a, b);                                                                   \
0194   return a;                                                                                     \
0195 }                                                                                               \
0196                                                                                                 \
0197 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0198 std::complex<T>& operator/=(std::complex<T>& a, const std::complex<T>& b) {                     \
0199   a = complex_divide(a, b);                                                                     \
0200   return  a;                                                                                    \
0201 }                                                                                               \
0202                                                                                                 \
0203 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0204 bool operator==(const std::complex<T>& a, const std::complex<T>& b) {                           \
0205   return numext::real(a) == numext::real(b) && numext::imag(a) == numext::imag(b);              \
0206 }                                                                                               \
0207                                                                                                 \
0208 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0209 bool operator==(const std::complex<T>& a, const T& b) {                                         \
0210   return numext::real(a) == b && numext::imag(a) == 0;                                          \
0211 }                                                                                               \
0212                                                                                                 \
0213 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0214 bool operator==(const T& a, const std::complex<T>& b) {                                         \
0215   return a  == numext::real(b) && 0 == numext::imag(b);                                         \
0216 }                                                                                               \
0217                                                                                                 \
0218 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0219 bool operator!=(const std::complex<T>& a, const std::complex<T>& b) {                           \
0220   return !(a == b);                                                                             \
0221 }                                                                                               \
0222                                                                                                 \
0223 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0224 bool operator!=(const std::complex<T>& a, const T& b) {                                         \
0225   return !(a == b);                                                                             \
0226 }                                                                                               \
0227                                                                                                 \
0228 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE                                                           \
0229 bool operator!=(const T& a, const std::complex<T>& b) {                                         \
0230   return !(a == b);                                                                             \
0231 }
0232 
0233 // Do not specialize for long double, since that reduces to double on device.
0234 EIGEN_CREATE_STD_COMPLEX_OPERATOR_SPECIALIZATIONS(float)
0235 EIGEN_CREATE_STD_COMPLEX_OPERATOR_SPECIALIZATIONS(double)
0236 
0237 #undef EIGEN_CREATE_STD_COMPLEX_OPERATOR_SPECIALIZATIONS
0238 
0239   
0240 }  // namespace complex_operator_detail
0241 
0242 EIGEN_USING_STD_COMPLEX_OPERATORS
0243 
0244 namespace numext {
0245 EIGEN_USING_STD_COMPLEX_OPERATORS
0246 }  // namespace numext
0247 
0248 namespace internal {
0249 EIGEN_USING_STD_COMPLEX_OPERATORS
0250 
0251 }  // namespace internal
0252 }  // namespace Eigen
0253 
0254 #endif  // !(EIGEN_COMP_ICC && _USE_COMPLEX_SPECIALIZATION_)
0255 
0256 #endif  // EIGEN_CUDACC && EIGEN_GPU_COMPILE_PHASE
0257 
0258 #endif  // EIGEN_COMPLEX_CUDA_H
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!