Algoritmo de Shor

Agora vamos voltar nossa atenção para o problema da fatoração de inteiros e ver como ele pode ser resolvido eficientemente em um computador quântico usando estimativa de fase. O algoritmo que vamos obter é o algoritmo de Shor para fatoração de inteiros. Shor não descreveu seu algoritmo especificamente em termos de estimativa de fase, mas é uma maneira natural e intuitiva de explicar como ele funciona.

Começaremos discutindo um problema intermediário conhecido como problema de busca de ordem e veremos como a estimativa de fase oferece uma solução para esse problema. Em seguida, veremos como uma solução eficiente para o problema de busca de ordem nos dá uma solução eficiente para o problema de fatoração de inteiros. (Quando a solução de um problema fornece uma solução para outro problema assim, dizemos que o segundo problema se reduz ao primeiro — portanto, neste caso, estamos reduzindo a fatoração de inteiros à busca de ordem.) Essa segunda parte do algoritmo de Shor não usa computação quântica de forma alguma; ela é completamente clássica. A computação quântica é necessária apenas para resolver a busca de ordem.

O problema de busca de ordem

Alguns conceitos básicos de teoria dos números

Para explicar o problema de busca de ordem e como ele pode ser resolvido usando estimativa de fase, é útil começar com alguns conceitos básicos de teoria dos números e apresentar uma notação prática ao longo do caminho.

Para começar, para qualquer inteiro positivo $N$, defina o conjunto $\mathbb{Z}_N$ da seguinte forma.

$$\mathbb{Z}_N = \{0,1,\ldots,N-1\}$$

Por exemplo, $\mathbb{Z}_1 = \{0\},\;$ $\mathbb{Z}_2 = \{0,1\},\;$ $\mathbb{Z}_3 = \{0,1,2\},\;$ e assim por diante.

Esses são conjuntos de números, mas podemos pensar neles como mais do que conjuntos. Em particular, podemos pensar em operações aritméticas em $\mathbb{Z}_N$ como adição e multiplicação — e se concordarmos em sempre tomar nossas respostas módulo $N$ (isto é, dividir por $N$ e tomar o resto como resultado), sempre permaneceremos dentro desse conjunto ao realizar essas operações. As duas operações específicas de adição e multiplicação, ambas tomadas módulo $N$, transformam $\mathbb{Z}_N$ em um anel, que é um tipo de objeto fundamentalmente importante em álgebra.

Por exemplo, $3$ e $5$ são elementos de $\mathbb{Z}_7,$ e se os multiplicarmos obtemos $3\cdot 5 = 15,$ que deixa um resto de $1$ quando dividido por $7.$ Às vezes expressamos isso da seguinte forma.

$$3 \cdot 5 \equiv 1 \; (\textrm{mod } 7)$$

Mas também podemos simplesmente escrever $3 \cdot 5 = 1$, desde que esteja claro que estamos trabalhando em $\mathbb{Z}_7$, apenas para manter nossa notação o mais simples possível.

Como exemplo, aqui estão as tabelas de adição e multiplicação para $\mathbb{Z}_6$.

$$\begin{array}{c|cccccc} + & 0 & 1 & 2 & 3 & 4 & 5 \\\hline 0 & 0 & 1 & 2 & 3 & 4 & 5 \\ 1 & 1 & 2 & 3 & 4 & 5 & 0 \\ 2 & 2 & 3 & 4 & 5 & 0 & 1 \\ 3 & 3 & 4 & 5 & 0 & 1 & 2 \\ 4 & 4 & 5 & 0 & 1 & 2 & 3 \\ 5 & 5 & 0 & 1 & 2 & 3 & 4 \\ \end{array} \qquad \begin{array}{c|cccccc} \cdot & 0 & 1 & 2 & 3 & 4 & 5 \\\hline 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & 1 & 2 & 3 & 4 & 5 \\ 2 & 0 & 2 & 4 & 0 & 2 & 4 \\ 3 & 0 & 3 & 0 & 3 & 0 & 3 \\ 4 & 0 & 4 & 2 & 0 & 4 & 2 \\ 5 & 0 & 5 & 4 & 3 & 2 & 1 \\ \end{array}$$

Entre os $N$ elementos de $\mathbb{Z}_N$, os elementos $a\in\mathbb{Z}_N$ que satisfazem $\gcd(a,N) = 1$ são especiais. Frequentemente, o conjunto que contém esses elementos é denotado com um asterisco assim.

$$\mathbb{Z}_N^{\ast} = \{a\in \mathbb{Z}_N : \gcd(a,N) = 1\}$$

Se focamos nossa atenção na operação de multiplicação, o conjunto $\mathbb{Z}_N^{\ast}$ forma um grupo — especificamente um grupo abeliano — que é outro tipo importante de objeto em álgebra. É um fato básico sobre esses conjuntos (e grupos finitos em geral) que, se escolhermos qualquer elemento $a\in\mathbb{Z}_N^{\ast}$ e multiplicarmos $a$ por si mesmo repetidamente, sempre chegaremos ao número $1$.

