Subversion Repositories QNX 8.QNX8 LLVM/Clang compiler suite

Rev

Blame | Last modification | View Log | Download | RSS feed

  1. //===------------------------- __complex_cmath.h --------------------------===//
  2. //
  3. // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
  4. // See https://llvm.org/LICENSE.txt for license information.
  5. // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
  6. //
  7. //===----------------------------------------------------------------------===//
  8. //
  9. // std::complex header copied from the libcxx source and simplified for use in
  10. // OpenMP target offload regions.
  11. //
  12. //===----------------------------------------------------------------------===//
  13.  
  14. #ifndef _OPENMP
  15. #error "This file is for OpenMP compilation only."
  16. #endif
  17.  
  18. #ifndef __cplusplus
  19. #error "This file is for C++ compilation only."
  20. #endif
  21.  
  22. #ifndef _LIBCPP_COMPLEX
  23. #define _LIBCPP_COMPLEX
  24.  
  25. #include <cmath>
  26. #include <type_traits>
  27.  
  28. #define __DEVICE__ static constexpr __attribute__((nothrow))
  29.  
  30. namespace std {
  31.  
  32. // abs
  33.  
  34. template <class _Tp> __DEVICE__ _Tp abs(const std::complex<_Tp> &__c) {
  35.   return hypot(__c.real(), __c.imag());
  36. }
  37.  
  38. // arg
  39.  
  40. template <class _Tp> __DEVICE__ _Tp arg(const std::complex<_Tp> &__c) {
  41.   return atan2(__c.imag(), __c.real());
  42. }
  43.  
  44. template <class _Tp>
  45. typename enable_if<is_integral<_Tp>::value || is_same<_Tp, double>::value,
  46.                    double>::type
  47. arg(_Tp __re) {
  48.   return atan2(0., __re);
  49. }
  50.  
  51. template <class _Tp>
  52. typename enable_if<is_same<_Tp, float>::value, float>::type arg(_Tp __re) {
  53.   return atan2f(0.F, __re);
  54. }
  55.  
  56. // norm
  57.  
  58. template <class _Tp> __DEVICE__ _Tp norm(const std::complex<_Tp> &__c) {
  59.   if (std::isinf(__c.real()))
  60.     return abs(__c.real());
  61.   if (std::isinf(__c.imag()))
  62.     return abs(__c.imag());
  63.   return __c.real() * __c.real() + __c.imag() * __c.imag();
  64. }
  65.  
  66. // conj
  67.  
  68. template <class _Tp> std::complex<_Tp> conj(const std::complex<_Tp> &__c) {
  69.   return std::complex<_Tp>(__c.real(), -__c.imag());
  70. }
  71.  
  72. // proj
  73.  
  74. template <class _Tp> std::complex<_Tp> proj(const std::complex<_Tp> &__c) {
  75.   std::complex<_Tp> __r = __c;
  76.   if (std::isinf(__c.real()) || std::isinf(__c.imag()))
  77.     __r = std::complex<_Tp>(INFINITY, copysign(_Tp(0), __c.imag()));
  78.   return __r;
  79. }
  80.  
  81. // polar
  82.  
  83. template <class _Tp>
  84. complex<_Tp> polar(const _Tp &__rho, const _Tp &__theta = _Tp()) {
  85.   if (std::isnan(__rho) || signbit(__rho))
  86.     return std::complex<_Tp>(_Tp(NAN), _Tp(NAN));
  87.   if (std::isnan(__theta)) {
  88.     if (std::isinf(__rho))
  89.       return std::complex<_Tp>(__rho, __theta);
  90.     return std::complex<_Tp>(__theta, __theta);
  91.   }
  92.   if (std::isinf(__theta)) {
  93.     if (std::isinf(__rho))
  94.       return std::complex<_Tp>(__rho, _Tp(NAN));
  95.     return std::complex<_Tp>(_Tp(NAN), _Tp(NAN));
  96.   }
  97.   _Tp __x = __rho * cos(__theta);
  98.   if (std::isnan(__x))
  99.     __x = 0;
  100.   _Tp __y = __rho * sin(__theta);
  101.   if (std::isnan(__y))
  102.     __y = 0;
  103.   return std::complex<_Tp>(__x, __y);
  104. }
  105.  
  106. // log
  107.  
  108. template <class _Tp> std::complex<_Tp> log(const std::complex<_Tp> &__x) {
  109.   return std::complex<_Tp>(log(abs(__x)), arg(__x));
  110. }
  111.  
  112. // log10
  113.  
  114. template <class _Tp> std::complex<_Tp> log10(const std::complex<_Tp> &__x) {
  115.   return log(__x) / log(_Tp(10));
  116. }
  117.  
  118. // sqrt
  119.  
  120. template <class _Tp>
  121. __DEVICE__ std::complex<_Tp> sqrt(const std::complex<_Tp> &__x) {
  122.   if (std::isinf(__x.imag()))
  123.     return std::complex<_Tp>(_Tp(INFINITY), __x.imag());
  124.   if (std::isinf(__x.real())) {
  125.     if (__x.real() > _Tp(0))
  126.       return std::complex<_Tp>(__x.real(), std::isnan(__x.imag())
  127.                                                ? __x.imag()
  128.                                                : copysign(_Tp(0), __x.imag()));
  129.     return std::complex<_Tp>(std::isnan(__x.imag()) ? __x.imag() : _Tp(0),
  130.                              copysign(__x.real(), __x.imag()));
  131.   }
  132.   return polar(sqrt(abs(__x)), arg(__x) / _Tp(2));
  133. }
  134.  
  135. // exp
  136.  
  137. template <class _Tp>
  138. __DEVICE__ std::complex<_Tp> exp(const std::complex<_Tp> &__x) {
  139.   _Tp __i = __x.imag();
  140.   if (std::isinf(__x.real())) {
  141.     if (__x.real() < _Tp(0)) {
  142.       if (!std::isfinite(__i))
  143.         __i = _Tp(1);
  144.     } else if (__i == 0 || !std::isfinite(__i)) {
  145.       if (std::isinf(__i))
  146.         __i = _Tp(NAN);
  147.       return std::complex<_Tp>(__x.real(), __i);
  148.     }
  149.   } else if (std::isnan(__x.real()) && __x.imag() == 0)
  150.     return __x;
  151.   _Tp __e = exp(__x.real());
  152.   return std::complex<_Tp>(__e * cos(__i), __e * sin(__i));
  153. }
  154.  
  155. // pow
  156.  
  157. template <class _Tp>
  158. std::complex<_Tp> pow(const std::complex<_Tp> &__x,
  159.                       const std::complex<_Tp> &__y) {
  160.   return exp(__y * log(__x));
  161. }
  162.  
  163. // __sqr, computes pow(x, 2)
  164.  
  165. template <class _Tp> std::complex<_Tp> __sqr(const std::complex<_Tp> &__x) {
  166.   return std::complex<_Tp>((__x.real() - __x.imag()) *
  167.                                (__x.real() + __x.imag()),
  168.                            _Tp(2) * __x.real() * __x.imag());
  169. }
  170.  
  171. // asinh
  172.  
  173. template <class _Tp>
  174. __DEVICE__ std::complex<_Tp> asinh(const std::complex<_Tp> &__x) {
  175.   const _Tp __pi(atan2(+0., -0.));
  176.   if (std::isinf(__x.real())) {
  177.     if (std::isnan(__x.imag()))
  178.       return __x;
  179.     if (std::isinf(__x.imag()))
  180.       return std::complex<_Tp>(__x.real(),
  181.                                copysign(__pi * _Tp(0.25), __x.imag()));
  182.     return std::complex<_Tp>(__x.real(), copysign(_Tp(0), __x.imag()));
  183.   }
  184.   if (std::isnan(__x.real())) {
  185.     if (std::isinf(__x.imag()))
  186.       return std::complex<_Tp>(__x.imag(), __x.real());
  187.     if (__x.imag() == 0)
  188.       return __x;
  189.     return std::complex<_Tp>(__x.real(), __x.real());
  190.   }
  191.   if (std::isinf(__x.imag()))
  192.     return std::complex<_Tp>(copysign(__x.imag(), __x.real()),
  193.                              copysign(__pi / _Tp(2), __x.imag()));
  194.   std::complex<_Tp> __z = log(__x + sqrt(__sqr(__x) + _Tp(1)));
  195.   return std::complex<_Tp>(copysign(__z.real(), __x.real()),
  196.                            copysign(__z.imag(), __x.imag()));
  197. }
  198.  
  199. // acosh
  200.  
  201. template <class _Tp>
  202. __DEVICE__ std::complex<_Tp> acosh(const std::complex<_Tp> &__x) {
  203.   const _Tp __pi(atan2(+0., -0.));
  204.   if (std::isinf(__x.real())) {
  205.     if (std::isnan(__x.imag()))
  206.       return std::complex<_Tp>(abs(__x.real()), __x.imag());
  207.     if (std::isinf(__x.imag())) {
  208.       if (__x.real() > 0)
  209.         return std::complex<_Tp>(__x.real(),
  210.                                  copysign(__pi * _Tp(0.25), __x.imag()));
  211.       else
  212.         return std::complex<_Tp>(-__x.real(),
  213.                                  copysign(__pi * _Tp(0.75), __x.imag()));
  214.     }
  215.     if (__x.real() < 0)
  216.       return std::complex<_Tp>(-__x.real(), copysign(__pi, __x.imag()));
  217.     return std::complex<_Tp>(__x.real(), copysign(_Tp(0), __x.imag()));
  218.   }
  219.   if (std::isnan(__x.real())) {
  220.     if (std::isinf(__x.imag()))
  221.       return std::complex<_Tp>(abs(__x.imag()), __x.real());
  222.     return std::complex<_Tp>(__x.real(), __x.real());
  223.   }
  224.   if (std::isinf(__x.imag()))
  225.     return std::complex<_Tp>(abs(__x.imag()),
  226.                              copysign(__pi / _Tp(2), __x.imag()));
  227.   std::complex<_Tp> __z = log(__x + sqrt(__sqr(__x) - _Tp(1)));
  228.   return std::complex<_Tp>(copysign(__z.real(), _Tp(0)),
  229.                            copysign(__z.imag(), __x.imag()));
  230. }
  231.  
  232. // atanh
  233.  
  234. template <class _Tp>
  235. __DEVICE__ std::complex<_Tp> atanh(const std::complex<_Tp> &__x) {
  236.   const _Tp __pi(atan2(+0., -0.));
  237.   if (std::isinf(__x.imag())) {
  238.     return std::complex<_Tp>(copysign(_Tp(0), __x.real()),
  239.                              copysign(__pi / _Tp(2), __x.imag()));
  240.   }
  241.   if (std::isnan(__x.imag())) {
  242.     if (std::isinf(__x.real()) || __x.real() == 0)
  243.       return std::complex<_Tp>(copysign(_Tp(0), __x.real()), __x.imag());
  244.     return std::complex<_Tp>(__x.imag(), __x.imag());
  245.   }
  246.   if (std::isnan(__x.real())) {
  247.     return std::complex<_Tp>(__x.real(), __x.real());
  248.   }
  249.   if (std::isinf(__x.real())) {
  250.     return std::complex<_Tp>(copysign(_Tp(0), __x.real()),
  251.                              copysign(__pi / _Tp(2), __x.imag()));
  252.   }
  253.   if (abs(__x.real()) == _Tp(1) && __x.imag() == _Tp(0)) {
  254.     return std::complex<_Tp>(copysign(_Tp(INFINITY), __x.real()),
  255.                              copysign(_Tp(0), __x.imag()));
  256.   }
  257.   std::complex<_Tp> __z = log((_Tp(1) + __x) / (_Tp(1) - __x)) / _Tp(2);
  258.   return std::complex<_Tp>(copysign(__z.real(), __x.real()),
  259.                            copysign(__z.imag(), __x.imag()));
  260. }
  261.  
  262. // sinh
  263.  
  264. template <class _Tp>
  265. __DEVICE__ std::complex<_Tp> sinh(const std::complex<_Tp> &__x) {
  266.   if (std::isinf(__x.real()) && !std::isfinite(__x.imag()))
  267.     return std::complex<_Tp>(__x.real(), _Tp(NAN));
  268.   if (__x.real() == 0 && !std::isfinite(__x.imag()))
  269.     return std::complex<_Tp>(__x.real(), _Tp(NAN));
  270.   if (__x.imag() == 0 && !std::isfinite(__x.real()))
  271.     return __x;
  272.   return std::complex<_Tp>(sinh(__x.real()) * cos(__x.imag()),
  273.                            cosh(__x.real()) * sin(__x.imag()));
  274. }
  275.  
  276. // cosh
  277.  
  278. template <class _Tp>
  279. __DEVICE__ std::complex<_Tp> cosh(const std::complex<_Tp> &__x) {
  280.   if (std::isinf(__x.real()) && !std::isfinite(__x.imag()))
  281.     return std::complex<_Tp>(abs(__x.real()), _Tp(NAN));
  282.   if (__x.real() == 0 && !std::isfinite(__x.imag()))
  283.     return std::complex<_Tp>(_Tp(NAN), __x.real());
  284.   if (__x.real() == 0 && __x.imag() == 0)
  285.     return std::complex<_Tp>(_Tp(1), __x.imag());
  286.   if (__x.imag() == 0 && !std::isfinite(__x.real()))
  287.     return std::complex<_Tp>(abs(__x.real()), __x.imag());
  288.   return std::complex<_Tp>(cosh(__x.real()) * cos(__x.imag()),
  289.                            sinh(__x.real()) * sin(__x.imag()));
  290. }
  291.  
  292. // tanh
  293.  
  294. template <class _Tp>
  295. __DEVICE__ std::complex<_Tp> tanh(const std::complex<_Tp> &__x) {
  296.   if (std::isinf(__x.real())) {
  297.     if (!std::isfinite(__x.imag()))
  298.       return std::complex<_Tp>(_Tp(1), _Tp(0));
  299.     return std::complex<_Tp>(_Tp(1),
  300.                              copysign(_Tp(0), sin(_Tp(2) * __x.imag())));
  301.   }
  302.   if (std::isnan(__x.real()) && __x.imag() == 0)
  303.     return __x;
  304.   _Tp __2r(_Tp(2) * __x.real());
  305.   _Tp __2i(_Tp(2) * __x.imag());
  306.   _Tp __d(cosh(__2r) + cos(__2i));
  307.   _Tp __2rsh(sinh(__2r));
  308.   if (std::isinf(__2rsh) && std::isinf(__d))
  309.     return std::complex<_Tp>(__2rsh > _Tp(0) ? _Tp(1) : _Tp(-1),
  310.                              __2i > _Tp(0) ? _Tp(0) : _Tp(-0.));
  311.   return std::complex<_Tp>(__2rsh / __d, sin(__2i) / __d);
  312. }
  313.  
  314. // asin
  315.  
  316. template <class _Tp>
  317. __DEVICE__ std::complex<_Tp> asin(const std::complex<_Tp> &__x) {
  318.   std::complex<_Tp> __z = asinh(complex<_Tp>(-__x.imag(), __x.real()));
  319.   return std::complex<_Tp>(__z.imag(), -__z.real());
  320. }
  321.  
  322. // acos
  323.  
  324. template <class _Tp>
  325. __DEVICE__ std::complex<_Tp> acos(const std::complex<_Tp> &__x) {
  326.   const _Tp __pi(atan2(+0., -0.));
  327.   if (std::isinf(__x.real())) {
  328.     if (std::isnan(__x.imag()))
  329.       return std::complex<_Tp>(__x.imag(), __x.real());
  330.     if (std::isinf(__x.imag())) {
  331.       if (__x.real() < _Tp(0))
  332.         return std::complex<_Tp>(_Tp(0.75) * __pi, -__x.imag());
  333.       return std::complex<_Tp>(_Tp(0.25) * __pi, -__x.imag());
  334.     }
  335.     if (__x.real() < _Tp(0))
  336.       return std::complex<_Tp>(__pi,
  337.                                signbit(__x.imag()) ? -__x.real() : __x.real());
  338.     return std::complex<_Tp>(_Tp(0),
  339.                              signbit(__x.imag()) ? __x.real() : -__x.real());
  340.   }
  341.   if (std::isnan(__x.real())) {
  342.     if (std::isinf(__x.imag()))
  343.       return std::complex<_Tp>(__x.real(), -__x.imag());
  344.     return std::complex<_Tp>(__x.real(), __x.real());
  345.   }
  346.   if (std::isinf(__x.imag()))
  347.     return std::complex<_Tp>(__pi / _Tp(2), -__x.imag());
  348.   if (__x.real() == 0 && (__x.imag() == 0 || isnan(__x.imag())))
  349.     return std::complex<_Tp>(__pi / _Tp(2), -__x.imag());
  350.   std::complex<_Tp> __z = log(__x + sqrt(__sqr(__x) - _Tp(1)));
  351.   if (signbit(__x.imag()))
  352.     return std::complex<_Tp>(abs(__z.imag()), abs(__z.real()));
  353.   return std::complex<_Tp>(abs(__z.imag()), -abs(__z.real()));
  354. }
  355.  
  356. // atan
  357.  
  358. template <class _Tp>
  359. __DEVICE__ std::complex<_Tp> atan(const std::complex<_Tp> &__x) {
  360.   std::complex<_Tp> __z = atanh(complex<_Tp>(-__x.imag(), __x.real()));
  361.   return std::complex<_Tp>(__z.imag(), -__z.real());
  362. }
  363.  
  364. // sin
  365.  
  366. template <class _Tp>
  367. __DEVICE__ std::complex<_Tp> sin(const std::complex<_Tp> &__x) {
  368.   std::complex<_Tp> __z = sinh(complex<_Tp>(-__x.imag(), __x.real()));
  369.   return std::complex<_Tp>(__z.imag(), -__z.real());
  370. }
  371.  
  372. // cos
  373.  
  374. template <class _Tp> std::complex<_Tp> cos(const std::complex<_Tp> &__x) {
  375.   return cosh(complex<_Tp>(-__x.imag(), __x.real()));
  376. }
  377.  
  378. // tan
  379.  
  380. template <class _Tp>
  381. __DEVICE__ std::complex<_Tp> tan(const std::complex<_Tp> &__x) {
  382.   std::complex<_Tp> __z = tanh(complex<_Tp>(-__x.imag(), __x.real()));
  383.   return std::complex<_Tp>(__z.imag(), -__z.real());
  384. }
  385.  
  386. } // namespace std
  387.  
  388. #endif
  389.