Implementation of Unitary in Shor's Algorithm in PennyLane

Question

I've been working on implementing Shor's Algorithm in PennyLane, but am struggling to understand how the circuit for 'U' has been constructed according to Qiskit. In the Qiskit textbook, they seek to factor 15, which, for my own sake, is what I decided to do as well (but in PennyLane). I used eight counting qubits, and four for the unitary to act upon. So far, my code looks like this:

n_count = 8

dev = qml.device('default.qubit', wires = n_count + 4, shots = 10)

def c_amod15(a, power):
    for iteration in range(power):
        if a in [2, 13]:
            qml.SWAP(wires = [8, 9])
            qml.SWAP(wires = [9, 10])
            qml.SWAP(wires = [10, 11])
        if a in [7, 8]:
            qml.SWAP(wires = [10, 11])
            qml.SWAP(wires = [9, 10])
            qml.SWAP(wires = [8, 9])
        if a == 11:
            qml.SWAP(wires = [9, 11])
            qml.SWAP(wires = [8, 10])
        if a in [7, 11, 13]:
            for k in range(8, 12):
                qml.PauliX(wires = k)
    
@qml.qnode(dev)
def circuit():
    N = 15
    np.random.seed(1)
    a = np.random.randint(2, 15)
    if math.gcd(a, N) != 1:
        raise ValueError("Non-trivial factor.")
    
    for k in range(0, 8):
        qml.Hadamard(wires = k)

    qml.PauliX(wires = n_count + 3)

    for k in range(0, 8):
        qml.ctrl(c_amod15, control = k)(a, 2**k)

    qml.adjoint(qml.QFT)(wires = [0, 1, 2, 3, 4, 5, 6, 7])

    return qml.sample(wires = range(0, 8))

circuit()
print(circuit())

fig, ax = qml.draw_mpl(circuit)()
fig.show()

The 'c_amod15' function was taken directly from Qiskit's hardcoded implementation. Am I applying this in the correct way? My results from this are:

[[1 1 1 1 1 1 1 1]
[1 1 1 1 1 1 1 1]
[0 0 0 0 0 0 0 0]
[1 1 1 1 1 1 1 1]
[0 0 0 0 0 0 0 0]
[0 0 0 0 0 0 0 1]
[0 0 0 0 0 0 0 1]
[1 1 1 1 1 0 1 1]
[1 1 1 1 1 1 1 1]
[1 1 1 1 1 1 1 1]]

I am still very new to this, and haven't done any of the classical post-processing yet, so any help would be much appreciated.

Answer 1

Yes, it seems to me you apply it in a right way, except some points I listed after. The circuit making $U$ you found is simply a way to check how Shor's algorithm works, and is simply built from the truth table of the operation $a*b\%15$, as shows for example this paper. As a result, it is diffult to generalize to other values of N, but it is another topic.

Maybe I didn't see everything, but here are some points I noticed:

• The if block if a == 11 should also be executed for $a=4$, and thus be if a in [4,11] (following Qiskit implementation of this circuit)
• I don't know what's wrong with my usage of pennylane QFT circuit, but when I implemented Shor's algorithm I got wrong results until I did my own QFT. Here it is:

def QFT(wires):
    for i in range(len(wires)):
        qml.Hadamard(wires=wires[i])
        for j in range(i):
            qml.CRZ(math.pi/2**(i - j),wires=[wires[i],wires[j]])

• I think (but it is quite personal) that for the sake of clarity it is better to separate quantum and classic part. So that would mean generate $a$ in a classic function, then call a qnode that only handles quantum operations, and finally post-process the result in a classic function again (maybe the same)

Keeping $n\_count$ as a variable holding the number of qubits doing the estimation, here is how to read the results:

1. Translate the binary to decimal numbers (e.g. your first result $11111111$->$255$)
2. These numbers are multiples of $\frac{2^{n\_count}}{r}$, were $r$ is the period we look for (i.e. it is the smallest integer such that $a^{r}\%N=1$). Divide them by $2^{n\_count}$, keeping the fraction (it gives $\frac{255}{256}$ in this case)
3. Keep the highest denominator of your 10 results, it is likely to be r (here it would be $256$, so obviously it is wrong, but I assume it is because $a=4$ and so $U$ was wrong or because the inverse QFT didn't do the job as told above)

Note: pennylane has a pre-implemented circuit for phase estimation.