Como primeiro exemplo, vamos tomar $N=6$. Temos que $5\in\mathbb{Z}_6^{\ast}$ porque $\gcd(5,6) = 1$, e se multiplicarmos $5$ por si mesmo obtemos $1$, como a tabela acima confirma.

$$5^2 = 1 \quad \text{(trabalhando dentro de $\mathbb{Z}_6$)}$$

Como segundo exemplo, vamos tomar $N = 21$. Se percorrermos os números de $0$ a $20$, os que têm MDC igual a $1$ com $21$ são os seguintes.

$$\mathbb{Z}_{21}^{\ast} = \{1,2,4,5,8,10,11,13,16,17,19,20\}$$

Para cada um desses elementos, é possível elevar esse número a uma potência inteira positiva para obter $1$. Aqui estão as menores potências para as quais isso funciona:

$$\begin{array}{ccc} 1^{1} = 1 \quad & 8^{2} = 1 \quad & 16^{3} = 1 \\[1mm] 2^{6} = 1 \quad & 10^{6} = 1 \quad & 17^{6} = 1 \\[1mm] 4^{3} = 1 \quad & 11^{6} = 1 \quad & 19^{6} = 1 \\[1mm] 5^{6} = 1 \quad & 13^{2} = 1 \quad & 20^{2} = 1 \end{array}$$

Naturalmente, estamos trabalhando dentro de $\mathbb{Z}_{21}$ em todas essas equações, o que não nos preocupamos em escrever — tomamos isso como implícito para evitar poluir a notação. Continuaremos fazendo isso ao longo do restante da aula.

Enunciado do problema e conexão com a estimativa de fase

Agora podemos enunciar o problema de busca de ordem.

Busca de ordem

Entrada: inteiros positivos $N$ e $a$ satisfazendo $\gcd(N,a) = 1$
Saída: o menor inteiro positivo $r$ tal que $a^r \equiv 1$ $(\textrm{mod } N)$

Alternativamente, em termos da notação que acabamos de introduzir acima, nos é dado $a \in \mathbb{Z}_N^{\ast}$, e estamos procurando o menor inteiro positivo $r$ tal que $a^r = 1$. Esse número $r$ é chamado de ordem de $a$ módulo $N$.

Para conectar o problema de busca de ordem à estimativa de fase, vamos pensar na operação definida em um sistema cujos estados clássicos correspondem a $\mathbb{Z}_N$, onde multiplicamos por um elemento fixo $a\in\mathbb{Z}_N^{\ast}$.

$$M_a \vert x\rangle = \vert ax \rangle \qquad \text{(para cada $x\in\mathbb{Z}_N$)}$$

Para ser claro, estamos fazendo a multiplicação em $\mathbb{Z}_N$, portanto está implícito que estamos tomando o produto módulo $N$ dentro do ket no lado direito da equação.

Por exemplo, se tomarmos $N = 15$ e $a=2$, então a ação de $M_2$ sobre a base padrão $\{\vert 0\rangle,\ldots,\vert 14\rangle\}$ é a seguinte.

$$\begin{array}{ccc} M_{2} \vert 0 \rangle = \vert 0\rangle \quad & M_{2} \vert 5 \rangle = \vert 10\rangle \quad & M_{2} \vert 10 \rangle = \vert 5\rangle \\[1mm] M_{2} \vert 1 \rangle = \vert 2\rangle \quad & M_{2} \vert 6 \rangle = \vert 12\rangle \quad & M_{2} \vert 11 \rangle = \vert 7\rangle \\[1mm] M_{2} \vert 2 \rangle = \vert 4\rangle \quad & M_{2} \vert 7 \rangle = \vert 14\rangle \quad & M_{2} \vert 12 \rangle = \vert 9\rangle \\[1mm] M_{2} \vert 3 \rangle = \vert 6\rangle \quad & M_{2} \vert 8 \rangle = \vert 1\rangle \quad & M_{2} \vert 13 \rangle = \vert 11\rangle \\[1mm] M_{2} \vert 4 \rangle = \vert 8\rangle \quad & M_{2} \vert 9 \rangle = \vert 3\rangle \quad & M_{2} \vert 14 \rangle = \vert 13\rangle \end{array}$$

Esta é uma operação unitária, desde que $\gcd(a,N)=1$; ela embaralha os elementos da base padrão $\{\vert 0\rangle,\ldots,\vert N-1\rangle\}$, portanto, como matriz, é uma matriz de permutação. É evidente pela sua definição que essa operação é determinística, e uma forma simples de ver que ela é invertível é pensar na ordem $r$ de $a$ módulo $N$, e reconhecer que o inverso de $M_a$ é $M_a^{r-1}$.

$$M_a^{r-1} M_a = M_a^r = M_{a^r} = M_1 = \mathbb{I}$$

