Commit fbee59ac authored by Martin Reinecke's avatar Martin Reinecke
Browse files

tweak tests

parent c2c6ace2
......@@ -2591,14 +2591,7 @@ py::array mycomplex2hartley(const py::array &in,
else throw runtime_error("unsupported data type");
}
py::array hartley2(const py::array &in, py::object axes_, double fct)
{
auto tmp = rfftn(in, axes_, fct);
if (in.dtype().is(f64))
return complex2hartley<double>(in, tmp, axes_);
else if (in.dtype().is(f32))
return complex2hartley<float>(in, tmp, axes_);
else throw runtime_error("unsupported data type");
}
{ return mycomplex2hartley(in, rfftn(in, axes_, fct), axes_); }
const char *pypocketfft_DS = R"DELIM(Fast Fourier and Hartley transforms.
......
......@@ -2,7 +2,6 @@ import pypocketfft
import pyfftw
import numpy as np
import pytest
from numpy.testing import assert_allclose
pmp = pytest.mark.parametrize
......@@ -10,58 +9,61 @@ shapes1D = ((10,), (127,))
shapes2D = ((128, 128), (128, 129), (1, 129), (129,1))
shapes3D = ((32,17,39),)
shapes = shapes1D+shapes2D+shapes3D
len1D=range(1,2048)
def test1D():
maxlen=8192
a=np.random.rand(maxlen)-0.5 + 1j*np.random.rand(maxlen)-0.5j
def _l2error(a,b):
return np.sqrt(np.sum(np.abs(a-b)**2)/np.sum(np.abs(a)**2))
@pmp("len", len1D)
def test1D(len):
a=np.random.rand(len)-0.5 + 1j*np.random.rand(len)-0.5j
b=a.astype(np.complex64)
for i in range(1,maxlen+1):
assert_allclose(pypocketfft.ifftn(pypocketfft.fftn(a[:i]),fct=1./i),a[:i], atol=2e-15, rtol=0)
assert_allclose(pypocketfft.irfftn(pypocketfft.rfftn(a[:i].real),fct=1./i,lastsize=i),a[:i].real, atol=2e-15, rtol=0)
assert_allclose(pypocketfft.ifftn(pypocketfft.fftn(b[:i]),fct=1./i),b[:i], atol=8e-7, rtol=0)
assert_allclose(pypocketfft.irfftn(pypocketfft.rfftn(b[:i].real),fct=1./i,lastsize=i),b[:i].real, atol=8e-7, rtol=0)
assert(_l2error(a, pypocketfft.ifftn(pypocketfft.fftn(a),fct=1./len))<1.5e-15)
assert(_l2error(a.real, pypocketfft.irfftn(pypocketfft.rfftn(a.real),fct=1./len,lastsize=len))<1e-15)
assert(_l2error(b, pypocketfft.ifftn(pypocketfft.fftn(b),fct=1./len))<5e-7)
assert(_l2error(b.real, pypocketfft.irfftn(pypocketfft.rfftn(b.real),fct=1./len,lastsize=len))<5e-7)
@pmp("shp", shapes)
def test_fftn(shp):
a=np.random.rand(*shp)-0.5 + 1j*np.random.rand(*shp)-0.5j
vmax = np.max(np.abs(pyfftw.interfaces.numpy_fft.fftn(a)))
assert_allclose(pyfftw.interfaces.numpy_fft.fftn(a), pypocketfft.fftn(a), atol=2e-15*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.fftn(a), pypocketfft.fftn(a))<1e-15*vmax)
a=a.astype(np.complex64)
assert_allclose(pyfftw.interfaces.numpy_fft.fftn(a), pypocketfft.fftn(a), atol=6e-7*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.fftn(a), pypocketfft.fftn(a))<5e-7*vmax)
@pmp("shp", shapes2D)
@pmp("axes", ((0,),(1,),(0,1),(1,0)))
def test_fftn2D(shp, axes):
a=np.random.rand(*shp)-0.5 + 1j*np.random.rand(*shp)-0.5j
vmax = np.max(np.abs(pyfftw.interfaces.numpy_fft.fftn(a, axes=axes)))
assert_allclose(pyfftw.interfaces.numpy_fft.fftn(a, axes=axes), pypocketfft.fftn(a, axes=axes), atol=2e-15*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.fftn(a, axes=axes), pypocketfft.fftn(a, axes=axes))<1e-15*vmax)
a=a.astype(np.complex64)
assert_allclose(pyfftw.interfaces.numpy_fft.fftn(a, axes=axes), pypocketfft.fftn(a, axes=axes), atol=6e-7*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.fftn(a, axes=axes), pypocketfft.fftn(a, axes=axes))<5e-7*vmax)
@pmp("shp", shapes)
def test_rfftn(shp):
a=np.random.rand(*shp)-0.5
vmax = np.max(np.abs(pyfftw.interfaces.numpy_fft.fftn(a)))
assert_allclose(pyfftw.interfaces.numpy_fft.rfftn(a), pypocketfft.rfftn(a), atol=2e-15*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.rfftn(a), pypocketfft.rfftn(a))<1e-15*vmax)
a=a.astype(np.float32)
assert_allclose(pyfftw.interfaces.numpy_fft.rfftn(a), pypocketfft.rfftn(a), atol=6e-7*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.rfftn(a), pypocketfft.rfftn(a))<5e-7*vmax)
@pmp("shp", shapes2D)
@pmp("axes", ((0,),(1,),(0,1),(1,0)))
def test_rfftn2D(shp, axes):
a=np.random.rand(*shp)-0.5
vmax = np.max(np.abs(pyfftw.interfaces.numpy_fft.rfftn(a, axes=axes)))
assert_allclose(pyfftw.interfaces.numpy_fft.rfftn(a, axes=axes), pypocketfft.rfftn(a, axes=axes), atol=2e-15*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.rfftn(a, axes=axes), pypocketfft.rfftn(a, axes=axes))<1e-15*vmax)
a=a.astype(np.float32)
assert_allclose(pyfftw.interfaces.numpy_fft.rfftn(a, axes=axes), pypocketfft.rfftn(a, axes=axes), atol=6e-7*vmax, rtol=0)
assert(_l2error(pyfftw.interfaces.numpy_fft.rfftn(a, axes=axes), pypocketfft.rfftn(a, axes=axes))<5e-7*vmax)
@pmp("shp", shapes)
def test_identity(shp):
a=np.random.rand(*shp)-0.5 + 1j*np.random.rand(*shp)-0.5j
vmax = np.max(np.abs(a))
assert_allclose(pypocketfft.ifftn(pypocketfft.fftn(a),fct=1./a.size), a, atol=2e-15*vmax, rtol=0)
assert(_l2error(pypocketfft.ifftn(pypocketfft.fftn(a),fct=1./a.size), a)<1e-15*vmax)
a=a.astype(np.complex64)
assert_allclose(pypocketfft.ifftn(pypocketfft.fftn(a),fct=1./a.size), a, atol=6e-7*vmax, rtol=0)
assert(_l2error(pypocketfft.ifftn(pypocketfft.fftn(a),fct=1./a.size), a)<5e-7*vmax)
@pmp("shp", shapes)
def test_identity_r(shp):
......@@ -70,8 +72,8 @@ def test_identity_r(shp):
vmax = np.max(np.abs(a))
for ax in range(a.ndim):
n = a.shape[ax]
assert_allclose(pypocketfft.irfftn(pypocketfft.rfftn(a,(ax,)),(ax,),lastsize=n,fct=1./n), a, atol=2e-15*vmax, rtol=0)
assert_allclose(pypocketfft.irfftn(pypocketfft.rfftn(b,(ax,)),(ax,),lastsize=n,fct=1./n), b, atol=6e-7*vmax, rtol=0)
assert(_l2error(pypocketfft.irfftn(pypocketfft.rfftn(a,(ax,)),(ax,),lastsize=n,fct=1./n), a)<1e-15*vmax)
assert(_l2error(pypocketfft.irfftn(pypocketfft.rfftn(b,(ax,)),(ax,),lastsize=n,fct=1./n), b)<5e-7*vmax)
@pmp("shp", shapes)
def test_identity_r2(shp):
......@@ -79,7 +81,7 @@ def test_identity_r2(shp):
a=pypocketfft.rfftn(pypocketfft.irfftn(a))
fct = 1./pypocketfft.irfftn(a).size
vmax = np.max(np.abs(a))
assert_allclose(pypocketfft.rfftn(pypocketfft.irfftn(a),fct=fct), a, atol=2e-15*vmax, rtol=0)
assert(_l2error(pypocketfft.rfftn(pypocketfft.irfftn(a),fct=fct), a)<1e-15*vmax)
@pmp("shp", shapes2D+shapes3D)
def test_hartley2(shp):
......@@ -88,21 +90,21 @@ def test_hartley2(shp):
v2 = pypocketfft.fftn(a.astype(np.complex128))
vmax = np.max(np.abs(v1))
v2 = v2.real+v2.imag
assert_allclose(v1, v2, atol=2e-15*vmax, rtol=0)
assert(_l2error(v1, v2)<1e-15*vmax)
@pmp("shp", shapes)
def test_hartley_identity(shp):
a=np.random.rand(*shp)-0.5
v1 = pypocketfft.hartley(pypocketfft.hartley(a))/a.size
vmax = np.max(np.abs(a))
assert_allclose(a, v1, atol=2e-15*vmax, rtol=0)
assert(_l2error(a, v1)<1e-15*vmax)
@pmp("shp", shapes)
def test_hartley2_identity(shp):
a=np.random.rand(*shp)-0.5
v1 = pypocketfft.hartley2(pypocketfft.hartley2(a))/a.size
vmax = np.max(np.abs(a))
assert_allclose(a, v1, atol=2e-15*vmax, rtol=0)
assert(_l2error(a, v1)<1e-15*vmax)
@pmp("shp", shapes2D)
@pmp("axes", ((0,),(1,),(0,1),(1,0)))
......@@ -110,4 +112,4 @@ def test_hartley2_2D(shp, axes):
a=np.random.rand(*shp)-0.5
fct = 1./np.prod(np.take(shp,axes))
vmax = np.max(np.abs(a))
assert_allclose(pypocketfft.hartley2(pypocketfft.hartley2(a, axes=axes), axes=axes, fct=fct),a , atol=2e-15*vmax, rtol=0)
assert(_l2error(pypocketfft.hartley2(pypocketfft.hartley2(a, axes=axes), axes=axes, fct=fct),a)<1e-15*vmax)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment