Reading Week 8 - More on Quantum Algorithms

Ch 6.8: RSA Security

Encryption is the science of converting information into a form that is readable only by the recipient. In modern times, the Data Encryption Standard (DES), meaning the same key is used to encrypt the plaintext and decrypt the cyphertext. DES was replaced later by Advanced Encryption Standard (AES) in 2001. It uses a substitution-permutation network which is a series of linked mathematical operations. Decryption occurs by reversing the order of the network.

We'll focus on the RSA algorithm. It is:

antisymmetric (has a private and public key)
popular, public-key based
crucial component of the public key infrastructure on which most internet transactions are done

Here the public key is published for all to see, while the private key is only known between the sender and recipient. To use it:

Alice encrypts her message, the plaintext $m$ using the public key
She sends the message, the encoded cyphertext $c$ to Bob.
Bob uses the private key to decipher the message.

It relies on the fact that large-enough prime numbers are difficult to factor, even for a computer. The General Number Field Sieve (GNFS) algorithm is the fastesst well-known classical integer factorization algorithm. It would take the following time to factor these sized-numbers:

768 bit number (232 digits) - 1500 years
1024 bit number - 1.5 Million years
DigiCert 2048 bit SSL certificate 617 digit number - 6.4 Quadrillion years

6.8.1: Modular Exponentiation

Assume a large composite integer number $N$ with two distinct primes $p, q$ where $N = p \cdot q$ . According to number theory, there exists a periodic function $F (x)$ with period $r$ such that:

F (x) = a^{x} \mod N

where $a$ is an integer, $a \neq N$ , coprime with $p$ or $q$ and coprime with $N$ . The period $r = lcm (p - 1, q - 1)$ and is called the Carmichael Totient Function.

Example

For instance, let $N = 21, p = 7, q = 3$ .
Then $N = p \cdot q = 7 \cdot 3 = 21$ .
Then $r = lcm (7 - 1, 3 - 1) = lcm (6, 2) = 6$
Assume $a = 5$ . Generating the table:
Pasted image 20240522211534.png
Notice in the table it's periodic every 6-th entry, which agrees with $r = 6$ .

Summary

To do the RSA process:

Select two large primes $p, q$ . For example, $p = 23, q = 31$ .
Calculate $N$ , the composite number, where $N = p \times q$ . Here $N = 23 \times 31 = 713$
Calculate the $r$ period via $r = lcm (p - 1, q - 1)$ . Here $r = lcm (22, 30) = 330$ .
Find $e$ public exponent where $e$ is coprime to $r$ ( $gcd (e, r) = 1$ )while $e < r$ . Here assume $e = 19$ since $gcd (19, 339) = 1$
Calculate $d$ , the private exponent using modular inverse, namely $d \times e \mod r = 1 \Leftrightarrow d = e^{- 1} \mod r$ . Here $d = 19^{- 1} \mod 330 = 139$ .

Now we have the public exponent $e = 19$ and private one $d = 139$ . With this, we can try to perform the encryption of the message "A"

Encryption: public key $(e, N) = (19, 713)$ and $m = 65$ which is the ASCII for "A" character.
- $c = m^{e} \mod N = 65^{19} \mod 713 = 198$
Decryption: private key $(d, N) = (139, 713)$ and have $c = 198$
- $m = c^{d} \mod N = 198^{139} \mod 713 = 65$

6.8.2: Period Estimation

If $p, q$ are known, it is easy to find $r$ . But only if $N$ is known, then it's hard, ensuring the security of RSA. Once $r$ is known, it's easy to calculate $d$ using $d = e^{- 1} \mod r$ and RSA fails to be secure. Finding the value of $r$ is called order finding or period estimation. Shor's algorithm is the quantum version of period estimation.

Ch 7.2: Quantum Fourier Transform

7.2.1: Classical Fourier Transformation