Há outra maneira de pensar no inverso que não requer nenhum conhecimento de $r$ (que, afinal, é o que estamos tentando calcular). Para todo elemento $a\in\mathbb{Z}_N^{\ast}$ existe sempre um único elemento $b\in\mathbb{Z}_N^{\ast}$ que satisfaz $ab=1$. Denotamos esse elemento $b$ por $a^{-1}$, e ele pode ser calculado eficientemente; uma extensão do algoritmo MDC de Euclides faz isso com custo quadrático em $\operatorname{lg}(N)$. E assim

$$M_{a^{-1}} M_a = M_{a^{-1}a} = M_1 = \mathbb{I}.$$

Portanto, a operação $M_a$ é tanto determinística quanto invertível. Isso implica que ela é descrita por uma matriz de permutação e, portanto, é unitária.

Agora vamos pensar nos autovetores e autovalores da operação $M_a$, assumindo que $a\in\mathbb{Z}_N^{\ast}$. Como acabamos de argumentar, essa suposição nos diz que $M_a$ é unitária.

Existem $N$ autovalores de $M_a$, possivelmente incluindo o mesmo autovalor repetido várias vezes, e em geral há alguma liberdade na seleção dos autovetores correspondentes — mas não precisamos nos preocupar com todas as possibilidades. Vamos começar de forma simples e identificar apenas um autovetor de $M_a$.

$$\vert \psi_0 \rangle = \frac{\vert 1 \rangle + \vert a \rangle + \cdots + \vert a^{r-1} \rangle}{\sqrt{r}}$$

O número $r$ é a ordem de $a$ módulo $N$, aqui e em todo o restante da aula. O autovalor associado a esse autovetor é $1$ porque ele não é alterado quando multiplicamos por $a$.

$$M_a \vert \psi_0 \rangle = \frac{\vert a \rangle + \cdots + \vert a^{r-1} \rangle + \vert a^r \rangle}{\sqrt{r}} = \frac{\vert a \rangle + \cdots + \vert a^{r-1} \rangle + \vert 1 \rangle}{\sqrt{r}} = \vert \psi_0 \rangle$$

Isso acontece porque $a^r = 1$, de modo que cada estado da base padrão $\vert a^k \rangle$ é deslocado para $\vert a^{k+1} \rangle$ para $k\leq r-1$, e $\vert a^{r-1} \rangle$ é deslocado de volta para $\vert 1\rangle$. Informalmente, é como se estivéssemos mexendo lentamente $\vert \psi_0 \rangle$, mas ele já está completamente misturado, então nada muda.

Aqui está outro exemplo de autovetor de $M_a$. Este é mais interessante no contexto de busca de ordem e estimativa de fase.

$$\vert \psi_1 \rangle = \frac{\vert 1 \rangle + \omega_r^{-1} \vert a \rangle + \cdots + \omega_r^{-(r-1)}\vert a^{r-1} \rangle}{\sqrt{r}}$$

Alternativamente, podemos escrever esse vetor usando uma somatória da seguinte forma.

$$\vert \psi_1 \rangle = \frac{1}{\sqrt{r}} \sum_{k = 0}^{r-1} \omega_r^{-k} \vert a^k \rangle$$

Aqui vemos o número complexo $\omega_r = e^{2\pi i/r}$ surgindo naturalmente, devido à forma como a multiplicação por $a$ funciona módulo $N$. Desta vez, o autovalor correspondente é $\omega_r$. Para ver isso, podemos primeiro calcular da seguinte forma.

$$M_a \vert \psi_1 \rangle = \frac{1}{\sqrt{r}}\sum_{k = 0}^{r-1} \omega_r^{-k} M_a\vert a^k \rangle = \frac{1}{\sqrt{r}}\sum_{k = 0}^{r-1} \omega_r^{-k} \vert a^{k+1} \rangle = \frac{1}{\sqrt{r}}\sum_{k = 1}^{r} \omega_r^{-(k - 1)} \vert a^{k} \rangle = \frac{1}{\sqrt{r}}\omega_r \sum_{k = 1}^{r} \omega_r^{-k} \vert a^{k} \rangle$$

Então, como $\omega_r^{-r} = 1 = \omega_r^0$ e $\vert a^r \rangle = \vert 1\rangle = \vert a^0\rangle$, vemos que

$$\frac{1}{\sqrt{r}}\sum_{k = 1}^{r} \omega_r^{-k} \vert a^{k} \rangle = \frac{1}{\sqrt{r}}\sum_{k = 0}^{r-1} \omega_r^{-k} \vert a^k \rangle = \vert\psi_1\rangle,$$

so $M_a \vert\psi_1\rangle = \omega_r \vert\psi_1\rangle$.

Usando o mesmo raciocínio, podemos identificar pares adicionais autovetor/autovalor para $M_a$. Para qualquer escolha de $j\in\{0,\ldots,r-1\}$ temos que

$$\vert \psi_j \rangle = \frac{1}{\sqrt{r}} \sum_{k = 0}^{r-1} \omega_r^{-jk} \vert a^k \rangle$$

é um autovetor de $M_a$ cujo autovalor correspondente é $\omega_r^j$.

