su(n)_deepth_new_circ_with_param_save.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-FDG6kZmtnfl"
      },
      "source": [
        "# Imports and Installation"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "kLpKQZYrqm6r",
        "outputId": "79cefc9e-5aa5-4a23-eaca-8c0d9d119f10",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Collecting torch\n",
            "  Using cached torch-2.3.1-cp311-none-macosx_11_0_arm64.whl.metadata (26 kB)\n",
            "Collecting filelock (from torch)\n",
            "  Downloading filelock-3.15.1-py3-none-any.whl.metadata (2.8 kB)\n",
            "Requirement already satisfied: typing-extensions>=4.8.0 in ./venv/lib/python3.11/site-packages (from torch) (4.12.2)\n",
            "Collecting sympy (from torch)\n",
            "  Using cached sympy-1.12.1-py3-none-any.whl.metadata (12 kB)\n",
            "Collecting networkx (from torch)\n",
            "  Using cached networkx-3.3-py3-none-any.whl.metadata (5.1 kB)\n",
            "Collecting jinja2 (from torch)\n",
            "  Using cached jinja2-3.1.4-py3-none-any.whl.metadata (2.6 kB)\n",
            "Collecting fsspec (from torch)\n",
            "  Using cached fsspec-2024.6.0-py3-none-any.whl.metadata (11 kB)\n",
            "Collecting MarkupSafe>=2.0 (from jinja2->torch)\n",
            "  Using cached MarkupSafe-2.1.5-cp311-cp311-macosx_10_9_universal2.whl.metadata (3.0 kB)\n",
            "Collecting mpmath<1.4.0,>=1.1.0 (from sympy->torch)\n",
            "  Using cached mpmath-1.3.0-py3-none-any.whl.metadata (8.6 kB)\n",
            "Using cached torch-2.3.1-cp311-none-macosx_11_0_arm64.whl (61.0 MB)\n",
            "Downloading filelock-3.15.1-py3-none-any.whl (15 kB)\n",
            "Using cached fsspec-2024.6.0-py3-none-any.whl (176 kB)\n",
            "Using cached jinja2-3.1.4-py3-none-any.whl (133 kB)\n",
            "Using cached networkx-3.3-py3-none-any.whl (1.7 MB)\n",
            "Using cached sympy-1.12.1-py3-none-any.whl (5.7 MB)\n",
            "Using cached MarkupSafe-2.1.5-cp311-cp311-macosx_10_9_universal2.whl (18 kB)\n",
            "Using cached mpmath-1.3.0-py3-none-any.whl (536 kB)\n",
            "Installing collected packages: mpmath, sympy, networkx, MarkupSafe, fsspec, filelock, jinja2, torch\n",
            "Successfully installed MarkupSafe-2.1.5 filelock-3.15.1 fsspec-2024.6.0 jinja2-3.1.4 mpmath-1.3.0 networkx-3.3 sympy-1.12.1 torch-2.3.1\n",
            "Note: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "pip install torch"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "hopcJg9Ev6ud",
        "outputId": "5a129f3e-ec94-4d19-a892-8efa1f8a9c9c",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Collecting qiskit_algorithms\n",
            "  Using cached qiskit_algorithms-0.3.0-py3-none-any.whl.metadata (4.2 kB)\n",
            "Collecting qiskit>=0.44 (from qiskit_algorithms)\n",
            "  Using cached qiskit-1.1.0-cp38-abi3-macosx_11_0_arm64.whl.metadata (12 kB)\n",
            "Collecting scipy>=1.4 (from qiskit_algorithms)\n",
            "  Using cached scipy-1.13.1-cp311-cp311-macosx_12_0_arm64.whl.metadata (60 kB)\n",
            "Collecting numpy>=1.17 (from qiskit_algorithms)\n",
            "  Downloading numpy-2.0.0-cp311-cp311-macosx_11_0_arm64.whl.metadata (60 kB)\n",
            "\u001b[2K     \u001b[90m━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━\u001b[0m \u001b[32m60.9/60.9 kB\u001b[0m \u001b[31m816.6 kB/s\u001b[0m eta \u001b[36m0:00:00\u001b[0ma \u001b[36m0:00:01\u001b[0m\n",
            "\u001b[?25hCollecting rustworkx>=0.14.0 (from qiskit>=0.44->qiskit_algorithms)\n",
            "  Using cached rustworkx-0.14.2-cp311-cp311-macosx_11_0_arm64.whl.metadata (10.0 kB)\n",
            "Requirement already satisfied: sympy>=1.3 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit_algorithms) (1.12.1)\n",
            "Collecting dill>=0.3 (from qiskit>=0.44->qiskit_algorithms)\n",
            "  Using cached dill-0.3.8-py3-none-any.whl.metadata (10 kB)\n",
            "Requirement already satisfied: python-dateutil>=2.8.0 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit_algorithms) (2.9.0.post0)\n",
            "Collecting stevedore>=3.0.0 (from qiskit>=0.44->qiskit_algorithms)\n",
            "  Using cached stevedore-5.2.0-py3-none-any.whl.metadata (2.3 kB)\n",
            "Requirement already satisfied: typing-extensions in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit_algorithms) (4.12.2)\n",
            "Collecting symengine>=0.11 (from qiskit>=0.44->qiskit_algorithms)\n",
            "  Using cached symengine-0.11.0-cp311-cp311-macosx_11_0_arm64.whl.metadata (1.2 kB)\n",
            "Requirement already satisfied: six>=1.5 in ./venv/lib/python3.11/site-packages (from python-dateutil>=2.8.0->qiskit>=0.44->qiskit_algorithms) (1.16.0)\n",
            "Collecting numpy>=1.17 (from qiskit_algorithms)\n",
            "  Using cached numpy-1.26.4-cp311-cp311-macosx_11_0_arm64.whl.metadata (114 kB)\n",
            "Collecting pbr!=2.1.0,>=2.0.0 (from stevedore>=3.0.0->qiskit>=0.44->qiskit_algorithms)\n",
            "  Using cached pbr-6.0.0-py2.py3-none-any.whl.metadata (1.3 kB)\n",
            "Requirement already satisfied: mpmath<1.4.0,>=1.1.0 in ./venv/lib/python3.11/site-packages (from sympy>=1.3->qiskit>=0.44->qiskit_algorithms) (1.3.0)\n",
            "Using cached qiskit_algorithms-0.3.0-py3-none-any.whl (308 kB)\n",
            "Using cached qiskit-1.1.0-cp38-abi3-macosx_11_0_arm64.whl (4.0 MB)\n",
            "Using cached scipy-1.13.1-cp311-cp311-macosx_12_0_arm64.whl (30.3 MB)\n",
            "Using cached dill-0.3.8-py3-none-any.whl (116 kB)\n",
            "Using cached rustworkx-0.14.2-cp311-cp311-macosx_11_0_arm64.whl (1.7 MB)\n",
            "Using cached numpy-1.26.4-cp311-cp311-macosx_11_0_arm64.whl (14.0 MB)\n",
            "Using cached stevedore-5.2.0-py3-none-any.whl (49 kB)\n",
            "Using cached symengine-0.11.0-cp311-cp311-macosx_11_0_arm64.whl (20.6 MB)\n",
            "Using cached pbr-6.0.0-py2.py3-none-any.whl (107 kB)\n",
            "Installing collected packages: symengine, pbr, numpy, dill, stevedore, scipy, rustworkx, qiskit, qiskit_algorithms\n",
            "Successfully installed dill-0.3.8 numpy-1.26.4 pbr-6.0.0 qiskit-1.1.0 qiskit_algorithms-0.3.0 rustworkx-0.14.2 scipy-1.13.1 stevedore-5.2.0 symengine-0.11.0\n",
            "Note: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "pip install qiskit_algorithms"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "reRQoomHhU4H",
        "outputId": "540b9734-f7b6-4d4c-a959-099624de9cbd",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "\u001b[33mWARNING: Skipping qiskit_nature as it is not installed.\u001b[0m\u001b[33m\n",
            "\u001b[0mNote: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "pip uninstall qiskit_nature --yes"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "nZsZ-6tvv6ud",
        "outputId": "0e735b3a-4fc6-4eb7-de6e-a55b4c1bbc0a",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "\u001b[33mWARNING: Skipping pyscf as it is not installed.\u001b[0m\u001b[33m\n",
            "\u001b[0mNote: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "pip uninstall pyscf --yes"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "SiIbXXQ3hYCN",
        "outputId": "f1554488-10d9-4633-e04a-a78a9ad25e68"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "\u001b[33mWARNING: Package(s) not found: qiskit_nature\u001b[0m\u001b[33m\n",
            "\u001b[0m\u001b[33mWARNING: Package(s) not found: pyscf\u001b[0m\u001b[33m\n",
            "\u001b[0m"
          ]
        }
      ],
      "source": [
        "!pip show qiskit_nature\n",
        "!pip show pyscf"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "05U9Y2AfhcU5",
        "outputId": "a2bfb73a-ff6c-4740-c5d4-eeb43023f8d1",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Collecting qiskit-nature\n",
            "  Using cached qiskit_nature-0.7.2-py3-none-any.whl.metadata (8.0 kB)\n",
            "Requirement already satisfied: qiskit>=0.44 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (1.1.0)\n",
            "Requirement already satisfied: qiskit-algorithms>=0.2.1 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (0.3.0)\n",
            "Requirement already satisfied: scipy>=1.4 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (1.13.1)\n",
            "Requirement already satisfied: numpy>=1.17 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (1.26.4)\n",
            "Requirement already satisfied: psutil>=5 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (6.0.0)\n",
            "Requirement already satisfied: setuptools>=40.1.0 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (70.0.0)\n",
            "Requirement already satisfied: typing-extensions in ./venv/lib/python3.11/site-packages (from qiskit-nature) (4.12.2)\n",
            "Collecting h5py (from qiskit-nature)\n",
            "  Using cached h5py-3.11.0-cp311-cp311-macosx_11_0_arm64.whl.metadata (2.5 kB)\n",
            "Requirement already satisfied: rustworkx>=0.12 in ./venv/lib/python3.11/site-packages (from qiskit-nature) (0.14.2)\n",
            "Requirement already satisfied: sympy>=1.3 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit-nature) (1.12.1)\n",
            "Requirement already satisfied: dill>=0.3 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit-nature) (0.3.8)\n",
            "Requirement already satisfied: python-dateutil>=2.8.0 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit-nature) (2.9.0.post0)\n",
            "Requirement already satisfied: stevedore>=3.0.0 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit-nature) (5.2.0)\n",
            "Requirement already satisfied: symengine>=0.11 in ./venv/lib/python3.11/site-packages (from qiskit>=0.44->qiskit-nature) (0.11.0)\n",
            "Requirement already satisfied: six>=1.5 in ./venv/lib/python3.11/site-packages (from python-dateutil>=2.8.0->qiskit>=0.44->qiskit-nature) (1.16.0)\n",
            "Requirement already satisfied: pbr!=2.1.0,>=2.0.0 in ./venv/lib/python3.11/site-packages (from stevedore>=3.0.0->qiskit>=0.44->qiskit-nature) (6.0.0)\n",
            "Requirement already satisfied: mpmath<1.4.0,>=1.1.0 in ./venv/lib/python3.11/site-packages (from sympy>=1.3->qiskit>=0.44->qiskit-nature) (1.3.0)\n",
            "Using cached qiskit_nature-0.7.2-py3-none-any.whl (2.2 MB)\n",
            "Using cached h5py-3.11.0-cp311-cp311-macosx_11_0_arm64.whl (2.9 MB)\n",
            "Installing collected packages: h5py, qiskit-nature\n",
            "Successfully installed h5py-3.11.0 qiskit-nature-0.7.2\n",
            "Note: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "pip install qiskit-nature"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "iAl7IvGOv6ue",
        "outputId": "01c0653c-9456-47e3-ae68-409115e72094",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Collecting pyscf\n",
            "  Using cached pyscf-2.6.1-py3-none-macosx_11_0_arm64.whl.metadata (6.3 kB)\n",
            "Requirement already satisfied: numpy!=1.16,!=1.17,>=1.13 in ./venv/lib/python3.11/site-packages (from pyscf) (1.26.4)\n",
            "Requirement already satisfied: scipy!=1.5.0,!=1.5.1 in ./venv/lib/python3.11/site-packages (from pyscf) (1.13.1)\n",
            "Requirement already satisfied: h5py>=2.7 in ./venv/lib/python3.11/site-packages (from pyscf) (3.11.0)\n",
            "Requirement already satisfied: setuptools in ./venv/lib/python3.11/site-packages (from pyscf) (70.0.0)\n",
            "Using cached pyscf-2.6.1-py3-none-macosx_11_0_arm64.whl (32.8 MB)\n",
            "Installing collected packages: pyscf\n",
            "Successfully installed pyscf-2.6.1\n",
            "Note: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "pip install pyscf"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "metadata": {
        "id": "0HIrwfxa5yNd",
        "tags": []
      },
      "outputs": [],
      "source": [
        "import qiskit\n",
        "from qiskit.quantum_info import SparsePauliOp, Operator\n",
        "from qiskit.circuit.library import UnitaryGate,PauliEvolutionGate\n",
        "import numpy as np\n",
        "from functools import lru_cache, reduce\n",
        "from itertools import product\n",
        "import torch\n",
        "from qiskit.circuit import Parameter, Gate\n",
        "from qiskit import QuantumCircuit\n",
        "from scipy.optimize import minimize\n",
        "from qiskit.circuit import ParameterVector\n",
        "from qiskit_algorithms.gradients import ParamShiftEstimatorGradient\n",
        "from scipy.linalg import expm\n",
        "from tqdm import tqdm\n",
        "from qiskit.primitives import Estimator\n",
        "from qiskit.circuit.library import GlobalPhaseGate, HamiltonianGate\n",
        "from qiskit_nature.second_q.operators import FermionicOp\n",
        "import marshal\n",
        "import csv"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XFh_PKHNtszy"
      },
      "source": [
        "# Utils"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "metadata": {
        "id": "VeQUiqiz6b20",
        "tags": []
      },
      "outputs": [],
      "source": [
        "_pauli_matrices = np.array(\n",
        "    [[[1, 0], [0, 1]], [[0, 1], [1, 0]], [[0, -1j], [1j, 0]], [[1, 0], [0, -1]]]\n",
        ")\n",
        "\"\"\"Single-qubit Paulis:    I                 X                   Y                  Z\"\"\"\n",
        "\n",
        "_pauli_letters = \"IXYZ\"\n",
        "\"\"\"Single-qubit Pauli letters that make up Pauli words.\"\"\"\n",
        "\n",
        "@lru_cache\n",
        "def pauli_basis_matrices(num_qubits):\n",
        "  if num_qubits < 1:\n",
        "        raise ValueError(f\"Require at least one wire, got {num_qubits}.\")\n",
        "  if num_qubits > 7:\n",
        "      raise ValueError(\n",
        "          f\"Creating the Pauli basis tensor for more than 7 wires (got {num_qubits}) is deactivated.\"\n",
        "      )\n",
        "  return reduce(np.kron, (_pauli_matrices for _ in range(num_qubits)))[1:]\n",
        "\n",
        "def detach_numpy(tensor):\n",
        "    tensor = tensor.detach().cpu()\n",
        "    if torch._C._functorch.is_gradtrackingtensor(tensor):\n",
        "        tensor = torch._C._functorch.get_unwrapped(tensor)\n",
        "        return np.array(tensor.storage().tolist()).reshape(tensor.shape)\n",
        "    return tensor.numpy()\n",
        "\n",
        "def pauli_decompose(matrix, num_qubits):\n",
        "  basis = pauli_basis_matrices(num_qubits)\n",
        "  coefficients = np.tensordot(torch.tensor(basis), matrix, [[1, 2], [-1, -2]]) * (1 / 2**num_qubits)\n",
        "  return coefficients\n",
        "\n",
        "@lru_cache\n",
        "def pauli_basis_strings(num_wires):\n",
        "  return [\"\".join(letters) for letters in product(_pauli_letters, repeat=num_wires)][1:]\n",
        "\n",
        "def get_file_path(file_name, data_directory):\n",
        "    if file_name:\n",
        "        if file_name[-5:] != '.data':\n",
        "            file_name = file_name + \".data\"\n",
        "    else:\n",
        "        raise (\"File name is not provided.\")\n",
        "\n",
        "    file_path = data_directory + '/' + file_name\n",
        "\n",
        "    return file_path\n",
        "sign_arr = ['+', '-']\n",
        "def load_operator(file_name=None, data_directory=None, plain_text=False, nqubits=4):\n",
        "    \"\"\"Load FermionOperator or QubitOperator from file.\n",
        "\n",
        "    Args:\n",
        "        file_name: The name of the saved file.\n",
        "        data_directory: Optional data directory to change from default data\n",
        "                        directory specified in config file.\n",
        "        plain_text: Whether the input file is plain text\n",
        "\n",
        "    Returns:\n",
        "        operator: The stored FermionOperator, BosonOperator,\n",
        "            QuadOperator, or QubitOperator\n",
        "\n",
        "    Raises:\n",
        "        TypeError: Operator of invalid type.\n",
        "    \"\"\"\n",
        "    file_path = get_file_path(file_name, data_directory)\n",
        "    if plain_text:\n",
        "        with open(file_path, 'r') as f:\n",
        "            data = f.read()\n",
        "            operator_type, operator_terms = data.split(\":\\n\")\n",
        "            operator = FermionicOp(operator_terms)\n",
        "    else:\n",
        "        with open(file_path, 'rb') as f:\n",
        "            data = marshal.load(f)\n",
        "            operator_type = data[0]\n",
        "            operator_terms = data[1]\n",
        "        if operator_type == 'FermionOperator':\n",
        "            max_term = 0\n",
        "            op_dict = {}\n",
        "            for term, coeff in operator_terms.items():\n",
        "                op = \"\"\n",
        "                for term_i in term:\n",
        "                  val, sign = term_i[0], term_i[1]\n",
        "                  op += sign_arr[sign] + \"_\" + str((nqubits-1-val)) + \" \"\n",
        "                  max_term = max(val, max_term)\n",
        "                op_dict[op] =  coeff\n",
        "\n",
        "            operator = FermionicOp(op_dict, max_term+1, copy=False)\n",
        "    return operator\n",
        "\n",
        "def problem_hamiltonian(n_qubits: int, seed: int, hamiltonian_directory: str):\n",
        "    print(f\"{str(n_qubits).zfill(2)}qubits_{str(seed).zfill(2)}\")\n",
        "    ham = load_operator(\n",
        "        file_name=f\"{str(n_qubits).zfill(2)}qubits_{str(seed).zfill(2)}\",\n",
        "        data_directory=hamiltonian_directory,\n",
        "        plain_text=False,\n",
        "        nqubits=n_qubits,\n",
        "    )\n",
        "    return ham"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OFMQ5nwhtwST"
      },
      "source": [
        "# SU Implementation"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "metadata": {
        "id": "itIydW0THC9w",
        "tags": []
      },
      "outputs": [],
      "source": [
        "class SU(Gate):\n",
        "    def __init__(self, theta:np.ndarray|list[Parameter], nqubits, label=None):\n",
        "        if isinstance(theta, np.ndarray):\n",
        "            self.theta = torch.from_numpy(np.array(theta))\n",
        "        super().__init__('SU', nqubits, theta, label=label)\n",
        "\n",
        "    def _define(self):\n",
        "        qc = QuantumCircuit(self.num_qubits)\n",
        "        qc.unitary(self.to_matrix(), range(self.num_qubits))\n",
        "\n",
        "        self.definition = qc\n",
        "\n",
        "    def to_matrix(self):\n",
        "        num_qubits = self.num_qubits\n",
        "        mat = self.compute_matrix(self.theta, num_qubits)\n",
        "        return detach_numpy(mat)\n",
        "\n",
        "    @staticmethod\n",
        "    def compute_matrix(theta, num_qubits):\n",
        "      if num_qubits > 5:\n",
        "        matrices = product(_pauli_matrices, repeat=num_qubits)\n",
        "        # Drop the identity from the generator of matrices\n",
        "        _ = next(matrices)\n",
        "        A = sum(t * reduce(np.kron, letters) for t, letters in zip(theta, matrices))\n",
        "      else:\n",
        "        paulis = torch.from_numpy(pauli_basis_matrices(num_qubits))\n",
        "        theta = theta.type(torch.complex128)\n",
        "        A = torch.tensordot(theta, paulis, dims=[[-1], [0]])\n",
        "\n",
        "      return torch.matrix_exp(torch.mul(A, 1j))\n",
        "\n",
        "    def get_one_parameter_generators(self):\n",
        "      theta = self.theta\n",
        "      num_qubits = self.num_qubits\n",
        "\n",
        "      def split_matrix(theta):\n",
        "        mat = self.compute_matrix(theta, num_qubits)\n",
        "        return torch.real(mat), torch.imag(mat)\n",
        "\n",
        "      rjac, ijac = torch.autograd.functional.jacobian(split_matrix, theta)\n",
        "      jac = rjac + 1j * ijac\n",
        "\n",
        "      U_dagger = self.compute_matrix(torch.mul(theta, -1), num_qubits)\n",
        "\n",
        "      # U_dagger = U_dagger.type(torch.complex128)\n",
        "      # jac = detach_numpy(jac)\n",
        "      jac = jac.type(torch.complex128)\n",
        "      U_dagger = U_dagger.type(torch.complex128)\n",
        "      ujac = torch.tensordot(U_dagger, jac, dims=([1], [0]))\n",
        "      ujac = ujac.permute(2, 0, 1)\n",
        "      return detach_numpy(ujac)\n",
        "\n",
        "      # ujac = np.tensordot(U_dagger, jac, axes=([1], [0]))\n",
        "      # ujac = np.transpose(ujac, [2, 0, 1])\n",
        "      # return ujac\n",
        "\n",
        "    def get_one_parameter_coeffs(self):\n",
        "      num_qubits = self.num_qubits\n",
        "      theta = self.theta\n",
        "      generators = self.get_one_parameter_generators()\n",
        "      Omega = pauli_decompose(generators, num_qubits)\n",
        "\n",
        "      # skew matrix condition on Omega\n",
        "      return Omega\n",
        "\n",
        "def su4(params, nqubits):\n",
        "    return SU(params, nqubits)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "metadata": {
        "id": "2dy0o0Ep6czj",
        "tags": []
      },
      "outputs": [],
      "source": [
        "def apply_su4_in_layers(theta, depth, nqubits, fix_pars_per_layer=False):\n",
        "    \"\"\"Apply a (callable) operation in layers.\n",
        "\n",
        "    Args:\n",
        "        theta (tensor_like): The arguments passed to the operations. The expected shape\n",
        "            is ``(depth, k, num_params_op)``, where ``k`` is determined by ``fix_pars_per_layer``\n",
        "            and ``num_params_op`` is the number of paramters each operation takes.\n",
        "        Op (callable): The operation to apply\n",
        "        depth (int): The number of layers to apply\n",
        "        fix_pars_per_layer (bool): Whether or not all operations applied in parallel share the\n",
        "            same parameters. If ``True``, the dimension ``k`` in the shape of ``theta`` is set to\n",
        "            two, otherwise it is set to ``nqubits``.\n",
        "\n",
        "    \"\"\"\n",
        "    qc = qiskit.QuantumCircuit(nqubits)\n",
        "    qubits = list(range(nqubits))\n",
        "\n",
        "    if fix_pars_per_layer:\n",
        "        size = 2\n",
        "    else:\n",
        "        size = nqubits\n",
        "\n",
        "    for d in range(depth):\n",
        "        # Even-odd qubit pairs\n",
        "        idx = 0\n",
        "        for q in qubits[0::2]:\n",
        "            su = su4(theta[(d*size + idx)*15:(d*size + idx +1)*15], 2)\n",
        "            qc = qc.compose(su, [q, (q+1+d)%nqubits])\n",
        "            if q == nqubits - 2:\n",
        "                idx += 1\n",
        "            elif not fix_pars_per_layer:\n",
        "                idx += 1\n",
        "\n",
        "        # Odd-even qubit pairs\n",
        "        for q in qubits[1::2]:\n",
        "            su = su4(theta[(d*size + idx)*15:(d*size + idx +1)*15],2)\n",
        "            qc = qc.compose(su, [q, (q+1+d)%nqubits])\n",
        "            if not fix_pars_per_layer:\n",
        "                idx += 1\n",
        "    return qc"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "g7UmA18Cv6ug"
      },
      "source": [
        "# Grad function"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OONj0ACAt05a"
      },
      "source": [
        "# Read observable"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 12,
      "metadata": {
        "id": "qkROisbTwaT1",
        "tags": []
      },
      "outputs": [],
      "source": [
        "from qiskit_nature.second_q.hamiltonians import FermiHubbardModel, HeisenbergModel\n",
        "from qiskit_nature.second_q.hamiltonians.lattices import BoundaryCondition,LineLattice\n",
        "from qiskit_nature.second_q.mappers import JordanWignerMapper\n",
        "from qiskit_nature.second_q.hamiltonians import QuadraticHamiltonian\n",
        "\n",
        "#read observable from file\n",
        "def make_observable(num_qubits, seed=0, mat=True):\n",
        "  hamiltonian_directory = \".\"\n",
        "  print(\"MN = \", num_qubits)\n",
        "  ham = problem_hamiltonian(num_qubits, seed, hamiltonian_directory)\n",
        "  H = JordanWignerMapper().map(ham)\n",
        "  if mat:\n",
        "    H_jw = H.to_matrix()\n",
        "    return H, H_jw\n",
        "  else:\n",
        "    return H"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "SJt0SicYt4y-"
      },
      "source": [
        "# VQE"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "metadata": {
        "id": "c6pkPAC5kdnZ",
        "tags": []
      },
      "outputs": [],
      "source": [
        "def print_cpu_usage():\n",
        "  import psutil\n",
        "  split_bar = '='*20\n",
        "  memory_info = psutil.virtual_memory()._asdict()\n",
        "  print(f\"{split_bar} Memory Usage {split_bar}\")\n",
        "  for k,v in memory_info.items():\n",
        "    print(k, v)\n",
        "  print(f\"{split_bar} CPU Usage {split_bar}\")\n",
        "  print(f\"CPU percent: {psutil.cpu_percent()}%\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 14,
      "metadata": {
        "id": "Z07U05xVt9lm",
        "tags": []
      },
      "outputs": [],
      "source": [
        "from qiskit_nature.second_q.circuit.library import HartreeFock, UCCSD\n",
        "from qiskit_algorithms.gradients import FiniteDiffEstimatorGradient\n",
        "from qiskit_algorithms import VQE\n",
        "import math\n",
        "\n",
        "def generate_init_params(num_qubits, num_layers):\n",
        "  \"\"\"\n",
        "  Generate initial parameters for the ansatz.\n",
        "  Returns:\n",
        "      list[float]: The initial parameters for the ansatz.\n",
        "  \"\"\"\n",
        "  NUM_WIRES = num_qubits\n",
        "  LAYERS_PER_BLOCK = 15*num_qubits\n",
        "\n",
        "  init_params = np.random.uniform(\n",
        "      0, 2 * np.pi, LAYERS_PER_BLOCK * num_layers\n",
        "  )\n",
        "  for i in range(num_layers):\n",
        "      layer_param_sum = np.sum(\n",
        "          init_params[i * LAYERS_PER_BLOCK : (i + 1) * LAYERS_PER_BLOCK]\n",
        "      )\n",
        "      init_params[i * LAYERS_PER_BLOCK : (i + 1) * LAYERS_PER_BLOCK] *= (\n",
        "          math.pi / (2 * NUM_WIRES) / layer_param_sum\n",
        "      )\n",
        "\n",
        "  return init_params\n",
        "\n",
        "\n",
        "def vqe(\n",
        "    num_qubits,\n",
        "    observable,\n",
        "    depth,\n",
        "    optimizer,\n",
        "    max_steps,\n",
        "    init_param=None):\n",
        "\n",
        "    # n = 4**num_SUgate_qubits - 1 = 15\n",
        "    param_len = depth*num_qubits*15\n",
        "\n",
        "    estimator = Estimator()\n",
        "\n",
        "    # List to store the minimized values after every iteration\n",
        "    minimized_values = []\n",
        "\n",
        "    def cost_function(params):\n",
        "        \"\"\"Cost function of the VQE.\"\"\"\n",
        "\n",
        "        # Initialize Hartree-Fock state\n",
        "        num_orbitals = num_qubits // 2\n",
        "        num_particles = num_qubits // 2\n",
        "        num_alpha = num_particles // 2\n",
        "        num_beta = num_particles // 2\n",
        "\n",
        "        mapper = JordanWignerMapper()\n",
        "        init_state = HartreeFock(num_orbitals, (num_alpha, num_beta) , qubit_mapper=mapper)\n",
        "\n",
        "        # Setup the variational form for UCCSD (CISD equivalent in quantum chemistry)\n",
        "        '''\n",
        "        var_form = UCCSD(\n",
        "            num_orbitals,\n",
        "            (num_alpha, num_beta),\n",
        "            mapper,\n",
        "            initial_state=init_state,\n",
        "        )\n",
        "        '''\n",
        "        circ = init_state.compose(apply_su4_in_layers(params, depth, num_qubits, fix_pars_per_layer=False))\n",
        "\n",
        "\n",
        "        estimator_expvals = estimator.run(circuits=circ, observables=observable).result().values[0]\n",
        "\n",
        "        # Store the minimized value\n",
        "        #minimized_values.append(estimator_expvals)\n",
        "        #print(\"Minimized Exp Value: \", estimator_expvals)\n",
        "        #print(\"-----------------------------------------\")\n",
        "\n",
        "        return estimator_expvals\n",
        "\n",
        "    def grad_function(params):\n",
        "        \"\"\" use FiniteDiffEstimatorGradient \"\"\"\n",
        "\n",
        "        # Initialize Hartree-Fock state\n",
        "        num_orbitals = num_qubits // 2\n",
        "        num_particles = num_qubits // 2\n",
        "        num_alpha = num_particles // 2\n",
        "        num_beta = num_particles // 2\n",
        "\n",
        "        mapper = JordanWignerMapper()\n",
        "        init_state = HartreeFock(num_orbitals, (num_alpha, num_beta) , qubit_mapper=mapper)\n",
        "\n",
        "        # Setup the variational form for UCCSD (CISD equivalent in quantum chemistry)\n",
        "        '''\n",
        "        var_form = UCCSD(\n",
        "            num_orbitals,\n",
        "            (num_alpha, num_beta),\n",
        "            mapper,\n",
        "            initial_state=init_state,\n",
        "        )\n",
        "        '''\n",
        "\n",
        "        # This could be optimized\n",
        "        indices = [*range(len(params))]\n",
        "\n",
        "        offset = np.identity(len(params))[indices, :]\n",
        "        plus = params + 0.001 * offset\n",
        "        minus = params - 0.001 * offset\n",
        "\n",
        "        # length of circ_list is param_len\n",
        "        circ_plus_list = [init_state.compose(apply_su4_in_layers(params_value, depth, num_qubits, fix_pars_per_layer=False)) for params_value in plus]\n",
        "        circ_minus_list = [init_state.compose(apply_su4_in_layers(params_value, depth, num_qubits, fix_pars_per_layer=False)) for params_value in minus]\n",
        "\n",
        "        job_plus = estimator.run(circ_plus_list, [observable] * param_len)\n",
        "        job_minus = estimator.run(circ_minus_list, [observable] * param_len)\n",
        "        results_plus = job_plus.result()\n",
        "        results_minus = job_minus.result()\n",
        "\n",
        "\n",
        "        return (results_plus.values-results_minus.values)/0.002\n",
        "        #return circ_plus_list\n",
        "\n",
        "\n",
        "    # Create initial parameters\n",
        "    if init_param is None:\n",
        "        init_params = np.array([0.1]*param_len)\n",
        "    else:\n",
        "        init_params = init_param\n",
        "\n",
        "    Optimizer = optimizer(maxiter=max_steps)\n",
        "    # Optimizer = optimizer(maxfun=max_steps)\n",
        "    opt_result = Optimizer.minimize(fun=cost_function, x0=init_params, jac=grad_function)\n",
        "\n",
        "    return opt_result"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "metadata": {
        "id": "2a1i4HFm16ue",
        "tags": []
      },
      "outputs": [],
      "source": [
        "def separat_vqe(\n",
        "    num_qubits,\n",
        "    observable,\n",
        "    depth,\n",
        "    optimizer,\n",
        "    max_steps,\n",
        "    separat_steps,\n",
        "    init_param=None):\n",
        "\n",
        "    with open(\"vqe_params_\" + datetime + \".csv\", \"a\", newline='') as fd:\n",
        "        statewriter = csv.writer(fd, delimiter=',', quoting=csv.QUOTE_MINIMAL)\n",
        "        if init_param is None:\n",
        "            result = vqe(num_qubits,\n",
        "                         observable,\n",
        "                         depth,\n",
        "                         optimizer,\n",
        "                         max_steps=separat_steps,\n",
        "                         init_param=None)\n",
        "        else:\n",
        "            result = vqe(num_qubits,\n",
        "                         observable,\n",
        "                         depth,\n",
        "                         optimizer,\n",
        "                         max_steps=separat_steps,\n",
        "                         init_param=init_param)\n",
        "        print(result.fun)\n",
        "        statewriter.writerow(result.x)\n",
        "\n",
        "\n",
        "        for i in range(separat_steps,max_steps,separat_steps):\n",
        "            result = vqe(num_qubits,\n",
        "                         observable,\n",
        "                         depth,\n",
        "                         optimizer,\n",
        "                         separat_steps,\n",
        "                         init_param=result.x)\n",
        "            print(result.fun)\n",
        "            statewriter.writerow(result.x)\n",
        "            print(\"completed steps:\",i)\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "metadata": {
        "id": "zkEQNh_06SNU"
      },
      "outputs": [],
      "source": [
        "def vqe_full(\n",
        "    num_qubits,\n",
        "    observable,\n",
        "    depth,\n",
        "    optimizer,\n",
        "    max_steps,\n",
        "    skip_steps,\n",
        "    starting_param=None):\n",
        "\n",
        "  for i in range(skip_steps,max_steps,skip_steps):\n",
        "    #run for 1 step\n",
        "    result = vqe(num_qubits,\n",
        "                  observable,\n",
        "                  depth,\n",
        "                  optimizer,\n",
        "                  skip_steps,\n",
        "                  init_param=starting_param)\n",
        "    print(result.fun)\n",
        "    # save init parms\n",
        "\n",
        "    with open(\"vqe_params_\" + str(i) + \".csv\", \"w\", newline='') as fd:\n",
        "      statewriter = csv.writer(fd, delimiter=',', quoting=csv.QUOTE_MINIMAL)\n",
        "      starting_params = result.x\n",
        "      statewriter.writerow(result.x)\n",
        "\n",
        "    print(\"completed steps:\",i)\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "DLl-TnzBwbXz",
        "outputId": "2867e176-fd4c-412e-d10a-8f78debd83a6",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "MN =  12\n",
            "12qubits_00\n"
          ]
        }
      ],
      "source": [
        "nqubits = 12  # Number of qubits\n",
        "fix_pars_per_layer = False\n",
        "observable = make_observable(nqubits, mat=False)\n",
        "# observable"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 23,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "zEyiniv6krbS",
        "outputId": "26bc9ee1-e108-4814-c611-560ed5cbf947"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "==================== Memory Usage ====================\n",
            "total 8589934592\n",
            "available 1320960000\n",
            "percent 84.6\n",
            "used 3265282048\n",
            "free 43958272\n",
            "active 1293746176\n",
            "inactive 1200652288\n",
            "wired 1971535872\n",
            "==================== CPU Usage ====================\n",
            "CPU percent: 78.5%\n"
          ]
        }
      ],
      "source": [
        "print_cpu_usage()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 24,
      "metadata": {
        "id": "iZTcGlpaGAyU",
        "tags": []
      },
      "outputs": [],
      "source": [
        "# from qiskit_algorithms.optimizers import SLSQP, P_BFGS, NELDER_MEAD\n",
        "# import time\n",
        "\n",
        "# start = time.time()\n",
        "# testt = vqe([0],\n",
        "#             nqubits,\n",
        "#             observable,\n",
        "#             [],\n",
        "#             depth=2,\n",
        "#             optimizer=NELDER_MEAD,\n",
        "#             max_steps=100\n",
        "#            )\n",
        "# # circuit.num_parameters\n",
        "# end = time.time()\n",
        "# print(\"Time Taken: \", end-start)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 25,
      "metadata": {
        "id": "vc0zR6KZv6uh",
        "tags": []
      },
      "outputs": [],
      "source": [
        "# SLSQP 100\n",
        "# testt.fun"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 26,
      "metadata": {
        "id": "mA_QMesa16ue",
        "tags": []
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "-8.284832195943544\n",
            "completed steps: 2\n",
            "-8.284832195943544\n",
            "completed steps: 4\n",
            "-8.284832195943544\n",
            "completed steps: 6\n",
            "-8.284832195943544\n",
            "completed steps: 8\n",
            "Time Taken:  146.6697609424591\n"
          ]
        }
      ],
      "source": [
        "from qiskit_algorithms.optimizers import SLSQP, P_BFGS, NELDER_MEAD\n",
        "import time\n",
        "\n",
        "start = time.time()\n",
        "datetime = time.strftime(\"%Y-%m-%d_%H-%M-%S\")\n",
        "\n",
        "testt = vqe_full(\n",
        "    nqubits,\n",
        "    observable,\n",
        "    depth=1,\n",
        "    optimizer=SLSQP,\n",
        "    max_steps=10,\n",
        "    skip_steps=2,\n",
        "    starting_param=None)\n",
        "\n",
        "# For running vqe without saving params \n",
        "'''\n",
        "test = vqe(nqubits,\n",
        "           observable,\n",
        "           depth=1,\n",
        "           optimizer=SLSQP,\n",
        "           max_steps=10,\n",
        "           init_param=None)\n",
        "'''\n",
        "end = time.time()\n",
        "print(\"Time Taken: \", end-start)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 141
        },
        "id": "by1La41516ue",
        "outputId": "2519a812-422f-4f9b-c80b-7dc906eac3dd",
        "tags": []
      },
      "outputs": [
        {
          "ename": "AttributeError",
          "evalue": "'NoneType' object has no attribute 'fun'",
          "output_type": "error",
          "traceback": [
            "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
            "\u001b[0;31mAttributeError\u001b[0m                            Traceback (most recent call last)",
            "\u001b[0;32m<ipython-input-22-ac9f2e130431>\u001b[0m in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mtestt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfun\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
            "\u001b[0;31mAttributeError\u001b[0m: 'NoneType' object has no attribute 'fun'"
          ]
        }
      ],
      "source": [
        "testt.fun"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "VkG73BJJv6uh"
      },
      "outputs": [],
      "source": [
        "energies = np.linalg.eigvalsh(observable_matrix)\n",
        "gs_energy = energies[0]\n",
        "gs_energy"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "k6j4FXyqMGVm"
      },
      "source": [
        "TODO:"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "DxQ2DlB9DRf1"
      },
      "source": [
        "1. Optimizer Table for complexity\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "sleOKCt_v6uh"
      },
      "outputs": [],
      "source": [
        "#2. check mps simulator"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "machine_shape": "hm",
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3 [Default]",
      "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.11.6"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}