|
| 1 | +{ |
| 2 | + "cells": [ |
| 3 | + { |
| 4 | + "cell_type": "markdown", |
| 5 | + "metadata": {}, |
| 6 | + "source": [ |
| 7 | + "# Tutorial for the Richardson-Gaudin model" |
| 8 | + ] |
| 9 | + }, |
| 10 | + { |
| 11 | + "cell_type": "markdown", |
| 12 | + "metadata": {}, |
| 13 | + "source": [ |
| 14 | + "This tutorial demonstrates how to use the Richardson-Gaudin model with the ModelHamiltonian package.\n", |
| 15 | + "Specifically, it shows how to set up the connectivity matrix and perform the FCI calculation using the pyscf library.\n", |
| 16 | + "\n", |
| 17 | + "The primary focus of this tutorial, however, it to demostrate the flexibility behind the Model Hamiltonin package. To do so, we will show how to use this library to build a Richardson Model that's given by:\n", |
| 18 | + "\n", |
| 19 | + "$$\n", |
| 20 | + "\\hat{H}_{B C S}=\\frac{1}{2} \\sum_{i=1}^N \\varepsilon_i\\left(a_{i \\uparrow}^{\\dagger} a_{1 \\uparrow}+a_{i \\downarrow}^{\\dagger} a_{i \\downarrow}\\right)-\\frac{g}{2} \\sum_{i, j=1}^N a_{i \\uparrow}^{\\dagger} a_{i \\downarrow}^{\\dagger} a_{j \\downarrow} a_{j \\uparrow}\n", |
| 21 | + "$$\n", |
| 22 | + "\n", |
| 23 | + "This definition is given in a [paper](https://arxiv.org/pdf/2506.09379) and can be used to validate numerical results that were kindly provided by [Paul Johnson](https://paul-a-johnson.fsg.ulaval.ca/en)." |
| 24 | + ] |
| 25 | + }, |
| 26 | + { |
| 27 | + "cell_type": "markdown", |
| 28 | + "metadata": {}, |
| 29 | + "source": [ |
| 30 | + "# 1. Match the definition\n", |
| 31 | + "\n", |
| 32 | + "In our case, we define the Richardson model as follows:\n", |
| 33 | + "$$\n", |
| 34 | + "\\hat{H}_{R G}=\\sum_p\\left(\\mu_p^Z-J_{p p}^{\\mathrm{eq}}\\right) S_p^Z+\\sum_{p q} J_{p q}^{\\mathrm{eq}} S_p^{+} S_q^{-}\n", |
| 35 | + "$$\n", |
| 36 | + "\n", |
| 37 | + "To match the above definition we need to move around some of the terms. Skipping the details, we can match the definition by setting the following parameters:\n", |
| 38 | + "\n", |
| 39 | + "$$\n", |
| 40 | + "\\begin{align*}\n", |
| 41 | + "\\mu_p &= \\varepsilon_p + J_{p p}^{\\mathrm{eq}} \\\\\n", |
| 42 | + "J_{p q}^{\\mathrm{eq}} &= -\\frac{g}{2} \\\\\n", |
| 43 | + "\\end{align*}\n", |
| 44 | + "$$\n", |
| 45 | + "\n", |
| 46 | + "# 2. Provided data\n", |
| 47 | + "For the values of $\\varepsilon_i$ = [0, 1, 2, 3, 4, 5] and ground state energy of 6 electron system, the following data was provided by the author:\n", |
| 48 | + "\n", |
| 49 | + "| g (unitless) | Energy (unitless) |\n", |
| 50 | + "| ----------------- | -------------- |\n", |
| 51 | + "| 3.536 | -14.20687553 |\n", |
| 52 | + "| 460.5454545454545 | -2755.77652711 |\n", |
| 53 | + "| 120.8888888888889 | -717.84780933 |\n", |
| 54 | + "| 8.182 | -41.80555486 |\n", |
| 55 | + "| 6.869 | -33.96821191 |\n", |
| 56 | + "| 480.5252525252525 | -2875.655157 |\n", |
| 57 | + "| 5.556 | -26.14992638 |\n", |
| 58 | + "| 1.011 | -0.14612126 |\n", |
| 59 | + "| 6.263 | -30.35668599 |\n", |
| 60 | + "| 4.849 | -21.95332479 |\n" |
| 61 | + ] |
| 62 | + }, |
| 63 | + { |
| 64 | + "cell_type": "markdown", |
| 65 | + "metadata": {}, |
| 66 | + "source": [ |
| 67 | + "# Calculation" |
| 68 | + ] |
| 69 | + }, |
| 70 | + { |
| 71 | + "cell_type": "code", |
| 72 | + "execution_count": 1, |
| 73 | + "metadata": {}, |
| 74 | + "outputs": [ |
| 75 | + { |
| 76 | + "name": "stdout", |
| 77 | + "output_type": "stream", |
| 78 | + "text": [ |
| 79 | + "| g | Energy | Computed | Error |\n", |
| 80 | + "|---------|--------------|--------------|---------|\n", |
| 81 | + "| 3.536 | -14.2069 | -14.2069 | 0 |\n", |
| 82 | + "| 460.545 | -2755.78 | -2755.78 | 0 |\n", |
| 83 | + "| 120.889 | -717.848 | -717.848 | 0 |\n", |
| 84 | + "| 8.182 | -41.8056 | -41.8056 | 0 |\n", |
| 85 | + "| 6.869 | -33.9682 | -33.9682 | 0 |\n", |
| 86 | + "| 480.525 | -2875.66 | -2875.66 | 0 |\n", |
| 87 | + "| 5.556 | -26.1499 | -26.1499 | 0 |\n", |
| 88 | + "| 1.011 | -0.146121 | -0.146121 | 0 |\n", |
| 89 | + "| 6.263 | -30.3567 | -30.3567 | 0 |\n", |
| 90 | + "| 4.849 | -21.9533 | -21.9533 | 0 |\n" |
| 91 | + ] |
| 92 | + } |
| 93 | + ], |
| 94 | + "source": [ |
| 95 | + "import numpy as np\n", |
| 96 | + "from moha.hamiltonians import HamRG\n", |
| 97 | + "from pyscf import fci\n", |
| 98 | + "from tabulate import tabulate\n", |
| 99 | + "\n", |
| 100 | + "# Define the connectivity matrix assuming a fully connected system\n", |
| 101 | + "connectivity = np.ones((6, 6))\n", |
| 102 | + "g_values = [\n", |
| 103 | + " 3.536,\n", |
| 104 | + " 460.5454545454545,\n", |
| 105 | + " 120.88888888888889,\n", |
| 106 | + " 8.182,\n", |
| 107 | + " 6.869000000000001,\n", |
| 108 | + " 480.52525252525254,\n", |
| 109 | + " 5.556000000000001,\n", |
| 110 | + " 1.011,\n", |
| 111 | + " 6.263000000000001,\n", |
| 112 | + " 4.849000000000001\n", |
| 113 | + "]\n", |
| 114 | + "\n", |
| 115 | + "energy_values = [\n", |
| 116 | + " -14.20687553,\n", |
| 117 | + " -2755.77652711,\n", |
| 118 | + " -717.84780933,\n", |
| 119 | + " -41.80555486,\n", |
| 120 | + " -33.96821191,\n", |
| 121 | + " -2875.655157,\n", |
| 122 | + " -26.14992638,\n", |
| 123 | + " -0.14612126,\n", |
| 124 | + " -30.35668599,\n", |
| 125 | + " -21.95332479\n", |
| 126 | + "]\n", |
| 127 | + "energy_computed = []\n", |
| 128 | + "eps = np.arange(6)\n", |
| 129 | + "# Create the Richardson-Gaudin Hamiltonian\n", |
| 130 | + "for g in g_values:\n", |
| 131 | + " # convert parameters to match the definition\n", |
| 132 | + " J = -g / 2\n", |
| 133 | + " mu = eps + J \n", |
| 134 | + "\n", |
| 135 | + " # Create the Hamiltonian\n", |
| 136 | + " ham = HamRG(J_eq = J*connectivity,\n", |
| 137 | + " mu = mu)\n", |
| 138 | + " e0 = ham.generate_zero_body_integral()\n", |
| 139 | + "\n", |
| 140 | + " h1 = ham.generate_one_body_integral(dense=True, basis='spatial basis')\n", |
| 141 | + " h2 = ham.generate_two_body_integral(dense=True, basis='spatial basis', sym=4)\n", |
| 142 | + "\n", |
| 143 | + " # convert to chemists notation\n", |
| 144 | + " h2_ch = np.einsum('ijkl->ikjl', h2)\n", |
| 145 | + "\n", |
| 146 | + "\n", |
| 147 | + " n_sites = 6\n", |
| 148 | + " occ = (n_sites//2, n_sites//2)\n", |
| 149 | + "\n", |
| 150 | + " # Solve the FCI problem\n", |
| 151 | + " e, ci = fci.direct_spin0.kernel(1*h1, 2*h2_ch, nelec=occ, norb=n_sites, nroots=6, ecore=0)\n", |
| 152 | + " energy_computed.append(e[0])\n", |
| 153 | + "\n", |
| 154 | + "# Print the results in a table format\n", |
| 155 | + "\n", |
| 156 | + "table_data = []\n", |
| 157 | + "for g, energy, computed in zip(g_values, energy_values, energy_computed):\n", |
| 158 | + " error = abs(computed - energy)\n", |
| 159 | + " table_data.append([f\"{g:.3f}\", f\"{energy:.6f}\", f\"{computed:.6f}\", f\"{error:.6f}\"])\n", |
| 160 | + "\n", |
| 161 | + "headers = [\"g\", \"Energy\", \"Computed\", \"Error\"]\n", |
| 162 | + "print(tabulate(table_data, headers=headers, tablefmt=\"github\"))" |
| 163 | + ] |
| 164 | + } |
| 165 | + ], |
| 166 | + "metadata": { |
| 167 | + "kernelspec": { |
| 168 | + "display_name": "base", |
| 169 | + "language": "python", |
| 170 | + "name": "python3" |
| 171 | + }, |
| 172 | + "language_info": { |
| 173 | + "codemirror_mode": { |
| 174 | + "name": "ipython", |
| 175 | + "version": 3 |
| 176 | + }, |
| 177 | + "file_extension": ".py", |
| 178 | + "mimetype": "text/x-python", |
| 179 | + "name": "python", |
| 180 | + "nbconvert_exporter": "python", |
| 181 | + "pygments_lexer": "ipython3", |
| 182 | + "version": "3.9.7" |
| 183 | + } |
| 184 | + }, |
| 185 | + "nbformat": 4, |
| 186 | + "nbformat_minor": 2 |
| 187 | +} |
0 commit comments