$$M_a \vert \psi_j \rangle = \omega_r^j \vert \psi_j \rangle$$

Existem outros autovetores de $M_a$, mas não precisamos nos preocupar com eles — vamos focar exclusivamente nos autovetores $\vert\psi_0\rangle,\ldots,\vert\psi_{r-1}\rangle$ que acabamos de identificar.

Encontrando a ordem por estimativa de fase

Para resolver o problema de encontrar a ordem de um dado $a\in\mathbb{Z}_N^{\ast}$, podemos aplicar o procedimento de estimativa de fase à operação $M_a$.

Para isso, precisamos implementar não apenas $M_a$ de forma eficiente com um circuito quântico, mas também $M_a^2,$ $M_a^4,$ $M_a^8,$ e assim por diante, indo tão longe quanto necessário para obter uma estimativa precisa o suficiente do procedimento de estimativa de fase. Aqui vamos explicar como isso pode ser feito, e determinaremos mais adiante exatamente quanta precisão é necessária.

Vamos começar com a operação $M_a$ por si só. Naturalmente, como estamos trabalhando com o modelo de circuito quântico, usaremos notação binária para codificar os números entre $0$ e $N-1$. O maior número que precisamos codificar é $N-1$, então o número de bits necessários é

$$n = \operatorname{lg}(N-1) = \lfloor \log(N-1) \rfloor + 1.$$

Por exemplo, se $N = 21$ temos $n = \operatorname{lg}(N-1) = 5$. Veja como fica a codificação dos elementos de $\mathbb{Z}_{21}$ como strings binárias de comprimento $5$.

$$\begin{gathered} 0 \mapsto 00000\\[1mm] 1 \mapsto 00001\\[1mm] \vdots\\[1mm] 20 \mapsto 10100 \end{gathered}$$

E agora, aqui está uma definição precisa de como $M_a$ é definida como uma operação de $n$ qubits.

$$M_a \vert x\rangle = \begin{cases} \vert ax \; (\textrm{mod}\;N)\rangle & 0\leq x < N\\[1mm] \vert x\rangle & N\leq x < 2^n \end{cases}$$

O ponto é que, embora só nos importe como $M_a$ funciona para $\vert 0\rangle,\ldots,\vert N-1\rangle$, precisamos especificar como ela funciona para os $2^n - N$ estados da base padrão restantes — e precisamos fazer isso de uma forma que ainda nos dê uma operação unitária. Definir $M_a$ de modo que ela não faça nada aos estados da base padrão restantes resolve isso.

Usando os algoritmos de multiplicação e divisão de inteiros discutidos na lição anterior, junto com a metodologia para implementações reversíveis e sem lixo computacional, podemos construir um circuito quântico que realiza $M_a$, para qualquer escolha de $a\in\mathbb{Z}_N^{\ast}$, com custo $O(n^2)$. Aqui está uma maneira de fazer isso.

1. Construa um circuito para realizar a operação

   $$\vert x \rangle \vert y \rangle \mapsto \vert x \rangle \vert y \oplus f_a(x)\rangle$$

   onde

   $$f_a(x) = \begin{cases} ax \; (\textrm{mod}\;N) & 0\leq x < N\\[1mm] x & N\leq x < 2^n \end{cases}$$

   usando o método descrito na lição anterior. Isso nos dá um circuito de tamanho $O(n^2)$.

2. Troque os dois sistemas de $n$ qubits usando $n$ portas de swap para trocar os qubits individualmente.

3. De forma semelhante ao primeiro passo, construa um circuito para a operação

   $$\vert x \rangle \vert y \rangle \mapsto \vert x \rangle \bigl\vert y \oplus f_{a^{-1}}(x)\bigr\rangle$$

   onde $a^{-1}$ é o inverso de $a$ em $\mathbb{Z}_N^{\ast}$.

Inicializando os $n$ qubits inferiores e compondo os três passos, obtemos esta transformação:

$$\vert x \rangle \vert 0^n \rangle \stackrel{\text{step 1}}{\mapsto} \vert x \rangle \vert f_a(x)\rangle \stackrel{\text{step 2}}{\mapsto} \vert f_a(x)\rangle \vert x \rangle \stackrel{\text{step 3}}{\mapsto} \vert f_a(x)\rangle \bigl\vert x \oplus f_{a^{-1}}(f_a(x)) \bigr\rangle = \vert f_a(x)\rangle\vert 0^n \rangle$$

O método requer qubits de espaço de trabalho, mas eles são devolvidos ao estado inicializado ao final, o que nos permite usar esses circuitos para estimativa de fase. O custo total do circuito obtido é $O(n^2)$.

Para realizar $M_a^2,$ $M_a^4,$ $M_a^8,$ e assim por diante, podemos usar exatamente o mesmo método, exceto que substituímos $a$ por $a^2,$ $a^4,$ $a^8,$ e assim por diante, como elementos de $\mathbb{Z}_N^{\ast}.$ Ou seja, para qualquer potência $k$ que escolhermos, podemos criar um circuito para $M_a^k$ não iterando $k$ vezes o circuito para $M_a,$ mas sim calculando $b = a^k \in \mathbb{Z}_N^{\ast}$ e então usando o circuito para $M_b.$