The Fourier Transformation (see Lecture Notes#page=113) is a reversible linear transformation. It converts signal $g (t)$ in the time domain into a signal $G (f)$ in the frequency domain. Hence, it has a lot of applications to image and signal processing. The equation (in case you don't use the link above) is:

X (ω) = F (x (t)) = \int_{- \infty}^{\infty} x (t) e^{- i ω t} d t

where $ω = 2 π f$ , and $f$ is the frequency and $t$ is the time. This equation is the forward transform, but we can do the inverse transform as well:

x (t) = F^{- 1} (X (ω)) = \int_{- \infty}^{\infty} X (ω) e^{i ω t} d ω

Note

We get, as an output of the Fourier transform, a period function of $x (t)$ , namely where $x (t) = x (t + T_{0})$ where $T_{0}$ is the period.

We can see that since we have $e^{j ω t} = \cos (ω t) + j \sin (ω t)$ we can conclude that complex sinusoidal functions are the basis functions of the signal. Think of them as unit vectors in a complex plane at a rate of $ω$ radians per second:

Pasted image 20240522213510.png

7.2.2: Discrete Fourier Transform (DFT)

The Discrete Fourier Transform is equivalent to the CFT (Continuous Fourier Transform); however it applies to $N$ samples taken at $T$ time intervals.

Let's assume $x (t)$ is the continuous signal under consideration. We take $N$ samples of this signal x[0], ..., x[N-1]. Each sample x[k] is taken after a time $T$ from the previous sample x[k-1]. the signal is periodic. Thus, we can assume the data in the range x[0], ..., x[N-1] are the same as the data x[N], ..., x[2N-1]. Hence, we can turn:

X (ω) = \int_{0}^{(N - 1) T} x (t) e^{- i ω t} d t = x [0] e^{- i 0} + x [1] e^{- i ω T} + \dots + x [k] e^{- i ω k T} + \dots + x [N - 1] e^{- i ω (N - 1) T}

Which becomes the summation:

X (ω) = \sum_{k = 0}^{N - 1} x [k] e^{- i ω k T}

The sample data is periodic (as we said above), so we should calculate the DFT for the fundamental frequency $\frac{2 π}{N T} rad/sec$ and its harmonics:

ω = 0, 1 \times \frac{2 π}{N T}, 2 \times \frac{2 π}{N T}, . . ., n \times \frac{2 π}{N T}, 1 \times \frac{2 π}{N T}, (N - 1) \times \frac{2 π}{N T}

Applying this equation to the summation above, we get (after normalization):

X [n] = \frac{1}{\sqrt{N}} \sum_{k = 0}^{N - 1} x [k] e^{- i \frac{2 π}{N} n \cdot k}

where $0 \leq n \leq N - 1$ . Here $X [n]$ is called the DFT of the sampled data response $x [k]$ . The sampling is done at uniform intervals of time $t = 0, . . ., N - 1$ . The output $X [n]$ is a vector of complex numbers that encodes the amplitude and phase of a sinusoidal wave of frequency $\frac{n}{N}$ :

[\begin{matrix} X [0] \\ X [1] \\ X [2] \\ ⋮ \\ X [N - 1] \end{matrix}] = \frac{1}{\sqrt{N}} [\begin{matrix} ω^{0} & ω^{0} & ω^{0} & \dots & ω^{0} \\ ω^{0} & ω^{1} & ω^{2} & \dots & ω^{N - 1} \\ ω^{0} & ⋮ & ⋮ & ⋱ & ⋮ \\ ω^{0} & ω^{(N - 1)} & ω^{2 (N - 1)} & \dots & ω^{(N - 1) (N - 1)} \end{matrix}] [\begin{matrix} x [0] \\ x [1] \\ x [2] \\ ⋮ \\ x [N - 1] \end{matrix}]

where $ω = e^{- i \frac{2 π}{N}}$ . The term $\frac{1}{\sqrt{N}}$ is a normalization constant and the DFT is unitary (just like a gate).

The Inverse DFT is given by:

x [k] = \frac{1}{N} \sum_{k = 0}^{N - 1} X [n] e^{i \frac{2 π}{N} n k}

where $0 \leq n \leq N - 1$ .

7.2.3: Quantum Fourier Transformation

Recall from the previous sections that the DFT takes in a vector $x [N] \in C^{N}$ and maps it to another vector $X [N] \in C^{N}$ . Similarly, the Quantum Fourier Transform acts on a quantum state $| x ⟩ = \sum_{i = 0}^{N - 1} x_{i} | i ⟩$ and maps it into a state $\sum_{i = 0}^{N - 1} X_{i} | i ⟩$ . The transformation is given by:

X_{k} = \frac{1}{\sqrt{N}} \sum_{n = 0}^{N - 1} x_{n} ω_{N}^{k n}

where $k = 0, . . ., N - 1$ and is called the primitive $N$ -th root of unity. Assuming $| x ⟩$ is a basis state, the QFT can be written as the mapping:

| x ⟩ = \frac{1}{\sqrt{N}} \sum_{y = 0}^{N - 1} ω_{N}^{x y} | y ⟩

In matrix form it's exactly as what we said the DFT matrix was:

U_{Q F T} = \frac{1}{\sqrt{N}} \sum_{y = 0}^{N - 1} ω_{N}^{x y} | y ⟩, U_{Q F T} = \frac{1}{\sqrt{N}} [\begin{matrix} ω^{0} & ω^{0} & ω^{0} & \dots & ω^{0} \\ ω^{0} & ω^{1} & ω^{2} & \dots & ω^{N - 1} \\ ω^{0} & ⋮ & ⋮ & ⋱ & ⋮ \\ ω^{0} & ω^{(N - 1)} & ω^{2 (N - 1)} & \dots & ω^{(N - 1) (N - 1)} \end{matrix}]

As an example, assume a qubit in state $| ψ ⟩ = a | 0 ⟩ + b | 1 ⟩$ . For this qubit $x_{0} = a$ and $x_{1} = b$ . We can calculate the QFT as:

\begin{aligned} X_{0} & = \frac{1}{\sqrt{2}} (a \cdot e^{2 π i \frac{0 \times 0}{2}} + a \cdot e^{2 π i \frac{0 \times 1}{2}}) \\ = \frac{1}{\sqrt{2}} (a \cdot e^{0} + b \cdot e^{0}) \\ = \frac{1}{\sqrt{2}} (a + b) \end{aligned}

Similarly, we'd get:

X_{1} = \dots = \frac{1}{\sqrt{2}} (a - b)

THus then:

U_{Q F T} | ψ ⟩ = \frac{1}{\sqrt{2}} (a + b) | 0 ⟩ + \frac{1}{\sqrt{2}} (a - b) | 1 ⟩

Which is equivalent to doing an $H$ gate!!!!

Recall that the $H$ gate transforms the $Z$ -basis (the computational basis $| 0 ⟩, | 1 ⟩$ ) into the $X$ -basis (the polar basis $| + ⟩, | - ⟩$ ). Similarly, we can define the QFT as the transformation from the computational basis to the Fourier basis. With this assumption, we can write our $U_{Q F T} | ψ ⟩$ for the 2-qubit case as:

U_{Q F T} | ψ ⟩ = \tilde{a} | 0 ⟩ + \tilde{b} | 1 ⟩

Let's calculate the QFT for large states $N = 2^{n}$ . Assume $| x ⟩ = | x_{1} . . . x_{n} ⟩$ . In this convention, $x_{1}$ is the most significant bit:

\begin{aligned} Q F T_{N} | x ⟩ & = \frac{1}{\sqrt{N}} \sum_{y = 0}^{N - 1} ω_{N}^{x y} | y ⟩ \\ = \frac{1}{\sqrt{2^{n}}} \sum_{y = 0}^{N - 1} e^{2 π i x y / 2^{n}} | y ⟩ \end{aligned}

We can expand $y$ in the computational basis ${0, 1}$ , the range of the sum becomes $y_{1}, y_{2}, . . ., y_{n} \in {0, 1}^{n}$ :

y = y_{1} 2^{n - 1} + y_{2} 2^{n - 2} + \dots + y_{n} 2^{0} = \sum_{i = 1}^{n} y_{i} 2^{n - i}

Again, here each $y_{i}$ is some toggle for a bit, either 0 or 1. Using this we get:

Q F T_{N} | x ⟩ = \frac{1}{\sqrt{2^{n}}} \sum_{y = 0}^{2^{n} - 1} \exp (2 π i \frac{x \sum_{k = 1}^{n} y_{k} 2^{n - k}}{2^{n}}) | y_{1} y_{2} . . . y_{n} ⟩

The exponential of a sum becomes a product of exponentials. With some cancellation:

Q F T_{N} | x ⟩ = \frac{1}{\sqrt{2^{n}}} \sum_{y = 0}^{2^{n} - 1} \prod_{k = 1}^{n} \exp (2 π i \frac{x y_{k}}{2^{k}}) | y_{1} y_{2} . . . y_{n} ⟩

When doing a product of vectors with constants, the constants combine and the vectors tensor up:

\prod_{i = 1}^{n} α_{i} | x_{1} . . . x_{n} ⟩ = α_{1} α_{2} \dots α_{n} | x_{1} ⟩ | x_{2} ⟩ \dots | x_{n} ⟩ = ⨂_{i = 1}^{n} α_{i} | x_{i} ⟩

Thus:

Q F T_{N} | x ⟩ = \frac{1}{\sqrt{2^{n}}} \sum_{y = 0}^{2^{n} - 1} ⨂_{k = 1}^{n} \exp (2 π i \frac{x y_{k}}{2^{k}}) | y_{k} ⟩

Then by rearranging the sum and product operators:

Q F T_{N} | x ⟩ = \frac{1}{\sqrt{2^{n}}} ⨂_{k = 1}^{n} (\sum_{y_{k} \in {0, 1}} \exp (2 π i \frac{x y_{k}}{2^{k}}) | y_{k} ⟩)

Expanding the sum part:

\begin{aligned} Q F T_{N} | x ⟩ & = \frac{1}{\sqrt{2^{n}}} ⨂_{k = 1}^{n} (\exp (2 π x \cdot \frac{0}{2^{k}}) | 0 ⟩ + \exp (2 π i x \cdot \frac{1}{2^{k}}) | 1 ⟩) \\ = \frac{1}{\sqrt{2^{n}}} ⨂_{k = 1}^{n} (| 0 ⟩ + \exp (2 π i \frac{x}{2^{k}}) | 1 ⟩) \end{aligned}

Note that we can expand the exponent part of the previous equation by expanding $x$ as a binary number:

\exp (2 π i \frac{x}{2^{k}}) = \exp (2 π i \sum_{l = 1}^{n} \frac{2^{n - l} x_{l}}{2^{k}}) = \exp (2 π i \sum_{l = 1}^{n} 2^{n - k - l} x_{l})

Like what we did in expanding the sum above, we can expand this into it's 0,1 parts:

= \exp (2 π i \sum_{l = 1}^{n} 2^{n - k - l} x_{l} + 2 π i \sum_{l = n - k + 1}^{n} 2^{n - k - l} x_{l}) = \exp (2 π i \sum_{l = 1}^{n} 2^{n - k - l} x_{l}) \times \exp (2 π i \sum_{l = n - k + 1}^{n} 2^{n - k - l} x_{l})

For the first exponent, notice there loop definition restricts $l \leq n - k$ so then $n - k - l \geq 0$ . So then the sum is just some integer, and because $\exp (2 π i m) = 1$ for any $m \in Z^{+}$ then the whole expression just becomes 1:

= \exp (2 π i \sum_{l = n - k + 1}^{n} 2^{n - k - l} x_{l})

Where the internal sum expands as:

\sum_{l = n - k + 1}^{n} 2^{n - k - l} x_{l} = \frac{x_{n - k + 1}}{2^{1}} + \frac{x_{n - k + 2}}{2^{2}} + \dots + \frac{x_{n}}{2^{n}}

Recall that a binary fraction is expanded as:

0. x_{1} x_{2} x_{3} \dots x_{n} = \sum_{i = 1}^{n} \frac{x_{i}}{2^{i}}

Thus we can put this into our final equation:

Q F T_{N} | x ⟩ = \frac{1}{\sqrt{2^{n}}} ⨂_{k = 1}^{n} (| 0 ⟩ + \exp (2 π i \cdot 0. x_{n - k + 1} \dots x_{n}) | 1 ⟩)

The inverse QFT can be implemented and derived by using $U_{Q F T}^{†}$ .

Ch 7.11: Quantum Phase Estimation (QPE)

Given an eigenvector $| ψ ⟩$ of unitary operator $U$ , such that $U | ψ ⟩ = e^{2 π i θ} | ψ ⟩$ , where $0 \leq θ \leq 1$ , the goal of this algorithm is to find the eigenvalue $e^{2 π i θ}$ for $| ψ ⟩$ . This is equivalent to finding the phase $θ$ to some accuracy. This phase estimation algorithm is used in Shor's Algorithm (see the next week for more on that).

Pasted image 20240523134702.png

The quantum phase estimation procedure is similar to the method used to perform arithmetic in the Fourier basis. Each of the successive blocks of the controlled $U$ gates kickbacks a phase $\frac{θ}{2^{n}}$ into the upper block, called a counting register. It sums up to a value proportional to the phase factor $e^{2 π i θ}$ . An inverse QFT is applied to the counting register to convert this to the computational basis again, which we measure and post-process.

Summary

The steps are:

Start $| ψ ⟩ \to {| 0 ⟩}^{\otimes n} | ψ ⟩$ .
Apply $H^{\otimes n}$ to the ${| 0 ⟩}^{\otimes n}$ part to put into superpositions. Now $\to \frac{1}{\sqrt{2^{n}}} (| 0 ⟩ + | 1 ⟩)^{\otimes n} | ψ ⟩$ .
Apply a cascade of $U^{2^{i} - 1}$ gates for each $q_{i}$ qubit.
Apply inverse QFT.

After step 2 we have:

\frac{1}{\sqrt{2^{n}}} (| 0 ⟩ + | 1 ⟩)^{\otimes n} | ψ ⟩

Since $U$ is unitary, then:

U^{2^{n}} | ψ ⟩ = U^{2^{n - 1}} e^{2 π i θ} | ψ ⟩ = U^{2^{n - 2}} e^{2 π i θ 2^{1}} e^{2 π i θ 2^{0}} | ψ ⟩

Applying this logic into our previous equation:

\to \frac{1}{\sqrt{2^{n}}} (| 0 ⟩ + | 1 ⟩) e^{2 π i θ 2^{n - 1}} | ψ ⟩ \otimes \dots \otimes (| 0 ⟩ + | 1 ⟩) e^{2 π i θ 2^{0}} | ψ ⟩

Rewriting:

\to \frac{1}{\sqrt{2^{n}}} [(| 0 ⟩ + e^{2 π i θ 2^{n - 1}} | 1 ⟩) | ψ ⟩ \otimes \dots \otimes (| 0 ⟩ + e^{2 π i θ 2^{0}} | 1 ⟩) | ψ ⟩]

\to \frac{1}{\sqrt{2^{n}}} \sum_{k = 0}^{2^{n} - 1} e^{2 π i θ k} | k ⟩ \otimes | ψ ⟩

If we do the $Q F T_{N}^{- 1}$ as the algorithm describes:

\begin{aligned} I Q F T_{N} | prev ⟩ & \to \frac{1}{\sqrt{2^{n}}} \sum_{k = 0}^{2^{n} - 1} e^{2 π i θ k} \frac{1}{\sqrt{2^{n}}} \sum_{y = 0}^{2^{n} - 1} e^{- 2 π \frac{i k y}{2^{n}}} | y ⟩ \otimes | ψ ⟩ \\ \to \frac{1}{2^{n}} \sum_{y = 0}^{2^{n} - 1} \sum_{k = 0}^{2^{n} - 1} e^{2 π i θ k} e^{\frac{- 2 π i k y}{2^{n}}} | y ⟩ \otimes | ψ ⟩ \\ \to \frac{1}{2^{n}} \sum_{y = 0}^{2^{n} - 1} \sum_{k = 0}^{2^{n} - 1} \exp (- 2 π \frac{i k}{2^{n}} (y - 2^{n} θ)) | y ⟩ \otimes | ψ ⟩ \end{aligned}

When $y = 2^{n} θ$ then the $\exp$ makes the whole thing 1, so then:

| 2^{n} θ ⟩ \otimes | ψ ⟩

Thus if we can measure this vector's probability, then $θ$ can be estimated as we know $2^{n}$ . The estimation error $Δ θ$ is namely $\frac{1}{2^{n}}$ since we only get a per-qubit estimation (we can only measure discrete state values). In more mathematical terms, we have the amplitudes $α_{i}$ where:

α_{i} \equiv \frac{1}{2^{n}} \sum_{k = 0}^{2^{n} - 1} \exp (2 π i (θ -))

(continue on pg. 348)

Ch 7.9: Shor's Algorithm

I would read Reading Week 8 - More on Quantum Algorithms#Ch 6.8 RSA Security if you haven't already. In there we learned that the period $r$ of a function $f (x) \equiv a^{x} \mod N$ could be used to factorize $N$ into two prime number $p, q$ such that $N = p \times q$ . Recall that $a$ was a number that was coprime to $N$ and that:

a^{r} \mod N = 1

Peter Shor proposed a novel quantum algorithm to factorize large integers in polynomial time. The largest currently factored integer using this method is 261980999226229. He solved the problem by shifting the problem domain to finding the period of the modular exponentiation function. The problem is as follows:

Problem

Given $N \in Z^{+}$ and two prime numbers $p, q$ where $p \times q = N$ , let $a \in Z_{N}$ . Let a moldular exponentiation function $f_{a} : Z_{+} \to Z_{N}$ where $f_{a} (x) = a^{x} \mod N$ , where $f_{a}$ satisfies the condition where $f (x) = f (x_{0}) \Leftrightarrow= x_{0} + k r$ where $r$ is the period of $f_{a}$ and $k$ is some integer.

Shor's algorithm is trying to determine the period $r$ by querying $f_{a}$ . The block diagram is:

Pasted image 20240530122853.png

An interesting feature of Shor's is the use of inverse QFT to find $r$ . The period $r$ can be used to find the factors for $N$ once known.

The algorithm uses two banks of qubit registers. The upper bank is called the source register while the lower one is the target register. The oracle (shown in green) performs modular exponentiation on the source register and stores the values in the target register. The process is done in one step using quantum parallelism.

The number of qubits $s$ in the source and target must be in the following range:

N^{2} \leq 2^{s} \leq 2 N^{2}

Using this range then $\frac{2^{s}}{r} > N$ and it guaruntees that when estimating $r$ that $N$ terms are contributing to the probability amplitude, even if $r \approx \frac{N}{2}$ . Using a large $s$ ensures a dense solution space so $r$ can be estimated accurately as possible.

Before starting the algorithm, we need to preprocess some thing. For one, we should check that $N$ is a prime number, using a classical computer (which can be done efficiently). The next thing is to identify $a$ and confirm that $N, a$ are coprime. Lastly, we should determine $s$ (by choosing any suitable value).

Running the Algorithm

We start with both registers in $| 0 ⟩$ for everything:

| ψ_{0} ⟩ = {| 0 ⟩}^{\otimes s} {| 0 ⟩}^{\otimes s}

We apply $N^{\otimes s}$ to the source to have an equal superposition:

| ψ_{1} ⟩ = H^{\otimes s} {| 0 ⟩}^{\otimes s} {| 0 ⟩}^{\otimes s} = \frac{1}{\sqrt{2^{s}}} \sum_{x = 0}^{2^{s} - 1} | x ⟩ \otimes {| 0 ⟩}^{\otimes s}

Next we apply modular exponentiation using the oracle. For all states of the source register we calculate $y \oplus f (x)$ and store that in the target $| y ⟩$ .

7.7 Review

We can calculate $a^{x} \mod N$ by:

Calculating $a^{x}$
Dividing by $N$ and giving the remainder

If we have $x$ in its binary form $x_{n - 1} x_{n - 2} . . . x_{0}$ we can expand the modular exponentiation to:

\begin{aligned} a^{x} \mod N & = (a^{2^{0}} \mod N)^{x_{0}} \cdot (a^{2^{1}} \mod N)^{x_{1}} \cdot \dots \cdot (a^{2^{n - 1}} \mod N)^{x_{n - 1}} \mod N \end{aligned}

We can do this idea in Quantum Circuits using the modular exponentiation block:

Pasted image 20240608181909.png

The idea here is that each of the controlled multiplier blocks performs a multiplication of the input register 2, by a constant function $a^{2^{i}} \mod N$ . Each controlled multiplier is controlled by a respective qubit $x_{i}$ in register 1. If we assume $| x ⟩$ is the control cubits and $| y ⟩$ as the target in register 2, then the operation produces:

C U_{a^{2^{i}}} | x ⟩ | y ⟩ = | x ⟩ (a^{2^{i}})^{x} y \mod N

For example, doing $2^{x} \mod 15$ is as follows. Consider the circuit:

Pasted image 20240608182146.png

Here q[0], q[1], q[2] forms $| x ⟩$ and the other qubits form $| y ⟩$ .

The state after $H^{\otimes 3}$ and the $X$ gate is:

| ψ_{1} ⟩ = \frac{1}{\sqrt{2^{3}}} (\sum_{x = 0}^{2^{3} - 1} | x ⟩) | 0001 ⟩

The application of the first block multiplier makes:

\begin{aligned} | ψ_{2} ⟩ & = C U_{0} \frac{1}{\sqrt{8}} (| 000 ⟩ | 0001 ⟩ + | 001 ⟩ | 0001 ⟩ + | 010 ⟩ | 0001 ⟩ \\ + | 011 ⟩ | 0001 ⟩ + | 100 ⟩ | 0001 ⟩ + | 101 ⟩ | 0001 ⟩ \\ + | 110 ⟩ | 0001 ⟩ + | 111 ⟩ | 0001 ⟩ \\ = \frac{1}{\sqrt{8}} (\begin{array}{c} | 000 ⟩ | 0001 ⟩ + \\ \underset{changed!}{\underset{⏟}{| 001 ⟩ | 0010 ⟩}} + \\ | 010 ⟩ | 0001 ⟩ + \\ | 011 ⟩ | 0010 ⟩ + \\ | 100 ⟩ | 0001 ⟩ + \\ | 101 ⟩ | 0010 ⟩ + \\ | 110 ⟩ | 0001 ⟩ + \\ | 111 ⟩ | 0010 ⟩ \end{array}) \end{aligned}

Applying $C U_{1}$ does something similar:

\begin{aligned} | ψ_{3} ⟩ & = \frac{1}{\sqrt{8}} ( \\ | 000 ⟩ | 0001 ⟩ + \\ | 001 ⟩ | 0010 ⟩ + \\ \underset{changed!}{\underset{⏟}{| 010 ⟩ | 0100 ⟩}} + \\ | 011 ⟩ | 1000 ⟩ + \\ | 100 ⟩ | 0001 ⟩ + \\ | 101 ⟩ | 0010 ⟩ + \\ | 110 ⟩ | 0100 ⟩ + \\ | 111 ⟩ | 1000 ⟩ \\ ) \end{aligned}

Applying the last $C U_{2}$ finishes the process:
Pasted image 20240608183526.png

We see if we try to measure $| y ⟩$ that we get only $| 1000 ⟩, | 0010 ⟩, | 0100 ⟩, | 0001 ⟩$ with equal probability. We can compare these values with the following table and see that the left ket $| x ⟩$ makes $| y ⟩ = | 2^{x} \mod 15 ⟩$ as we wanted (keeping in mind the flippded order):

Pasted image 20240608183708.png

Back to the Algorithm

We can used this generalized exponentiation process as follows:

\begin{aligned} | ψ^{'} ⟩ & = \frac{1}{\sqrt{2^{s}}} \sum_{x = 0}^{2^{s} - 1} | x ⟩ f (x) {| 0 ⟩}^{\otimes s} \\ = \frac{1}{\sqrt{2^{s}}} \sum_{x = 0}^{2^{s} - 1} | y \oplus f (x) ⟩ \\ = \frac{1}{\sqrt{2^{s}}} \sum_{x = 0}^{2^{s} - 1} | f (x) ⟩ \end{aligned}

as we desired. After doing this, we measure the target register. This is because we want to see the point where we get a return back to $c = 1$ , if at all. It has the side effect of collapsing our superposition of the source register to the states for which $f (x_{0}) = a^{x_{0}} \mod N$ . Since $f (x)$ is periodic via $r$ , then the possible values are of the form $(x_{0} + k r)$ . Note that $x_{0}$ is the least $x$ such that $f (x)$ is what was measured in the target register. This ensures $k$ is strictly nonnegative. Thus our new state is:

| ψ^{″} ⟩ = \frac{1}{\sqrt{\frac{2^{s}}{r}}} \sum_{k = 0}^{\frac{2^{s}}{r} - 1} | x_{0} + k r ⟩ | f (x_{0}) ⟩

Since the period of $f (x)$ is $r$ , the measurement reduces the range of the source register to $\frac{2^{s}}{r}$ . There are $2^{s}$ states, so the maximum number of periods we can have is $\frac{2^{s}}{r}$ .

Pasted image 20240608184417.png

But this is why we want to apply the inverse QFT! Because we are guarunteed to get a period for some value $k$ , then we can apply inverse AFT to identify this period $r$ , which is what we are after. This gate does:

| x ⟩ = \frac{1}{\sqrt{N}} \sum_{y = 0}^{N - 1} e^{- 2 π i \frac{x y}{N}} | y ⟩

for $x = 0, . . ., N - 1$ . Applying the $i Q F T$ gives:

| ψ^{‴} ⟩ = \frac{1}{\sqrt{2^{s}}} \frac{1}{\sqrt{Q}} \sum_{k = 0}^{Q - 1} \sum_{y = 0}^{2^{s} - 1} (e^{- 2 π i \frac{x_{0} y}{2^{s}}} \cdot e^{- 2 π i \frac{k r y}{2^{s}}}) | y ⟩

Which (saving you the trouble) reduces down to:

= \frac{1}{\sqrt{2^{s}Q}}\sum\limits_{y=0}^{2^{s}-1} e^{\frac{-2\pi i x_{0} y}{2^{s}}} \sum\limits_{k=0}^{Q-1}\left(\exp\left(-2\pi i \frac{ry}{2^{s}}\right)\right)^{k}\ket{y} $${#1fdfasdf} The right side summation is a geometric series which is non-zero iff $\exp\frac{-2\pi i ry}{2^{s}}= 1$. This is possible if $\frac{ry}{2^{s}}$ is an integer, so values of $y$ which are multiples of $\frac{2^{s}}{r}$ can have non-zero amplitudes. In other words, if $\frac{2^{s}}{r}$ is an integer, by constructive wave interference, there is a peak values in the value of $y$, which will thus become measured more often. We can try to calculate this probability: {#7dead2}

P(y) = \frac{1}{2^{s}Q}\left| \sum\limits_{k=0}^{Q-1} \exp\left(-2\pi i \frac{kry}{2^{s}}\right)\right|^

W h e n $ \frac{r y}{2^{s}} $ i s a n i n t e g e r, t h e p h a s e s a d d u p t o 1 :

P(y) = \frac{r}{2^{s}2^{s}}\left[\frac{2^{s}}{r}\frac{2^{s}}{r}\right] = \frac{1}

We can plug this back into ourequationwe had for $\ket{\psi'''}$:

\ket{\psi'''} = \frac{1}{\sqrt{r}} \sum\limits_{z=0}^{r-1} \exp\left(-2\pi i \frac{x_{0}z}{r}\right)\left|z \frac{2^{s}}{r}\right\rangle

T h e v a l u e f o r $ c $, t h e $ a^{x} \mod N $ v a l u e, i s :

c = z \frac{2^{s}}

W e k n o w $ c $ a n d $ 2^{s} $ t h e n w e c a n d e t e r m i n e t h a t :

\frac{c}{2^{s}} = \frac{z}

W e c a n t r y t o o b t a i n $ r $ f r o m h e r e u s i n g c o n t i n u e d f r a c t i o n s . I f $ z, r $ d o n^{'} t s h a r e c o m m o n d i v i s o r s, n a m e l y $ gcd (z, r) = 1 $, t h e n :

\frac{c}{2^{s}} = a_{0} + \frac{1}{a_{1} + \frac{1}{a_{2} + \frac{1}{a_{3} + \dots}}}

T h e a p p r o x i m a t i o n t o t h e r a t i o $ \frac{z}{r} $ c a n b e c a l c u l a t e d u s i n g t h i s c o n t i n u e d f r a c t i o n :

\begin{align}
z_{0} &= a_{0} \
r_{0} &= 1 \
z_{1} &= a_{0}a_{1} + 1 \
r_{1} &= a_{1} \
\vdots \
z_{n} &= a_{n}z_{n-1} + z_{n-2} \
r_{n} &= a_{n}r_{n-1} + r_{n-2}
\end