QuState / PhastFT

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

Mention output scaling and normalization in docs #13

Closed astral4 closed 4 months ago

astral4 commented 5 months ago

It seems that phastft and rustfft use the same "normalization factor" for outputs: $\frac{1}{\sqrt{n}}$, where $n$ is the sequence length. The rustfft docs mention normalizing outputs, and it would be nice if this were covered in the readme/docs for phastft too.

smu160 commented 5 months ago

Hi @astral4,

Thank you for bringing this up. Can you please provide a code snippet that I can run to make sure I understand what you are referring to exactly? Going off of the link you provided, it says "RustFFT does not normalize outputs.". I know that as of right now, the inverse FFT in our implementation does not scale all the outputs by $\frac{1}{N}$.

Thank you!

Best, Saveliy

astral4 commented 5 months ago

Here's a demonstration:

use phastft::planner::Direction;
use rustfft::num_complex::Complex;
use std::iter::zip;

fn main() {
    // example data
    let mut real_data = [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8];
    let mut imag_data = [-0.1, -0.2, -0.3, -0.4, -0.5, -0.6, -0.7, -0.8];

   // same data, but in the format that rustfft wants
    let mut complex_data = [
        Complex::new(0.1, -0.1),
        Complex::new(0.2, -0.2),
        Complex::new(0.3, -0.3),
        Complex::new(0.4, -0.4),
        Complex::new(0.5, -0.5),
        Complex::new(0.6, -0.6),
        Complex::new(0.7, -0.7),
        Complex::new(0.8, -0.8),
    ];

    // use phastft to calculate the Fourier transform
    phastft::fft(&mut real_data, &mut imag_data, Direction::Forward);

    // convert result to facilitate side-by-side comparison with rustfft
    let phastft_results: Vec<_> = zip(real_data, imag_data)
        .map(|(real, imag)| Complex::new(real, imag))
        .collect();

    // use rustfft to calculate the Fourier transform
    rustfft::FftPlanner::new()
        .plan_fft_forward(8)
        .process(&mut complex_data);

    // log results
    println!("{phastft_results:#?}");
    println!("{complex_data:#?}");
}

I got the same answer with phastft and rustfft, so I assume that both of them have the same scaling factor. The rustfft docs that I linked states that this scaling factor is $\frac{1}{\sqrt{n}}$.

Also, I don't think either phastft or rustfft scales the output by $\frac{1}{n}$ when computing the inverse Fourier transform—the rustfft docs imply that the scaling factor applies twice when you first compute the Fourier transform and then the inverse, so the final scaling factor is $\frac{1}{\sqrt{n}}\cdot\frac{1}{\sqrt{n}}=\boxed{\frac{1}{n}}$. I hope this clears up any confusion!

smu160 commented 4 months ago

RustFFT does not normalize outputs. Callers must manually normalize the results by scaling each element by 1/len().sqrt(). Multiple normalization steps can be merged into one via pairwise multiplication, so when doing a forward FFT followed by an inverse callers can normalize once by scaling each element by 1/len()

Hi @astral4,

The rustfft docs seem to imply there is no scaling at all. Everything must be done manually. I assume you got the same results using phastft and rustfft because phastft does not do any scaling -- all scaling is left up to the user.

Let me know if this sounds right to you. Thanks!

Best, Saveliy

astral4 commented 4 months ago

That makes sense. Thank you for the docs update!