O cálculo de potências $a^k \in \mathbb{Z}_N$ é o problema de exponenciação modular mencionado na lição anterior. Esse cálculo pode ser feito classicamente, usando o algoritmo de exponenciação modular mencionado na lição anterior (frequentemente chamado de algoritmo de potência na teoria computacional dos números). Na verdade, precisamos apenas das potências de $a$ cujos expoentes são potências de $2$, especificamente $a^2, a^4, \ldots a^{2^{m-1}} \in \mathbb{Z}_N^{\ast}$, e podemos obter essas potências elevando ao quadrado iterativamente $m-1$ vezes. Cada elevação ao quadrado pode ser realizada por um circuito booleano de tamanho $O(n^2)$.

Em essência, o que estamos efetivamente fazendo aqui é transferir o problema de iterar $M_a$ até $2^{m-1}$ vezes para um cálculo clássico eficiente. E é uma sorte que isso seja possível! Para uma escolha arbitrária de circuito quântico no problema de estimativa de fase, isso provavelmente não seria possível — e nesse caso o custo resultante para a estimativa de fase cresce exponencialmente com o número de qubits de controle $m$.

Solução dado um autovetor conveniente

Para entender como podemos resolver o problema de encontrar a ordem usando estimativa de fase, vamos começar supondo que executamos o procedimento de estimativa de fase na operação $M_a$ usando o autovetor $\vert\psi_1\rangle$. Obter esse autovetor não é tarefa fácil, como veremos, então essa não será a história completa — mas é útil começar por aqui.

O autovalor de $M_a$ correspondente ao autovetor $\vert \psi_1\rangle$ é

$$\omega_r = e^{2\pi i \frac{1}{r}}.$$

Ou seja, $\omega_r = e^{2\pi i \theta}$ com $\theta = 1/r.$ Então, se executarmos o procedimento de estimativa de fase em $M_a$ usando o autovetor $\vert\psi_1\rangle,$ obteremos uma aproximação de $1/r.$ Calculando o recíproco, conseguiremos descobrir $r$ — desde que nossa aproximação seja boa o suficiente.

Em mais detalhes, quando executamos o procedimento de estimativa de fase usando $m$ qubits de controle, o que obtemos é um número $y\in\{0,\ldots,2^m-1\}$. Tomamos então $y/2^m$ como estimativa para $\theta$, que é $1/r$ no caso em questão. Para descobrir $r$ a partir dessa aproximação, o natural é calcular o recíproco da nossa aproximação e arredondar para o inteiro mais próximo.

$$\left\lfloor \frac{2^m}{y} + \frac{1}{2} \right\rfloor$$

Por exemplo, suponha que $r = 6$ e realizamos a estimativa de fase em $M_a$ com o autovetor $\vert\psi_1\rangle$ usando $m = 5$ bits de controle. A melhor aproximação de $5$ bits para $1/r = 1/6$ é $5/32$, e temos uma chance razoável (cerca de $68\%$ nesse caso) de obter o resultado $y=5$ da estimativa de fase. Temos

$$\frac{2^m}{y} = \frac{32}{5} = 6.4,$$

e arredondando para o inteiro mais próximo obtemos $6$, que é a resposta correta.

Por outro lado, se não usarmos precisão suficiente, podemos não obter a resposta certa. Por exemplo, se tomarmos $m = 4$ qubits de controle na estimativa de fase, podemos obter a melhor aproximação de $4$ bits para $1/r = 1/6$, que é $3/16$. Calculando o recíproco obtemos

$$\frac{2^m}{y} = \frac{16}{3} = 5.333 \cdots$$

e arredondando para o inteiro mais próximo obtemos a resposta incorreta $5$.

Então, quanta precisão precisamos para obter a resposta certa? Sabemos que a ordem $r$ é um inteiro, e intuitivamente o que precisamos é de precisão suficiente para distinguir $1/r$ de possibilidades próximas, incluindo $1/(r+1)$ e $1/(r-1)$. O número mais próximo de $1/r$ com o qual precisamos nos preocupar é $1/(r+1)$, e a distância entre esses dois números é

$$\frac{1}{r} - \frac{1}{r+1} = \frac{1}{r(r+1)}.$$

Então, se quisermos garantir que não confundimos $1/r$ com $1/(r+1)$, é suficiente usar precisão o bastante para garantir que a melhor aproximação $y/2^m$ para $1/r$ seja mais próxima de $1/r$ do que de $1/(r+1)$. Se usarmos precisão suficiente para garantir que

$$\left\vert \frac{y}{2^m} - \frac{1}{r} \right\vert < \frac{1}{2 r (r+1)},$$

de forma que o erro seja menor que a metade da distância entre $1/r$ e $1/(r+1)$, então $y/2^m$ estará mais próximo de $1/r$ do que de qualquer outra possibilidade, incluindo $1/(r+1)$ e $1/(r-1)$.

