After performing Quantum Phase Estimation (QPE) and learning its role in factoring numbers (breaking the number 15, for example, down to 3 x 5 — reversing multiplication), I wanted to use it to do my first practical application with quantum computing. The practical use of factoring is in cryptography; with sufficient quantum computing power (which doesn’t — yet — exist), a factoring algorithm could break the strongest RSA encryption.

The problem is, I can’t find any information on how to actually apply QPE. There are many equations and complex explanations out there, but there’s nothing that seems to define, at least at my current level of understanding, how to actually configure a circuit to perform the quantum component of the algorithm.

Consequently, much of this experiment is guesswork. I tried to initialize the quantum bits (qubits) in such a way as to produce the desired outcome, and then I filtered the output to try to find the answer I was looking for.

The two big questions are:

Where do you put the number that you want factored? How do you find the period in the histogram?

The four blocks, from left to right, are the initialization, the controlled unitary, an inverse Quantum Fourier Transform (QFT), and then the measurements. Block 1 contains my first educated guess, that the small register is the number to be factored.

This guess is based on one reference I found somewhere that the n bits of the small register should be large enough to hold N, the number being worked with. Therefore, I set the five qubits of the small register to 10101, which is binary for 21.

The reason for my lack of confidence is that I also saw a reference somewhere to the small register needing to be all ones. But, if the small register is all ones and the large register is all superpositions, where else could N possibly go?

The histogram with 10 classical registers is completely illegible. I did quite a bit of experimentation with 3-qubit small registers, and there are always spikes in the middle of and on the right side of the histogram. For this experiment, the values of the spike don’t seem to be helpful at all. But, since I couldn’t zoom in and determine precise values, I decided to try filtering the results.

This histogram only contains measurements of the first five qubits, the right half of the previous histogram. The spike is near and on 11111, which is binary for 31. That’s not helpful at all.

But, this histogram is the left half of that first histogram. That middle spike looks like it contains 00110, which is binary for 6. That is actually the result that I am looking for! Let’s zoom in a little more, and then I will explain why I am looking for the number 6.

This histogram measures only the last 3 qubits of the large register. I thought I finally had it, but a closer inspection reveals 011 (3), not 110 (6). There are so many digits, that I seem to have misread the previous histogram.

So, this histogram was my final attempt. It measures only the very first 3 qubits. Ignoring the highest spike, which is all ones, the second highest spike is finally 110, the 6 that I am looking for.

Does my search for the number 6 feel like desparation, almost cheating? It feels like it. Until I find definitive answers anywhere, all I can do is test this out on larger numbers and see what happens. For consistency, I will have to read only the least significant qubits of the large register.

On the other hand, I have seen the inverse QFT set up differently. It is possible that the histogram is backward. If this experiment doesn’t work on a larger number, that’s one area to look at. The most significant qubits backward give me the same value as the the least significant qubits forward.

Before getting to the code, allow me to explain that the number 6 is the period p of the number 21. The simplified factoring algorithm is that you guess a random number g and then use a quantum computer to calculate the period p.

Then, gp/2+1 and gp/2-1 are your guesses for finding the factors of N, which is 21 for this experiment. If we guess 2, mostly because I already know it works, we get 2 raised to that elusive number 6, but divided by 2. That’s 26/2 = 23 = 8. We then find the greatest common divisor (gcd) of 21 and 8+1, as well as 21 and 8-1. That’s the bolded formula at the start of this paragraph.

The gcd of 21 and 9 is 3. And, the gcd of 21 and 7 is 7. We are left with 3 and 7, the factors of 21.

The full algorithm is quite a bit more complex. For starters, your random guess g must be smaller than the number to be factored N. And the resultant period p cannot be odd, because you would have a non-integer exponent. Therefore, the full algorithm includes tests to help you find the guess g and period p that ultimately output the factors of N. Because of the inverse QFT (Block 3), this is exponentially faster than any classical factoring algorithm.

