#pragma once
#include <cmath>
#define INFINITY (__builtin_inff())
#define NAN (__builtin_nanf (""))
namespace std {
// Taken from libc++
template <class _Tp>
class complex {
public:
typedef _Tp value_type;
private:
value_type __re_;
value_type __im_;
public:
complex(const value_type &__re = value_type(), const value_type &__im = value_type())
: __re_(__re), __im_(__im) {}
template <class _Xp>
complex(const complex<_Xp> &__c)
: __re_(__c.real()), __im_(__c.imag()) {}
value_type real() const { return __re_; }
value_type imag() const { return __im_; }
void real(value_type __re) { __re_ = __re; }
void imag(value_type __im) { __im_ = __im; }
complex &operator=(const value_type &__re) {
__re_ = __re;
__im_ = value_type();
return *this;
}
complex &operator+=(const value_type &__re) {
__re_ += __re;
return *this;
}
complex &operator-=(const value_type &__re) {
__re_ -= __re;
return *this;
}
complex &operator*=(const value_type &__re) {
__re_ *= __re;
__im_ *= __re;
return *this;
}
complex &operator/=(const value_type &__re) {
__re_ /= __re;
__im_ /= __re;
return *this;
}
template <class _Xp>
complex &operator=(const complex<_Xp> &__c) {
__re_ = __c.real();
__im_ = __c.imag();
return *this;
}
template <class _Xp>
complex &operator+=(const complex<_Xp> &__c) {
__re_ += __c.real();
__im_ += __c.imag();
return *this;
}
template <class _Xp>
complex &operator-=(const complex<_Xp> &__c) {
__re_ -= __c.real();
__im_ -= __c.imag();
return *this;
}
template <class _Xp>
complex &operator*=(const complex<_Xp> &__c) {
*this = *this * complex(__c.real(), __c.imag());
return *this;
}
template <class _Xp>
complex &operator/=(const complex<_Xp> &__c) {
*this = *this / complex(__c.real(), __c.imag());
return *this;
}
};
template <class _Tp>
inline complex<_Tp>
operator+(const complex<_Tp> &__x, const complex<_Tp> &__y) {
complex<_Tp> __t(__x);
__t += __y;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator+(const complex<_Tp> &__x, const _Tp &__y) {
complex<_Tp> __t(__x);
__t += __y;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator+(const _Tp &__x, const complex<_Tp> &__y) {
complex<_Tp> __t(__y);
__t += __x;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator-(const complex<_Tp> &__x, const complex<_Tp> &__y) {
complex<_Tp> __t(__x);
__t -= __y;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator-(const complex<_Tp> &__x, const _Tp &__y) {
complex<_Tp> __t(__x);
__t -= __y;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator-(const _Tp &__x, const complex<_Tp> &__y) {
complex<_Tp> __t(-__y);
__t += __x;
return __t;
}
template <class _Tp>
complex<_Tp>
operator*(const complex<_Tp> &__z, const complex<_Tp> &__w) {
_Tp __a = __z.real();
_Tp __b = __z.imag();
_Tp __c = __w.real();
_Tp __d = __w.imag();
_Tp __ac = __a * __c;
_Tp __bd = __b * __d;
_Tp __ad = __a * __d;
_Tp __bc = __b * __c;
_Tp __x = __ac - __bd;
_Tp __y = __ad + __bc;
if (std::isnan(__x) && std::isnan(__y)) {
bool __recalc = false;
if (std::isinf(__a) || std::isinf(__b)) {
__a = copysign(std::isinf(__a) ? _Tp(1) : _Tp(0), __a);
__b = copysign(std::isinf(__b) ? _Tp(1) : _Tp(0), __b);
if (std::isnan(__c))
__c = copysign(_Tp(0), __c);
if (std::isnan(__d))
__d = copysign(_Tp(0), __d);
__recalc = true;
}
if (std::isinf(__c) || std::isinf(__d)) {
__c = copysign(std::isinf(__c) ? _Tp(1) : _Tp(0), __c);
__d = copysign(std::isinf(__d) ? _Tp(1) : _Tp(0), __d);
if (std::isnan(__a))
__a = copysign(_Tp(0), __a);
if (std::isnan(__b))
__b = copysign(_Tp(0), __b);
__recalc = true;
}
if (!__recalc && (std::isinf(__ac) || std::isinf(__bd) ||
std::isinf(__ad) || std::isinf(__bc))) {
if (std::isnan(__a))
__a = copysign(_Tp(0), __a);
if (std::isnan(__b))
__b = copysign(_Tp(0), __b);
if (std::isnan(__c))
__c = copysign(_Tp(0), __c);
if (std::isnan(__d))
__d = copysign(_Tp(0), __d);
__recalc = true;
}
if (__recalc) {
__x = _Tp(INFINITY) * (__a * __c - __b * __d);
__y = _Tp(INFINITY) * (__a * __d + __b * __c);
}
}
return complex<_Tp>(__x, __y);
}
template <class _Tp>
inline complex<_Tp>
operator*(const complex<_Tp> &__x, const _Tp &__y) {
complex<_Tp> __t(__x);
__t *= __y;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator*(const _Tp &__x, const complex<_Tp> &__y) {
complex<_Tp> __t(__y);
__t *= __x;
return __t;
}
template <class _Tp>
complex<_Tp>
operator/(const complex<_Tp> &__z, const complex<_Tp> &__w) {
int __ilogbw = 0;
_Tp __a = __z.real();
_Tp __b = __z.imag();
_Tp __c = __w.real();
_Tp __d = __w.imag();
_Tp __logbw = logb(fmax(fabs(__c), fabs(__d)));
if (std::isfinite(__logbw)) {
__ilogbw = static_cast<int>(__logbw);
__c = scalbn(__c, -__ilogbw);
__d = scalbn(__d, -__ilogbw);
}
_Tp __denom = __c * __c + __d * __d;
_Tp __x = scalbn((__a * __c + __b * __d) / __denom, -__ilogbw);
_Tp __y = scalbn((__b * __c - __a * __d) / __denom, -__ilogbw);
if (std::isnan(__x) && std::isnan(__y)) {
if ((__denom == _Tp(0)) && (!std::isnan(__a) || !std::isnan(__b))) {
__x = copysign(_Tp(INFINITY), __c) * __a;
__y = copysign(_Tp(INFINITY), __c) * __b;
} else if ((std::isinf(__a) || std::isinf(__b)) && std::isfinite(__c) && std::isfinite(__d)) {
__a = copysign(std::isinf(__a) ? _Tp(1) : _Tp(0), __a);
__b = copysign(std::isinf(__b) ? _Tp(1) : _Tp(0), __b);
__x = _Tp(INFINITY) * (__a * __c + __b * __d);
__y = _Tp(INFINITY) * (__b * __c - __a * __d);
} else if (std::isinf(__logbw) && __logbw > _Tp(0) && std::isfinite(__a) && std::isfinite(__b)) {
__c = copysign(std::isinf(__c) ? _Tp(1) : _Tp(0), __c);
__d = copysign(std::isinf(__d) ? _Tp(1) : _Tp(0), __d);
__x = _Tp(0) * (__a * __c + __b * __d);
__y = _Tp(0) * (__b * __c - __a * __d);
}
}
return complex<_Tp>(__x, __y);
}
template <class _Tp>
inline complex<_Tp>
operator/(const complex<_Tp> &__x, const _Tp &__y) {
return complex<_Tp>(__x.real() / __y, __x.imag() / __y);
}
template <class _Tp>
inline complex<_Tp>
operator/(const _Tp &__x, const complex<_Tp> &__y) {
complex<_Tp> __t(__x);
__t /= __y;
return __t;
}
template <class _Tp>
inline complex<_Tp>
operator+(const complex<_Tp> &__x) {
return __x;
}
template <class _Tp>
inline complex<_Tp>
operator-(const complex<_Tp> &__x) {
return complex<_Tp>(-__x.real(), -__x.imag());
}
template <class _Tp>
inline bool
operator==(const complex<_Tp> &__x, const complex<_Tp> &__y) {
return __x.real() == __y.real() && __x.imag() == __y.imag();
}
template <class _Tp>
inline bool
operator==(const complex<_Tp> &__x, const _Tp &__y) {
return __x.real() == __y && __x.imag() == 0;
}
template <class _Tp>
inline bool
operator==(const _Tp &__x, const complex<_Tp> &__y) {
return __x == __y.real() && 0 == __y.imag();
}
template <class _Tp>
inline bool
operator!=(const complex<_Tp> &__x, const complex<_Tp> &__y) {
return !(__x == __y);
}
template <class _Tp>
inline bool
operator!=(const complex<_Tp> &__x, const _Tp &__y) {
return !(__x == __y);
}
template <class _Tp>
inline bool
operator!=(const _Tp &__x, const complex<_Tp> &__y) {
return !(__x == __y);
}
template <class _Tp> _Tp abs(const std::complex<_Tp> &__c);
// arg
template <class _Tp> _Tp arg(const std::complex<_Tp> &__c);
// norm
template <class _Tp> _Tp norm(const std::complex<_Tp> &__c);
// conj
template <class _Tp> std::complex<_Tp> conj(const std::complex<_Tp> &__c);
// proj
template <class _Tp> std::complex<_Tp> proj(const std::complex<_Tp> &__c);
// polar
template <class _Tp>
complex<_Tp> polar(const _Tp &__rho, const _Tp &__theta = _Tp());
// log
template <class _Tp> std::complex<_Tp> log(const std::complex<_Tp> &__x);
// log10
template <class _Tp> std::complex<_Tp> log10(const std::complex<_Tp> &__x);
// sqrt
template <class _Tp>
std::complex<_Tp> sqrt(const std::complex<_Tp> &__x);
// exp
template <class _Tp>
std::complex<_Tp> exp(const std::complex<_Tp> &__x);
// pow
template <class _Tp>
std::complex<_Tp> pow(const std::complex<_Tp> &__x,
const std::complex<_Tp> &__y);
// __sqr, computes pow(x, 2)
template <class _Tp> std::complex<_Tp> __sqr(const std::complex<_Tp> &__x);
// asinh
template <class _Tp>
std::complex<_Tp> asinh(const std::complex<_Tp> &__x);
// acosh
template <class _Tp>
std::complex<_Tp> acosh(const std::complex<_Tp> &__x);
// atanh
template <class _Tp>
std::complex<_Tp> atanh(const std::complex<_Tp> &__x);
// sinh
template <class _Tp>
std::complex<_Tp> sinh(const std::complex<_Tp> &__x);
// cosh
template <class _Tp>
std::complex<_Tp> cosh(const std::complex<_Tp> &__x);
// tanh
template <class _Tp>
std::complex<_Tp> tanh(const std::complex<_Tp> &__x);
// asin
template <class _Tp>
std::complex<_Tp> asin(const std::complex<_Tp> &__x);
// acos
template <class _Tp>
std::complex<_Tp> acos(const std::complex<_Tp> &__x);
// atan
template <class _Tp>
std::complex<_Tp> atan(const std::complex<_Tp> &__x);
// sin
template <class _Tp>
std::complex<_Tp> sin(const std::complex<_Tp> &__x);
// cos
template <class _Tp> std::complex<_Tp> cos(const std::complex<_Tp> &__x);
// tan
template <class _Tp>
std::complex<_Tp> tan(const std::complex<_Tp> &__x);
} // namespace std