Podemos verificar isso da seguinte forma. Suponha que

$$\frac{y}{2^m} = \frac{1}{r} + \varepsilon$$

para $\varepsilon$ satisfazendo

$$\vert\varepsilon\vert < \frac{1}{2 r (r+1)}.$$

Ao calcular o recíproco obtemos

$$\frac{2^m}{y} = \frac{1}{\frac{1}{r} + \varepsilon} = \frac{r}{1+\varepsilon r} = r - \frac{\varepsilon r^2}{1+\varepsilon r}.$$

Maximizando no numerador e minimizando no denominador, podemos limitar o quanto estamos distantes de $r$ da seguinte forma.

$$\left\vert \frac{\varepsilon r^2}{1+\varepsilon r} \right\vert \leq \frac{ \frac{r^2}{2 r(r+1)}}{1 - \frac{r}{2r(r+1)}} = \frac{r}{2 r + 1} < \frac{1}{2}$$

Estamos a menos de $1/2$ de $r$, então, como esperado, obteremos $r$ ao arredondar.

Infelizmente, como ainda não sabemos o que é $r,$ não podemos usá-lo para nos dizer quanta precisão precisamos. O que podemos fazer em vez disso é usar o fato de que $r$ deve ser menor que $N$ para garantir que usamos precisão suficiente. Em particular, se usarmos precisão suficiente para garantir que a melhor aproximação $y/2^m$ para $1/r$ satisfaça

$$\left\vert \frac{y}{2^m} - \frac{1}{r} \right\vert \leq \frac{1}{2N^2},$$

então teremos precisão suficiente para determinar $r$ corretamente ao calcular o recíproco. Tomar $m = 2\operatorname{lg}(N)+1$ garante que temos uma boa chance de obter uma estimativa com essa precisão usando o método descrito anteriormente. (Tomar $m = 2\operatorname{lg}(N)$ é suficiente se você estiver confortável com um limite inferior de 40% na probabilidade de sucesso.)

Solução geral

Como acabamos de ver, se temos o autovetor $\vert \psi_1 \rangle$ de $M_a$, podemos aprender $r$ por meio da estimativa de fase, desde que usemos qubits de controle suficientes para fazer isso com precisão adequada. Infelizmente, não é fácil obter o autovetor $\vert\psi_1\rangle$, então precisamos descobrir como prosseguir.

Vamos supor momentaneamente que procedemos como acima, exceto com o autovetor $\vert\psi_k\rangle$ no lugar de $\vert\psi_1\rangle$, para qualquer escolha de $k\in\{0,\ldots,r-1\}$ que queiramos considerar. O resultado que obtemos do procedimento de estimativa de fase será uma aproximação

$$\frac{y}{2^m} \approx \frac{k}{r}.$$

Trabalhando com a suposição de que não conhecemos nem $k$ nem $r$, isso pode ou não nos permitir identificar $r$. Por exemplo, se $k = 0$ obteremos uma aproximação $y/2^m$ de $0$, o que infelizmente não nos diz nada. Esse, no entanto, é um caso incomum; para outros valores de $k$, pelo menos conseguiremos aprender algo sobre $r$.

Podemos usar um algoritmo conhecido como algoritmo de frações contínuas para converter nossa aproximação $y/2^m$ em frações próximas — incluindo $k/r$ se a aproximação for boa o suficiente. Não vamos explicar o algoritmo de frações contínuas aqui. Em vez disso, aqui está uma declaração de um fato conhecido sobre esse algoritmo.

Fato

Dados um inteiro $N\geq 2$ e um número real $\alpha\in(0,1)$, existe no máximo uma escolha de inteiros $u,v\in\{0,\ldots,N-1\}$ com $v\neq 0$ e $\gcd(u,v)=1$ satisfazendo $\vert \alpha - u/v\vert < \frac{1}{2N^2}$. Dados $\alpha$ e $N$, o algoritmo de frações contínuas encontra $u$ e $v$, ou informa que eles não existem. Esse algoritmo pode ser implementado como um circuito booleano de tamanho $O((\operatorname{lg}(N))^3)$.

Se temos uma aproximação muito próxima $y/2^m$ de $k/r$, e executamos o algoritmo de frações contínuas para $N$ e $\alpha = y/2^m$, obteremos $u$ e $v$, como descritos no fato. Uma análise do fato nos permite concluir que

$$\frac{u}{v} = \frac{k}{r}.$$

Observe em particular que não necessariamente aprendemos $k$ e $r$ separadamente — aprendemos apenas $k/r$ em sua forma irredutível.

Por exemplo, como já notamos, não vamos aprender nada com $k=0$. Mas esse é o único valor de $k$ em que isso acontece. Quando $k$ é diferente de zero, ele pode ter fatores comuns com $r$, mas o número $v$ que obtemos do algoritmo de frações contínuas deve pelo menos dividir $r$.

