From 0566909e936ead73551c76c0e7445b89e1cd3cf5 Mon Sep 17 00:00:00 2001 From: Python Infra CI Date: Mon, 5 Oct 2026 16:48:21 +0000 Subject: [PATCH 1/2] fix: cast unsupported-dtype input to complex128 even when `out=` is give --- mkl_fft/_pydfti.pyx | 17 +++-- .../tests/test_c2c_out_unsupported_dtype.py | 72 +++++++++++++++++++ 2 files changed, 83 insertions(+), 6 deletions(-) create mode 100644 mkl_fft/tests/test_c2c_out_unsupported_dtype.py diff --git a/mkl_fft/_pydfti.pyx b/mkl_fft/_pydfti.pyx index 36e152fb..a3dc5183 100644 --- a/mkl_fft/_pydfti.pyx +++ b/mkl_fft/_pydfti.pyx @@ -425,11 +425,11 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): x_arr = _process_arguments(x, n, axis, &axis_, &n_, &in_place, &xnd, 0) x_type = cnp.PyArray_TYPE(x_arr) - if out is not None: - in_place = 0 - elif x_type is cnp.NPY_CFLOAT or x_type is cnp.NPY_CDOUBLE: + if x_type is cnp.NPY_CFLOAT or x_type is cnp.NPY_CDOUBLE: + if out is not None: + in_place = 0 # we can operate in place if requested. - if in_place: + elif in_place: if not cnp.PyArray_ISONESEGMENT(x_arr): in_place = 0 if internal_overlap(x_arr) else 1 elif x_type is cnp.NPY_FLOAT or x_type is cnp.NPY_DOUBLE: @@ -438,7 +438,12 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): in_place = 0 else: # we must cast the input and allocate the output, - # so we cast to complex double and operate in place + # so we cast to complex double and operate in place. + # This cast must happen even when `out` is given: `x_type` is + # otherwise left as this unsupported type, and the dispatch below + # (which only recognizes float/double/cfloat/cdouble) would then + # silently skip calling into MKL, leaving `out`/`f_arr` filled with + # whatever uninitialized memory `_allocate_result` handed back. try: x_arr = cnp.PyArray_FROM_OTF( x_arr, cnp.NPY_CDOUBLE, @@ -450,7 +455,7 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): "or real sequence of single or double precision" ) x_type = cnp.PyArray_TYPE(x_arr) - in_place = 1 + in_place = 0 if out is not None else 1 if in_place: _cache_capsule = _tls_dfti_cache_capsule() diff --git a/mkl_fft/tests/test_c2c_out_unsupported_dtype.py b/mkl_fft/tests/test_c2c_out_unsupported_dtype.py new file mode 100644 index 00000000..1b80e6c7 --- /dev/null +++ b/mkl_fft/tests/test_c2c_out_unsupported_dtype.py @@ -0,0 +1,72 @@ +# Copyright (c) 2025, 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. + +""" +Regression test for a silent-garbage-output bug in ``_c2c_fft1d_impl`` +(mkl_fft/_pydfti.pyx). + +``mkl_fft.fft``/``mkl_fft.ifft`` accept input of any dtype and, for dtypes +other than float32/float64/complex64/complex128, cast the input to +complex128 before computing the transform (cf. test_vector12 and +test_vector7/8 in test_fft1d.py, which exercise this cast when ``out`` is +not given). + +Before the fix, that cast only happened when ``out`` was *not* passed: the +original code checked ``if out is not None: in_place = 0`` *before* looking +at ``x_type`` at all, so the branch that performs +``PyArray_FROM_OTF(..., NPY_CDOUBLE, ...)`` was skipped whenever the caller +supplied ``out``. ``x_type`` was then left as the unsupported dtype's type +code, which the dispatch further down (``if x_type is NPY_DOUBLE: ... elif +x_type is NPY_CDOUBLE: ...``) does not recognize, so none of its branches +ran and the ``status`` variable kept its initial value of 0. Because a +"success" status was reported without MKL ever having been called, the +function returned ``out`` unmodified -- whatever uninitialized values +``np.empty`` happened to produce -- instead of raising an error or +computing the transform. +""" + +import numpy as np +import pytest + +import mkl_fft + + +@pytest.mark.parametrize("dt", ["i4", "i8", "f2"]) +@pytest.mark.parametrize("func, npfunc", [("fft", "fft"), ("ifft", "ifft")]) +def test_c2c_out_with_unsupported_input_dtype(dt, func, npfunc): + """fft/ifft with `out=` must compute the real transform for dtypes + that require an internal cast to complex128 (e.g. integer, float16), + not silently return the uninitialized `out` buffer.""" + x = np.arange(1, 17, dtype=dt) + expected = getattr(np.fft, npfunc)(x.astype(np.float64)) + + # fill `out` with a sentinel value that could never, by chance, match + # the expected FFT output, so a correctness regression is detected + # reliably rather than only "most of the time". + out = np.full(x.shape, -1 - 1j, dtype=np.complex128) + result = getattr(mkl_fft, func)(x, out=out) + + assert result is out + np.testing.assert_allclose(result, expected, rtol=1e-10, atol=1e-10) From 2417932979e953c49fb3d2b590adfa91747229ef Mon Sep 17 00:00:00 2001 From: Anton Volkov Date: Tue, 6 Oct 2026 14:24:00 +0200 Subject: [PATCH 2/2] fix: check `out` strides against the array passed to MKL in c2c FFT `_c2c_fft1d_impl` decided whether `out` can be written directly by MKL by comparing its strides with the original `x`, but MKL reads `x_arr`, which may be a cast or zero-padded contiguous copy of `x`. When `x` and `out` shared a non-contiguous or F-ordered layout, rank >= 3 transforms gave wrong results and wrote outside of `out`. With the cast now also done when `out` is given, this affected integer/bool/float16 input as well as the pre-existing case of padded complex input. Also restore the original dtype-dispatch structure and apply the `out` override after the cast, as `_direct_fftnd` does; raise instead of silently returning when the out-of-place dispatch sees an unsupported type; move the regression tests into test_fft1d.py/test_fftnd.py with coverage for strided, F-ordered, truncated and padded layouts and for fftn/ifftn; and add CHANGELOG entries. --- CHANGELOG.md | 2 + mkl_fft/_pydfti.pyx | 24 ++++--- .../tests/test_c2c_out_unsupported_dtype.py | 72 ------------------- mkl_fft/tests/test_fft1d.py | 67 +++++++++++++++++ mkl_fft/tests/test_fftnd.py | 22 ++++++ 5 files changed, 104 insertions(+), 83 deletions(-) delete mode 100644 mkl_fft/tests/test_c2c_out_unsupported_dtype.py diff --git a/CHANGELOG.md b/CHANGELOG.md index 3d1795ff..1a843956 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -22,6 +22,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 * Fixed possible memory leaks when `PyMem_Malloc` fails, and raise `MemoryError` [gh-373](https://github.com/IntelPython/mkl_fft/pull/373) * Fixed possible memory leaks when multi-iterator constructors fail [gh-373](https://github.com/IntelPython/mkl_fft/pull/373) * Fixed N-D transforms returning a success status when a scratch allocation fails [gh-373](https://github.com/IntelPython/mkl_fft/pull/373) +* Fixed `fft`, `ifft`, and the `fftn`/`fft2` family returning `out` unmodified when `out` is given and the input dtype is not `float32`, `float64`, `complex64` or `complex128` (e.g. integer, `bool`, `float16`) [gh-387](https://github.com/IntelPython/mkl_fft/pull/387) +* Fixed wrong results and writes outside of `out` in `fft`/`ifft` when the input is cast or zero-padded into a contiguous copy while `out` has the same non-contiguous layout as the original input [gh-387](https://github.com/IntelPython/mkl_fft/pull/387) ## [2.3.2] - 2026-08-04 diff --git a/mkl_fft/_pydfti.pyx b/mkl_fft/_pydfti.pyx index a3dc5183..9671bbfd 100644 --- a/mkl_fft/_pydfti.pyx +++ b/mkl_fft/_pydfti.pyx @@ -426,10 +426,8 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): x_type = cnp.PyArray_TYPE(x_arr) if x_type is cnp.NPY_CFLOAT or x_type is cnp.NPY_CDOUBLE: - if out is not None: - in_place = 0 # we can operate in place if requested. - elif in_place: + if in_place: if not cnp.PyArray_ISONESEGMENT(x_arr): in_place = 0 if internal_overlap(x_arr) else 1 elif x_type is cnp.NPY_FLOAT or x_type is cnp.NPY_DOUBLE: @@ -438,12 +436,7 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): in_place = 0 else: # we must cast the input and allocate the output, - # so we cast to complex double and operate in place. - # This cast must happen even when `out` is given: `x_type` is - # otherwise left as this unsupported type, and the dispatch below - # (which only recognizes float/double/cfloat/cdouble) would then - # silently skip calling into MKL, leaving `out`/`f_arr` filled with - # whatever uninitialized memory `_allocate_result` handed back. + # so we cast to complex double and operate in place try: x_arr = cnp.PyArray_FROM_OTF( x_arr, cnp.NPY_CDOUBLE, @@ -455,7 +448,11 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): "or real sequence of single or double precision" ) x_type = cnp.PyArray_TYPE(x_arr) - in_place = 0 if out is not None else 1 + in_place = 1 + + # checked only after the cast above, which is needed even if out is given + if out is not None: + in_place = 0 if in_place: _cache_capsule = _tls_dfti_cache_capsule() @@ -506,8 +503,9 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): _validate_out_array(out, x, out_dtype, axis=axis_, n=n_) # out array that is used in OneMKL c2c FFT must have the exact same # stride as input array. If not, we need to allocate a new array. + # Compare with x_arr, which may be a cast or padded copy of x. # TODO: check to see if this condition can be relaxed - if _get_element_strides(x) == _get_element_strides(out): + if _get_element_strides(x_arr) == _get_element_strides(out): f_arr = out else: f_arr = _allocate_result(x_arr, n_, axis_, f_type) @@ -548,6 +546,10 @@ def _c2c_fft1d_impl(x, n=None, axis=-1, direction=+1, double fsc=1.0, out=None): status = cdouble_cdouble_mkl_fft1d_out( x_arr, n_, axis_, f_arr, fsc, _cache ) + else: + raise ValueError( + "An input argument x is not of a supported type" + ) else: if x_type is cnp.NPY_FLOAT: if direction < 0: diff --git a/mkl_fft/tests/test_c2c_out_unsupported_dtype.py b/mkl_fft/tests/test_c2c_out_unsupported_dtype.py deleted file mode 100644 index 1b80e6c7..00000000 --- a/mkl_fft/tests/test_c2c_out_unsupported_dtype.py +++ /dev/null @@ -1,72 +0,0 @@ -# Copyright (c) 2025, 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. - -""" -Regression test for a silent-garbage-output bug in ``_c2c_fft1d_impl`` -(mkl_fft/_pydfti.pyx). - -``mkl_fft.fft``/``mkl_fft.ifft`` accept input of any dtype and, for dtypes -other than float32/float64/complex64/complex128, cast the input to -complex128 before computing the transform (cf. test_vector12 and -test_vector7/8 in test_fft1d.py, which exercise this cast when ``out`` is -not given). - -Before the fix, that cast only happened when ``out`` was *not* passed: the -original code checked ``if out is not None: in_place = 0`` *before* looking -at ``x_type`` at all, so the branch that performs -``PyArray_FROM_OTF(..., NPY_CDOUBLE, ...)`` was skipped whenever the caller -supplied ``out``. ``x_type`` was then left as the unsupported dtype's type -code, which the dispatch further down (``if x_type is NPY_DOUBLE: ... elif -x_type is NPY_CDOUBLE: ...``) does not recognize, so none of its branches -ran and the ``status`` variable kept its initial value of 0. Because a -"success" status was reported without MKL ever having been called, the -function returned ``out`` unmodified -- whatever uninitialized values -``np.empty`` happened to produce -- instead of raising an error or -computing the transform. -""" - -import numpy as np -import pytest - -import mkl_fft - - -@pytest.mark.parametrize("dt", ["i4", "i8", "f2"]) -@pytest.mark.parametrize("func, npfunc", [("fft", "fft"), ("ifft", "ifft")]) -def test_c2c_out_with_unsupported_input_dtype(dt, func, npfunc): - """fft/ifft with `out=` must compute the real transform for dtypes - that require an internal cast to complex128 (e.g. integer, float16), - not silently return the uninitialized `out` buffer.""" - x = np.arange(1, 17, dtype=dt) - expected = getattr(np.fft, npfunc)(x.astype(np.float64)) - - # fill `out` with a sentinel value that could never, by chance, match - # the expected FFT output, so a correctness regression is detected - # reliably rather than only "most of the time". - out = np.full(x.shape, -1 - 1j, dtype=np.complex128) - result = getattr(mkl_fft, func)(x, out=out) - - assert result is out - np.testing.assert_allclose(result, expected, rtol=1e-10, atol=1e-10) diff --git a/mkl_fft/tests/test_fft1d.py b/mkl_fft/tests/test_fft1d.py index 16464e75..ccf4332c 100644 --- a/mkl_fft/tests/test_fft1d.py +++ b/mkl_fft/tests/test_fft1d.py @@ -367,6 +367,73 @@ def test_fft_out_strided(axis, func): assert_allclose(result, expected) +@pytest.mark.parametrize("n", [None, 8, 24]) +@pytest.mark.parametrize("dt", ["?", "i1", "u1", "i4", "i8", "f2"]) +@pytest.mark.parametrize("func", ["fft", "ifft"]) +def test_fft_out_cast_input(func, dt, n): + # dtypes other than f4/f8/c8/c16 are cast to complex128 + x = np.arange(1, 17).astype(dt) + # sentinel: catches out being returned untouched + out = np.full(16 if n is None else n, -1 - 1j) + result = getattr(mkl_fft, func)(x, n=n, out=out) + expected = getattr(np.fft, func)(x.astype(np.float64), n=n) + + assert result is out + assert_allclose(result, expected, atol=1e-12) + + +@pytest.mark.parametrize("axis", [0, 1, 2]) +@pytest.mark.parametrize("dt", ["i8", "c16"]) +@pytest.mark.parametrize("func", ["fft", "ifft"]) +def test_fft_out_strided_input_copied(func, dt, axis): + # x and out have the same strides, but x is cast (i8) or padded (c16) + # into a contiguous copy, so out cannot be handed to MKL as is + shape = (20, 33, 54) + base = np.full(shape, -1 - 1j) + out = base[::2, ::3, ::4] + x = rnd.randint(-50, 50, size=shape).astype(dt)[::2, ::3, ::4] + n = None + if dt == "c16": + ind = [slice(None)] * x.ndim + ind[axis] = slice(0, x.shape[axis] - 3) + x = x[tuple(ind)] + n = out.shape[axis] + + result = getattr(mkl_fft, func)(x, n=n, axis=axis, out=out) + expected = getattr(np.fft, func)(x.astype(np.complex128), n=n, axis=axis) + + assert result is out + assert_allclose(result, expected, atol=1e-10) + # nothing outside of out was written to + out[...] = -1 - 1j + assert np.all(base == -1 - 1j) + + +@pytest.mark.parametrize("axis", [0, 1, 2]) +@pytest.mark.parametrize("func", ["fft", "ifft"]) +def test_fft_out_fortran_cast_input(func, axis): + x = np.asfortranarray(rnd.randint(-50, 50, size=(4, 5, 6))) + out = np.full(x.shape, -1 - 1j, order="F") + result = getattr(mkl_fft, func)(x, axis=axis, out=out) + expected = getattr(np.fft, func)(x.astype(np.float64), axis=axis) + + assert result is out + assert_allclose(result, expected, atol=1e-10) + + +@pytest.mark.skipif( + np.can_cast(np.longdouble, np.complex128), + reason="long double is the same as double", +) +@pytest.mark.parametrize("dt", [np.longdouble, np.clongdouble]) +@pytest.mark.parametrize("use_out", [False, True]) +def test_fft_longdouble_unsupported(dt, use_out): + x = np.ones(8, dtype=dt) + out = np.empty(8, dtype=np.complex128) if use_out else None + with pytest.raises(ValueError, match="single or double precision"): + mkl_fft.fft(x, out=out) + + @requires_numpy_2 @pytest.mark.parametrize("axis", [0, 1, 2]) def test_rfft_out_strided(axis): diff --git a/mkl_fft/tests/test_fftnd.py b/mkl_fft/tests/test_fftnd.py index 4aa6f159..0b2f9b0c 100644 --- a/mkl_fft/tests/test_fftnd.py +++ b/mkl_fft/tests/test_fftnd.py @@ -309,6 +309,28 @@ def test_out_strided(axes, func): assert_allclose(result, expected, strict=True) +@pytest.mark.parametrize("strided", [False, True]) +@pytest.mark.parametrize("axes", [None, (0, 1), (0, 2), (1, 2)]) +@pytest.mark.parametrize("func", ["fftn", "ifftn"]) +def test_out_int_input(func, axes, strided): + # integer input is cast to complex128 before being transformed into out + shape = (20, 30, 40) + x = rnd.randint(-50, 50, size=shape) + base = np.full(shape, -1 - 1j) + out = base + if strided: + x = x[::2, ::3, ::4] + out = base[::2, ::3, ::4] + result = getattr(mkl_fft, func)(x, axes=axes, out=out) + expected = getattr(np.fft, func)(x.astype(np.float64), axes=axes) + + assert result is out + assert_allclose(result, expected, rtol=reps_64, atol=1e-9) + # nothing outside of out was written to + out[...] = -1 - 1j + assert np.all(base == -1 - 1j) + + @pytest.mark.parametrize( "dtype", [np.float32, np.float64, np.complex64, np.complex128] )