QuState / PhastFT

A high-performance, "quantum-inspired" Fast Fourier Transform (FFT) library written in pure and safe Rust.
Apache License 2.0
203 stars 8 forks source link

Discrepancy in Inverse FFT output order #24

Closed smu160 closed 5 months ago

smu160 commented 5 months ago

The output order of running FFT followed by Inverse FFT seems to be scrambled. As an example, see the following outputs:

Original signal

Signal Reals: [1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0, 15.0, 16.0]
Signal Imags: [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

Using PhastFT

FFT Output Reals: [136.0, -7.999999999999999, -7.999999999999999, -7.999999999999999, -8.0, -7.999999999999999, -8.0, -8.000000000000002, -8.0, -7.999999999999999, -7.999999999999999, -7.999999999999999, -8.0, -7.999999999999999, -8.0, -8.000000000000002]
FFT Output Imags: [0.0, 40.218715937006785, 19.31370849898476, 11.972846101323913, 8.0, 5.345429103354393, 3.3137084989847594, 1.5912989390372587, 0.0, -1.5912989390372658, -3.313708498984761, -5.345429103354395, -8.0, -11.972846101323915, -19.31370849898476, -40.21871593700678]

IFFT Output Reals: [1.0, 2.0, 3.0, 4.000000000000001, 13.0, 14.0, 15.0, 16.0, 9.0, 10.0, 11.0, 12.0, 5.0, 6.0, 7.000000000000001, 8.0]
IFFT Output Imags: [-3.3306690738754696e-16, 3.3306690738754696e-16, 1.1102230246251565e-16, 1.110223024625157e-16, -1.1102230246251565e-16, -1.1102230246251565e-16, -3.3306690738754696e-16, 7.771561172376096e-16, 1.1102230246251565e-16, -5.551115123125783e-16, -3.3306690738754696e-16, -3.330669073875469e-16, 3.3306690738754696e-16, -1.1102230246251565e-16, 5.551115123125783e-16, -1.110223024625157e-16]

Using RustFFT

RustFFT FFT Reals: [136.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0, -8.0]
RustFFT FFT Imags: [0.0, 40.218715937006785, 19.31370849898476, 11.972846101323912, 8.0, 5.345429103354389, 3.313708498984761, 1.5912989390372658, 0.0, -1.5912989390372658, -3.3137084989847594, -5.345429103354389, -8.0, -11.972846101323912, -19.31370849898476, -40.218715937006785]

RustFFT IFFT Reals: [1.0, 2.0000000000000004, 3.000000000000001, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0, 15.0, 16.0]
RustFFT IFFT Imags: [1.1102230246251565e-16, 0.0, 0.0, 0.0, -1.1102230246251565e-16, 0.0, 0.0, 0.0, 1.1102230246251565e-16, 0.0, 0.0, 0.0, -1.1102230246251565e-16, 0.0, 0.0, 0.0]

Expected Behavior The Inverse FFT output should match the original input signal, as RustFFT produces.

Actual Behavior The Inverse FFT output that PhastFT produces does not match the expected output; the real components appears to be out of order.

Steps to Reproduce Please run the following test

    #[test]
    fn ifft_with_manual_scaling() {
        let n = 4;
        let big_n = 1 << n;
        let mut reals: Vec<_> = (1..=big_n).map(|i| i as f64).collect();
        let mut imags: Vec<_> = (1..=big_n).map(|_| 0.0f64).collect();

        println!("Signal Reals: {:?}", reals);
        println!("Signal Imags: {:?}\n", imags);

        // Perform FFT
        fft_64(&mut reals, &mut imags, Direction::Forward);

        println!("FFT Output Reals: {:?}", reals);
        println!("FFT Output Imags: {:?}", imags);

        // Perform IFFT
        fft_64(&mut reals, &mut imags, Direction::Reverse);

        let scale_factor = reals.len() as f64;

        // Manually scale the output by the dataset length
        for i in 0..reals.len() {
            reals[i] /= scale_factor;
            imags[i] /= scale_factor;
        }

        println!("IFFT Output Reals: {:?}", reals);
        println!("IFFT Output Imags: {:?}\n", imags);

        let mut signal: Vec<_> = (1..=big_n).map(|i| Complex::new(i as f64, 0.0)).collect();

        let mut planner = FftPlanner::new();
        let fft = planner.plan_fft_forward(signal.len());
        fft.process(&mut signal);
        println!(
            "RustFFT FFT Reals: {:?}",
            signal.iter().map(|z| z.re).collect::<Vec<_>>()
        );
        println!(
            "RustFFT FFT Imags: {:?}",
            signal.iter().map(|z| z.im).collect::<Vec<_>>()
        );

        let fft = planner.plan_fft_inverse(signal.len());
        fft.process(&mut signal);
        for z in signal.iter_mut() {
            z.re /= scale_factor;
            z.im /= scale_factor;
        }
        println!(
            "RustFFT IFFT Reals: {:?}",
            signal.iter().map(|z| z.re).collect::<Vec<_>>()
        );
        println!(
            "RustFFT IFFT Imags: {:?}",
            signal.iter().map(|z| z.im).collect::<Vec<_>>()
        );
    }

Versions PhastFT: v0.2.0 RustFFT: v6.2.0