diff --git a/dpnp/CMakeLists.txt b/dpnp/CMakeLists.txt index 46939410e73d..97ea9759236b 100644 --- a/dpnp/CMakeLists.txt +++ b/dpnp/CMakeLists.txt @@ -62,5 +62,4 @@ add_subdirectory(backend/extensions/ufunc) add_subdirectory(dpnp_algo) add_subdirectory(dpnp_utils) -add_subdirectory(fft) add_subdirectory(random) diff --git a/dpnp/backend/CMakeLists.txt b/dpnp/backend/CMakeLists.txt index cef60c7632a4..64d1e81c9c28 100644 --- a/dpnp/backend/CMakeLists.txt +++ b/dpnp/backend/CMakeLists.txt @@ -27,7 +27,6 @@ set(DPNP_SRC kernels/dpnp_krnl_arraycreation.cpp kernels/dpnp_krnl_common.cpp kernels/dpnp_krnl_elemwise.cpp - kernels/dpnp_krnl_fft.cpp kernels/dpnp_krnl_indexing.cpp kernels/dpnp_krnl_mathematical.cpp kernels/dpnp_krnl_random.cpp diff --git a/dpnp/backend/include/dpnp_gen_1arg_2type_tbl.hpp b/dpnp/backend/include/dpnp_gen_1arg_2type_tbl.hpp index 5efb489ba366..018a625a606a 100644 --- a/dpnp/backend/include/dpnp_gen_1arg_2type_tbl.hpp +++ b/dpnp/backend/include/dpnp_gen_1arg_2type_tbl.hpp @@ -85,7 +85,6 @@ #endif -MACRO_1ARG_2TYPES_OP(dpnp_copyto_c, input_elem, q.submit(kernel_func)) MACRO_1ARG_2TYPES_OP(dpnp_sqrt_c, sycl::sqrt(input_elem), oneapi::mkl::vm::sqrt(q, input1_size, input1_data, result)) diff --git a/dpnp/backend/include/dpnp_iface.hpp b/dpnp/backend/include/dpnp_iface.hpp index d3b342bee303..e99921baa2da 100644 --- a/dpnp/backend/include/dpnp_iface.hpp +++ b/dpnp/backend/include/dpnp_iface.hpp @@ -58,7 +58,6 @@ typedef ssize_t shape_elem_type; #include -#include "dpnp_iface_fft.hpp" #include "dpnp_iface_random.hpp" /** diff --git a/dpnp/backend/include/dpnp_iface_fft.hpp b/dpnp/backend/include/dpnp_iface_fft.hpp deleted file mode 100644 index 308ec7897f27..000000000000 --- a/dpnp/backend/include/dpnp_iface_fft.hpp +++ /dev/null @@ -1,141 +0,0 @@ -//***************************************************************************** -// Copyright (c) 2016-2024, Intel Corporation -// All rights reserved. -// -// 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. -// -// 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 HOLDER 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. -//***************************************************************************** - -/* - * This header file is for interface Cython with C++. - * It should not contains any backend specific headers (like SYCL or math - * library) because all included headers will be exposed in Cython compilation - * procedure - * - * We would like to avoid backend specific things in higher level Cython - * modules. Any backend interface functions and types should be defined here. - * - * Also, this file should contains documentation on functions and types - * which are used in the interface - */ - -#pragma once -#ifndef BACKEND_IFACE_FFT_H // Cython compatibility -#define BACKEND_IFACE_FFT_H - -/** - * @defgroup BACKEND_FFT_API Backend C++ library interface FFT API - * @{ - * This section describes Backend API for Discrete Fourier Transform (DFT) part. - * @} - */ - -/** - * @ingroup BACKEND_FFT_API - * @brief 1D discrete Fourier Transform. - * - * Compute the one-dimensional discrete Fourier Transform. - * - * @param[in] q_ref Reference to SYCL queue. - * @param[in] array1_in Input array. - * @param[out] result_out Output array. - * @param[in] input_shape Array with shape information for input array. - * @param[in] result_shape Array with shape information for result - * array. - * @param[in] shape_size Number of elements in @ref input_shape or - * @ref result_shape arrays. - * @param[in] axis Axis ID to compute by. - * @param[in] input_boundarie Limit number of elements for @ref axis. - * @param[in] inverse Using inverse algorithm. - * @param[in] norm Normalization mode. 0 - backward, 1 - - * forward, 2 - ortho. - * @param[in] dep_event_vec_ref Reference to vector of SYCL events. - */ -template -INP_DLLEXPORT DPCTLSyclEventRef - dpnp_fft_fft_c(DPCTLSyclQueueRef q_ref, - const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *result_shape, - size_t shape_size, - long axis, - long input_boundarie, - size_t inverse, - const size_t norm, - const DPCTLEventVectorRef dep_event_vec_ref); - -template -INP_DLLEXPORT void dpnp_fft_fft_c(const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *output_shape, - size_t shape_size, - long axis, - long input_boundarie, - size_t inverse, - const size_t norm); - -/** - * @ingroup BACKEND_FFT_API - * @brief 1D discrete Fourier Transform. - * - * Compute the one-dimensional discrete Fourier Transform for real input. - * - * @param[in] q_ref Reference to SYCL queue. - * @param[in] array1_in Input array. - * @param[out] result_out Output array. - * @param[in] input_shape Array with shape information for input array. - * @param[in] result_shape Array with shape information for result - * array. - * @param[in] shape_size Number of elements in @ref input_shape or - * @ref result_shape arrays. - * @param[in] axis Axis ID to compute by. - * @param[in] input_boundarie Limit number of elements for @ref axis. - * @param[in] inverse Using inverse algorithm. - * @param[in] norm Normalization mode. 0 - backward, 1 - - * forward, 2 - ortho. - * @param[in] dep_event_vec_ref Reference to vector of SYCL events. - */ -template -INP_DLLEXPORT DPCTLSyclEventRef - dpnp_fft_rfft_c(DPCTLSyclQueueRef q_ref, - const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *result_shape, - size_t shape_size, - long axis, - long input_boundarie, - size_t inverse, - const size_t norm, - const DPCTLEventVectorRef dep_event_vec_ref); - -template -INP_DLLEXPORT void dpnp_fft_fft_c(const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *output_shape, - size_t shape_size, - long axis, - long input_boundarie, - size_t inverse, - const size_t norm); -#endif // BACKEND_IFACE_FFT_H diff --git a/dpnp/backend/include/dpnp_iface_fptr.hpp b/dpnp/backend/include/dpnp_iface_fptr.hpp index f0b6cc67399d..bc5d0042a35a 100644 --- a/dpnp/backend/include/dpnp_iface_fptr.hpp +++ b/dpnp/backend/include/dpnp_iface_fptr.hpp @@ -65,9 +65,6 @@ enum class DPNPFuncName : size_t DPNP_FN_CHOOSE, /**< Used in numpy.choose() impl */ DPNP_FN_CHOOSE_EXT, /**< Used in numpy.choose() impl, requires extra parameters */ - DPNP_FN_COPYTO, /**< Used in numpy.copyto() impl */ - DPNP_FN_COPYTO_EXT, /**< Used in numpy.copyto() impl, requires extra - parameters */ DPNP_FN_CORRELATE, /**< Used in numpy.correlate() impl */ DPNP_FN_CORRELATE_EXT, /**< Used in numpy.correlate() impl, requires extra parameters */ @@ -81,9 +78,6 @@ enum class DPNPFuncName : size_t DPNP_FN_ERF, /**< Used in scipy.special.erf impl */ DPNP_FN_ERF_EXT, /**< Used in scipy.special.erf impl, requires extra parameters */ - DPNP_FN_FFT_FFT, /**< Used in numpy.fft.fft() impl */ - DPNP_FN_FFT_FFT_EXT, /**< Used in numpy.fft.fft() impl, requires extra - parameters */ DPNP_FN_INITVAL, /**< Used in numpy ones, ones_like, zeros, zeros_like impls */ DPNP_FN_INITVAL_EXT, /**< Used in numpy ones, ones_like, zeros, zeros_like diff --git a/dpnp/backend/kernels/dpnp_krnl_elemwise.cpp b/dpnp/backend/kernels/dpnp_krnl_elemwise.cpp index 4a02c1dfe43c..7227b0dcd99d 100644 --- a/dpnp/backend/kernels/dpnp_krnl_elemwise.cpp +++ b/dpnp/backend/kernels/dpnp_krnl_elemwise.cpp @@ -231,82 +231,6 @@ using dpctl::tensor::kernels::alignment_utils::required_alignment; static void func_map_init_elemwise_1arg_2type(func_map_t &fmap) { - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_BLN][eft_BLN] = { - eft_BLN, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_BLN][eft_INT] = { - eft_INT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_BLN][eft_LNG] = { - eft_LNG, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_BLN][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_BLN][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_INT][eft_BLN] = { - eft_BLN, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_INT][eft_INT] = { - eft_INT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_INT][eft_LNG] = { - eft_LNG, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_INT][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_INT][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_LNG][eft_BLN] = { - eft_BLN, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_LNG][eft_INT] = { - eft_INT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_LNG][eft_LNG] = { - eft_LNG, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_LNG][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_LNG][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_FLT][eft_BLN] = { - eft_BLN, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_FLT][eft_INT] = { - eft_INT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_FLT][eft_LNG] = { - eft_LNG, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_FLT][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_FLT][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_DBL][eft_BLN] = { - eft_BLN, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_DBL][eft_INT] = { - eft_INT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_DBL][eft_LNG] = { - eft_LNG, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_DBL][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_DBL][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_default}; - fmap[DPNPFuncName::DPNP_FN_COPYTO][eft_C128][eft_C128] = { - eft_C128, - (void *) - dpnp_copyto_c_default, std::complex>}; - - // dpnp_copyto_c is required by dpnp_fft_fft_c and dpnp_fft_rfft_c - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_BLN][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_INT][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_LNG][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_FLT][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_DBL][eft_DBL] = { - eft_DBL, (void *)dpnp_copyto_c_ext}; - - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_BLN][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_INT][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_LNG][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_COPYTO_EXT][eft_FLT][eft_FLT] = { - eft_FLT, (void *)dpnp_copyto_c_ext}; - fmap[DPNPFuncName::DPNP_FN_SQRT][eft_INT][eft_INT] = { eft_DBL, (void *)dpnp_sqrt_c_default}; fmap[DPNPFuncName::DPNP_FN_SQRT][eft_LNG][eft_LNG] = { diff --git a/dpnp/backend/kernels/dpnp_krnl_fft.cpp b/dpnp/backend/kernels/dpnp_krnl_fft.cpp deleted file mode 100644 index 3d81dde39250..000000000000 --- a/dpnp/backend/kernels/dpnp_krnl_fft.cpp +++ /dev/null @@ -1,613 +0,0 @@ -//***************************************************************************** -// Copyright (c) 2016-2024, Intel Corporation -// All rights reserved. -// -// 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. -// -// 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 HOLDER 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. -//***************************************************************************** - -#include - -#include "dpnp_fptr.hpp" -#include "dpnp_utils.hpp" -#include "dpnpc_memory_adapter.hpp" -#include "queue_sycl.hpp" -#include - -namespace mkl_dft = oneapi::mkl::dft; - -typedef mkl_dft::descriptor - desc_dp_cmplx_t; -typedef mkl_dft::descriptor - desc_sp_cmplx_t; -typedef mkl_dft::descriptor - desc_dp_real_t; -typedef mkl_dft::descriptor - desc_sp_real_t; - -#ifdef _WIN32 -#ifndef M_PI // Windows compatibility -#define M_PI 3.14159265358979323846 -#endif -#endif - -template -class dpnp_fft_fft_c_kernel; - -template -static void dpnp_fft_fft_sycl_c(DPCTLSyclQueueRef q_ref, - const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *output_shape, - size_t shape_size, - const size_t result_size, - const size_t input_size, - long axis, - long input_boundarie, - size_t inverse) -{ - if (!(input_size && result_size && shape_size)) { - return; - } - - sycl::event event; - - const double kernel_pi = inverse ? -M_PI : M_PI; - - sycl::queue queue = *(reinterpret_cast(q_ref)); - - _DataType_input *array_1 = - static_cast<_DataType_input *>(const_cast(array1_in)); - _DataType_output *result = static_cast<_DataType_output *>(result_out); - - // kernel specific temporal data - shape_elem_type *output_shape_offsets = reinterpret_cast( - sycl::malloc_shared(shape_size * sizeof(shape_elem_type), queue)); - shape_elem_type *input_shape_offsets = reinterpret_cast( - sycl::malloc_shared(shape_size * sizeof(shape_elem_type), queue)); - // must be a thread local storage. - shape_elem_type *axis_iterator = - reinterpret_cast(sycl::malloc_shared( - result_size * shape_size * sizeof(shape_elem_type), queue)); - - get_shape_offsets_inkernel(output_shape, shape_size, output_shape_offsets); - get_shape_offsets_inkernel(input_shape, shape_size, input_shape_offsets); - - sycl::range<1> gws(result_size); - auto kernel_parallel_for_func = [=](sycl::id<1> global_id) { - size_t output_id = global_id[0]; - - double sum_real = 0.0; - double sum_imag = 0.0; - // need to replace this array by thread local storage - shape_elem_type *axis_iterator_thread = - axis_iterator + (output_id * shape_size); - - size_t xyz_id; - for (size_t i = 0; i < shape_size; ++i) { - xyz_id = get_xyz_id_by_id_inkernel(output_id, output_shape_offsets, - shape_size, i); - axis_iterator_thread[i] = xyz_id; - } - - const long axis_length = input_boundarie; - for (long it = 0; it < axis_length; ++it) { - double in_real = 0.0; - double in_imag = 0.0; - - axis_iterator_thread[axis] = it; - - const size_t input_it = get_id_by_xyz_inkernel( - axis_iterator_thread, shape_size, input_shape_offsets); - - if (it < input_shape[axis]) { - if constexpr (std::is_same<_DataType_input, - std::complex>::value) { - const _DataType_input *cmplx_ptr = array_1 + input_it; - const double *dbl_ptr = - reinterpret_cast(cmplx_ptr); - in_real = *dbl_ptr; - in_imag = *(dbl_ptr + 1); - } - else if constexpr (std::is_same<_DataType_input, - std::complex>::value) { - const _DataType_input *cmplx_ptr = array_1 + input_it; - const float *dbl_ptr = - reinterpret_cast(cmplx_ptr); - in_real = *dbl_ptr; - in_imag = *(dbl_ptr + 1); - } - else { - in_real = array_1[input_it]; - } - } - - xyz_id = get_xyz_id_by_id_inkernel(output_id, output_shape_offsets, - shape_size, axis); - const size_t output_local_id = xyz_id; - const double angle = - 2.0 * kernel_pi * it * output_local_id / axis_length; - - const double angle_cos = sycl::cos(angle); - const double angle_sin = sycl::sin(angle); - - sum_real += in_real * angle_cos + in_imag * angle_sin; - sum_imag += -in_real * angle_sin + in_imag * angle_cos; - } - - if (inverse) { - sum_real = sum_real / input_boundarie; - sum_imag = sum_imag / input_boundarie; - } - - result[output_id] = _DataType_output(sum_real, sum_imag); - }; - - auto kernel_func = [&](sycl::handler &cgh) { - cgh.parallel_for< - class dpnp_fft_fft_c_kernel<_DataType_input, _DataType_output>>( - gws, kernel_parallel_for_func); - }; - - event = queue.submit(kernel_func); - event.wait(); - - sycl::free(input_shape_offsets, queue); - sycl::free(output_shape_offsets, queue); - sycl::free(axis_iterator, queue); - - return; -} - -template -static void - dpnp_fft_fft_mathlib_cmplx_to_cmplx_c(DPCTLSyclQueueRef q_ref, - const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *result_shape, - const size_t shape_size, - const size_t input_size, - const size_t result_size, - size_t inverse, - const size_t norm) -{ - // avoid warning unused variable - (void)result_shape; - (void)input_size; - (void)result_size; - - if (!shape_size) { - return; - } - - sycl::queue queue = *(reinterpret_cast(q_ref)); - - _DataType_input *array_1 = - static_cast<_DataType_input *>(const_cast(array1_in)); - _DataType_output *result = static_cast<_DataType_output *>(result_out); - - const size_t n_iter = - std::accumulate(input_shape, input_shape + shape_size - 1, 1, - std::multiplies()); - - const size_t shift = input_shape[shape_size - 1]; - - double backward_scale = 1.; - double forward_scale = 1.; - - if (norm == 0) // norm = "backward" - { - backward_scale = 1. / shift; - } - else if (norm == 1) // norm = "forward" - { - forward_scale = 1. / shift; - } - else // norm = "ortho" - { - if (inverse) { - backward_scale = 1. / sqrt(shift); - } - else { - forward_scale = 1. / sqrt(shift); - } - } - - std::vector fft_events(n_iter); - - for (size_t i = 0; i < n_iter; ++i) { - std::unique_ptr<_Descriptor_type> desc = - std::make_unique<_Descriptor_type>(shift); - desc->set_value(mkl_dft::config_param::BACKWARD_SCALE, backward_scale); - desc->set_value(mkl_dft::config_param::FORWARD_SCALE, forward_scale); - desc->set_value(mkl_dft::config_param::PLACEMENT, DFTI_NOT_INPLACE); - desc->commit(queue); - - if (inverse) { - fft_events[i] = - mkl_dft::compute_backward<_Descriptor_type, _DataType_input, - _DataType_output>( - *desc, array_1 + i * shift, result + i * shift); - } - else { - fft_events[i] = - mkl_dft::compute_forward<_Descriptor_type, _DataType_input, - _DataType_output>( - *desc, array_1 + i * shift, result + i * shift); - } - } - - sycl::event::wait(fft_events); -} - -template -class dpnp_fft_fft_mathlib_real_to_cmplx_c_kernel; - -template -static DPCTLSyclEventRef - dpnp_fft_fft_mathlib_real_to_cmplx_c(DPCTLSyclQueueRef q_ref, - const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *result_shape, - const size_t shape_size, - const size_t input_size, - const size_t result_size, - size_t inverse, - const size_t norm, - const size_t real) -{ - // avoid warning unused variable - (void)input_size; - - DPCTLSyclEventRef event_ref = nullptr; - if (!shape_size) { - return event_ref; - } - - sycl::queue queue = *(reinterpret_cast(q_ref)); - - _DataType_input *array_1 = - static_cast<_DataType_input *>(const_cast(array1_in)); - _DataType_output *result = static_cast<_DataType_output *>(result_out); - - const size_t n_iter = - std::accumulate(input_shape, input_shape + shape_size - 1, 1, - std::multiplies()); - - const size_t input_shift = input_shape[shape_size - 1]; - const size_t result_shift = result_shape[shape_size - 1]; - - double backward_scale = 1.; - double forward_scale = 1.; - - if (norm == 0) // norm = "backward" - { - if (inverse) { - forward_scale = 1. / result_shift; - } - else { - backward_scale = 1. / result_shift; - } - } - else if (norm == 1) // norm = "forward" - { - if (inverse) { - backward_scale = 1. / result_shift; - } - else { - forward_scale = 1. / result_shift; - } - } - else // norm = "ortho" - { - forward_scale = 1. / sqrt(result_shift); - } - - std::vector fft_events(n_iter); - - for (size_t i = 0; i < n_iter; ++i) { - std::unique_ptr<_Descriptor_type> desc = - std::make_unique<_Descriptor_type>(input_shift); - desc->set_value(mkl_dft::config_param::BACKWARD_SCALE, backward_scale); - desc->set_value(mkl_dft::config_param::FORWARD_SCALE, forward_scale); - desc->set_value(mkl_dft::config_param::PLACEMENT, DFTI_NOT_INPLACE); - desc->commit(queue); - - // real result_size = 2 * result_size, because real type of "result" is - // twice wider than '_DataType_output' - fft_events[i] = - mkl_dft::compute_forward<_Descriptor_type, _DataType_input, - _DataType_output>( - *desc, array_1 + i * input_shift, - result + i * result_shift * 2); - } - - sycl::event::wait(fft_events); - - if (real) // the output size of the rfft function is input_size/2 + 1 so we - // don't need to fill the second half of the output - { - return event_ref; - } - - size_t n_conj = - result_shift % 2 == 0 ? result_shift / 2 - 1 : result_shift / 2; - - sycl::event event; - - sycl::range<2> gws(n_iter, n_conj); - - auto kernel_parallel_for_func = [=](sycl::id<2> global_id) { - size_t i = global_id[0]; - { - size_t j = global_id[1]; - { - *(reinterpret_cast *>(result) + - result_shift * (i + 1) - (j + 1)) = - std::conj( - *(reinterpret_cast *>( - result) + - result_shift * i + (j + 1))); - } - } - }; - - auto kernel_func = [&](sycl::handler &cgh) { - cgh.parallel_for>( - gws, kernel_parallel_for_func); - }; - - event = queue.submit(kernel_func); - - if (inverse) { - event.wait(); - event = oneapi::mkl::vm::conj( - queue, result_size, - reinterpret_cast *>(result), - reinterpret_cast *>(result)); - } - - event_ref = reinterpret_cast(&event); - return DPCTLEvent_Copy(event_ref); -} - -template -DPCTLSyclEventRef dpnp_fft_fft_c(DPCTLSyclQueueRef q_ref, - const void *array1_in, - void *result_out, - const shape_elem_type *input_shape, - const shape_elem_type *result_shape, - size_t shape_size, - long axis, - long input_boundarie, - size_t inverse, - const size_t norm, - const DPCTLEventVectorRef dep_event_vec_ref) -{ - static_assert(is_complex<_DataType_output>::value, - "Output data type must be a complex type."); - - DPCTLSyclEventRef event_ref = nullptr; - - if (!shape_size || !array1_in || !result_out) { - return event_ref; - } - - const size_t result_size = - std::accumulate(result_shape, result_shape + shape_size, 1, - std::multiplies()); - const size_t input_size = - std::accumulate(input_shape, input_shape + shape_size, 1, - std::multiplies()); - - if constexpr (std::is_same<_DataType_output, std::complex>::value || - std::is_same<_DataType_output, std::complex>::value) - { - if constexpr (std::is_same<_DataType_input, - std::complex>::value && - std::is_same<_DataType_output, - std::complex>::value) - { - dpnp_fft_fft_mathlib_cmplx_to_cmplx_c< - _DataType_input, _DataType_output, desc_dp_cmplx_t>( - q_ref, array1_in, result_out, input_shape, result_shape, - shape_size, input_size, result_size, inverse, norm); - } - /* complex-to-complex, single precision */ - else if constexpr (std::is_same<_DataType_input, - std::complex>::value && - std::is_same<_DataType_output, - std::complex>::value) - { - dpnp_fft_fft_mathlib_cmplx_to_cmplx_c< - _DataType_input, _DataType_output, desc_sp_cmplx_t>( - q_ref, array1_in, result_out, input_shape, result_shape, - shape_size, input_size, result_size, inverse, norm); - } - /* real-to-complex, double precision */ - else if constexpr (std::is_same<_DataType_input, double>::value && - std::is_same<_DataType_output, - std::complex>::value) - { - event_ref = - dpnp_fft_fft_mathlib_real_to_cmplx_c<_DataType_input, double, - desc_dp_real_t>( - q_ref, array1_in, result_out, input_shape, result_shape, - shape_size, input_size, result_size, inverse, norm, 0); - } - /* real-to-complex, single precision */ - else if constexpr (std::is_same<_DataType_input, float>::value && - std::is_same<_DataType_output, - std::complex>::value) - { - event_ref = - dpnp_fft_fft_mathlib_real_to_cmplx_c<_DataType_input, float, - desc_sp_real_t>( - q_ref, array1_in, result_out, input_shape, result_shape, - shape_size, input_size, result_size, inverse, norm, 0); - } - else if constexpr (std::is_same<_DataType_input, int32_t>::value || - std::is_same<_DataType_input, int64_t>::value) - { - using CastType = typename _DataType_output::value_type; - - CastType *array1_copy = reinterpret_cast( - dpnp_memory_alloc_c(q_ref, input_size * sizeof(CastType))); - - shape_elem_type *copy_strides = reinterpret_cast( - dpnp_memory_alloc_c(q_ref, sizeof(shape_elem_type))); - *copy_strides = 1; - shape_elem_type *copy_shape = reinterpret_cast( - dpnp_memory_alloc_c(q_ref, sizeof(shape_elem_type))); - *copy_shape = input_size; - shape_elem_type copy_shape_size = 1; - event_ref = dpnp_copyto_c<_DataType_input, CastType>( - q_ref, array1_copy, input_size, copy_shape_size, copy_shape, - copy_strides, array1_in, input_size, copy_shape_size, - copy_shape, copy_strides, NULL, dep_event_vec_ref); - DPCTLEvent_WaitAndThrow(event_ref); - DPCTLEvent_Delete(event_ref); - - event_ref = dpnp_fft_fft_mathlib_real_to_cmplx_c< - CastType, CastType, - std::conditional_t::value, - desc_dp_real_t, desc_sp_real_t>>( - q_ref, array1_copy, result_out, input_shape, result_shape, - shape_size, input_size, result_size, inverse, norm, 0); - - DPCTLEvent_WaitAndThrow(event_ref); - DPCTLEvent_Delete(event_ref); - event_ref = nullptr; - - dpnp_memory_free_c(q_ref, array1_copy); - dpnp_memory_free_c(q_ref, copy_strides); - dpnp_memory_free_c(q_ref, copy_shape); - } - else { - dpnp_fft_fft_sycl_c<_DataType_input, _DataType_output>( - q_ref, array1_in, result_out, input_shape, result_shape, - shape_size, result_size, input_size, axis, input_boundarie, - inverse); - } - } - - return event_ref; -} - -template -void dpnp_fft_fft_c(const void *array1_in, - void *result1, - const shape_elem_type *input_shape, - const shape_elem_type *output_shape, - size_t shape_size, - long axis, - long input_boundarie, - size_t inverse, - const size_t norm) -{ - DPCTLSyclQueueRef q_ref = reinterpret_cast(&DPNP_QUEUE); - DPCTLEventVectorRef dep_event_vec_ref = nullptr; - DPCTLSyclEventRef event_ref = - dpnp_fft_fft_c<_DataType_input, _DataType_output>( - q_ref, array1_in, result1, input_shape, output_shape, shape_size, - axis, input_boundarie, inverse, norm, dep_event_vec_ref); - DPCTLEvent_WaitAndThrow(event_ref); - DPCTLEvent_Delete(event_ref); -} - -template -void (*dpnp_fft_fft_default_c)(const void *, - void *, - const shape_elem_type *, - const shape_elem_type *, - size_t, - long, - long, - size_t, - const size_t) = - dpnp_fft_fft_c<_DataType_input, _DataType_output>; - -template -DPCTLSyclEventRef (*dpnp_fft_fft_ext_c)(DPCTLSyclQueueRef, - const void *, - void *, - const shape_elem_type *, - const shape_elem_type *, - size_t, - long, - long, - size_t, - const size_t, - const DPCTLEventVectorRef) = - dpnp_fft_fft_c<_DataType_input, _DataType_output>; - -void func_map_init_fft_func(func_map_t &fmap) -{ - fmap[DPNPFuncName::DPNP_FN_FFT_FFT][eft_INT][eft_INT] = { - eft_C128, - (void *)dpnp_fft_fft_default_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT][eft_LNG][eft_LNG] = { - eft_C128, - (void *)dpnp_fft_fft_default_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT][eft_FLT][eft_FLT] = { - eft_C64, (void *)dpnp_fft_fft_default_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT][eft_DBL][eft_DBL] = { - eft_C128, (void *)dpnp_fft_fft_default_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT][eft_C64][eft_C64] = { - eft_C64, - (void *) - dpnp_fft_fft_default_c, std::complex>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT][eft_C128][eft_C128] = { - eft_C128, - (void *) - dpnp_fft_fft_default_c, std::complex>}; - - fmap[DPNPFuncName::DPNP_FN_FFT_FFT_EXT][eft_INT][eft_INT] = { - eft_C128, (void *)dpnp_fft_fft_ext_c>, - eft_C64, (void *)dpnp_fft_fft_ext_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT_EXT][eft_LNG][eft_LNG] = { - eft_C128, (void *)dpnp_fft_fft_ext_c>, - eft_C64, (void *)dpnp_fft_fft_ext_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT_EXT][eft_FLT][eft_FLT] = { - eft_C64, (void *)dpnp_fft_fft_ext_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT_EXT][eft_DBL][eft_DBL] = { - eft_C128, (void *)dpnp_fft_fft_ext_c>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT_EXT][eft_C64][eft_C64] = { - eft_C64, - (void *)dpnp_fft_fft_ext_c, std::complex>}; - fmap[DPNPFuncName::DPNP_FN_FFT_FFT_EXT][eft_C128][eft_C128] = { - eft_C128, - (void *)dpnp_fft_fft_ext_c, std::complex>}; - - return; -} diff --git a/dpnp/backend/src/dpnp_fptr.hpp b/dpnp/backend/src/dpnp_fptr.hpp index 212ce848a5c4..f57f9cd449b5 100644 --- a/dpnp/backend/src/dpnp_fptr.hpp +++ b/dpnp/backend/src/dpnp_fptr.hpp @@ -327,7 +327,6 @@ static constexpr DPNPFuncType get_floating_res_type() */ void func_map_init_arraycreation(func_map_t &fmap); void func_map_init_elemwise(func_map_t &fmap); -void func_map_init_fft_func(func_map_t &fmap); void func_map_init_indexing_func(func_map_t &fmap); void func_map_init_linalg(func_map_t &fmap); void func_map_init_mathematical(func_map_t &fmap); diff --git a/dpnp/backend/src/dpnp_iface_fptr.cpp b/dpnp/backend/src/dpnp_iface_fptr.cpp index e2659511d9c4..2e691064e23f 100644 --- a/dpnp/backend/src/dpnp_iface_fptr.cpp +++ b/dpnp/backend/src/dpnp_iface_fptr.cpp @@ -129,7 +129,6 @@ static func_map_t func_map_init() func_map_init_arraycreation(fmap); func_map_init_elemwise(fmap); - func_map_init_fft_func(fmap); func_map_init_indexing_func(fmap); func_map_init_linalg(fmap); func_map_init_mathematical(fmap); diff --git a/dpnp/dpnp_algo/dpnp_algo.pxd b/dpnp/dpnp_algo/dpnp_algo.pxd index afee043421bc..8fb06cc58e44 100644 --- a/dpnp/dpnp_algo/dpnp_algo.pxd +++ b/dpnp/dpnp_algo/dpnp_algo.pxd @@ -37,7 +37,6 @@ cdef extern from "dpnp_iface_fptr.hpp" namespace "DPNPFuncName": # need this na DPNP_FN_CORRELATE_EXT DPNP_FN_EDIFF1D_EXT DPNP_FN_ERF_EXT - DPNP_FN_FFT_FFT_EXT DPNP_FN_MEDIAN_EXT DPNP_FN_MODF_EXT DPNP_FN_PARTITION_EXT diff --git a/dpnp/fft/CMakeLists.txt b/dpnp/fft/CMakeLists.txt deleted file mode 100644 index 3b6146a2a854..000000000000 --- a/dpnp/fft/CMakeLists.txt +++ /dev/null @@ -1,7 +0,0 @@ -# Building dpnp_algo_fft Cython extension - -build_dpnp_cython_ext_with_backend( - dpnp_algo_fft - ${CMAKE_CURRENT_SOURCE_DIR}/dpnp_algo_fft.pyx - dpnp/fft - ) diff --git a/dpnp/fft/dpnp_algo_fft.pyx b/dpnp/fft/dpnp_algo_fft.pyx deleted file mode 100644 index ae1294591343..000000000000 --- a/dpnp/fft/dpnp_algo_fft.pyx +++ /dev/null @@ -1,105 +0,0 @@ -# cython: language_level=3 -# cython: linetrace=True -# -*- coding: utf-8 -*- -# ***************************************************************************** -# Copyright (c) 2016-2024, Intel Corporation -# All rights reserved. -# -# 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. -# -# 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 HOLDER 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. -# ***************************************************************************** - -"""Module Backend (FFT part) - -This module contains interface functions between C backend layer -and the rest of the library - -""" - - -cimport dpnp.dpnp_utils as utils -from dpnp.dpnp_algo cimport * - -__all__ = [ - "dpnp_fft_deprecated", -] - -ctypedef c_dpctl.DPCTLSyclEventRef(*fptr_dpnp_fft_fft_t)(c_dpctl.DPCTLSyclQueueRef, void *, void * , - shape_elem_type * , shape_elem_type * , size_t, long, - long, size_t, size_t, const c_dpctl.DPCTLEventVectorRef) - - -cpdef utils.dpnp_descriptor dpnp_fft_deprecated(utils.dpnp_descriptor input, - size_t input_boundarie, - size_t output_boundarie, - long axis, - size_t inverse, - size_t norm): - - cdef shape_type_c input_shape = input.shape - cdef shape_type_c output_shape = input_shape - - cdef long axis_norm = utils.normalize_axis((axis,), input_shape.size())[0] - output_shape[axis_norm] = output_boundarie - - # convert string type names (dtype) to C enum DPNPFuncType - cdef DPNPFuncType param1_type = dpnp_dtype_to_DPNPFuncType(input.dtype) - - # get the FPTR data structure - cdef DPNPFuncData kernel_data = get_dpnp_function_ptr(DPNP_FN_FFT_FFT_EXT, param1_type, param1_type) - - input_obj = input.get_array() - - # get FPTR function and return type - cdef (DPNPFuncType, void *) ret_type_and_func = utils.get_ret_type_and_func(kernel_data, - input_obj.sycl_device.has_aspect_fp64) - cdef DPNPFuncType return_type = ret_type_and_func[0] - cdef fptr_dpnp_fft_fft_t func = < fptr_dpnp_fft_fft_t > ret_type_and_func[1] - - # create result array with type given by FPTR data - cdef utils.dpnp_descriptor result = utils.create_output_descriptor(output_shape, - return_type, - None, - device=input_obj.sycl_device, - usm_type=input_obj.usm_type, - sycl_queue=input_obj.sycl_queue) - - result_sycl_queue = result.get_array().sycl_queue - - cdef c_dpctl.SyclQueue q = result_sycl_queue - cdef c_dpctl.DPCTLSyclQueueRef q_ref = q.get_queue_ref() - - # call FPTR function - cdef c_dpctl.DPCTLSyclEventRef event_ref = func(q_ref, - input.get_data(), - result.get_data(), - input_shape.data(), - output_shape.data(), - input_shape.size(), - axis_norm, - input_boundarie, - inverse, - norm, - NULL) # dep_events_ref - - with nogil: c_dpctl.DPCTLEvent_WaitAndThrow(event_ref) - c_dpctl.DPCTLEvent_Delete(event_ref) - - return result diff --git a/dpnp/fft/dpnp_iface_fft.py b/dpnp/fft/dpnp_iface_fft.py index c117e2b19770..08f5245ae021 100644 --- a/dpnp/fft/dpnp_iface_fft.py +++ b/dpnp/fft/dpnp_iface_fft.py @@ -33,14 +33,11 @@ it contains: - Interface functions - documentation for the functions - - The functions parameters check """ # pylint: disable=invalid-name -from enum import Enum - import numpy import dpnp @@ -50,9 +47,6 @@ call_origin, checker_throw_axis_error, ) -from dpnp.fft.dpnp_algo_fft import ( - dpnp_fft_deprecated, -) from .dpnp_utils_fft import ( dpnp_fft, @@ -80,24 +74,22 @@ ] -# TODO: remove pylint disable, once new implementation is ready -# pylint: disable=missing-class-docstring -class Norm(Enum): - backward = 0 - forward = 1 - ortho = 2 +_SWAP_DIRECTION_MAP = { + "backward": "forward", + None: "forward", + "ortho": "ortho", + "forward": "backward", +} -# TODO: remove pylint disable, once new implementation is ready -# pylint: disable=missing-function-docstring -def get_validated_norm(norm): - if norm is None or norm == "backward": - return Norm.backward - if norm == "forward": - return Norm.forward - if norm == "ortho": - return Norm.ortho - raise ValueError("Unknown norm value.") +def _swap_direction(norm): + try: + return _SWAP_DIRECTION_MAP[norm] + except KeyError: + raise ValueError( + f'Invalid norm value {norm}; should be None, "backward", ' + '"ortho" or "forward".' + ) from None def fft(a, n=None, axis=-1, norm=None, out=None): @@ -425,58 +417,98 @@ def fftshift(x, axes=None): return dpnp.roll(x, shift, axes) -def hfft(x, n=None, axis=-1, norm=None): +def hfft(a, n=None, axis=-1, norm=None, out=None): """ - Compute the one-dimensional discrete Fourier Transform of a signal that has - Hermitian symmetry. + Compute the FFT of a signal that has Hermitian symmetry, i.e., + a real spectrum. For full documentation refer to :obj:`numpy.fft.hfft`. - Limitations - ----------- - Parameter `x` is supported either as :class:`dpnp.ndarray`. - Parameter `norm` is unsupported. - Only `dpnp.float64`, `dpnp.float32`, `dpnp.int64`, `dpnp.int32`, - `dpnp.complex128` data types are supported. - Otherwise the function will be executed sequentially on CPU. + Parameters + ---------- + a : {dpnp.ndarray, usm_ndarray} + Input array. + n : {None, int}, optional + Length of the transformed axis of the output. + For `n` output points, ``n//2+1`` input points are necessary. If the + input is longer than this, it is cropped. If it is shorter than this, + it is padded with zeros. If `n` is not given, it is taken to be + ``2*(m-1)`` where ``m`` is the length of the input along the axis + specified by `axis`. Default: ``None``. + axis : int, optional + Axis over which to compute the FFT. If not given, the last axis is + used. Default: ``-1``. + norm : {None, "backward", "ortho", "forward"}, optional + Normalization mode (see :obj:`dpnp.fft`). + Indicates which direction of the forward/backward pair of transforms + is scaled and with what normalization factor. ``None`` is an alias of + the default option ``"backward"``. + Default: ``"backward"``. + out : {None, dpnp.ndarray, usm_ndarray}, optional + If provided, the result will be placed in this array. It should be + of the appropriate shape and dtype. + Default: ``None``. - """ + Returns + ------- + out : dpnp.ndarray + The truncated or zero-padded input, transformed along the axis + indicated by `axis`, or the last one if `axis` is not specified. + The length of the transformed axis is `n`, or, if `n` is not given, + ``2*(m-1)`` where ``m`` is the length of the transformed axis of the + input. To get an odd number of output points, `n` must be specified, + for instance as ``2*m - 1`` in the typical case. - x_desc = dpnp.get_dpnp_descriptor(x, copy_when_nondefault_queue=False) - # TODO: enable implementation - # pylint: disable=condition-evals-to-constant - if x_desc and 0: - norm_ = get_validated_norm(norm) + See Also + -------- + :obj:`dpnp.fft` : For definition of the DFT and conventions used. + :obj:`dpnp.fft.rfft` : The one-dimensional FFT of real input. + :obj:`dpnp.fft.ihfft` :The inverse of :obj:`dpnp.fft.hfft`. - if axis is None: - axis_param = -1 # the most right dimension (default value) - else: - axis_param = axis - if n is None: - input_boundarie = x_desc.shape[axis_param] - else: - input_boundarie = n + Notes + ----- + :obj:`dpnp.fft.hfft`/:obj:`dpnp.fft.ihfft` are a pair analogous to + :obj:`dpnp.fft.rfft`/:obj:`dpnp.fft.irfft`, but for the opposite case: + here the signal has Hermitian symmetry in the time domain and is real in + the frequency domain. So here it's :obj:`dpnp.fft.hfft` for which you must + supply the length of the result if it is to be odd. - if x.size < 1: - pass # let fallback to handle exception - elif input_boundarie < 1: - pass # let fallback to handle exception - elif norm is not None: - pass - else: - output_boundarie = input_boundarie + * even: ``ihfft(hfft(a, 2*len(a) - 2)) == a``, within roundoff error, + * odd: ``ihfft(hfft(a, 2*len(a) - 1)) == a``, within roundoff error. + + The correct interpretation of the Hermitian input depends on the length of + the original data, as given by `n`. This is because each input shape could + correspond to either an odd or even length signal. By default, + :obj:`dpnp.fft.hfft` assumes an even output length which puts the last + entry at the Nyquist frequency; aliasing with its symmetric counterpart. + By Hermitian symmetry, the value is thus treated as purely real. To avoid + losing information, the correct length of the real input **must** be given. + + Examples + -------- + >>> import dpnp as np + >>> signal = np.array([1, 2, 3, 4, 3, 2]) + >>> np.fft.fft(signal) + array([15.+0.j, -4.+0.j, 0.+0.j, -1.-0.j, 0.+0.j, -4.+0.j]) + >>> np.fft.hfft(signal[:4]) # Input first half of signal + array([15., -4., 0., -1., 0., -4.]) + >>> np.fft.hfft(signal, 6) # Input entire signal and truncate + array([15., -4., 0., -1., 0., -4.]) + + >>> signal = np.array([[1, 1.j], [-1.j, 2]]) + >>> np.conj(signal.T) - signal # check Hermitian symmetry + array([[ 0.-0.j, -0.+0.j], # may vary + [ 0.+0.j, 0.-0.j]]) + >>> freq_spectrum = np.fft.hfft(signal) + >>> freq_spectrum + array([[ 1., 1.], + [ 2., -2.]]) - return dpnp_fft_deprecated( - x_desc, - input_boundarie, - output_boundarie, - axis_param, - False, - norm_.value, - ).get_pyobj() + """ - return call_origin(numpy.fft.hfft, x, n, axis, norm) + new_norm = _swap_direction(norm) + return irfft(dpnp.conjugate(a), n=n, axis=axis, norm=new_norm, out=out) def ifft(a, n=None, axis=-1, norm=None, out=None): @@ -691,60 +723,76 @@ def ifftshift(x, axes=None): return dpnp.roll(x, shift, axes) -def ihfft(x, n=None, axis=-1, norm=None): +def ihfft(a, n=None, axis=-1, norm=None, out=None): """ - Compute inverse one-dimensional discrete Fourier Transform of a signal that - has Hermitian symmetry. + Compute the inverse FFT of a signal that has Hermitian symmetry. For full documentation refer to :obj:`numpy.fft.ihfft`. - Limitations - ----------- - Parameter `x` is supported either as :class:`dpnp.ndarray`. - Parameter `norm` is unsupported. - Only `dpnp.float64`, `dpnp.float32`, `dpnp.int64`, `dpnp.int32`, - `dpnp.complex128` data types are supported. - Otherwise the function will be executed sequentially on CPU. + Parameters + ---------- + a : {dpnp.ndarray, usm_ndarray} + Input array. + n : {None, int}, optional + Length of the inverse FFT, the number of points along + transformation axis in the input to use. If `n` is smaller than + the length of the input, the input is cropped. If it is larger, + the input is padded with zeros. If `n` is not given, the length of + the input along the axis specified by `axis` is used. + Default: ``None``. + axis : int, optional + Axis over which to compute the FFT. If not given, the last axis is + used. Default: ``-1``. + norm : {None, "backward", "ortho", "forward"}, optional + Normalization mode (see :obj:`dpnp.fft`). + Indicates which direction of the forward/backward pair of transforms + is scaled and with what normalization factor. ``None`` is an alias of + the default option ``"backward"``. + Default: ``"backward"``. + out : {None, dpnp.ndarray or usm_ndarray of complex dtype}, optional + If provided, the result will be placed in this array. It should be + of the appropriate shape and dtype. + Default: ``None``. - """ + Returns + ------- + out : dpnp.ndarray of complex dtype + The truncated or zero-padded input, transformed along the axis + indicated by `axis`, or the last one if `axis` is not specified. + The length of the transformed axis is ``n//2 + 1``. - x_desc = dpnp.get_dpnp_descriptor(x, copy_when_nondefault_queue=False) - # TODO: enable implementation - # pylint: disable=condition-evals-to-constant - if x_desc and 0: - norm_ = get_validated_norm(norm) + See Also + -------- + :obj:`dpnp.fft` : For definition of the DFT and conventions used. + :obj:`dpnp.fft.hfft` : Compute the FFT of a signal that has + Hermitian symmetry. + :obj:`dpnp.fft.irfft` : The inverse of :obj:`dpnp.fft.rfft`. - if axis is None: - axis_param = -1 # the most right dimension (default value) - else: - axis_param = axis + Notes + ----- + :obj:`dpnp.fft.hfft`/:obj:`dpnp.fft.ihfft` are a pair analogous to + :obj:`dpnp.fft.rfft`/:obj:`dpnp.fft.irfft`, but for the opposite case: + here the signal has Hermitian symmetry in the time domain and is real in + the frequency domain. So here it's :obj:`dpnp.fft.hfft` for which you must + supply the length of the result if it is to be odd. - if n is None: - input_boundarie = x_desc.shape[axis_param] - else: - input_boundarie = n + * even: ``ihfft(hfft(a, 2*len(a) - 2)) == a``, within roundoff error, + * odd: ``ihfft(hfft(a, 2*len(a) - 1)) == a``, within roundoff error. - if x_desc.size < 1: - pass # let fallback to handle exception - elif input_boundarie < 1: - pass # let fallback to handle exception - elif norm is not None: - pass - elif n is not None: - pass - else: - output_boundarie = input_boundarie + Examples + -------- + >>> import dpnp as np + >>> spectrum = np.array([ 15, -4, 0, -1, 0, -4]) + >>> np.fft.ifft(spectrum) + array([1.+0.j, 2.+0.j, 3.+0.j, 4.+0.j, 3.+0.j, 2.+0.j]) # may vary + >>> np.fft.ihfft(spectrum) + array([1.-0.j, 2.-0.j, 3.-0.j, 4.-0.j]) # may vary - return dpnp_fft_deprecated( - x_desc, - input_boundarie, - output_boundarie, - axis_param, - True, - norm_.value, - ).get_pyobj() + """ - return call_origin(numpy.fft.ihfft, x, n, axis, norm) + new_norm = _swap_direction(norm) + res = rfft(a, n=n, axis=axis, norm=new_norm, out=out) + return dpnp.conjugate(res, out=out) def irfft(a, n=None, axis=-1, norm=None, out=None): @@ -771,9 +819,9 @@ def irfft(a, n=None, axis=-1, norm=None, out=None): Input array. n : {None, int}, optional Length of the transformed axis of the output. - For `n` output points, ``n//2+1`` input points are necessary. If the + For `n` output points, ``n//2+1`` input points are necessary. If the input is longer than this, it is cropped. If it is shorter than this, - it is padded with zeros. If `n` is not given, it is taken to be + it is padded with zeros. If `n` is not given, it is taken to be ``2*(m-1)`` where ``m`` is the length of the input along the axis specified by `axis`. Default: ``None``. axis : int, optional @@ -822,7 +870,7 @@ def irfft(a, n=None, axis=-1, norm=None, out=None): thus resample a series to `m` points via Fourier interpolation by: ``a_resamp = irfft(rfft(a), m)``. - The correct interpretation of the hermitian input depends on the length of + The correct interpretation of the Hermitian input depends on the length of the original data, as given by `n`. This is because each input shape could correspond to either an odd or even length signal. By default, :obj:`dpnp.fft.irfft` assumes an even output length which puts the last diff --git a/tests/test_fft.py b/tests/test_fft.py index 1f7b389ed75c..5d142f7e31a8 100644 --- a/tests/test_fft.py +++ b/tests/test_fft.py @@ -343,141 +343,89 @@ def test_fft_validate_out(self): assert_raises(TypeError, dpnp.fft.fft, a, out=out) -class TestRfft: - def setup_method(self): - numpy.random.seed(42) +class TestFftfreq: + @pytest.mark.parametrize("func", ["fftfreq", "rfftfreq"]) + @pytest.mark.parametrize("n", [10, 20]) + @pytest.mark.parametrize("d", [0.5, 2]) + def test_fftfreq(self, func, n, d): + expected = getattr(dpnp.fft, func)(n, d) + result = getattr(numpy.fft, func)(n, d) + assert_dtype_allclose(expected, result) - @pytest.mark.parametrize( - "dtype", get_all_dtypes(no_bool=True, no_complex=True) - ) - @pytest.mark.parametrize( - "shape", [(64,), (8, 8), (4, 16), (4, 4, 4), (2, 4, 4, 2)] - ) - def test_fft_rfft(self, dtype, shape): - np_data = numpy.arange(64, dtype=dtype).reshape(shape) - dpnp_data = dpnp.arange(64, dtype=dtype).reshape(shape) + @pytest.mark.parametrize("func", ["fftfreq", "rfftfreq"]) + def test_error(self, func): + # n should be an integer + assert_raises(ValueError, getattr(dpnp.fft, func), 10.0) - np_res = numpy.fft.rfft(np_data) - dpnp_res = dpnp.fft.rfft(dpnp_data) + # d should be an scalar + assert_raises(ValueError, getattr(dpnp.fft, func), 10, (2,)) - assert_dtype_allclose(dpnp_res, np_res, check_only_type_kind=True) - @pytest.mark.parametrize( - "dtype", get_all_dtypes(no_bool=True, no_complex=True) - ) - @pytest.mark.parametrize("n", [None, 5, 20]) - @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) - def test_fft_1D(self, dtype, n, norm): - x = dpnp.linspace(-1, 1, 11, dtype=dtype) - a = dpnp.sin(x) - a_np = dpnp.asnumpy(a) +class TestFftshift: + @pytest.mark.parametrize("func", ["fftshift", "ifftshift"]) + @pytest.mark.parametrize("axes", [None, 1, (0, 1)]) + def test_fftshift(self, func, axes): + x = dpnp.arange(12).reshape(3, 4) + x_np = x.asnumpy() + expected = getattr(dpnp.fft, func)(x, axes=axes) + result = getattr(numpy.fft, func)(x_np, axes=axes) + assert_dtype_allclose(expected, result) - result = dpnp.fft.rfft(a, n=n, norm=norm) - expected = numpy.fft.rfft(a_np, n=n, norm=norm) - assert_dtype_allclose(result, expected, check_only_type_kind=True) +class TestHfft: + def setup_method(self): + numpy.random.seed(42) + + @pytest.mark.parametrize("dtype", get_all_dtypes(no_complex=True)) @pytest.mark.parametrize("n", [None, 5, 20]) @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) - def test_fft_bool(self, n, norm): - a = dpnp.ones(11, dtype=dpnp.bool) + def test_hfft_1D(self, dtype, n, norm): + x = dpnp.linspace(-1, 1, 11, dtype=dtype) + a = dpnp.sin(x) a_np = dpnp.asnumpy(a) - result = dpnp.fft.rfft(a, n=n, norm=norm) - expected = numpy.fft.rfft(a_np, n=n, norm=norm) - assert_dtype_allclose(result, expected, check_only_type_kind=True) - - @pytest.mark.parametrize("dtype", get_float_dtypes()) - @pytest.mark.parametrize("n", [None, 5, 8]) - @pytest.mark.parametrize("axis", [-1, 1, 0]) - @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) - @pytest.mark.parametrize("order", ["C", "F"]) - def test_fft_1D_on_2D_array(self, dtype, n, axis, norm, order): - a_np = numpy.arange(12, dtype=dtype).reshape(3, 4, order=order) - a = dpnp.asarray(a_np) - - result = dpnp.fft.rfft(a, n=n, axis=axis, norm=norm) - expected = numpy.fft.rfft(a_np, n=n, axis=axis, norm=norm) - assert_dtype_allclose(result, expected, check_only_type_kind=True) - - @pytest.mark.parametrize("dtype", get_float_dtypes()) - @pytest.mark.parametrize("n", [None, 5, 8]) - @pytest.mark.parametrize("axis", [0, 1, 2]) - @pytest.mark.parametrize("norm", ["forward", "backward", "ortho"]) - @pytest.mark.parametrize("order", ["C", "F"]) - def test_fft_1D_on_3D_array(self, dtype, n, axis, norm, order): - a_np = numpy.arange(24, dtype=dtype).reshape(2, 3, 4, order=order) - a = dpnp.asarray(a_np) - - result = dpnp.fft.rfft(a, n=n, axis=axis, norm=norm) - expected = numpy.fft.rfft(a_np, n=n, axis=axis, norm=norm) - assert_dtype_allclose(result, expected, check_only_type_kind=True) - - @pytest.mark.parametrize("n", [None, 5, 20]) - def test_fft_usm_ndarray(self, n): - x = dpt.linspace(-1, 1, 11) - a_usm = dpt.asarray(dpt.sin(x)) - a_np = dpt.asnumpy(a_usm) - out_shape = a_usm.shape[0] // 2 + 1 if n is None else n // 2 + 1 - out_dtype = map_dtype_to_device(dpnp.complex128, a_usm.sycl_device) - out = dpt.empty(out_shape, dtype=out_dtype) - - result = dpnp.fft.rfft(a_usm, n=n, out=out) - assert out is result.get_array() - expected = numpy.fft.rfft(a_np, n=n) + result = dpnp.fft.hfft(a, n=n, norm=norm) + expected = numpy.fft.hfft(a_np, n=n, norm=norm) + # check_only_type_kind=True since numpy always returns float64 + # but dpnp return float32 if input is float32 assert_dtype_allclose(result, expected, check_only_type_kind=True) - @pytest.mark.parametrize("dtype", get_float_dtypes()) + @pytest.mark.parametrize("dtype", get_complex_dtypes()) @pytest.mark.parametrize("n", [None, 5, 20]) @pytest.mark.parametrize("norm", ["forward", "backward", "ortho"]) - def test_fft_1D_out(self, dtype, n, norm): + def test_hfft_1D_complex(self, dtype, n, norm): x = dpnp.linspace(-1, 1, 11) a = dpnp.sin(x) + 1j * dpnp.cos(x) a = dpnp.asarray(a, dtype=dtype) a_np = dpnp.asnumpy(a) - out_shape = a.shape[0] // 2 + 1 if n is None else n // 2 + 1 - out_dtype = dpnp.complex64 if dtype == dpnp.float32 else dpnp.complex128 - out = dpnp.empty(out_shape, dtype=out_dtype) - - result = dpnp.fft.rfft(a, n=n, norm=norm, out=out) - assert out is result - expected = numpy.fft.rfft(a_np, n=n, norm=norm) + result = dpnp.fft.hfft(a, n=n, norm=norm) + expected = numpy.fft.hfft(a_np, n=n, norm=norm) assert_dtype_allclose(result, expected, check_only_type_kind=True) - @pytest.mark.parametrize("dtype", get_float_dtypes()) - @pytest.mark.parametrize("n", [None, 5, 8]) - @pytest.mark.parametrize("axis", [-1, 0]) + @pytest.mark.parametrize( + "dtype", get_all_dtypes(no_bool=True, no_complex=True) + ) + @pytest.mark.parametrize("n", [None, 5, 20]) @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) - @pytest.mark.parametrize("order", ["C", "F"]) - def test_fft_1D_on_2D_array_out(self, dtype, n, axis, norm, order): - a_np = numpy.arange(12, dtype=dtype).reshape(3, 4, order=order) - a = dpnp.asarray(a_np) - - out_shape = list(a.shape) - out_shape[axis] = a.shape[axis] // 2 + 1 if n is None else n // 2 + 1 - out_shape = tuple(out_shape) - out_dtype = dpnp.complex64 if dtype == dpnp.float32 else dpnp.complex128 - out = dpnp.empty(out_shape, dtype=out_dtype) + def test_ihfft_1D(self, dtype, n, norm): + x = dpnp.linspace(-1, 1, 11, dtype=dtype) + a = dpnp.sin(x) + a_np = dpnp.asnumpy(a) - result = dpnp.fft.rfft(a, n=n, axis=axis, norm=norm, out=out) - assert out is result - expected = numpy.fft.rfft(a_np, n=n, axis=axis, norm=norm) + result = dpnp.fft.ihfft(a, n=n, norm=norm) + expected = numpy.fft.ihfft(a_np, n=n, norm=norm) assert_dtype_allclose(result, expected, check_only_type_kind=True) - @pytest.mark.parametrize("xp", [numpy, dpnp]) - def test_fft_error(self, xp): - a = xp.ones((4, 3), dtype=xp.complex64) - # invalid dtype of input array for r2c FFT - if xp == dpnp: - # stock NumPy-1.26 ignores imaginary part - # IntelĀ® NumPy, dpnp, stock NumPy-2.0 return TypeError - assert_raises(TypeError, xp.fft.rfft, a) + @pytest.mark.parametrize("n", [None, 5, 20]) + @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) + def test_ihfft_bool(self, n, norm): + a = dpnp.ones(11, dtype=dpnp.bool) + a_np = dpnp.asnumpy(a) - def test_fft_validate_out(self): - # Invalid shape for r2c FFT - a = dpnp.ones((10,), dtype=dpnp.float32) - out = dpnp.empty((10,), dtype=dpnp.complex64) - assert_raises(ValueError, dpnp.fft.rfft, a, out=out) + result = dpnp.fft.ihfft(a, n=n, norm=norm) + expected = numpy.fft.ihfft(a_np, n=n, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) class TestIrfft: @@ -602,30 +550,138 @@ def test_fft_validate_out(self): assert_raises(TypeError, dpnp.fft.irfft, a, out=out) -class TestFftfreq: - @pytest.mark.parametrize("func", ["fftfreq", "rfftfreq"]) - @pytest.mark.parametrize("n", [10, 20]) - @pytest.mark.parametrize("d", [0.5, 2]) - def test_fftfreq(self, func, n, d): - expected = getattr(dpnp.fft, func)(n, d) - result = getattr(numpy.fft, func)(n, d) - assert_dtype_allclose(expected, result) +class TestRfft: + def setup_method(self): + numpy.random.seed(42) - @pytest.mark.parametrize("func", ["fftfreq", "rfftfreq"]) - def test_error(self, func): - # n should be an integer - assert_raises(ValueError, getattr(dpnp.fft, func), 10.0) + @pytest.mark.parametrize( + "dtype", get_all_dtypes(no_bool=True, no_complex=True) + ) + @pytest.mark.parametrize( + "shape", [(64,), (8, 8), (4, 16), (4, 4, 4), (2, 4, 4, 2)] + ) + def test_fft_rfft(self, dtype, shape): + np_data = numpy.arange(64, dtype=dtype).reshape(shape) + dpnp_data = dpnp.arange(64, dtype=dtype).reshape(shape) - # d should be an scalar - assert_raises(ValueError, getattr(dpnp.fft, func), 10, (2,)) + np_res = numpy.fft.rfft(np_data) + dpnp_res = dpnp.fft.rfft(dpnp_data) + assert_dtype_allclose(dpnp_res, np_res, check_only_type_kind=True) -class TestFftshift: - @pytest.mark.parametrize("func", ["fftshift", "ifftshift"]) - @pytest.mark.parametrize("axes", [None, 1, (0, 1)]) - def test_fftshift(self, func, axes): - x = dpnp.arange(12).reshape(3, 4) - x_np = x.asnumpy() - expected = getattr(dpnp.fft, func)(x, axes=axes) - result = getattr(numpy.fft, func)(x_np, axes=axes) - assert_dtype_allclose(expected, result) + @pytest.mark.parametrize( + "dtype", get_all_dtypes(no_bool=True, no_complex=True) + ) + @pytest.mark.parametrize("n", [None, 5, 20]) + @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) + def test_fft_1D(self, dtype, n, norm): + x = dpnp.linspace(-1, 1, 11, dtype=dtype) + a = dpnp.sin(x) + a_np = dpnp.asnumpy(a) + + result = dpnp.fft.rfft(a, n=n, norm=norm) + expected = numpy.fft.rfft(a_np, n=n, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("n", [None, 5, 20]) + @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) + def test_fft_bool(self, n, norm): + a = dpnp.ones(11, dtype=dpnp.bool) + a_np = dpnp.asnumpy(a) + + result = dpnp.fft.rfft(a, n=n, norm=norm) + expected = numpy.fft.rfft(a_np, n=n, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("dtype", get_float_dtypes()) + @pytest.mark.parametrize("n", [None, 5, 8]) + @pytest.mark.parametrize("axis", [-1, 1, 0]) + @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) + @pytest.mark.parametrize("order", ["C", "F"]) + def test_fft_1D_on_2D_array(self, dtype, n, axis, norm, order): + a_np = numpy.arange(12, dtype=dtype).reshape(3, 4, order=order) + a = dpnp.asarray(a_np) + + result = dpnp.fft.rfft(a, n=n, axis=axis, norm=norm) + expected = numpy.fft.rfft(a_np, n=n, axis=axis, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("dtype", get_float_dtypes()) + @pytest.mark.parametrize("n", [None, 5, 8]) + @pytest.mark.parametrize("axis", [0, 1, 2]) + @pytest.mark.parametrize("norm", ["forward", "backward", "ortho"]) + @pytest.mark.parametrize("order", ["C", "F"]) + def test_fft_1D_on_3D_array(self, dtype, n, axis, norm, order): + a_np = numpy.arange(24, dtype=dtype).reshape(2, 3, 4, order=order) + a = dpnp.asarray(a_np) + + result = dpnp.fft.rfft(a, n=n, axis=axis, norm=norm) + expected = numpy.fft.rfft(a_np, n=n, axis=axis, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("n", [None, 5, 20]) + def test_fft_usm_ndarray(self, n): + x = dpt.linspace(-1, 1, 11) + a_usm = dpt.asarray(dpt.sin(x)) + a_np = dpt.asnumpy(a_usm) + out_shape = a_usm.shape[0] // 2 + 1 if n is None else n // 2 + 1 + out_dtype = map_dtype_to_device(dpnp.complex128, a_usm.sycl_device) + out = dpt.empty(out_shape, dtype=out_dtype) + + result = dpnp.fft.rfft(a_usm, n=n, out=out) + assert out is result.get_array() + expected = numpy.fft.rfft(a_np, n=n) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("dtype", get_float_dtypes()) + @pytest.mark.parametrize("n", [None, 5, 20]) + @pytest.mark.parametrize("norm", ["forward", "backward", "ortho"]) + def test_fft_1D_out(self, dtype, n, norm): + x = dpnp.linspace(-1, 1, 11) + a = dpnp.sin(x) + 1j * dpnp.cos(x) + a = dpnp.asarray(a, dtype=dtype) + a_np = dpnp.asnumpy(a) + + out_shape = a.shape[0] // 2 + 1 if n is None else n // 2 + 1 + out_dtype = dpnp.complex64 if dtype == dpnp.float32 else dpnp.complex128 + out = dpnp.empty(out_shape, dtype=out_dtype) + + result = dpnp.fft.rfft(a, n=n, norm=norm, out=out) + assert out is result + expected = numpy.fft.rfft(a_np, n=n, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("dtype", get_float_dtypes()) + @pytest.mark.parametrize("n", [None, 5, 8]) + @pytest.mark.parametrize("axis", [-1, 0]) + @pytest.mark.parametrize("norm", [None, "forward", "ortho"]) + @pytest.mark.parametrize("order", ["C", "F"]) + def test_fft_1D_on_2D_array_out(self, dtype, n, axis, norm, order): + a_np = numpy.arange(12, dtype=dtype).reshape(3, 4, order=order) + a = dpnp.asarray(a_np) + + out_shape = list(a.shape) + out_shape[axis] = a.shape[axis] // 2 + 1 if n is None else n // 2 + 1 + out_shape = tuple(out_shape) + out_dtype = dpnp.complex64 if dtype == dpnp.float32 else dpnp.complex128 + out = dpnp.empty(out_shape, dtype=out_dtype) + + result = dpnp.fft.rfft(a, n=n, axis=axis, norm=norm, out=out) + assert out is result + expected = numpy.fft.rfft(a_np, n=n, axis=axis, norm=norm) + assert_dtype_allclose(result, expected, check_only_type_kind=True) + + @pytest.mark.parametrize("xp", [numpy, dpnp]) + def test_fft_error(self, xp): + a = xp.ones((4, 3), dtype=xp.complex64) + # invalid dtype of input array for r2c FFT + if xp == dpnp: + # stock NumPy-1.26 ignores imaginary part + # IntelĀ® NumPy, dpnp, stock NumPy-2.0 return TypeError + assert_raises(TypeError, xp.fft.rfft, a) + + def test_fft_validate_out(self): + # Invalid shape for r2c FFT + a = dpnp.ones((10,), dtype=dpnp.float32) + out = dpnp.empty((10,), dtype=dpnp.complex64) + assert_raises(ValueError, dpnp.fft.rfft, a, out=out) diff --git a/tests/test_sycl_queue.py b/tests/test_sycl_queue.py index 6dcdaff2da81..57e2e6dc4c9f 100644 --- a/tests/test_sycl_queue.py +++ b/tests/test_sycl_queue.py @@ -1223,15 +1223,17 @@ def test_out_multi_dot(device): assert_sycl_queue_equal(result.sycl_queue, exec_q) -@pytest.mark.parametrize("func", ["fft", "ifft", "rfft", "irfft"]) +@pytest.mark.parametrize( + "func", ["fft", "ifft", "rfft", "irfft", "hfft", "ihfft"] +) @pytest.mark.parametrize( "device", valid_devices, ids=[device.filter_string for device in valid_devices], ) def test_fft(func, device): - dtype = numpy.float64 if func == "rfft" else numpy.complex128 - data = numpy.arange(100, dtype=dtype) + dtype = numpy.float64 if func in ["rfft", "ihfft"] else numpy.complex128 + data = numpy.arange(20, dtype=dtype) dpnp_data = dpnp.array(data, device=device) expected = getattr(numpy.fft, func)(data) diff --git a/tests/test_usm_type.py b/tests/test_usm_type.py index dbd97dd8896c..b81b9f25edda 100644 --- a/tests/test_usm_type.py +++ b/tests/test_usm_type.py @@ -938,10 +938,12 @@ def test_eigenvalue(func, shape, usm_type): assert a.usm_type == dp_val.usm_type -@pytest.mark.parametrize("func", ["fft", "ifft", "rfft", "irfft"]) +@pytest.mark.parametrize( + "func", ["fft", "ifft", "rfft", "irfft", "hfft", "ihfft"] +) @pytest.mark.parametrize("usm_type", list_of_usm_types, ids=list_of_usm_types) def test_fft(func, usm_type): - dtype = dp.float32 if func == "rfft" else dp.complex64 + dtype = dp.float32 if func in ["rfft", "ihfft"] else dp.complex64 dpnp_data = dp.arange(100, usm_type=usm_type, dtype=dtype) result = getattr(dp.fft, func)(dpnp_data) diff --git a/tests/third_party/cupy/fft_tests/test_fft.py b/tests/third_party/cupy/fft_tests/test_fft.py index 52a2da0a0033..4a25f738dcc7 100644 --- a/tests/third_party/cupy/fft_tests/test_fft.py +++ b/tests/third_party/cupy/fft_tests/test_fft.py @@ -317,39 +317,49 @@ def test_irfftn(self, dtype): xp.fft.irfftn(a, s=self.s, axes=self.axes, norm=self.norm) +@pytest.mark.usefixtures("skip_forward_backward") @testing.parameterize( *testing.product( { "n": [None, 5, 10, 15], "shape": [(10,), (10, 10)], - "norm": [None, "ortho"], + "norm": [None, "backward", "ortho", "forward", ""], } ) ) -@pytest.mark.usefixtures("allow_fall_back_on_numpy") class TestHfft: @testing.for_all_dtypes() @testing.numpy_cupy_allclose( rtol=1e-4, - atol=1e-7, + atol=2e-6, + accept_error=ValueError, + contiguous_check=False, type_check=has_support_aspect64(), ) def test_hfft(self, xp, dtype): a = testing.shaped_random(self.shape, xp, dtype) out = xp.fft.hfft(a, n=self.n, norm=self.norm) + if xp is np and dtype in [np.float16, np.float32, np.complex64]: + out = out.astype(np.float32) + return out @testing.for_all_dtypes(no_complex=True) @testing.numpy_cupy_allclose( rtol=1e-4, atol=1e-7, + accept_error=ValueError, + contiguous_check=False, type_check=has_support_aspect64(), ) def test_ihfft(self, xp, dtype): a = testing.shaped_random(self.shape, xp, dtype) out = xp.fft.ihfft(a, n=self.n, norm=self.norm) + if xp is np and dtype in [np.float16, np.float32, np.complex64]: + out = out.astype(np.complex64) + return out