import matplotlib.pyplot as plt
import numpy as np
from typing import *

from qiskit import QuantumCircuit
from qiskit.quantum_info import Statevector, Operator
from qiskit_aer import AerSimulator
sim = AerSimulator()

from util import zero, one, qft

Shor’s Part 3/5: Modular Exponentiation#

References

  1. Introduction to Classical and Quantum Computing, Chapter 7.9.3

  2. Qiskit notebook on Shor’s Algorithm

  3. Introduction to Quantum Information Science: Lecture 19 by Scott Aaronson

  4. Quantum Computation and Quantum Information: Chapter 5, Nielsen and Chuang

Quantum Oracle for Modular Multiplication#

  1. Suppose we have a unitary

Uf|y=|ay(modN)

so that Uf encodes modular multiplication by a. This oracle will be useful for quantum order finding.

  1. Question: why don’t we need an xor oracle?

Example#

Consider

f(x)=ax(modN)

where a=7 and N=15.

a = 7; N = 15
xs = [x for x in range(N)]
ys = [(a ** x) % N for x in xs]
plt.plot(xs, ys, marker='x'); plt.title("a=7 and N=15")
Text(0.5, 1.0, 'a=7 and N=15')
../_images/8bbe650d2d30cb29fda1d7c34e1eb3e206229814001bdf049f553bdbcc8aa191.png

Quantum Oracle#

  1. Let y=y4y3y2y1 be the binary expansion of y in 4 bits where y4 is the most significant bit.

  2. Our goal is to construct the unitary

Uf|y=|0111y(mod1111)

where 7=0111 and 15=1111.

def amod15(a, power):
    # Adapted From: https://qiskit.org/textbook/ch-algorithms/shor.html
    if a not in [2,4,7,8,11,13]:
        raise ValueError("'a' must be 2,4,7,8,11 or 13")
    U = QuantumCircuit(4)        
    for iteration in range(power):
        if a in [2,13]:
            U.swap(2,3); U.swap(1,2); U.swap(0,1)
        if a in [7,8]:
            U.swap(0,1); U.swap(1,2); U.swap(2,3)
        if a in [4, 11]:
            U.swap(1,3); U.swap(0,2)
        if a in [7,11,13]:
            for q in range(4):
                U.x(q)
    U = U.to_gate(); U.name = "%i^%i mod 15" % (a, power); 
    return U
qc_amod15 = QuantumCircuit(4)
qc_amod15.append(amod15(7, 2**0), [0, 1, 2, 3])
U_f = Operator(qc_amod15)
qc_amod15.draw(output="mpl", style="iqp")
../_images/952fba8cab3ff4d926eed06007bb8531ad9ae2c9560e9404639e7d8e8a5a6c79.png
print((1 * 7) % 15) # Check that we get binary 7 (little endian)
(zero ^ zero ^ zero ^ one).evolve(U_f).draw("latex")
7
|0111
print((3 * 7) % 15) # Check that we get binary 6 (little endian)
(zero ^ zero ^ one ^ one).evolve(U_f).draw("latex")
6
|0110

Iterating Applications of Uf#

  1. Observe that if we iteratively apply Uf to itself

Uf|0001=|a(modN)Uf2|0001=|a2(modN)=Ufs1|0001=|as1(modN)Ufs|0001=|as(modN)=|0001.
  1. Thus we can encode the original function as

f(x)=Ufx|0001

so that we can sequence Uf x times to simulate f(x).

# Starting vector  (reminder: little endian)
(zero ^ zero ^ zero ^ one).draw("latex")
|0001
# U_f |0001>   (reminder: little endian)
(zero ^ zero ^ zero ^ one).evolve(U_f).draw("latex")
|0111
U_f2 = U_f.compose(U_f)
# U_f^2 |0001>   (reminder: little endian)
(zero ^ zero ^ zero ^ one).evolve(U_f2).draw("latex")
|0100
U_f3 = U_f.compose(U_f2)
# U_f^3 |0001>   (reminder: little endian)
(zero ^ zero ^ zero ^ one).evolve(U_f3).draw("latex")
|1101
U_f4 = U_f.compose(U_f3)
# U_f^4 |0011> = |0001>   (reminder: little endian)
# Recall that the order was 4
(zero ^ zero ^ zero ^ one).evolve(U_f4).draw("latex")
|0001

Aside: Eigenvectors and Eigenvalues#

We say that x is an eigenvector and λ is an eigenvalue of a matrix A if

Ax=λx.

Eigenvectors and eigenvalues of Uf#

  1. |0001 is an eigenvector of Uf with eigenvalue 1.

  2. Are there other eigenvectors/eigenvalues of Uf?

  3. Let

|u=1sk=0s1Ufk|1

be an equal superposition of the powers of Uf. Then

(definition)Uf|u=Uf(1sk=0s1Ufk|1)(linearity)=(1sk=0s1Ufk+1|1)(Uf0=Ufs)=(1sk=0s1Ufk|1)(definition)=|u
(1/2.*((zero ^ zero ^ zero ^ one).evolve(U_f) +
(zero ^ zero ^ zero ^ one).evolve(U_f2) +
(zero ^ zero ^ zero ^ one).evolve(U_f3) +
(zero ^ zero ^ zero ^ one))).draw("latex")
12|0001+12|0100+12|0111+12|1101
# Checking eigenvector
(1/2.*((zero ^ zero ^ zero ^ one).evolve(U_f) +
(zero ^ zero ^ zero ^ one).evolve(U_f2) +
(zero ^ zero ^ zero ^ one).evolve(U_f3) +
(zero ^ zero ^ zero ^ one))).evolve(U_f).draw("latex")
12|0001+12|0100+12|0111+12|1101

More useful eigenvectors and eigenvalues of Uf?#

  1. It is a fact that unitary matrices have unit norm eigenvalues.

  2. Consequently, we might wonder if it possible to have phases e2πik as eigenvalues.

  3. Where might we get such eigenvalues?

  4. Let

|u=1sk=0s1e2πiksUfk|1

be an equal superposition of the powers of Uf with phase shifts e2πiks .

  1. Then

Uf|u=e2πis|u

which means we leak out information about the period s in the global phase.

  1. This is the change-of-basis given by the QFT.

qft4 = QuantumCircuit(4)
qft4 = qft(qft4, 4)
qft4.draw(output="mpl", style="iqp")
../_images/14dffb2b84b2c51d5a9cbafe75a2dccd8d5454ddf78433c23363a50c19d12646.png
(zero ^ zero ^ zero ^ one).evolve(U_f).evolve(qft4).draw("latex")
14|0000+(0.2309698831+0.0956708581i)|0001+(0.17677669530.1767766953i)|0010+(0.0956708581+0.2309698831i)|0011i4|0100+(0.0956708581+0.2309698831i)|0101++(0.09567085810.2309698831i)|1011+i4|1100+(0.09567085810.2309698831i)|1101+(0.1767766953+0.1767766953i)|1110+(0.23096988310.0956708581i)|1111

Final fact: Sum of |u#

One final observation is that

|1=1s=0s1|u.

Summary#

We reviewed the quantum modular exponentiation algorithm.