Está longe de ser óbvio, mas é verdade que se temos a capacidade de aprender $u$ e $v$ com $u/v = k/r$ para $k\in\{0,\ldots,r-1\}$ escolhido uniformemente ao acaso, então é muito provável que consigamos recuperar $r$ após apenas algumas amostras. Em particular, se nosso palpite para $r$ for o mínimo múltiplo comum de todos os valores do denominador $v$ que observamos, estaremos certos com alta probabilidade. Intuitivamente falando, alguns valores de $k$ não são bons porque compartilham fatores comuns com $r$, e esses fatores comuns ficam ocultos quando aprendemos $u$ e $v$. Mas escolhas aleatórias de $k$ provavelmente não vão ocultar fatores de $r$ por muito tempo, e a probabilidade de não adivinharmos $r$ corretamente ao tomar o mínimo múltiplo comum dos denominadores observados cai exponencialmente com o número de amostras.

Resta abordar a questão de como obtemos um autovetor $\vert\psi_k\rangle$ de $M_a$ para executar o procedimento de estimativa de fase. Como se vê, na prática não precisamos criá-los!

O que faremos em vez disso é executar o procedimento de estimativa de fase no estado $\vert 1\rangle$, com o qual queremos dizer a codificação binária de $n$ bits do número $1$, no lugar de um autovetor $\vert\psi\rangle$ de $M_a$. Até agora, falamos apenas em executar o procedimento de estimativa de fase em um autovetor específico, mas nada nos impede de executar o procedimento em um estado de entrada que não seja autovetor de $M_a$, e é isso que estamos fazendo aqui com o estado $\vert 1\rangle$. (Esse não é um autovetor de $M_a$ a menos que $a=1$, o que não é uma escolha de nosso interesse.)

A justificativa para escolher o estado $\vert 1\rangle$ no lugar de um autovetor de $M_a$ é que a seguinte equação é verdadeira.

$$\vert 1\rangle = \frac{1}{\sqrt{r}} \sum_{k = 0}^{r-1} \vert \psi_k\rangle$$

Uma maneira de verificar essa equação é comparar os produtos internos dos dois lados com cada estado da base padrão, usando fórmulas mencionadas anteriormente na lição para ajudar a avaliar os resultados do lado direito. Como consequência, obteremos exatamente os mesmos resultados de medição que se tivéssemos escolhido $k\in\{0,\ldots,r-1\}$ uniformemente ao acaso e usado $\vert\psi_k\rangle$ como autovetor.

Em mais detalhes, vamos imaginar que executamos o procedimento de estimativa de fase com o estado $\vert 1\rangle$ no lugar de um dos autovetores $\vert\psi_k\rangle$. Após a transformada de Fourier quântica inversa ser realizada, isso nos deixa com o estado

$$\frac{1}{\sqrt{r}} \sum_{k = 0}^{r-1} \vert \psi_k\rangle \vert \gamma_k\rangle,$$

onde

$$\vert\gamma_k\rangle = \frac{1}{2^m} \sum_{y=0}^{2^m - 1} \sum_{x=0}^{2^m-1} e^{2\pi i x (k/r - y/2^m)} \vert y\rangle.$$

O vetor $\vert\gamma_k\rangle$ representa o estado dos $m$ qubits superiores após a transformada de Fourier quântica inversa ter sido realizada sobre eles.

Portanto, em virtude do fato de que $\{\vert\psi_0\rangle,\ldots,\vert\psi_{r-1}\rangle\}$ é um conjunto ortonormal, encontramos que uma medição dos $m$ qubits superiores produz uma aproximação $y/2^m$ do valor $k/r$ onde $k\in\{0,\ldots,r-1\}$ é escolhido uniformemente ao acaso. Como já discutimos, isso nos permite aprender $r$ com alto grau de confiança após várias execuções independentes, que era nosso objetivo.

Custo total

O custo para implementar cada unitário controlado $M_a^k$ é $O(n^2)$. Há $m$ operações unitárias controladas, e temos $m = O(n)$, então o custo total para as operações unitárias controladas é $O(n^3)$. Além disso, temos $m$ portas Hadamard (que contribuem $O(n)$ para o custo), e a transformada de Fourier quântica inversa contribui $O(n^2)$ para o custo. Assim, o custo das operações unitárias controladas domina o custo de todo o procedimento — que é portanto $O(n^3)$.

Além do próprio circuito quântico, há alguns cálculos clássicos que precisam ser realizados ao longo do caminho. Isso inclui o cálculo das potências $a^k$ em $\mathbb{Z}_N$ para $k = 2, 4, 8, \ldots, 2^{m-1}$, que são necessárias para criar as portas unitárias controladas, bem como o algoritmo de frações contínuas que converte aproximações de $\theta$ em frações. Esses cálculos podem ser realizados por circuitos booleanos com custo total de $O(n^3)$.

Como é típico, todos esses limites podem ser melhorados usando algoritmos assintoticamente mais rápidos; esses limites assumem que estamos usando algoritmos padrão para operações aritméticas básicas.

Fatoração por busca de ordem

