{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "e72c9e73-98e1-49ae-bc3c-511954eeb15c",
      "metadata": {
        "tags": [
          "remove_cell"
        ]
      },
      "source": [
        "---\n",
        "title: \"Algoritmo de Shor\"\n",
        "description: \"Este tutorial se centra en demostrar el algoritmo de Shor mediante la factorización del número 15 en un ordenador cuántico.\"\n",
        "---\n",
        "\n",
        "{/* cspell:ignore textrm */}\n",
        "\n",
        "<span id=\"shors-algorithm\" />\n",
        "\n",
        "# Algoritmo de Shor\n",
        "\n",
        "*Estimación de uso: Tres segundos en un procesador Eagle r3 (NOTA: Esto es sólo una estimación. Su tiempo de ejecución puede variar)*\n",
        "\n",
        "<span id=\"learning-outcomes\" />\n",
        "\n",
        "## Resultados del aprendizaje\n",
        "\n",
        "Una vez completado este tutorial, los usuarios deberían comprender:\n",
        "\n",
        "* Los fundamentos matemáticos del algoritmo de Shor para la factorización de números enteros\n",
        "* Cómo ejecutar una instancia de ejemplo de este algoritmo en un equipo\n",
        "\n",
        "<span id=\"prerequisites\" />\n",
        "\n",
        "## Requisitos previos\n",
        "\n",
        "Recomendamos a los usuarios que se familiaricen con los siguientes temas antes de seguir este tutorial:\n",
        "\n",
        "* [Fundamentos de los algoritmos cuánticos](/learning/courses/fundamentals-of-quantum-algorithms).\n",
        "* [Estimación de fase y factorización](/learning/courses/fundamentals-of-quantum-algorithms/phase-estimation-and-factoring/introduction). En este tutorial tratamos parte de este tema.\n",
        "\n",
        "<span id=\"background\" />\n",
        "\n",
        "## En segundo plano\n",
        "\n",
        "[El algoritmo de Shor](https://epubs.siam.org/doi/abs/10.1137/S0036144598347011), desarrollado por Peter Shor en 1994, es un algoritmo cuántico revolucionario que permite factorizar números enteros en tiempo polinomial. Su importancia radica en su capacidad para factorizar números enteros grandes de forma exponencialmente más rápida que cualquier algoritmo clásico conocido, lo que pone en peligro la seguridad de sistemas criptográficos de uso generalizado como el RSA, que se basan en la dificultad de factorizar números grandes. Si se resolviera este problema de manera eficiente en un ordenador cuántico lo suficientemente potente, el algoritmo de Shor podría revolucionar campos como la criptografía, la ciberseguridad y las matemáticas computacionales, lo que pondría de relieve el poder transformador de la computación cuántica.\n",
        "\n",
        "Este tutorial se centra en la demostración del algoritmo de Shor mediante la factorización de 15 en un ordenador cuántico.\n",
        "\n",
        "En primer lugar, definimos el problema de búsqueda de orden y construimos los circuitos correspondientes a partir del protocolo cuántico de estimación de fase. A continuación, ejecutamos los circuitos de búsqueda de órdenes en hardware real utilizando circuitos de menor profundidad que podamos transpilar. La última sección completa el algoritmo de Shor conectando el problema de búsqueda de órdenes con la factorización entera.\n",
        "\n",
        "Terminamos el tutorial con una discusión sobre otras demostraciones del algoritmo de Shor en hardware real, centrándonos tanto en implementaciones genéricas como en aquellas adaptadas a la factorización de enteros específicos como 15 y 21.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "22f3a6f3-0ba5-4826-a4f3-9dcc62f51c70",
      "metadata": {},
      "source": [
        "Nota: Este tutorial se centra más en la implementación y demostración de los circuitos relativos al algoritmo de Shor. Para profundizar en el material, consulte el curso [Fundamentos de los algoritmos cuánticos](/learning/courses/fundamentals-of-quantum-algorithms/phase-estimation-and-factoring/introduction), del Dr. John Watrous, y los artículos de la sección [Referencias](#references).\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "b41b8639-ff72-4cd5-8bd3-3c16c01bedc0",
      "metadata": {},
      "source": [
        "<span id=\"requirements\" />\n",
        "\n",
        "### Requisitos\n",
        "\n",
        "Antes de empezar este tutorial, asegúrate de que tienes instalado lo siguiente:\n",
        "\n",
        "* Qiskit SDK v2.0 o posterior, con soporte [de visualización](/docs/api/qiskit/visualization)\n",
        "* Qiskit Runtime v0.40 o posterior (`pip install qiskit-ibm-runtime`)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "93dcb213-ddcd-4329-8bcb-7c61fa298ee5",
      "metadata": {},
      "source": [
        "<span id=\"setup\" />\n",
        "\n",
        "### Configuración\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "0860914d-cf1f-4bee-907d-3c442cd92cc6",
      "metadata": {},
      "outputs": [],
      "source": [
        "import numpy as np\n",
        "import pandas as pd\n",
        "from fractions import Fraction\n",
        "from math import floor, gcd, log\n",
        "\n",
        "from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister\n",
        "from qiskit.circuit.library import QFT, UnitaryGate\n",
        "from qiskit.transpiler import CouplingMap, generate_preset_pass_manager\n",
        "from qiskit.visualization import plot_histogram\n",
        "\n",
        "from qiskit_ibm_runtime import QiskitRuntimeService\n",
        "from qiskit_ibm_runtime import SamplerV2 as Sampler"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a28eeaa5-2da0-4b6b-8877-547a86fbcd76",
      "metadata": {},
      "source": [
        "<span id=\"step-1-map-classical-inputs-to-a-quantum-problem\" />\n",
        "\n",
        "## Paso 1: Asignar entradas clásicas a un problema cuántico\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c09ec908-51ba-4ea3-bc84-142e8f8002bf",
      "metadata": {},
      "source": [
        "El algoritmo de Shor para la factorización de enteros utiliza un problema intermedio conocido como problema de *búsqueda de órdenes*. En esta sección, demostramos cómo resolver el problema de búsqueda de orden utilizando *la estimación cuántica de fase*.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e1aa8e70-a760-43ca-bcf2-def23d49ba36",
      "metadata": {},
      "source": [
        "<span id=\"phase-estimation-problem\" />\n",
        "\n",
        "### Problema de estimación de fase\n",
        "\n",
        "En el problema de estimación de fase, se nos da un estado cuántico $\\ket{\\psi}$ de $n$ qubits, junto con un circuito cuántico unitario que actúa sobre $n$ qubits. Nos han prometido que $\\ket{\\psi}$ es un vector propio de la matriz unitaria $U$ que describe la acción del circuito, y nuestro objetivo es calcular o aproximar el valor propio $\\lambda = e^{2 \\pi i \\theta}$ al que corresponde $\\ket{\\psi}$. En otras palabras, el circuito debe dar como salida una aproximación al número $\\theta \\in [0, 1)$ que satisfaga $U \\ket{\\psi}= e^{2 \\pi i \\theta} \\ket{\\psi}.$ El objetivo del circuito de estimación de fase es aproximar $\\theta$ en $m$ bits. Matemáticamente hablando, nos gustaría encontrar $y$ tal que $\\theta \\approx y / 2^m$, donde $y \\in {0, 1, 2, \\dots, 2^{m-1}}$. La siguiente imagen muestra el circuito cuántico que estima $y$ en $m$ bits realizando una medida en $m$ qubits.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7ae92a5f-55b4-4597-be58-fe47b450dc39",
      "metadata": {},
      "source": [
        "![Circuito cuántico de estimación de fase](https://quantum.cloud.ibm.com/learning/images/courses/fundamentals-of-quantum-algorithms/phase-estimation-and-factoring/phase-estimation-procedure.svg)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c527ff28-8038-44fc-8d62-1c9dda2b600e",
      "metadata": {},
      "source": [
        "En el circuito anterior, los qubits superiores $m$ se inician en el estado $\\ket{0^m}$, y los qubits inferiores $n$ se inician en $\\ket{\\psi}$, que se promete que es un eigenvector de $U$. El primer ingrediente en el circuito de estimación de fase son las operaciones unitarias controladas que son responsables de realizar un *retroceso de fase* a su qubit de control correspondiente. Estos unitarios controlados se exponencian en función de la posición del qubit de control, desde el bit menos significativo hasta el bit más significativo. Dado que $\\ket{\\psi}$ es un eigenvector de $U$, el estado de los qubits inferiores $n$ no se ve afectado por esta operación, pero la información de fase del eigenvalor se propaga a los qubits superiores $m$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "40b4a6e3-6f46-4971-8054-4a7ae9dd535e",
      "metadata": {},
      "source": [
        "Resulta que después de la operación de retroceso de fase mediante unidades controladas, todos los estados posibles de los qubits superiores $m$ son ortonormales entre sí para cada vector propio $\\ket{\\psi}$ del unitario $U$. Por lo tanto, estos estados son perfectamente distinguibles, y podemos rotar la base que forman de nuevo a la base computacional para hacer una medición. Un análisis matemático muestra que esta matriz de rotación corresponde a la transformada cuántica de Fourier (QFT) inversa en el espacio de Hilbert $2^m$ -dimensional. La intuición detrás de esto es que la estructura periódica de los operadores de exponenciación modular está codificada en el estado cuántico, y la QFT convierte esta periodicidad en picos medibles en el dominio de la frecuencia.\n",
        "\n",
        "Para una comprensión más profunda de por qué se emplea el circuito QFT en el algoritmo de Shor, remitimos al lector al curso [Fundamentos de algoritmos cuánticos](/learning/courses/fundamentals-of-quantum-algorithms/phase-estimation-and-factoring/introduction).\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "acfa31d3-7980-42ad-b418-77a24d16513a",
      "metadata": {},
      "source": [
        "Ahora estamos preparados para utilizar el circuito de estimación de fase para la búsqueda de órdenes.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "535674d6-cc8a-45f5-9fc7-f5fd3053f4c4",
      "metadata": {},
      "source": [
        "<span id=\"order-finding-problem\" />\n",
        "\n",
        "### Problema con el pedido\n",
        "\n",
        "Para definir el problema de búsqueda de órdenes, empezamos con algunos conceptos de teoría de números. En primer lugar, para cualquier número entero positivo dado $N$, defina el conjunto $\\mathbb{Z}_N$ como $\\mathbb{Z}_N = \\{0, 1, 2, \\dots, N-1\\}.$ Todas las operaciones aritméticas en $\\mathbb{Z}_N$ se realizan modulo $N$. En particular, todos los elementos $a \\in \\mathbb{Z}_n$ que son coprimos con $N$ son especiales y constituyen $\\mathbb{Z}^*_N$ como $\\mathbb{Z}^*_N = \\{ a \\in \\mathbb{Z}_N : \\mathrm{gcd}(a, N)=1 \\}.$ Para un elemento $a \\in \\mathbb{Z}^*_N$, el menor número entero positivo $r$ tal que $a^r \\equiv 1 \\; (\\mathrm{mod} \\; N)$ se define como el *orden* de $a$ módulo $N$. Como veremos más adelante, encontrar el orden de un $a \\in \\mathbb{Z}^*_N$ nos permitirá factorizar $N$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6de4cbed-6591-40da-a994-2ff139c9f064",
      "metadata": {},
      "source": [
        "Para construir el circuito de búsqueda de orden a partir del circuito de estimación de fase, necesitamos dos consideraciones. En primer lugar, necesitamos definir el unitario $U$ que nos permitirá encontrar el orden $r$, y en segundo lugar, necesitamos definir un vector propio $\\ket{\\psi}$ de $U$ para preparar el estado inicial del circuito de estimación de fase.\n",
        "\n",
        "Para conectar el problema de búsqueda de orden con la estimación de fase, consideramos la operación definida sobre un sistema cuyos estados clásicos corresponden a $\\mathbb{Z}_N$, donde multiplicamos por un elemento fijo $a \\in \\mathbb{Z}^*_N$. En particular, definimos este operador de multiplicación $M_a$ tal que $M_a \\ket{x} = \\ket{ax \\; (\\mathrm{mod} \\; N)}$ para cada $x \\in \\mathbb{Z}_N$. Nótese que está implícito que estamos tomando el producto modulo $N$ dentro del ket en el lado derecho de la ecuación. Un análisis matemático muestra que $M_a$ es un operador unitario. Además, resulta que $M_a$ tiene pares de vectores y valores propios que nos permiten conectar el orden $r$ de $a$ con el problema de estimación de fase. En concreto, para cualquier elección de $j \\in \\{0, \\dots, r-1\\}$, tenemos que $\\ket{\\psi_j} = \\frac{1}{\\sqrt{r}} \\sum^{r-1}_{k=0} \\omega^{-jk}_{r} \\ket{a^k}$ es un vector propio de $M_a$ cuyo valor propio correspondiente es $\\omega^{j}_{r}$, donde $\\omega^{j}_{r} = e^{2 \\pi i \\frac{j}{r}}.$\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "57942d8d-b2fc-4863-bfdd-f99aa941e0c0",
      "metadata": {},
      "source": [
        "Por observación, vemos que un par eigenvector/eigenvalor conveniente es el estado $\\ket{\\psi_1}$ con $\\omega^{1}_{r} = e^{2 \\pi i \\frac{1}{r}}$. Por lo tanto, si pudiéramos encontrar el vector propio $\\ket{\\psi_1}$, podríamos estimar la fase $\\theta=1/r$ con nuestro circuito cuántico y, por lo tanto, obtener una estimación del orden $r$. Sin embargo, no es fácil hacerlo, y tenemos que considerar una alternativa.\n",
        "\n",
        "Consideremos cuál sería el resultado del circuito si preparamos el estado computacional $\\ket{1}$ como estado inicial. No se trata de un estado propio de $M_a$, sino de la superposición uniforme de los estados propios que acabamos de describir. En otras palabras, se cumple la siguiente relación. $\\ket{1} = \\frac{1}{\\sqrt{r}} \\sum^{r-1}_{k=0} \\ket{\\psi_k}$ La implicación de la ecuación anterior es que si fijamos el estado inicial en $\\ket{1}$, obtendremos precisamente el mismo resultado de medición que si hubiéramos elegido $k \\in \\{ 0, \\dots, r-1\\}$ uniformemente al azar y utilizado $\\ket{\\psi_k}$ como un vector propio en el circuito de estimación de fase. En otras palabras, una medición de los $m$ qubits superiores produce una aproximación $y / 2^m$ al valor $k / r$ donde $k \\in \\{ 0, \\dots, r-1\\}$ se elige uniformemente al azar. Esto nos permite aprender $r$ con un alto grado de confianza tras varias ejecuciones independientes, que era nuestro objetivo.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c037fd9d-2ad6-4258-a11d-bd91944e7f01",
      "metadata": {},
      "source": [
        "<span id=\"modular-exponentiation-operators\" />\n",
        "\n",
        "### Operadores de exponenciación modular\n",
        "\n",
        "Hasta ahora, hemos vinculado el problema de estimación de fase al problema de búsqueda de orden definiendo $U = M_a$ y $\\ket{\\psi} = \\ket{1}$ en nuestro circuito cuántico. Por lo tanto, el último ingrediente que queda es encontrar una manera eficiente de definir exponenciales modulares de $M_a$ como $M_a^k$ para $k = 1, 2, 4, \\dots, 2^{m-1}$. Para realizar este cálculo, encontramos que para cualquier potencia $k$ que elijamos, podemos crear un circuito para $M_a^k$ no iterando $k$ veces el circuito para $M_a$, sino calculando $b = a^k \\; \\mathrm{mod} \\; N$ y luego usando el circuito para $M_b$. Dado que sólo necesitamos las potencias que son potencias de 2, podemos hacer esto de forma clásica y eficiente utilizando el cuadrado iterativo.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "09d93ee7-ab42-43a8-a44f-ad303c9407c7",
      "metadata": {},
      "source": [
        "<span id=\"step-2-optimize-problem-for-quantum-hardware-execution\" />\n",
        "\n",
        "## Paso 2: Optimizar el problema para la ejecución en hardware cuántico\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "80d39ec2-715f-4014-9dbc-37db001528b1",
      "metadata": {},
      "source": [
        "<span id=\"specific-example-with-$n-=-15$-and-$a=2$\" />\n",
        "\n",
        "### Ejemplo específico con $N = 15$ y $a=2$\n",
        "\n",
        "Podemos detenernos aquí para discutir un ejemplo específico y construir el circuito de búsqueda de orden para $N=15$. Obsérvese que los posibles $a \\in \\mathbb{Z}_N^*$ no triviales para $N=15$ son $a \\in \\{2, 4, 7, 8, 11, 13, 14 \\}$. Para este ejemplo, elegimos $a=2$. Construiremos el operador $M_2$ y los operadores de exponenciación modular $M_2^k$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "5fcd0c79-0c62-42f1-b6fa-c0f4c58e1994",
      "metadata": {},
      "source": [
        "La acción de $M_2$ sobre los estados de la base de cálculo es la siguiente. $M_2 \\ket{0} = \\ket{0} \\quad M_2 \\ket{5} = \\ket{10} \\quad M_2 \\ket{10} = \\ket{5}$ $M_2 \\ket{1} = \\ket{2} \\quad M_2 \\ket{6} = \\ket{12} \\quad M_2 \\ket{11} = \\ket{7}$ $M_2 \\ket{2} = \\ket{4} \\quad M_2 \\ket{7} = \\ket{14} \\quad M_2 \\ket{12} = \\ket{9}$ $M_2 \\ket{3} = \\ket{6} \\quad M_2 \\ket{8} = \\ket{1} \\quad M_2 \\ket{13} = \\ket{11}$ $M_2 \\ket{4} = \\ket{8} \\quad M_2 \\ket{9} = \\ket{3} \\quad M_2 \\ket{14} = \\ket{13}$ Por observación, podemos ver que los estados base se barajan, por lo que tenemos una matriz de permutación. Podemos construir esta operación en cuatro qubits con puertas de intercambio. A continuación, construimos las operaciones $M_2$ y $M_2$ controlada.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "id": "ec5d416b-7fc9-43f3-b291-e0edc8ad195f",
      "metadata": {},
      "outputs": [],
      "source": [
        "def M2mod15():\n",
        "    \"\"\"\n",
        "    M2 (mod 15)\n",
        "    \"\"\"\n",
        "    b = 2\n",
        "    U = QuantumCircuit(4)\n",
        "\n",
        "    U.swap(2, 3)\n",
        "    U.swap(1, 2)\n",
        "    U.swap(0, 1)\n",
        "\n",
        "    U = U.to_gate()\n",
        "    U.name = f\"M_{b}\"\n",
        "\n",
        "    return U"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "id": "0a8885f1-91d4-40bd-912d-dc5eea05f5bd",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/0a8885f1-91d4-40bd-912d-dc5eea05f5bd-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 3,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Get the M2 operator\n",
        "M2 = M2mod15()\n",
        "\n",
        "# Add it to a circuit and plot\n",
        "circ = QuantumCircuit(4)\n",
        "circ.compose(M2, inplace=True)\n",
        "circ.decompose(reps=2).draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "id": "d0ae7456-053a-4389-8653-a1c9c8ff757c",
      "metadata": {},
      "outputs": [],
      "source": [
        "def controlled_M2mod15():\n",
        "    \"\"\"\n",
        "    Controlled M2 (mod 15)\n",
        "    \"\"\"\n",
        "    b = 2\n",
        "    U = QuantumCircuit(4)\n",
        "\n",
        "    U.swap(2, 3)\n",
        "    U.swap(1, 2)\n",
        "    U.swap(0, 1)\n",
        "\n",
        "    U = U.to_gate()\n",
        "    U.name = f\"M_{b}\"\n",
        "    c_U = U.control()\n",
        "\n",
        "    return c_U"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "id": "ab7fe331-2f9e-47ca-ba3b-f5d67992062a",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/ab7fe331-2f9e-47ca-ba3b-f5d67992062a-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 5,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Get the controlled-M2 operator\n",
        "controlled_M2 = controlled_M2mod15()\n",
        "\n",
        "# Add it to a circuit and plot\n",
        "circ = QuantumCircuit(5)\n",
        "circ.compose(controlled_M2, inplace=True)\n",
        "circ.decompose(reps=1).draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "139b6d31-33d0-49ea-82b1-5bb50b5b6c9e",
      "metadata": {},
      "source": [
        "Las puertas que actúan sobre más de dos qubits se descompondrán en puertas de dos qubits.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "id": "13b4841d-a4ac-46bd-b4d0-d111b3017189",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/13b4841d-a4ac-46bd-b4d0-d111b3017189-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 6,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "circ.decompose(reps=2).draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e361cc09-e357-4ba7-950b-9fa83898fa88",
      "metadata": {},
      "source": [
        "Ahora necesitamos construir los operadores de exponenciación modular. Para obtener suficiente precisión en la estimación de fase, utilizaremos ocho qubits para la medición de la estimación. Por lo tanto, necesitamos construir $M_b$ con $b = a^{2^k} \\; (\\mathrm{mod} \\; N)$ para cada $k = 0, 1, \\dots, 7$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "id": "72d989c8-ef62-44d7-887c-5fa13db818e9",
      "metadata": {},
      "outputs": [],
      "source": [
        "def a2kmodN(a, k, N):\n",
        "    \"\"\"Compute a^{2^k} (mod N) by repeated squaring\"\"\"\n",
        "    for _ in range(k):\n",
        "        a = int(np.mod(a**2, N))\n",
        "    return a"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "id": "69fa8c9f-4107-4339-96db-c6f950a71261",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "[2, 4, 1, 1, 1, 1, 1, 1]\n"
          ]
        }
      ],
      "source": [
        "k_list = range(8)\n",
        "b_list = [a2kmodN(2, k, 15) for k in k_list]\n",
        "\n",
        "print(b_list)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e28ed7e5-f7ca-4b7b-b875-f14c5229588c",
      "metadata": {},
      "source": [
        "Como podemos ver en la lista de valores $b$, además de $M_2$ que construimos previamente, también necesitamos construir $M_4$ y $M_1$. Nótese que $M_1$ actúa trivialmente sobre los estados de la base computacional, por lo que es simplemente el operador identidad.\n",
        "\n",
        "$M_4$ actúa sobre los estados base computacionales de la siguiente manera. $M_4 \\ket{0} = \\ket{0} \\quad M_4 \\ket{5} = \\ket{5} \\quad M_4 \\ket{10} = \\ket{10}$ $M_4 \\ket{1} = \\ket{4} \\quad M_4 \\ket{6} = \\ket{9} \\quad M_4 \\ket{11} = \\ket{14}$ $M_4 \\ket{2} = \\ket{8} \\quad M_4 \\ket{7} = \\ket{13} \\quad M_4 \\ket{12} = \\ket{3}$ $M_4 \\ket{3} = \\ket{12} \\quad M_4 \\ket{8} = \\ket{2} \\quad M_4 \\ket{13} = \\ket{7}$ $M_4 \\ket{4} = \\ket{1} \\quad M_4 \\ket{9} = \\ket{6} \\quad M_4 \\ket{14} = \\ket{11}$\n",
        "\n",
        "Por lo tanto, esta permutación se puede construir con la siguiente operación de intercambio.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "id": "4f3a5fe4-5449-4869-ae94-6fdefd6765f0",
      "metadata": {},
      "outputs": [],
      "source": [
        "def M4mod15():\n",
        "    \"\"\"\n",
        "    M4 (mod 15)\n",
        "    \"\"\"\n",
        "    b = 4\n",
        "    U = QuantumCircuit(4)\n",
        "\n",
        "    U.swap(1, 3)\n",
        "    U.swap(0, 2)\n",
        "\n",
        "    U = U.to_gate()\n",
        "    U.name = f\"M_{b}\"\n",
        "\n",
        "    return U"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "id": "be041e3d-28b1-453e-983e-184c2366aeb9",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/be041e3d-28b1-453e-983e-184c2366aeb9-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 10,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Get the M4 operator\n",
        "M4 = M4mod15()\n",
        "\n",
        "# Add it to a circuit and plot\n",
        "circ = QuantumCircuit(4)\n",
        "circ.compose(M4, inplace=True)\n",
        "circ.decompose(reps=2).draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "id": "0efb7000-7d13-4b73-95d4-f9f747ec5119",
      "metadata": {},
      "outputs": [],
      "source": [
        "def controlled_M4mod15():\n",
        "    \"\"\"\n",
        "    Controlled M4 (mod 15)\n",
        "    \"\"\"\n",
        "    b = 4\n",
        "    U = QuantumCircuit(4)\n",
        "\n",
        "    U.swap(1, 3)\n",
        "    U.swap(0, 2)\n",
        "\n",
        "    U = U.to_gate()\n",
        "    U.name = f\"M_{b}\"\n",
        "    c_U = U.control()\n",
        "\n",
        "    return c_U"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 12,
      "id": "8d943b00-a502-4157-8a0d-13fb1f55e705",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/8d943b00-a502-4157-8a0d-13fb1f55e705-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 12,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Get the controlled-M4 operator\n",
        "controlled_M4 = controlled_M4mod15()\n",
        "\n",
        "# Add it to a circuit and plot\n",
        "circ = QuantumCircuit(5)\n",
        "circ.compose(controlled_M4, inplace=True)\n",
        "circ.decompose(reps=1).draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7d2a4f30-afb1-49fc-abbb-b9070b58fe8d",
      "metadata": {},
      "source": [
        "Las puertas que actúan sobre más de dos qubits se descompondrán en puertas de dos qubits.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "id": "68399eef-5e55-4c95-a8a4-c8efaebd34b9",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/68399eef-5e55-4c95-a8a4-c8efaebd34b9-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 13,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "circ.decompose(reps=2).draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "1f12fb24-257b-4b12-b7c5-e01b975e3216",
      "metadata": {},
      "source": [
        "Hemos visto que los operadores $M_b$ para un $b \\in \\mathbb{Z}^*_N$ dado son operaciones de permutación. Debido al tamaño relativamente pequeño del problema de permutación que tenemos aquí, ya que $N=15$ requiere sólo cuatro qubits, pudimos sintetizar estas operaciones directamente con `SWAP` gates por inspección. En general, éste podría no ser un enfoque escalable. En su lugar, podríamos necesitar construir la matriz de permutación explícitamente, y utilizar la clase `UnitaryGate` de Qiskit y los métodos de transpilación para sintetizar esta matriz de permutación. Sin embargo, esto puede dar lugar a circuitos mucho más profundos. A continuación\n",
        "se proporciona un ejemplo.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 14,
      "id": "33a328b8-11d8-4b0c-a277-5d76b2c07c5e",
      "metadata": {},
      "outputs": [],
      "source": [
        "def mod_mult_gate(b, N):\n",
        "    \"\"\"\n",
        "    Modular multiplication gate from permutation matrix.\n",
        "    \"\"\"\n",
        "    if gcd(b, N) > 1:\n",
        "        print(f\"Error: gcd({b},{N}) > 1\")\n",
        "    else:\n",
        "        n = floor(log(N - 1, 2)) + 1\n",
        "        U = np.full((2**n, 2**n), 0)\n",
        "        for x in range(N):\n",
        "            U[b * x % N][x] = 1\n",
        "        for x in range(N, 2**n):\n",
        "            U[x][x] = 1\n",
        "        G = UnitaryGate(U)\n",
        "        G.name = f\"M_{b}\"\n",
        "        return G"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "id": "c184f6dd-9f80-4487-ac0b-0dd94170b0f0",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "qubits: 4\n",
            "2q-depth: 94\n",
            "2q-size: 96\n",
            "Operator counts: OrderedDict({'cx': 45, 'swap': 32, 'u': 24, 'u1': 7, 'u3': 4, 'unitary': 3, 'circuit-335': 1, 'circuit-338': 1, 'circuit-341': 1, 'circuit-344': 1, 'circuit-347': 1, 'circuit-350': 1, 'circuit-353': 1, 'circuit-356': 1, 'circuit-359': 1, 'circuit-362': 1, 'circuit-365': 1, 'circuit-368': 1, 'circuit-371': 1, 'circuit-374': 1, 'circuit-377': 1, 'circuit-380': 1})\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/c184f6dd-9f80-4487-ac0b-0dd94170b0f0-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 15,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Let's build M2 using the permutation matrix definition\n",
        "M2_other = mod_mult_gate(2, 15)\n",
        "\n",
        "# Add it to a circuit\n",
        "circ = QuantumCircuit(4)\n",
        "circ.compose(M2_other, inplace=True)\n",
        "circ = circ.decompose()\n",
        "\n",
        "# Transpile the circuit and get the depth\n",
        "coupling_map = CouplingMap.from_line(4)\n",
        "pm = generate_preset_pass_manager(coupling_map=coupling_map)\n",
        "transpiled_circ = pm.run(circ)\n",
        "\n",
        "print(f\"qubits: {circ.num_qubits}\")\n",
        "print(\n",
        "    f\"2q-depth: {transpiled_circ.depth(lambda x: x.operation.num_qubits==2)}\"\n",
        ")\n",
        "print(f\"2q-size: {transpiled_circ.size(lambda x: x.operation.num_qubits==2)}\")\n",
        "print(f\"Operator counts: {transpiled_circ.count_ops()}\")\n",
        "transpiled_circ.decompose().draw(\n",
        "    output=\"mpl\", fold=-1, style=\"clifford\", idle_wires=False\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "fe8caaf6-4bb6-4d27-8177-71f1b530556d",
      "metadata": {},
      "source": [
        "Comparemos estos recuentos con la profundidad del circuito compilado de nuestra implementación manual de la puerta $M_2$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "id": "0235c931-0adb-4972-9fce-32a0341822bf",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "qubits: 4\n",
            "2q-depth: 9\n",
            "2q-size: 9\n",
            "Operator counts: OrderedDict({'cx': 9})\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/0235c931-0adb-4972-9fce-32a0341822bf-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 16,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Get the M2 operator from our manual construction\n",
        "M2 = M2mod15()\n",
        "\n",
        "# Add it to a circuit\n",
        "circ = QuantumCircuit(4)\n",
        "circ.compose(M2, inplace=True)\n",
        "circ = circ.decompose(reps=3)\n",
        "\n",
        "# Transpile the circuit and get the depth\n",
        "coupling_map = CouplingMap.from_line(4)\n",
        "pm = generate_preset_pass_manager(coupling_map=coupling_map)\n",
        "transpiled_circ = pm.run(circ)\n",
        "\n",
        "print(f\"qubits: {circ.num_qubits}\")\n",
        "print(\n",
        "    f\"2q-depth: {transpiled_circ.depth(lambda x: x.operation.num_qubits==2)}\"\n",
        ")\n",
        "print(f\"2q-size: {transpiled_circ.size(lambda x: x.operation.num_qubits==2)}\")\n",
        "print(f\"Operator counts: {transpiled_circ.count_ops()}\")\n",
        "transpiled_circ.draw(\n",
        "    output=\"mpl\", fold=-1, style=\"clifford\", idle_wires=False\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c3f0f349-a26c-4176-aa9d-407e17184b91",
      "metadata": {},
      "source": [
        "Como podemos ver, el enfoque de la matriz de permutación resultó en un circuito significativamente profundo incluso para una única puerta $M_2$ en comparación con nuestra implementación manual del mismo. Por lo tanto, continuaremos con nuestra aplicación anterior de las operaciones $M_b$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0dafd632-797b-402c-a88a-821a62c7265a",
      "metadata": {},
      "source": [
        "Ahora, estamos listos para construir el circuito de búsqueda de orden completo utilizando nuestros operadores de exponenciación modular controlada previamente definidos. En el siguiente código, también importamos el [circuito QFT](/docs/api/qiskit/qiskit.circuit.library.QFT) de la librería Qiskit Circuit, que utiliza puertas Hadamard en cada qubit, una serie de puertas controlled-U1 (o Z, dependiendo de la fase), y una capa de puertas swap.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 17,
      "id": "0e854aed-c11b-494c-8c80-adeb8eb0e8fe",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/0e854aed-c11b-494c-8c80-adeb8eb0e8fe-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 17,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Order finding problem for N = 15 with a = 2\n",
        "N = 15\n",
        "a = 2\n",
        "\n",
        "# Number of qubits\n",
        "num_target = floor(log(N - 1, 2)) + 1  # for modular exponentiation operators\n",
        "num_control = 2 * num_target  # for enough precision of estimation\n",
        "\n",
        "# List of M_b operators in order\n",
        "k_list = range(num_control)\n",
        "b_list = [a2kmodN(2, k, 15) for k in k_list]\n",
        "\n",
        "# Initialize the circuit\n",
        "control = QuantumRegister(num_control, name=\"C\")\n",
        "target = QuantumRegister(num_target, name=\"T\")\n",
        "output = ClassicalRegister(num_control, name=\"out\")\n",
        "circuit = QuantumCircuit(control, target, output)\n",
        "\n",
        "# Initialize the target register to the state |1>\n",
        "circuit.x(num_control)\n",
        "\n",
        "# Add the Hadamard gates and controlled versions of the\n",
        "# multiplication gates\n",
        "for k, qubit in enumerate(control):\n",
        "    circuit.h(k)\n",
        "    b = b_list[k]\n",
        "    if b == 2:\n",
        "        circuit.compose(\n",
        "            M2mod15().control(), qubits=[qubit] + list(target), inplace=True\n",
        "        )\n",
        "    elif b == 4:\n",
        "        circuit.compose(\n",
        "            M4mod15().control(), qubits=[qubit] + list(target), inplace=True\n",
        "        )\n",
        "    else:\n",
        "        continue  # M1 is the identity operator\n",
        "\n",
        "# Apply the inverse QFT to the control register\n",
        "circuit.compose(QFT(num_control, inverse=True), qubits=control, inplace=True)\n",
        "\n",
        "# Measure the control register\n",
        "circuit.measure(control, output)\n",
        "\n",
        "circuit.draw(\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "78cf31ef-b2e7-4256-996a-5f89500efb52",
      "metadata": {},
      "source": [
        "Nótese que omitimos las operaciones de exponenciación modular controlada de los qubits de control restantes porque $M_1$ es el operador identidad.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "792a5cb7-2ac0-4964-b747-2780193a0401",
      "metadata": {},
      "source": [
        "Tenga en cuenta que más adelante en este tutorial, ejecutaremos este circuito en el backend `ibm_marrakesh` . Para ello, transpilamos el circuito según este backend específico e informamos de la profundidad del circuito y el recuento de puertas.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "95925dd5-7ba9-4746-b96e-ba50400fa5ac",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "2q-depth: 187\n",
            "2q-size: 260\n",
            "Operator counts: OrderedDict({'sx': 521, 'rz': 354, 'cz': 260, 'measure': 8, 'x': 4})\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/95925dd5-7ba9-4746-b96e-ba50400fa5ac-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 18,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "service = QiskitRuntimeService()\n",
        "backend = service.backend(\"ibm_marrakesh\")\n",
        "pm = generate_preset_pass_manager(optimization_level=2, backend=backend)\n",
        "\n",
        "transpiled_circuit = pm.run(circuit)\n",
        "\n",
        "print(\n",
        "    f\"2q-depth: {transpiled_circuit.depth(lambda x: x.operation.num_qubits==2)}\"\n",
        ")\n",
        "print(\n",
        "    f\"2q-size: {transpiled_circuit.size(lambda x: x.operation.num_qubits==2)}\"\n",
        ")\n",
        "print(f\"Operator counts: {transpiled_circuit.count_ops()}\")\n",
        "transpiled_circuit.draw(\n",
        "    output=\"mpl\", fold=-1, style=\"clifford\", idle_wires=False\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "f96a2a7e-0363-49d6-bc21-f0f207a84610",
      "metadata": {},
      "source": [
        "<span id=\"step-3-execute-using-qiskit-primitives\" />\n",
        "\n",
        "## Paso 3: Ejecutar utilizando Qiskit primitives\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a3fd5248-11f2-40b6-857a-6efc25e31bdb",
      "metadata": {},
      "source": [
        "En primer lugar, analizamos lo que obtendríamos teóricamente si ejecutáramos este circuito en un simulador ideal. A continuación, tenemos una serie de resultados de simulación del circuito anterior utilizando 1024 tomas. Como vemos, obtenemos una distribución aproximadamente uniforme en cuatro cadenas de bits sobre los qubits de control.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 19,
      "id": "c8720362-684f-4114-8abb-83d71d21c0df",
      "metadata": {},
      "outputs": [],
      "source": [
        "# Obtained from the simulator\n",
        "counts = {\"00000000\": 264, \"01000000\": 268, \"10000000\": 249, \"11000000\": 243}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 20,
      "id": "0d6d2702-02e4-47de-8f7e-0b256657ef0f",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/0d6d2702-02e4-47de-8f7e-0b256657ef0f-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 20,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "plot_histogram(counts)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "5ddd6952-fb4d-447c-8781-968ca5bb88a3",
      "metadata": {},
      "source": [
        "Midiendo los qubits de control, obtenemos una estimación de fase de ocho bits del operador $M_a$. Podemos convertir esta representación binaria a decimal para hallar la fase medida. Como podemos ver en el histograma anterior, se midieron cuatro cadenas de bits diferentes, y cada una de ellas corresponde a un valor de fase, como se indica a continuación.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 21,
      "id": "59ae21fc-e3d0-48e3-8cc6-dbce59c747ca",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "            Register Output           Phase\n",
            "0  00000000(bin) =   0(dec)    0/256 = 0.00\n",
            "1  01000000(bin) =  64(dec)   64/256 = 0.25\n",
            "2  10000000(bin) = 128(dec)  128/256 = 0.50\n",
            "3  11000000(bin) = 192(dec)  192/256 = 0.75\n"
          ]
        }
      ],
      "source": [
        "# Rows to be displayed in table\n",
        "rows = []\n",
        "# Corresponding phase of each bitstring\n",
        "measured_phases = []\n",
        "\n",
        "for output in counts:\n",
        "    decimal = int(output, 2)  # Convert bitstring to decimal\n",
        "    phase = decimal / (2**num_control)  # Find corresponding eigenvalue\n",
        "    measured_phases.append(phase)\n",
        "    # Add these values to the rows in our table:\n",
        "    rows.append(\n",
        "        [\n",
        "            f\"{output}(bin) = {decimal:>3}(dec)\",\n",
        "            f\"{decimal}/{2 ** num_control} = {phase:.2f}\",\n",
        "        ]\n",
        "    )\n",
        "\n",
        "# Print the rows in a table\n",
        "headers = [\"Register Output\", \"Phase\"]\n",
        "df = pd.DataFrame(rows, columns=headers)\n",
        "print(df)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "f0ca65bb-abb4-4f4b-bc32-4749f36e5780",
      "metadata": {},
      "source": [
        "Recordemos que la fase medida cualquiera corresponde a $\\theta = k / r$ donde $k$ se muestrea uniformemente al azar a partir de $\\{0, 1, \\dots, r-1 \\}$. Por lo tanto, podemos utilizar el algoritmo de fracciones continuas para intentar encontrar $k$ y el orden $r$. Python tiene esta funcionalidad incorporada. Podemos utilizar el módulo `fractions` para convertir un flotador en un objeto `Fraction` , por ejemplo:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "id": "999f04bc-f912-465f-9cba-90d08bda3758",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "Fraction(5998794703657501, 9007199254740992)"
            ]
          },
          "execution_count": 22,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "Fraction(0.666)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c5452b79-967f-4fe9-86d7-f9797987a8b5",
      "metadata": {},
      "source": [
        "Debido a que esto da fracciones que devuelven el resultado exactamente (en este caso, `0.6660000...`), esto puede dar resultados retorcidos como el de arriba. Podemos utilizar el método `.limit_denominator()` para obtener la fracción que más se parezca a nuestro flotante, con un denominador inferior a un determinado valor:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 23,
      "id": "1352aa9e-7c98-4862-8ac6-8d17fce61f3d",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "Fraction(2, 3)"
            ]
          },
          "execution_count": 23,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Get fraction that most closely resembles 0.666\n",
        "# with denominator < 15\n",
        "Fraction(0.666).limit_denominator(15)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ab1ce9d5-7a24-4717-a94b-576f8c940ec2",
      "metadata": {},
      "source": [
        "Esto es mucho más bonito. El orden (r) debe ser menor que N, por lo que fijaremos el denominador máximo en `15`:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 24,
      "id": "6c20bff2-29b7-45ea-b1ed-e802e0d1a0f9",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "   Phase Fraction  Guess for r\n",
            "0   0.00      0/1            1\n",
            "1   0.25      1/4            4\n",
            "2   0.50      1/2            2\n",
            "3   0.75      3/4            4\n"
          ]
        }
      ],
      "source": [
        "# Rows to be displayed in a table\n",
        "rows = []\n",
        "\n",
        "for phase in measured_phases:\n",
        "    frac = Fraction(phase).limit_denominator(15)\n",
        "    rows.append(\n",
        "        [phase, f\"{frac.numerator}/{frac.denominator}\", frac.denominator]\n",
        "    )\n",
        "\n",
        "# Print the rows in a table\n",
        "headers = [\"Phase\", \"Fraction\", \"Guess for r\"]\n",
        "df = pd.DataFrame(rows, columns=headers)\n",
        "print(df)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7590599a-bb21-4d9c-8a78-2ff406acc4af",
      "metadata": {},
      "source": [
        "Podemos ver que dos de los valores propios medidos nos dieron el resultado correcto: $r=4$, y podemos ver que el algoritmo de Shor para encontrar el orden tiene una posibilidad de fallar. Estos malos resultados se deben a que $k = 0$, o porque $k$ y $r$ no son coprimos - y en lugar de $r$, se nos da un factor de $r$. La solución más fácil para esto es simplemente repetir el experimento hasta que obtengamos un resultado satisfactorio para $r$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "174b446c-ca09-4816-8983-33f61eaecfa6",
      "metadata": {},
      "source": [
        "Hasta ahora, hemos implementado el problema de búsqueda de orden para $N=15$ con $a=2$ utilizando el circuito de estimación de fase en un simulador. El último paso del algoritmo de Shor consistirá en relacionar el problema de búsqueda de órdenes con el problema de factorización de enteros. Esta última parte del algoritmo es puramente clásica y puede resolverse en un ordenador clásico una vez obtenidas las medidas de fase en un ordenador cuántico. Por lo tanto, aplazamos la última parte del algoritmo hasta después de demostrar cómo podemos ejecutar el circuito de búsqueda de órdenes en hardware real.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7c826de7-0caa-4337-826e-3e81024ce125",
      "metadata": {},
      "source": [
        "<span id=\"hardware-runs\" />\n",
        "\n",
        "### Ejecución de hardware\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a8d6d6a9-9494-4733-89a2-790ac8adf929",
      "metadata": {},
      "source": [
        "Ahora podemos ejecutar el circuito de búsqueda de órdenes que transpilamos anteriormente para `ibm_marrakesh`. Aquí nos centramos en [el desacoplamiento dinámico](/docs/guides/error-mitigation-and-suppression-techniques#dynamical-decoupling) (DD) para la supresión de errores, y en [el giro de compuertas](/docs/guides/error-mitigation-and-suppression-techniques#pauli-twirling) para la mitigación de errores. La DD consiste en aplicar secuencias de pulsos de control sincronizados con precisión a un dispositivo cuántico, lo que elimina las interacciones ambientales no deseadas y la decoherencia. El giro de puertas, por su parte, aleatoriza puertas cuánticas específicas para transformar los errores coherentes en errores Pauli, que se acumulan linealmente en lugar de cuadráticamente. Ambas técnicas se combinan a menudo para mejorar la coherencia y fidelidad de los cálculos cuánticos.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "af61f1de-be24-4fc6-b2bb-af6845163c3c",
      "metadata": {},
      "outputs": [],
      "source": [
        "# Sampler primitive to obtain the probability distribution\n",
        "sampler = Sampler(backend)\n",
        "\n",
        "# Turn on dynamical decoupling with sequence XpXm\n",
        "sampler.options.dynamical_decoupling.enable = True\n",
        "sampler.options.dynamical_decoupling.sequence_type = \"XpXm\"\n",
        "# Enable gate twirling\n",
        "sampler.options.twirling.enable_gates = True\n",
        "\n",
        "# Assign tags before executing\n",
        "sampler.options.environment.job_tags = [\"TUT_SA\"]\n",
        "\n",
        "pub = transpiled_circuit\n",
        "job = sampler.run([pub], shots=1024)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 25,
      "id": "e5f53a20-40a4-440f-a6ea-349c20efbc95",
      "metadata": {},
      "outputs": [],
      "source": [
        "result = job.result()[0]\n",
        "counts = result.data[\"out\"].get_counts()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 26,
      "id": "559d7030-1f67-44e8-afa7-6afc7a334677",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/shors-algorithm/extracted-outputs/559d7030-1f67-44e8-afa7-6afc7a334677-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 26,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "plot_histogram(counts, figsize=(35, 5))"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "3614aa7e-0410-4e9e-b33a-3caf407e086a",
      "metadata": {},
      "source": [
        "Como podemos ver, obtuvimos las mismas cadenas de bits con los recuentos más altos. Como el hardware cuántico tiene ruido, hay algunas fugas a otras cadenas de bits, que podemos filtrar estadísticamente.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 27,
      "id": "f9cf9a07-5251-47bc-9713-c802f8f1a37c",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "{'00000000': 58, '01000000': 41, '11000000': 42, '10000000': 40}\n"
          ]
        }
      ],
      "source": [
        "# Dictionary of bitstrings and their counts to keep\n",
        "counts_keep = {}\n",
        "# Threshold to filter\n",
        "threshold = np.max(list(counts.values())) / 2\n",
        "\n",
        "for key, value in counts.items():\n",
        "    if value > threshold:\n",
        "        counts_keep[key] = value\n",
        "\n",
        "print(counts_keep)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c32cca67-9c9d-4373-9a0a-3174e884dd91",
      "metadata": {},
      "source": [
        "<span id=\"step-4-post-process-and-return-result-in-desired-classical-format\" />\n",
        "\n",
        "## Paso 4: Procesamiento posterior y devolución del resultado en el formato clásico deseado\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "5a7b116b-02d0-47dc-8cab-e792ba36d868",
      "metadata": {},
      "source": [
        "<span id=\"integer-factorization\" />\n",
        "\n",
        "### Factorización de números enteros\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "691df180-6ac3-49ac-ad4d-2d78348ceea7",
      "metadata": {},
      "source": [
        "Hasta ahora, hemos discutido cómo podemos implementar el problema de búsqueda de orden utilizando un circuito de estimación de fase. Ahora, conectamos el problema de búsqueda de orden a la factorización entera, lo que completa el algoritmo de Shor. Tenga en cuenta que esta parte del algoritmo es clásica.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "fe3472a3-f8cd-4db9-8178-c025ccf2ec55",
      "metadata": {},
      "source": [
        "Ahora lo demostraremos con nuestro ejemplo de $N = 15$ y $a = 2$. Recordemos que la fase que medimos es $k / r$, donde $a^r \\; (\\textrm{mod} \\; N) = 1$ y $k$ es un número entero aleatorio entre $0$ y $r - 1$. A partir de esta ecuación, tenemos $(a^r - 1) \\; (\\textrm{mod} \\; N) = 0,$, lo que significa que $N$ debe dividir a $a^r-1$. Si $r$ también es par, entonces podemos escribir $a^r -1 = (a^{r/2}-1)(a^{r/2}+1).$ Si $r$ no es par, no podemos ir más allá y debemos intentarlo de nuevo con un valor diferente para $a$; de lo contrario, hay una alta probabilidad de que el máximo común divisor de $N$ y $a^{r/2}-1$, o $a^{r/2}+1$ sea un factor propio de $N$.\n",
        "\n",
        "Dado que algunas ejecuciones del algoritmo fallarán estadísticamente, repetiremos este algoritmo hasta que se encuentre al menos un factor de $N$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "dda4ad11-ef2b-4bfb-ba3c-5cb1f5810871",
      "metadata": {},
      "source": [
        "La celda siguiente repite el algoritmo hasta que se encuentra al menos un factor de $N=15$. Utilizaremos los resultados de la ejecución de hardware anterior para adivinar la fase y el factor correspondiente en cada iteración.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 28,
      "id": "5d67cd71-8651-4a10-a913-a5bdaa6d6b38",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "\n",
            "ATTEMPT 0:\n",
            "Phase: theta = 0.0\n",
            "Order of 2 modulo 15 estimated as: r = 1\n",
            "\n",
            "ATTEMPT 1:\n",
            "Phase: theta = 0.25\n",
            "Order of 2 modulo 15 estimated as: r = 4\n",
            "*** Non-trivial factor found: 3 ***\n"
          ]
        }
      ],
      "source": [
        "a = 2\n",
        "N = 15\n",
        "\n",
        "FACTOR_FOUND = False\n",
        "num_attempt = 0\n",
        "\n",
        "while not FACTOR_FOUND:\n",
        "    print(f\"\\nATTEMPT {num_attempt}:\")\n",
        "    # Here, we get the bitstring by iterating over outcomes\n",
        "    # of a previous hardware run with multiple shots.\n",
        "    # Instead, we can also perform a single-shot measurement\n",
        "    # here in the loop.\n",
        "    bitstring = list(counts_keep.keys())[num_attempt]\n",
        "    num_attempt += 1\n",
        "    # Find the phase from measurement\n",
        "    decimal = int(bitstring, 2)\n",
        "    phase = decimal / (2**num_control)  # phase = k / r\n",
        "    print(f\"Phase: theta = {phase}\")\n",
        "\n",
        "    # Guess the order from phase\n",
        "    frac = Fraction(phase).limit_denominator(N)\n",
        "    r = frac.denominator  # order = r\n",
        "    print(f\"Order of {a} modulo {N} estimated as: r = {r}\")\n",
        "\n",
        "    if phase != 0:\n",
        "        # Guesses for factors are gcd(a^{r / 2} ± 1, 15)\n",
        "        if r % 2 == 0:\n",
        "            x = pow(a, r // 2, N) - 1\n",
        "            d = gcd(x, N)\n",
        "            if d > 1:\n",
        "                FACTOR_FOUND = True\n",
        "                print(f\"*** Non-trivial factor found: {x} ***\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "00f31b56-e1a9-479e-8a82-bb708263f1a3",
      "metadata": {},
      "source": [
        "<span id=\"discussion\" />\n",
        "\n",
        "## Debate\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "69189c94-6e78-410d-a365-8649d4c69163",
      "metadata": {},
      "source": [
        "<span id=\"related-work\" />\n",
        "\n",
        "### Trabajos relacionados\n",
        "\n",
        "En esta sección, analizamos otros trabajos que han demostrado el algoritmo de Shor en hardware real.\n",
        "\n",
        "El trabajo seminal [\\[3\\]](#references) de IBM® demostró por primera vez el algoritmo de Shor, factorizando el número 15 en sus factores primos 3 y 5 utilizando un ordenador cuántico de resonancia magnética nuclear (RMN) de siete qubits. Otro experimento [\\[4\\]](#references) factorizó 15 utilizando qubits fotónicos. Empleando un único qubit reciclado varias veces y codificando el registro de trabajo en estados de mayor dimensión, los investigadores redujeron el número necesario de qubits a un tercio del del protocolo estándar, utilizando un algoritmo compilado de dos fotones. Un artículo significativo en la demostración del algoritmo de Shor es [\\[5\\]](#references), que utiliza la técnica de estimación de fase iterativa de Kitaev [\\[8\\]](#references) para reducir el requisito de qubits del algoritmo. Los autores utilizaron siete qubits de control y cuatro qubits de caché, junto con la implementación de multiplicadores modulares. Esta implementación, sin embargo, requiere mediciones en mitad del circuito con operaciones feed-forward y reciclado de qubits con operaciones de reinicio. Esta demostración se realizó en un ordenador cuántico de trampa de iones.\n",
        "\n",
        "Un trabajo más reciente [\\[6\\]](#references) se centró en la factorización de 15, 21 y 35 en el hardware IBM Quantum®. Al igual que en trabajos anteriores, los investigadores utilizaron una versión compilada del algoritmo que empleaba una transformada cuántica de Fourier semiclásica como la propuesta por Kitaev para minimizar el número de qubits físicos y puertas. Un trabajo más reciente [\\[7\\]](#references) también realizó una demostración de prueba de concepto para factorizar el entero 21. Esta demostración también implicó el uso de una versión compilada de la rutina de estimación de fase cuántica, y se basó en la demostración anterior de [\\[4\\]](#references). Los autores fueron más allá de este trabajo utilizando una configuración de puertas Toffoli aproximadas con desplazamientos de fase residuales. El algoritmo se implementó en procesadores cuánticos IBM utilizando sólo cinco qubits, y se verificó con éxito la presencia de entrelazamiento entre los qubits de control y registro.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2dab9ac1-d3e4-4b4e-babe-472d69b026ab",
      "metadata": {},
      "source": [
        "<span id=\"scaling-of-the-algorithm\" />\n",
        "\n",
        "### Escalado del algoritmo\n",
        "\n",
        "El cifrado RSA suele requerir claves del orden de 2048 a 4096 bits. Intentar factorizar un número de 2048 bits con el algoritmo de Shor dará como resultado un circuito cuántico con millones de qubits, incluida la sobrecarga de corrección de errores y una profundidad de circuito del orden de mil millones, que está más allá de los límites de ejecución del hardware cuántico actual. Por lo tanto, el algoritmo de Shor necesitará métodos de construcción de circuitos optimizados o una corrección de errores cuántica robusta para ser viable en la práctica a la hora de romper sistemas criptográficos modernos. Le remitimos a [\\[9\\]](#references) para una discusión más detallada sobre la estimación de recursos para el algoritmo de Shor.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "cb53f62c-e560-43cf-a614-d31f266a569c",
      "metadata": {},
      "source": [
        "<span id=\"challenge\" />\n",
        "\n",
        "## Reto\n",
        "\n",
        "¡Enhorabuena por terminar el tutorial! Ahora es un buen momento para poner a prueba tus conocimientos. ¿Podrías intentar construir el circuito para factorizar 21? Puede seleccionar un $a$ de su elección. Tendrás que decidir la precisión de bits del algoritmo para elegir el número de qubits, y tendrás que diseñar los operadores de exponenciación modular $M_a$. Te animamos a que pruebes esto por ti mismo, y luego leas sobre las metodologías mostradas en la Fig. 9 de [\\[6\\]](#references) y la Fig. 2 de [\\[7\\]](#references).\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "09d347b8-51a4-4eb5-8422-76508f7ba1ad",
      "metadata": {},
      "outputs": [],
      "source": [
        "def M_a_mod21():\n",
        "    \"\"\"\n",
        "    M_a (mod 21)\n",
        "    \"\"\"\n",
        "\n",
        "    # Your code here\n",
        "    pass"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "dc16f509-7441-4221-8113-e09b118397a0",
      "metadata": {},
      "source": [
        "<span id=\"references\" />\n",
        "\n",
        "## Referencias\n",
        "\n",
        "1. Shor, Peter W. \"[Algoritmos de tiempo polinómico para la factorización de primos y logaritmos discretos en un ordenador cuántico](https://epubs.siam.org/doi/abs/10.1137/S0036144598347011) \" Revisión SIAM 41.2 (1999): 303-332.\n",
        "2. IBM Quantum Curso [«Fundamentos de los algoritmos cuánticos»,](/learning/courses/fundamentals-of-quantum-algorithms/phase-estimation-and-factoring/introduction) impartido por el Dr. John Watrous.\n",
        "3. Vandersypen, Lieven MK, et al. \"[Experimental realization of Shor's quantum factoring algorithm using nuclear magnetic resonance](https://www.nature.com/articles/414883a) \" Nature 414.6866 (2001): 883-887.\n",
        "4. Martín-López, Enrique, et al. \"[Experimental realization of Shor's quantum factoring algorithm using qubit](https://www.nature.com/articles/nphoton.2012.259) recycling\" Nature photonics 6.11 (2012): 773-776.\n",
        "5. Monz, Thomas, et al. \"[Realización de un algoritmo Shor escalable](https://www.science.org/doi/full/10.1126/science.aad9480) \" Science 351.6277 (2016): 1068-1070.\n",
        "6. Amico, Mirko, Zain H. Saleem y Muir Kumph. \"[Estudio experimental del algoritmo de factorización de Shor utilizando la experiencia IBM Q](https://journals.aps.org/pra/abstract/10.1103/PhysRevA.100.012305) \" Physical Review A 100.1 (2019): 012305.\n",
        "7. Skosana, Unathi y Mark Tame. \"[Demostración del algoritmo de factorización de Shor para N=21 en procesadores cuánticos IBM](https://www.nature.com/articles/s41598-021-95973-w) \" Informes científicos 11.1 (2021): 16599.\n",
        "8. Kitaev, A. Yu. \"[Medidas cuánticas y el problema del estabilizador abeliano](https://arxiv.org/abs/quant-ph/9511026) \" arXiv preprint quant-ph/9511026 (1995).\n",
        "9. Gidney, Craig, y Martin Ekerå. \"[Cómo factorizar enteros RSA de 2048 bits en 8 horas usando 20 millones de qubits ruidosos](https://doi.org/10.22331/q-2021-04-15-433) \" Quantum 5 (2021): 433.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "id": "a1b8767d",
      "source": "© IBM Corp., 2017-2026"
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3"
    },
    "hours": 1,
    "qpuSeconds": 3
  },
  "nbformat": 4,
  "nbformat_minor": 4
}