I just scrolled up to make sure that I never stated that this was going to be easy!

I’ll wrap up this already-lengthy blog post with the code. Again, I have questions to which I cannot find answers. Therefore, this is my first attempt at figuring it out. I will have to try this again on a larger number N, but I am ultimately limited by the number of qubits available, even on a simulator.

OPENQASM 2.0;

include "qelib1.inc";

qreg q[15];

creg c[10];

// initialize ancilla qubits

h q[0];

h q[1];

h q[2];

h q[3];

h q[4];

h q[5];

h q[6];

h q[7];

h q[8];

h q[9];

// eigenstate of the unitary operator: 21 (10101)

x q[10];

x q[12];

x q[14];

barrier q; // unitary operator

ccx q[9], q[10], q[12];

ccx q[11], q[12], q[13];

cz q[13], q[14];

ccx q[11], q[12], q[13];

ccx q[9], q[10], q[12];

barrier q; // inverse Quantum Fourier Transform (QFT)

cu1(-pi/512) q[9], q[0];

cu1(-pi/256) q[9], q[1];

cu1(-pi/128) q[9], q[2];

cu1(-pi/64) q[9], q[3];

cu1(-pi/32) q[9], q[4];

cu1(-pi/16) q[9], q[5];

cu1(-pi/8) q[9], q[6];

cu1(-pi/4) q[9], q[7];

cu1(-pi/2) q[9], q[8];

cu1(-pi/256) q[8], q[0];

cu1(-pi/128) q[8], q[1];

cu1(-pi/64) q[8], q[2];

cu1(-pi/32) q[8], q[3];

cu1(-pi/16) q[8], q[4];

cu1(-pi/8) q[8], q[5];

cu1(-pi/4) q[8], q[6];

cu1(-pi/2) q[8], q[7];

cu1(-pi/128) q[7], q[0];

cu1(-pi/64) q[7], q[1];

cu1(-pi/32) q[7], q[2];

cu1(-pi/16) q[7], q[3];

cu1(-pi/8) q[7], q[4];

cu1(-pi/4) q[7], q[5];

cu1(-pi/2) q[7], q[6];

cu1(-pi/64) q[6], q[0];

cu1(-pi/32) q[6], q[1];

cu1(-pi/16) q[6], q[2];

cu1(-pi/8) q[6], q[3];

cu1(-pi/4) q[6], q[4];

cu1(-pi/2) q[6], q[5];

cu1(-pi/32) q[5], q[0];

cu1(-pi/16) q[5], q[1];

cu1(-pi/8) q[5], q[2];

cu1(-pi/4) q[5], q[3];

cu1(-pi/2) q[5], q[4];

cu1(-pi/16) q[4], q[0];

cu1(-pi/8) q[4], q[1];

cu1(-pi/4) q[4], q[2];

cu1(-pi/2) q[4], q[3];

cu1(-pi/8) q[3], q[0];

cu1(-pi/4) q[3], q[1];

cu1(-pi/2) q[3], q[2];

cu1(-pi/4) q[2], q[0];

cu1(-pi/2) q[2], q[1];

cu1(-pi/2) q[1], q[0];

barrier q; // output should be 6 (binary 110)

h q[0];

//measure q[0] -> c[0];

h q[1];

//measure q[1] -> c[1];

h q[2];

//measure q[2] -> c[2];

h q[3];

//measure q[3] -> c[3];

h q[4];

//measure q[4] -> c[4];

h q[5];

//measure q[5] -> c[5];

h q[6];

//measure q[6] -> c[6];

h q[7];

measure q[7] -> c[7];

h q[8];

measure q[8] -> c[8];

h q[9];

measure q[9] -> c[9];

Yes, if you want to factor larger numbers, and assuming you have the available qubits to do that, you will want to do this in Python. I use OpenQASM because you can copy-and-paste the Python for everything except the controlled unitary (Block 2 of 4), and there’s no fun in that.