A última coisa que precisamos discutir é como resolver o problema de busca de ordem nos ajuda a fatorar números. Essa parte é completamente clássica — não tem nada de específico com computação quântica.

A ideia básica é a seguinte. Queremos fatorar o número $N$, e podemos fazer isso recursivamente. Especificamente, podemos focar na tarefa de dividir $N$, o que significa encontrar dois inteiros $b,c\geq 2$ tais que $N = bc$. Isso não é possível se $N$ for um número primo, mas podemos testar de forma eficiente se $N$ é primo usando um algoritmo de teste de primalidade antes, e se $N$ não for primo tentaremos dividi-lo. Assim que dividirmos $N$, basta aplicar a recursão em $b$ e $c$ até que todos os fatores sejam primos e obtenhamos a fatoração prima de $N$.

Dividir números pares é fácil: basta retornar $2$ e $N/2$.

Também é fácil dividir potências perfeitas, ou seja, números da forma $N = s^j$ para inteiros $s,j\geq 2,$ aproximando as raízes $N^{1/2},$ $N^{1/3},$ $N^{1/4},$ e assim por diante, e verificando os inteiros próximos como candidatos para $s.$ Não precisamos ir além de $\log(N)$ passos nessa sequência, pois nesse ponto a raiz cai abaixo de $2$ e não revela candidatos adicionais.

É bom que consigamos fazer ambas essas coisas, pois a busca de ordem não nos ajudará a fatorar números pares nem potências de primos, onde o número $s$ é primo. Se $N$ for ímpar e não for uma potência de primo, porém, a busca de ordem nos permite dividir $N$.

Algoritmo probabilístico para dividir um inteiro ímpar e composto N que não é potência de primo

1. Escolha aleatoriamente $a\in\{2,\ldots,N-1\}$.

2. Calcule $d=\gcd(a,N)$.

3. Se $d > 1$, retorne $b = d$ e $c = N/d$ e pare. Caso contrário, continue para o próximo passo sabendo que $a\in\mathbb{Z}_N^{\ast}$.

4. Seja $r$ a ordem de $a$ módulo $N$. (É aqui que precisamos da busca de ordem.)

5. Se $r$ for par:

   5.1 Calcule $x = a^{r/2} - 1$ módulo $N$
   5.2 Calcule $d = \gcd(x,N)$.
   5.3 Se $d>1$, retorne $b=d$ e $c = N/d$ e pare.

6. Se esse ponto for atingido, o algoritmo falhou em encontrar um fator de $N$.

Uma execução desse algoritmo pode falhar em encontrar um fator de $N$. Especificamente, isso ocorre em duas situações:

• A ordem de $a$ módulo $N$ é ímpar.
• A ordem de $a$ módulo $N$ é par e $\gcd\bigl(a^{r/2} - 1, N\bigr) = 1$.

Usando teoria dos números elementar, pode-se provar que, para uma escolha aleatória de $a$, com probabilidade de pelo menos $1/2$ nenhum desses eventos ocorre. Na verdade, a probabilidade de que algum desses eventos ocorra é no máximo $2^{-(m-1)}$, onde $m$ é o número de fatores primos distintos de $N$ — e é por isso que a hipótese de que $N$ não é uma potência de primo é necessária. (A hipótese de que $N$ é ímpar também é necessária para que esse resultado seja válido.)

Isso significa que cada execução tem pelo menos 50% de chance de dividir $N$. Portanto, se rodarmos o algoritmo $t$ vezes, escolhendo $a$ aleatoriamente a cada vez, conseguiremos dividir $N$ com probabilidade de pelo menos $1 - 2^{-t}$.

A ideia fundamental por trás do algoritmo é a seguinte. Se temos uma escolha de $a$ para a qual a ordem $r$ de $a$ módulo $N$ é par, então $r/2$ é um inteiro e podemos considerar os números

$$a^{r/2} - 1\; (\textrm{mod}\; N) \quad \text{e} \quad a^{r/2} + 1\; (\textrm{mod}\; N).$$

Usando a fórmula $Z^2 - 1 = (Z+1)(Z-1)$, concluímos que

$$\bigl(a^{r/2} - 1\bigr) \bigl(a^{r/2} + 1\bigr) = a^r - 1.$$

Agora, sabemos que $a^r \; (\textrm{mod}\; N) = 1$ pela definição da ordem — o que é equivalente a dizer que $N$ divide $a^r - 1$ exatamente. Isso significa que $N$ divide o produto

$$\bigl(a^{r/2} - 1\bigr) \bigl(a^{r/2} + 1\bigr).$$

Para que isso seja verdade, todos os fatores primos de $N$ devem ser também fatores primos de $a^{r/2} - 1$ ou de $a^{r/2} + 1$ (ou de ambos) — e para uma escolha aleatória de $a$ é improvável que todos os fatores primos de $N$ dividam um mesmo termo e nenhum divida o outro. Caso contrário, desde que alguns fatores primos de $N$ dividam o primeiro termo e outros dividam o segundo, será possível encontrar um fator não trivial de $N$ calculando o MDC com o primeiro termo.