Incorrect results when using complex dtype in numpy

Some of the calculations that use complex numbers in numpy produce results that are incorrect, which makes the corresponding numpy tests fail.

```
_________________________ TestComplexFunctions.test_it _________________________

self =

def test_it(self):
for f in self.funcs:
if f is np.arccosh:
x = 1.5
else:
x = .5
fr = f(x)
fz = f(complex(x))
> assert_almost_equal(fz.real, fr, err_msg='real part %s' % f)
E AssertionError:
E Arrays are not almost equal to 7 decimals real part
E ACTUAL: -1.1102230246251565e-16
E DESIRED: 0.5235987755982989

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2424: AssertionError
____________________ TestComplexFunctions.test_branch_cuts _____________________

self =

def test_branch_cuts(self):
# check branch cuts and continuity on them
_check_branch_cut(np.log, -0.5, 1j, 1, -1, True)
_check_branch_cut(np.log2, -0.5, 1j, 1, -1, True)
_check_branch_cut(np.log10, -0.5, 1j, 1, -1, True)
_check_branch_cut(np.log1p, -1.5, 1j, 1, -1, True)
_check_branch_cut(np.sqrt, -0.5, 1j, 1, -1, True)

> _check_branch_cut(np.arcsin, [ -2, 2], [1j, 1j], 1, -1, True)

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2444:
f = , x0 = array([-2.+0.j, 2.+0.j])
dx = array([0.+1.j, 0.+1.j]), re_sign = 1, im_sign = -1, sig_zero_ok = True
dtype =

def _check_branch_cut(f, x0, dx, re_sign=1, im_sign=-1, sig_zero_ok=False,
dtype=complex):
"""
Check for a branch cut in a function.

Assert that `x0` lies on a branch cut of function `f` and `f` is
continuous from the direction `dx`.

Parameters
----------
f : func
Function to check
x0 : array-like
Point on branch cut
dx : array-like
Direction to check continuity in
re_sign, im_sign : {1, -1}
Change of sign of the real or imaginary part expected
sig_zero_ok : bool
Whether to check if the branch cut respects signed zero (if applicable)
dtype : dtype
Dtype to check (should be complex)

"""
x0 = np.atleast_1d(x0).astype(dtype)
dx = np.atleast_1d(dx).astype(dtype)
if np.dtype(dtype).char == 'F':
scale = np.finfo(dtype).eps * 1e2
atol = np.float32(1e-2)
else:
scale = np.finfo(dtype).eps * 1e3
atol = 1e-4

y0 = f(x0)
yp = f(x0 + dx*scale*np.absolute(x0)/np.absolute(dx))
ym = f(x0 - dx*scale*np.absolute(x0)/np.absolute(dx))

assert_(np.all(np.absolute(y0.real - yp.real) < atol), (y0, yp))
assert_(np.all(np.absolute(y0.imag - yp.imag) < atol), (y0, yp))
> assert_(np.all(np.absolute(y0.real - ym.real*re_sign) < atol), (y0, ym))
E AssertionError: (array([-1.3169579-1.57079633j, -1.3169579+1.57079633j]), array([1.3169579-1.57079633j, 1.3169579+1.57079633j]))

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2672: AssertionError
_______________ TestComplexFunctions.test_branch_cuts_complex64 ________________

self =

def test_branch_cuts_complex64(self):
# check branch cuts and continuity on them
_check_branch_cut(np.log, -0.5, 1j, 1, -1, True, np.complex64)
_check_branch_cut(np.log2, -0.5, 1j, 1, -1, True, np.complex64)
_check_branch_cut(np.log10, -0.5, 1j, 1, -1, True, np.complex64)
_check_branch_cut(np.log1p, -1.5, 1j, 1, -1, True, np.complex64)
_check_branch_cut(np.sqrt, -0.5, 1j, 1, -1, True, np.complex64)
> _check_branch_cut(np.arcsin, [ -2, 2], [1j, 1j], 1, -1, True, np.complex64)

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2469:
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

f = , x0 = array([-2.+0.j, 2.+0.j], dtype=complex64)
dx = array([0.+1.j, 0.+1.j], dtype=complex64), re_sign = 1, im_sign = -1
sig_zero_ok = True, dtype =

def _check_branch_cut(f, x0, dx, re_sign=1, im_sign=-1, sig_zero_ok=False,
dtype=complex):
"""
Check for a branch cut in a function.

Assert that `x0` lies on a branch cut of function `f` and `f` is
continuous from the direction `dx`.

Parameters
----------
f : func
Function to check
x0 : array-like
Point on branch cut
dx : array-like
Direction to check continuity in
re_sign, im_sign : {1, -1}
Change of sign of the real or imaginary part expected
sig_zero_ok : bool
Whether to check if the branch cut respects signed zero (if applicable)
dtype : dtype
Dtype to check (should be complex)

"""
x0 = np.atleast_1d(x0).astype(dtype)
dx = np.atleast_1d(dx).astype(dtype)
dx = np.atleast_1d(dx).astype(dtype)

if np.dtype(dtype).char == 'F':
scale = np.finfo(dtype).eps * 1e2
atol = np.float32(1e-2)
else:
scale = np.finfo(dtype).eps * 1e3
atol = 1e-4

y0 = f(x0)
yp = f(x0 + dx*scale*np.absolute(x0)/np.absolute(dx))
ym = f(x0 - dx*scale*np.absolute(x0)/np.absolute(dx))

assert_(np.all(np.absolute(y0.real - yp.real) < atol), (y0, yp))
assert_(np.all(np.absolute(y0.imag - yp.imag) < atol), (y0, yp))
> assert_(np.all(np.absolute(y0.real - ym.real*re_sign) < atol), (y0, ym))
E AssertionError: (array([-1.3169578-1.5707964j, -1.3169578+1.5707964j], dtype=complex64), array([1.316958-1.5707825j, 1.316958+1.5707825j], dtype=comp
lex64))
___________________ TestComplexFunctions.test_against_cmath ____________________

self =

def test_against_cmath(self):
import cmath

points = [-1-1j, -1+1j, +1-1j, +1+1j]
name_map = {'arcsin': 'asin', 'arccos': 'acos', 'arctan': 'atan',
'arcsinh': 'asinh', 'arccosh': 'acosh', 'arctanh': 'atanh'}
atol = 4*np.finfo(complex).eps
for func in self.funcs:
fname = func.__name__.split('.')[-1]
cname = name_map.get(fname, fname)
try:
cfunc = getattr(cmath, cname)
except AttributeError:
continue
for p in points:
a = complex(func(np.complex_(p)))
b = cfunc(p)
> assert_(abs(a - b) < atol, "%s %s: %s; cmath: %s" % (fname, p, a, b))
E AssertionError: arcsin (-1-1j): (1.0612750619050357-0.6662394324925153j); cmath: (-0.6662394324925153-1.0612750619050357j)

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2503: AssertionError
_________________ TestComplexFunctions.test_loss_of_precision __________________

self =

def test_loss_of_precision(self):
for dtype in [np.complex64, np.complex_]:
> self.check_loss_of_precision(dtype)

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2592:
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2549: in check_loss_of_precision
check(x_series, 2.1*eps)
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

x = array([9.99999968e-21, 1.21736865e-20, 1.48198648e-20, 1.80412373e-20,
2.19628376e-20, 2.67368701e-20, 3.254862...55170e-04, 4.54267400e-04,
5.53010846e-04, 6.73218106e-04, 8.19554611e-04, 9.97700030e-04],
dtype=float32)
rtol = 2.5033950805664064e-07

def check(x, rtol):
x = x.astype(real_dtype)

z = x.astype(dtype)
d = np.absolute(np.arcsinh(x)/np.arcsinh(z).real - 1)
assert_(np.all(d < rtol), (np.argmax(d), x[np.argmax(d)], d.max(),
> 'arcsinh'))
E AssertionError: (0, 1e-20, inf, 'arcsinh')

lib/python3.7/site-packages/numpy/core/tests/test_umath.py:2520: AssertionError
```

A number of other tests that use complex dtypes pass however.

Full test logs can be found in https://github.com/iodide-project/pyodide/issues/69#issuecomment-420616635
Tested with Python 3.7 and numpy 1.15.1.