diff --git a/LICENSE.txt b/LICENSE.txt index f0eb281..85fa35f 100644 --- a/LICENSE.txt +++ b/LICENSE.txt @@ -1,4 +1,4 @@ -Copyright (c) 2017, Intel Corporation +Copyright (c) 2017-2019, Intel Corporation Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: diff --git a/conda-recipe/meta.yaml b/conda-recipe/meta.yaml index 42a8814..b025578 100644 --- a/conda-recipe/meta.yaml +++ b/conda-recipe/meta.yaml @@ -1,4 +1,4 @@ -{% set version = "1.0dev" %} +{% set version = "1.1.0" %} {% set buildnumber = 0 %} @@ -19,6 +19,7 @@ build: - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/_pydfti.* - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/_numpy_fft.py - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/_scipy_fft.py + - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/_scipy_fft_backend.py - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/setup.py - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/tests/test_fft1d.py - {{ SP_DIR.replace('\\', '/') if win else SP_DIR }}/mkl_fft/__init__.pyc [py27] @@ -36,12 +37,13 @@ requirements: - python - setuptools - intelpython - - mkl-devel [not nomkl] + - mkl-devel - cython - numpy x.x run: - python - - mkl [not nomkl] + - mkl + - mkl-service - intelpython - numpy x.x diff --git a/mkl_fft/__init__.py b/mkl_fft/__init__.py index a23246d..0f3109a 100644 --- a/mkl_fft/__init__.py +++ b/mkl_fft/__init__.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/_float_utils.py b/mkl_fft/_float_utils.py index 7f50a37..7f58b1b 100644 --- a/mkl_fft/_float_utils.py +++ b/mkl_fft/_float_utils.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/_numpy_fft.py b/mkl_fft/_numpy_fft.py index 7f222eb..15ae26c 100644 --- a/mkl_fft/_numpy_fft.py +++ b/mkl_fft/_numpy_fft.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/_pydfti.pyx b/mkl_fft/_pydfti.pyx index f2826b2..1b818b7 100644 --- a/mkl_fft/_pydfti.pyx +++ b/mkl_fft/_pydfti.pyx @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/_scipy_fft.py b/mkl_fft/_scipy_fft.py index 8f9f825..97a7ff8 100644 --- a/mkl_fft/_scipy_fft.py +++ b/mkl_fft/_scipy_fft.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/_scipy_fft_backend.py b/mkl_fft/_scipy_fft_backend.py new file mode 100644 index 0000000..102a577 --- /dev/null +++ b/mkl_fft/_scipy_fft_backend.py @@ -0,0 +1,265 @@ +#!/usr/bin/env python +# Copyright (c) 2019-2020, Intel Corporation +# +# Redistribution and use in source and binary forms, with or without +# modification, are permitted provided that the following conditions are met: +# +# * Redistributions of source code must retain the above copyright notice, +# this list of conditions and the following disclaimer. +# * Redistributions in binary form must reproduce the above copyright +# notice, this list of conditions and the following disclaimer in the +# documentation and/or other materials provided with the distribution. +# * Neither the name of Intel Corporation nor the names of its contributors +# may be used to endorse or promote products derived from this software +# without specific prior written permission. +# +# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE +# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE +# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL +# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR +# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER +# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + +from . import _pydfti +from . import _float_utils +import mkl + +import scipy.fft as _fft + +# Complete the namespace (these are not actually used in this module) +from scipy.fft import ( + dct, idct, dst, idst, dctn, idctn, dstn, idstn, + hfft2, ihfft2, hfftn, ihfftn, + fftshift, ifftshift, fftfreq, rfftfreq, + get_workers, set_workers +) + +from numpy.core import (array, asarray, shape, conjugate, take, sqrt, prod) + +__all__ = ['fft', 'ifft', 'fft2', 'ifft2', 'fftn', 'ifftn', + 'rfft', 'irfft', 'rfft2', 'irfft2', 'rfftn', 'irfftn', + 'hfft', 'ihfft', 'hfft2', 'ihfft2', 'hfftn', 'ihfftn', + 'dct', 'idct', 'dst', 'idst', 'dctn', 'idctn', 'dstn', 'idstn', + 'fftshift', 'ifftshift', 'fftfreq', 'rfftfreq', 'get_workers', + 'set_workers', 'next_fast_len'] + +__ua_domain__ = 'numpy.scipy.fft' +__implemented = dict() + +def __ua_function__(method, args, kwargs): + """Fetch registered UA function.""" + fn = __implemented.get(method, None) + if fn is None: + return NotImplemented + return fn(*args, **kwargs) + + +def _implements(scipy_func): + """Decorator adds function to the dictionary of implemented UA functions""" + def inner(func): + __implemented[scipy_func] = func + return func + + return inner + + +def _unitary(norm): + if norm not in (None, "ortho"): + raise ValueError("Invalid norm value %s, should be None or \"ortho\"." + % norm) + return norm is not None + + +def _cook_nd_args(a, s=None, axes=None, invreal=0): + if s is None: + shapeless = 1 + if axes is None: + s = list(a.shape) + else: + s = take(a.shape, axes) + else: + shapeless = 0 + s = list(s) + if axes is None: + axes = list(range(-len(s), 0)) + if len(s) != len(axes): + raise ValueError("Shape and axes have different lengths.") + if invreal and shapeless: + s[-1] = (a.shape[axes[-1]] - 1) * 2 + return s, axes + + +def _tot_size(x, axes): + s = x.shape + if axes is None: + return x.size + return prod([s[ai] for ai in axes]) + + +def _workers_to_num_threads(w): + if w is None: + return mkl.domain_get_max_threads(domain='fft') + return int(w) + + +class Workers: + def __init__(self, workers): + self.workers = workers + self.n_threads = _workers_to_num_threads(workers) + + def __enter__(self): + try: + mkl.domain_set_num_threads(self.n_threads, domain='fft') + except: + raise ValueError("Class argument {} result in invalid number of threads {}".format(self.workers, self.n_threads)) + + def __exit__(self, *args): + # restore default + max_num_threads = mkl.domain_get_max_threads(domain='fft') + mkl.domain_set_num_threads(max_num_threads, domain='fft') + + +@_implements(_fft.fft) +def fft(a, n=None, axis=-1, norm=None, overwrite_x=False, workers=None): + x = _float_utils.__upcast_float16_array(a) + with Workers(workers): + output = _pydfti.fft(x, n=n, axis=axis, overwrite_x=overwrite_x) + if _unitary(norm): + output *= 1 / sqrt(output.shape[axis]) + return output + + +@_implements(_fft.ifft) +def ifft(a, n=None, axis=-1, norm=None, overwrite_x=False, workers=None): + x = _float_utils.__upcast_float16_array(a) + with Workers(workers): + output = _pydfti.ifft(x, n=n, axis=axis, overwrite_x=overwrite_x) + if _unitary(norm): + output *= sqrt(output.shape[axis]) + return output + + +@_implements(_fft.fft2) +def fft2(a, s=None, axes=(-2,-1), norm=None, overwrite_x=False, workers=None): + x = _float_utils.__upcast_float16_array(a) + with Workers(workers): + output = _pydfti.fftn(x, shape=s, axes=axes, overwrite_x=overwrite_x) + if _unitary(norm): + factor = 1 + for axis in axes: + factor *= 1 / sqrt(output.shape[axis]) + output *= factor + return output + + +@_implements(_fft.ifft2) +def ifft2(a, s=None, axes=(-2,-1), norm=None, overwrite_x=False, workers=None): + x = _float_utils.__upcast_float16_array(a) + with Workers(workers): + output = _pydfti.ifftn(x, shape=s, axes=axes, overwrite_x=overwrite_x) + if _unitary(norm): + factor = 1 + _axes = range(output.ndim) if axes is None else axes + for axis in _axes: + factor *= sqrt(output.shape[axis]) + output *= factor + return output + + +@_implements(_fft.fftn) +def fftn(a, s=None, axes=None, norm=None, overwrite_x=False, workers=None): + x = _float_utils.__upcast_float16_array(a) + with Workers(workers): + output = _pydfti.fftn(x, shape=s, axes=axes, overwrite_x=overwrite_x) + if _unitary(norm): + factor = 1 + _axes = range(output.ndim) if axes is None else axes + for axis in _axes: + factor *= 1 / sqrt(output.shape[axis]) + output *= factor + return output + + +@_implements(_fft.ifftn) +def ifftn(a, s=None, axes=None, norm=None, overwrite_x=False, workers=None): + x = _float_utils.__upcast_float16_array(a) + with Workers(workers): + output = _pydfti.ifftn(x, shape=s, axes=axes, overwrite_x=overwrite_x) + if _unitary(norm): + factor = 1 + _axes = range(output.ndim) if axes is None else axes + for axis in _axes: + factor *= sqrt(output.shape[axis]) + output *= factor + return output + + +@_implements(_fft.rfft) +def rfft(a, n=None, axis=-1, norm=None, workers=None): + x = _float_utils.__upcast_float16_array(a) + unitary = _unitary(norm) + x = _float_utils.__downcast_float128_array(x) + if unitary and n is None: + x = asarray(x) + n = x.shape[axis] + with Workers(workers): + output = _pydfti.rfft_numpy(x, n=n, axis=axis) + if unitary: + output *= 1 / sqrt(n) + return output + + +@_implements(_fft.irfft) +def irfft(a, n=None, axis=-1, norm=None, workers=None): + x = _float_utils.__upcast_float16_array(a) + x = _float_utils.__downcast_float128_array(x) + with Workers(workers): + output = _pydfti.irfft_numpy(x, n=n, axis=axis) + if _unitary(norm): + output *= sqrt(output.shape[axis]) + return output + + +@_implements(_fft.rfft2) +def rfft2(a, s=None, axes=(-2, -1), norm=None, workers=None): + x = _float_utils.__upcast_float16_array(a) + x = _float_utils.__downcast_float128_array(a) + return rfftn(x, s, axes, norm, workers) + + +@_implements(_fft.irfft2) +def irfft2(a, s=None, axes=(-2, -1), norm=None, workers=None): + x = _float_utils.__upcast_float16_array(a) + x = _float_utils.__downcast_float128_array(x) + return irfftn(x, s, axes, norm, workers) + + +@_implements(_fft.rfftn) +def rfftn(a, s=None, axes=None, norm=None, workers=None): + unitary = _unitary(norm) + x = _float_utils.__upcast_float16_array(a) + x = _float_utils.__downcast_float128_array(x) + if unitary: + x = asarray(x) + s, axes = _cook_nd_args(x, s, axes) + with Workers(workers): + output = _pydfti.rfftn_numpy(x, s, axes) + if unitary: + n_tot = prod(asarray(s, dtype=output.dtype)) + output *= 1 / sqrt(n_tot) + return output + + +@_implements(_fft.irfftn) +def irfftn(a, s=None, axes=None, norm=None, workers=None): + x = _float_utils.__upcast_float16_array(a) + x = _float_utils.__downcast_float128_array(x) + with Workers(workers): + output = _pydfti.irfftn_numpy(x, s, axes) + if _unitary(norm): + output *= sqrt(_tot_size(output, axes)) + return output diff --git a/mkl_fft/setup.py b/mkl_fft/setup.py index f6f6a14..8be6f39 100644 --- a/mkl_fft/setup.py +++ b/mkl_fft/setup.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/src/mklfft.c.src b/mkl_fft/src/mklfft.c.src index 4e97a72..94244bc 100644 --- a/mkl_fft/src/mklfft.c.src +++ b/mkl_fft/src/mklfft.c.src @@ -1,5 +1,5 @@ /* - Copyright (c) 2017, Intel Corporation + Copyright (c) 2017-2020, Intel Corporation Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: @@ -35,6 +35,7 @@ #include #ifdef DEBUG +#include #include #define _debug_print(...) printf(__VA_ARGS__) #else @@ -471,6 +472,9 @@ int @name@_mkl_@mode@_in(PyArrayObject* x_inout, npy_intp n, int axis, DftiCache get_basic_array_data(x_inout, &x_rank, &x_shape, &x_strides, &x_itemsize, &x_size); + if (x_size == 0) /* nothing to do */ + return status; + assert( x_size > 0 ); /* assert that x is non-empty */ assert( 0 <= axis && axis < x_rank ); assert( x_itemsize == sizeof(@MKL_TYPE@) ); @@ -625,6 +629,11 @@ int @REALIN@_@COMPLEXOUT@_mkl_@mode@_out( assert(xout_rank == xin_rank); + if (xin_size == 0) { + /* nothing to do */ + assert(xin_size == xout_size); + return status; + } assert( xin_size > 0 ); /* assert that array is non-empty */ assert(0 <= axis && axis < xin_rank); assert( 0 < n && n <= xin_shape[axis] ); @@ -886,6 +895,12 @@ int @COMPLEXIN@_@COMPLEXOUT@_mkl_@mode@_out( assert(xout_rank == xin_rank); + if (xin_size == 0) { + /* nothing to do */ + assert(xout_size == 0); + return status; + } + assert( xin_size > 0 ); /* assert that array is non-empty */ assert(0 <= axis && axis < xin_rank); assert( 0 < n && n <= xin_shape[axis] ); @@ -1059,6 +1074,9 @@ int @name@_mkl_@mode@_in(PyArrayObject* x_inout, npy_intp n, int axis, DftiCache get_basic_array_data(x_inout, &x_rank, &x_shape, &x_strides, &x_itemsize, &x_size); + if (x_size == 0) /* nothing to do */ + return status; + assert( x_size > 0 ); /* assert that */ assert(0 <= axis && axis < x_rank); assert(x_itemsize == sizeof(@MKL_TYPE@)); @@ -1195,6 +1213,12 @@ int assert(xout_rank == xin_rank); + if (xin_size == 0) { + /* nothing to do */ + assert(xout_size == 0); + return status; + } + assert( xin_size > 0 ); /* assert that array is non-empty */ assert(0 <= axis && axis < xin_rank); assert( 0 < n && n <= xout_shape[axis] ); @@ -1369,6 +1393,12 @@ int @name@_@name@_mkl_@mode@_out( assert(xout_rank == xin_rank); + if(xin_size == 0) { + /* nothing to do */ + assert(xout_size == 0); + return status; + } + assert( xin_size > 0 ); /* assert that array is non-empty */ assert(0 <= axis && axis < xin_rank); assert( 0 < n && n <= xin_shape[axis] ); @@ -1543,6 +1573,9 @@ int &x_strides, &x_itemsize, &x_size); + if (x_size == 0) /* nothing to do */ + return status; + assert(x_size > 0); assert(x_itemsize == sizeof(@MKL_IN_TYPE@)); @@ -1649,6 +1682,11 @@ int get_basic_array_data(x_out, &xout_rank, &xout_shape, &xout_strides, &xout_itemsize, &xout_size); + if(xin_size == 0) { + /* nothing to do */ + assert(xout_size == 0); + return status; + } assert(xin_size > 0); assert(xin_itemsize == sizeof(@MKL_IN_TYPE@)); @@ -1778,6 +1816,11 @@ int get_basic_array_data(x_out, &xout_rank, &xout_shape, &xout_strides, &xout_itemsize, &xout_size); + if(xin_size == 0) { + /* nothing to do */ + assert(xout_size == 0); + return status; + } assert(xin_size > 0); assert(xin_itemsize == sizeof(@MKL_IN_TYPE@)); diff --git a/mkl_fft/src/mklfft.h b/mkl_fft/src/mklfft.h index e53f66d..1ed3073 100644 --- a/mkl_fft/src/mklfft.h +++ b/mkl_fft/src/mklfft.h @@ -1,5 +1,5 @@ /* - Copyright (c) 2017, Intel Corporation + Copyright (c) 2017-2019, Intel Corporation Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/src/multi_iter.h b/mkl_fft/src/multi_iter.h index 2940b68..8cbde22 100644 --- a/mkl_fft/src/multi_iter.h +++ b/mkl_fft/src/multi_iter.h @@ -1,5 +1,5 @@ /* - Copyright (c) 2017, Intel Corporation + Copyright (c) 2017-2019, Intel Corporation Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/tests/test_fft1d.py b/mkl_fft/tests/test_fft1d.py index de94a98..72aae6b 100644 --- a/mkl_fft/tests/test_fft1d.py +++ b/mkl_fft/tests/test_fft1d.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/mkl_fft/tests/test_fftnd.py b/mkl_fft/tests/test_fftnd.py index 089dff5..7bd644d 100644 --- a/mkl_fft/tests/test_fftnd.py +++ b/mkl_fft/tests/test_fftnd.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2019, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: diff --git a/setup.py b/setup.py index 7d0a910..ae9a06a 100644 --- a/setup.py +++ b/setup.py @@ -1,5 +1,5 @@ #!/usr/bin/env python -# Copyright (c) 2017, Intel Corporation +# Copyright (c) 2017-2020, Intel Corporation # # Redistribution and use in source and binary forms, with or without # modification, are permitted provided that the following conditions are met: @@ -24,7 +24,6 @@ # OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE # OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. -from __future__ import division, print_function, absolute_import import io import re @@ -86,7 +85,7 @@ def setup_package(): classifiers = [_f for _f in CLASSIFIERS.split('\n') if _f], platforms = ["Windows", "Linux", "Mac OS-X"], test_suite = 'nose.collector', - python_requires = '>=2.7,!=3.0.*,!=3.1.*,!=3.2.*,!=3.3.*', + python_requires = '>=3.5', install_requires = ['numpy'], configuration = configuration )