Construct and implement generalized FFFT

[kevinsung]

I don't see a reason not to implement an FFFT that works for any number of modes (not just powers of 2). Has anyone considered this? I think the T-count would also be much better than using Givens rotations.

[Strilanc]

Do you have a construction that works for numbers that aren't a power of 2?

In the classical FFT world, generalizing to values that aren't a power of 2 involves extending to a power of 2 and putting in chirps. https://math.stackexchange.com/questions/77118/non-power-of-2-ffts I tried using the same techniques for a modular QFT once, but ran into the roadblock that the chirps aren't unitary. I just assumed the same would happen with the FFFT.

[kevinsung]

I haven't looked very deeply into it, but I just assumed that we'd be able to adapt a classical FFT algorithm. The QFT is pretty different from the FFFT (with n qubits, the QFT performs the DFT on 2^n coefficients, while the FFFT performs the DFT on n coefficients), so I wouldn't be surprised if we could make it work.

[babbush]

Implementations of the FFT that are not based on radix-2 decimation can get pretty messy. @idk3 and I thought about this at length last summer but we didn't come up with any terribly compelling solutions. If you do, please let us know!

[Strilanc]

@babbush Do you have a short mathematical description of what the non-power-of-2 operation is supposed to do? I could look at it.

[babbush]

The FFT definition has nothing to do with powers of 2 (although most of the algorithms to implement do). For instance, Equation C3 in https://journals.aps.org/prx/pdf/10.1103/PhysRevX.8.011044 is an example. The main definition here is fine too: https://en.wikipedia.org/wiki/Fast_Fourier_transform

I need to go - I can answer more questions tomorrow.

[Strilanc]

I still don't really get creation/annihilation operators and exactly how they map back into unitary operations on qubits. But my guess based on this:

that we want to map the basis state where the k'th qubit is set to a superposition of each possible qubit being set with equal magnitude, but with each case having a phase that's been rotated by an amount proportional to k as you go from case to case. Is that right?

What I don't quite get is how to then extend that to a computational basis state where 2 qubits are set instead of 1.

[idk3]

@Strilanc you can get a radix-3 FFT from the regular radix-2 definition - instead of splitting in two parts (with x(2n) and x(2n+1) getting the phase) you instead divide into three parts, x(3n), x(3n+1), and x(3n+2). You have more phase factors, which are no longer nice phases pi/2^k, and you also need to do a 3-way recombination in the recursive step. Both these parts seem a bit troubling, but we could implement them using Givens rotations (in blocks of 3 at a time). This should allow you to develop an arbitrary-radix fermionic Fourier transform with O(n log n) T cost. (The T cost comes in when you do the recombination and phasing using Givens rotations.)

http://cecas.clemson.edu/~janoski/reu/2008/FFT-book.pdf is a pretty detailed reference that gives the butterfly diagrams etc for radices 2 through 6 and 8.

[kevinsung]

@Strilanc to answer your question about the creation/annihilation operators...

What you said about mapping a basis state with one 1 and the rest 0s is correct, assuming the Jordan-Wigner Transform. Under the JWT, the j-th qubit holds the occupancy of the j-th fermionic mode. That is not true about other transforms; you can see the OpenFermion demo The Jordan-Wigner and Bravyi-Kitaev Transforms for an explanation.

If you start with a computational basis state where 2 qubits are set to 1, then you will end up with a superposition of computational basis states where 2 qubits are set to 1. You can work out the coefficients by multiplying the transformed creation operators together. That is, if qubits i and j are set, then multiply c^\dagger_i and c^\dagger_j to get an expression in terms of the a^\dagger's (you'll get a superposition of terms with two a^\daggers), and then apply those a^\daggers to the vacuum state. Under the JWT, the vacuum state is the all zeros state, and applying a^\dagger_i to it means setting the i-th qubit to 1 (potentially introducing a minus sign). I think the result can be described as a convolution or something like that. Again, the demo might be helpful for understanding this.

[kevinsung]

The Cooley-Tukey algorithm actually works by recursively breaking N down into its factors N = N_1 * N_2 and performing smaller FFT's of size N_1 and N_2: https://en.wikipedia.org/wiki/Cooley%E2%80%93Tukey_FFT_algorithm#General_factorizations. The base case is when N is prime. I'm pretty sure this recursive procedure translates to the FFFT. Classically, there are FFTs for the prime case with running time N log N. I'm not sure if that's possible for the FFFT, but at worst we can use the Givens rotation procedure for this base case.

[kevinsung]

The task here is to perform the transformation in Eq. (11) of https://arxiv.org/pdf/0804.1888.pdf. The right boxed part of Figure 1 shows the structure of the circuit when the number of modes is a power of 2. As you can see, it's possible to perform this case with depth n log n and linear connectivity, where n is the total number of modes. The circuit structure is exactly the same as the classical Cooley-Tukey algorithm. I think the relation with the Cooley-Tukey circuit is a bit clearer in Figure 1 of https://arxiv.org/pdf/1310.7605.pdf, although that figure doesn't show how to do it with linear connectivity. The term "FFFT" comes from https://arxiv.org/pdf/1706.00023.pdf, which also contains some useful information.

As I said above, I believe it's possible to get it to work for numbers other than power of 2. The idea is that the FFFT on n = n1 * n2 modes can be done by performing two smaller FFFT's of sizes n1 and n2. I think the first step is to work out what the circuit for this decomposition looks like. If it works out, then implement the recursive structure of the circuit. The base case should be handled in a general way, i.e., perhaps by defining a "prime number FFFT gate" which we can later define custom decompositions for.

The discussion so far is just on the 1-d FFFT under the Jordan-Wigner Transform. This is clearly the first step. We can worry about the 2-d and 3-d versions later (and possibly other fermionic transforms eventually).

I think the function signature for the 1-d FFFT under the Jordan-Wigner Transform should just be ffft(qubits: Sequence[cirq.QubitId])

Also, it would be good to keep https://github.com/quantumlib/OpenFermion-Cirq/issues/302 in mind while implementing.

[kevinsung]

Here is a short note I wrote on the operation of rotating together adjacent fermionic modes under the JWT. It might be useful for understanding what's going on. givens_rotations.pdf