[libcxx] Improve accuracy of complex asinh and acosh

Summary:
Currently std::asinh and std::acosh use std::pow to compute x^2. This
results in a significant error when computing e.g. asinh(i) or
acosh(-1).

This patch expresses x^2 directly via x.real() and x.imag(), like it
is done in libstdc++/glibc, and adds tests that checks the accuracy.

Reviewers: EricWF, mclow.lists

Reviewed By: mclow.lists

Subscribers: christof, cfe-commits

Differential Revision: https://reviews.llvm.org/D41629

git-svn-id: https://llvm.org/svn/llvm-project/libcxx/trunk@325510 91177308-0d34-0410-b5e6-96231b3b80d8
This commit is contained in:
Mikhail Maltsev
2018-02-19 15:41:36 +00:00
parent 6878e852d1
commit c4658abe60
4 changed files with 113 additions and 3 deletions

View File

@@ -1125,6 +1125,17 @@ pow(const _Tp& __x, const complex<_Up>& __y)
return _VSTD::pow(result_type(__x), result_type(__y));
}
// __sqr, computes pow(x, 2)
template<class _Tp>
inline _LIBCPP_INLINE_VISIBILITY
complex<_Tp>
__sqr(const complex<_Tp>& __x)
{
return complex<_Tp>((__x.real() - __x.imag()) * (__x.real() + __x.imag()),
_Tp(2) * __x.real() * __x.imag());
}
// asinh
template<class _Tp>
@@ -1150,7 +1161,7 @@ asinh(const complex<_Tp>& __x)
}
if (__libcpp_isinf_or_builtin(__x.imag()))
return complex<_Tp>(copysign(__x.imag(), __x.real()), copysign(__pi/_Tp(2), __x.imag()));
complex<_Tp> __z = log(__x + sqrt(pow(__x, _Tp(2)) + _Tp(1)));
complex<_Tp> __z = log(__x + sqrt(__sqr(__x) + _Tp(1)));
return complex<_Tp>(copysign(__z.real(), __x.real()), copysign(__z.imag(), __x.imag()));
}
@@ -1184,7 +1195,7 @@ acosh(const complex<_Tp>& __x)
}
if (__libcpp_isinf_or_builtin(__x.imag()))
return complex<_Tp>(abs(__x.imag()), copysign(__pi/_Tp(2), __x.imag()));
complex<_Tp> __z = log(__x + sqrt(pow(__x, _Tp(2)) - _Tp(1)));
complex<_Tp> __z = log(__x + sqrt(__sqr(__x) - _Tp(1)));
return complex<_Tp>(copysign(__z.real(), _Tp(0)), copysign(__z.imag(), __x.imag()));
}
@@ -1318,7 +1329,7 @@ acos(const complex<_Tp>& __x)
return complex<_Tp>(__pi/_Tp(2), -__x.imag());
if (__x.real() == 0 && (__x.imag() == 0 || isnan(__x.imag())))
return complex<_Tp>(__pi/_Tp(2), -__x.imag());
complex<_Tp> __z = log(__x + sqrt(pow(__x, _Tp(2)) - _Tp(1)));
complex<_Tp> __z = log(__x + sqrt(__sqr(__x) - _Tp(1)));
if (signbit(__x.imag()))
return complex<_Tp>(abs(__z.imag()), abs(__z.real()));
return complex<_Tp>(abs(__z.imag()), -abs(__z.real()));