Complete dbi_strategies_compare notebook

examples/dbi/dbi_strategies_compare.ipynb

@@ -0,0 +1,254 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# DBI strategies comparison\n",
+    "\n",
+    "This notebook is a comparison of the so-far developed diagonalization strategies for DBI, including the canonical, Pauli-Z, and magnetic field strategies. On top of these, we also show case the use of invariant DBI generators such as 'BHMM'."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from copy import deepcopy\n",
+    "\n",
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
+    "\n",
+    "from qibo import hamiltonians, set_backend\n",
+    "from qibo.hamiltonians import Hamiltonian, SymbolicHamiltonian\n",
+    "from qibo.quantum_info import random_hermitian\n",
+    "from qibo.models.dbi.double_bracket import DoubleBracketGeneratorType, DoubleBracketScheduling, DoubleBracketIteration\n",
+    "from qibo.models.dbi.utils import *"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def visualize_matrix(matrix, title=\"\"):\n",
+    "    \"\"\"Visualize hamiltonian in a heatmap form.\"\"\"\n",
+    "    fig, ax = plt.subplots(figsize=(5,5))\n",
+    "    ax.set_title(title)\n",
+    "    try:\n",
+    "        im = ax.imshow(np.absolute(matrix), cmap=\"inferno\")\n",
+    "    except TypeError:\n",
+    "        im = ax.imshow(np.absolute(matrix.get()), cmap=\"inferno\")\n",
+    "    fig.colorbar(im, ax=ax)\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Test on random Hamiltonian\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# backend\n",
+    "set_backend(\"qibojit\", \"numba\")\n",
+    "# initialize dbi object\n",
+    "nqubits = 5\n",
+    "h0 = random_hermitian(2**nqubits, seed=2)\n",
+    "dbi = DoubleBracketIteration(Hamiltonian(nqubits=nqubits, matrix=h0))\n",
+    "print(\"Initial off diagonal norm\", dbi.off_diagonal_norm)\n",
+    "visualize_matrix(dbi.h.matrix, title=f'Random hamiltonian with L={nqubits}')"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# iterations steps\n",
+    "NSTEPS = 15\n",
+    "# choose polynomial scheduling\n",
+    "scheduling = DoubleBracketScheduling.use_hyperopt"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Canonical"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# initialize DBI class for the canonical case\n",
+    "dbi_canonical = DoubleBracketIteration(Hamiltonian(nqubits=nqubits, matrix=h0), mode=DoubleBracketGeneratorType.canonical, scheduling=scheduling)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Canonical\n",
+    "off_diagonal_norm_history_canonical = [dbi_canonical.off_diagonal_norm]\n",
+    "steps_canonical_plot = [0]\n",
+    "for s in range(NSTEPS):\n",
+    "    # same settings as iteration from list\n",
+    "    step = dbi_canonical.choose_step(d=dbi.diagonal_h_matrix, max_evals=50)\n",
+    "    dbi_canonical(step=step)\n",
+    "    print(f\"New optimized step at iteration {s+1}/{NSTEPS}: {step}, loss {dbi_canonical.off_diagonal_norm}\")\n",
+    "    off_diagonal_norm_history_canonical.append(dbi_canonical.off_diagonal_norm)\n",
+    "    steps_canonical_plot.append(steps_canonical_plot[-1]+step)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Pauli-Z"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# initialize DBI class for the Pauli-Z strategy\n",
+    "dbi_pauli = DoubleBracketIteration(Hamiltonian(nqubits=nqubits, matrix=h0), mode=DoubleBracketGeneratorType.single_commutator, scheduling=scheduling)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "generate_local_Z = generate_Z_operators(nqubits)\n",
+    "Z_ops = list(generate_local_Z.values())\n",
+    "Z_names = list(generate_local_Z.keys())\n",
+    "Z_optimal = []\n",
+    "# add in initial values for plotting\n",
+    "off_diagonal_norm_history_pauli = [dbi_pauli.off_diagonal_norm]\n",
+    "steps_pauli_plot = [0]\n",
+    "scheduling = DoubleBracketScheduling.use_hyperopt\n",
+    "for _ in range(NSTEPS):\n",
+    "    dbi_pauli, idx, step, flip_sign = select_best_dbr_generator(dbi_pauli, Z_ops, scheduling=scheduling, compare_canonical=False, max_evals=50)\n",
+    "    off_diagonal_norm_history_pauli.append(dbi_pauli.off_diagonal_norm)\n",
+    "    steps_pauli_plot.append(steps_pauli_plot[-1]+step)\n",
+    "    if flip_sign < 0:\n",
+    "        Z_optimal.append('-' + Z_names[idx])\n",
+    "    else:\n",
+    "        Z_optimal.append(Z_names[idx])\n",
+    "    print(f\"New optimized step at iteration {_+1}/{NSTEPS}: {step} with operator {Z_optimal[-1]}, loss {dbi_pauli.off_diagonal_norm}\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Magnetic field"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# initialize DBI class for the canonical case\n",
+    "dbi_gradient = DoubleBracketIteration(Hamiltonian(nqubits=nqubits, matrix=h0), mode=DoubleBracketGeneratorType.single_commutator, scheduling=scheduling)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "onsite_Z_ops = generate_onsite_Z_ops(nqubits)\n",
+    "d_coef = onsite_Z_decomposition(dbi_gradient.h.matrix, onsite_Z_ops)\n",
+    "d = sum([d_coef[i] * onsite_Z_ops[i] for i in range(nqubits)])"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "off_diagonal_norm_history_gradient = [dbi_gradient.off_diagonal_norm]\n",
+    "steps_gradient_plot= [0]\n",
+    "for _ in range(NSTEPS):\n",
+    "    step, d_coef, d = gradient_descent_onsite_Z(dbi_gradient, d_coef, d, onsite_Z_ops=onsite_Z_ops, max_evals=50, n_taylor=5)\n",
+    "    off_diagonal_norm_history_gradient.append(dbi_gradient.off_diagonal_norm)\n",
+    "    print(f\"New optimized step at iteration {_+1}/{NSTEPS}: {step} with d_coef {d_coef}, loss {dbi_gradient.off_diagonal_norm}\")\n",
+    "    steps_gradient_plot.append(steps_gradient_plot[-1]+step)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "plt.title(str(nqubits) + ' random Hamiltonian diagonalization')\n",
+    "plt.plot(off_diagonal_norm_history_canonical, label='canonical')\n",
+    "plt.plot(off_diagonal_norm_history_pauli, label='Pauli-Z')\n",
+    "plt.plot(off_diagonal_norm_history_gradient, label='gradient')\n",
+    "plt.legend()\n",
+    "plt.xlabel('Iteration')\n",
+    "plt.ylabel(r'$|| \\sigma(e^{sW}He^{-sW}) || $')"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "plt.title(str(nqubits) + ' random Hamiltonian diagonalization')\n",
+    "plt.plot(steps_canonical_plot, off_diagonal_norm_history_canonical, marker='o', label='canonical')\n",
+    "plt.plot(steps_pauli_plot, off_diagonal_norm_history_pauli, marker='o', label='Pauli-Z')\n",
+    "plt.plot(steps_gradient_plot,off_diagonal_norm_history_gradient, marker='o', label='gradient')\n",
+    "plt.legend()\n",
+    "plt.xlabel('Iteration')\n",
+    "plt.ylabel(r'$|| \\sigma(e^{sW}He^{-sW}) || $')"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "DBF_qibo",
+   "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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}