{
"cells": [
{
"cell_type": "markdown",
"id": "de193554-b271-4295-95e4-8904f0f6ee8a",
"metadata": {},
"source": [
"---\n",
"title: \"Geometría molecular\"\n",
"description: \"En esta lección variamos la geometría de una molécula simple, minimizando la energía en cada paso.\"\n",
"---\n",
"\n",
"{/* cspell:ignore pxxr prqs nelecas mcscf chmax Dmax vmax ecore ncas Excp disp */}\n",
"\n",
"\n",
"\n",
"# Determinación de la geometría molecular\n",
"\n",
"En la sección anterior, implementamos VQE para determinar la energía de estado básico de una molécula. Ese es un uso válido de la computación cuántica, pero aún más útil sería determinar la estructura de una molécula.\n",
"\n",
"\n",
"\n",
"## Paso 1: Asignar entradas clásicas a un problema cuántico\n",
"\n",
"Siguiendo con nuestro ejemplo básico del hidrógeno diatómico, el único parámetro geométrico que varía es la longitud del enlace. Para ello, procedemos como antes, pero utilizando una variable en nuestra construcción inicial de la molécula (una longitud de enlace, *x*, en el argumento). Se trata de un cambio bastante sencillo, pero requiere que la variable se incluya en funciones a lo largo de todo el proceso, ya que comienza en la construcción del hamiltoniano fermiónico y se propaga a través del mapeo y, finalmente, a la función de coste.\n",
"\n",
"En primer lugar, cargamos algunos de los paquetes que hemos utilizado antes y definimos la función Cholesky.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 30,
"id": "0a5d39bd-8c26-404f-8057-c29e3af70df4",
"metadata": {},
"outputs": [],
"source": [
"from qiskit.quantum_info import SparsePauliOp\n",
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"\n",
"#!pip install pyscf==2.4.0\n",
"from pyscf import ao2mo, gto, mcscf, scf\n",
"\n",
"\n",
"def cholesky(V, eps):\n",
" # see https://arxiv.org/pdf/1711.02242.pdf section B2\n",
" # see https://arxiv.org/abs/1808.02625\n",
" # see https://arxiv.org/abs/2104.08957\n",
" no = V.shape[0]\n",
" chmax, ng = 20 * no, 0\n",
" W = V.reshape(no**2, no**2)\n",
" L = np.zeros((no**2, chmax))\n",
" Dmax = np.diagonal(W).copy()\n",
" nu_max = np.argmax(Dmax)\n",
" vmax = Dmax[nu_max]\n",
" while vmax > eps:\n",
" L[:, ng] = W[:, nu_max]\n",
" if ng > 0:\n",
" L[:, ng] -= np.dot(L[:, 0:ng], (L.T)[0:ng, nu_max])\n",
" L[:, ng] /= np.sqrt(vmax)\n",
" Dmax[: no**2] -= L[: no**2, ng] ** 2\n",
" ng += 1\n",
" nu_max = np.argmax(Dmax)\n",
" vmax = Dmax[nu_max]\n",
" L = L[:, :ng].reshape((no, no, ng))\n",
" print(\n",
" \"accuracy of Cholesky decomposition \",\n",
" np.abs(np.einsum(\"prg,qsg->prqs\", L, L) - V).max(),\n",
" )\n",
" return L, ng\n",
"\n",
"\n",
"def identity(n):\n",
" return SparsePauliOp.from_list([(\"I\" * n, 1)])\n",
"\n",
"\n",
"def creators_destructors(n, mapping=\"jordan_wigner\"):\n",
" c_list = []\n",
" if mapping == \"jordan_wigner\":\n",
" for p in range(n):\n",
" if p == 0:\n",
" ell, r = \"I\" * (n - 1), \"\"\n",
" elif p == n - 1:\n",
" ell, r = \"\", \"Z\" * (n - 1)\n",
" else:\n",
" ell, r = \"I\" * (n - p - 1), \"Z\" * p\n",
" cp = SparsePauliOp.from_list([(ell + \"X\" + r, 0.5), (ell + \"Y\" + r, -0.5j)])\n",
" c_list.append(cp)\n",
" else:\n",
" raise ValueError(\"Unsupported mapping.\")\n",
" d_list = [cp.adjoint() for cp in c_list]\n",
" return c_list, d_list"
]
},
{
"cell_type": "markdown",
"id": "40e25b54-7927-4dbb-a26f-1c6b33f7f349",
"metadata": {},
"source": [
"Ahora para definir nuestro Hamiltoniano, usaremos PySCF exactamente como en el ejemplo anterior, pero ahora incluiremos una variable, `x`, para que juegue el papel de nuestra distancia interatómica. Esto devolverá la energía del núcleo, la energía de un solo electrón y las energías de dos electrones como antes.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "dbd10d0c-feb1-4a86-9bd9-b61101a08b95",
"metadata": {},
"outputs": [],
"source": [
"def ham_terms(x: float):\n",
" distance = x\n",
" a = distance / 2\n",
" mol = gto.Mole()\n",
" mol.build(\n",
" verbose=0,\n",
" atom=[\n",
" [\"H\", (0, 0, -a)],\n",
" [\"H\", (0, 0, a)],\n",
" ],\n",
" basis=\"sto-6g\",\n",
" spin=0,\n",
" charge=0,\n",
" symmetry=\"Dooh\",\n",
" )\n",
"\n",
" # mf = scf.RHF(mol)\n",
" # mx = mcscf.CASCI(mf, ncas=2, nelecas=(1, 1))\n",
" # mx.kernel()\n",
"\n",
" mf = scf.RHF(mol)\n",
" mf.kernel()\n",
" if not mf.converged:\n",
" raise RuntimeError(f\"SCF did not converge for distance {x}\")\n",
"\n",
" mx = mcscf.CASCI(mf, ncas=2, nelecas=(1, 1))\n",
" casci_energy = mx.kernel()\n",
" if casci_energy is None:\n",
" raise RuntimeError(f\"CASCI failed for distance {x}\")\n",
"\n",
" # Other variables that might come in handy:\n",
" # active_space = range(mol.nelectron // 2 - 1, mol.nelectron // 2 + 1)\n",
" # E1 = mf.kernel()\n",
" # mo = mx.sort_mo(active_space, base=0)\n",
" # E2 = mx.kernel(mo)[:2]\n",
"\n",
" h1e, ecore = mx.get_h1eff()\n",
" h2e = ao2mo.restore(1, mx.get_h2eff(), mx.ncas)\n",
" return ecore, h1e, h2e"
]
},
{
"cell_type": "markdown",
"id": "db49a702-e60e-4c9b-8ff2-94aa2ada022c",
"metadata": {},
"source": [
"Recordemos que la construcción anterior está haciendo un Hamiltoniano fermiónico basado en la especie atómica, la geometría y los orbitales electrónicos. A continuación, mapeamos este Hamiltoniano fermiónico en operadores de Pauli. Esta función `build_hamiltonian` también incluirá una variable geométrica como argumento.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 32,
"id": "84e6a56b-eaea-4c0a-8502-567cfc5140a2",
"metadata": {},
"outputs": [],
"source": [
"def build_hamiltonian(distx: float) -> SparsePauliOp:\n",
" ecore = ham_terms(distx)[0]\n",
" h1e = ham_terms(distx)[1]\n",
" h2e = ham_terms(distx)[2]\n",
"\n",
" ncas, _ = h1e.shape\n",
"\n",
" C, D = creators_destructors(2 * ncas, mapping=\"jordan_wigner\")\n",
" Exc = []\n",
" for p in range(ncas):\n",
" Excp = [C[p] @ D[p] + C[ncas + p] @ D[ncas + p]]\n",
" for r in range(p + 1, ncas):\n",
" Excp.append(\n",
" C[p] @ D[r]\n",
" + C[ncas + p] @ D[ncas + r]\n",
" + C[r] @ D[p]\n",
" + C[ncas + r] @ D[ncas + p]\n",
" )\n",
" Exc.append(Excp)\n",
"\n",
" # low-rank decomposition of the Hamiltonian\n",
" Lop, ng = cholesky(h2e, 1e-6)\n",
" t1e = h1e - 0.5 * np.einsum(\"pxxr->pr\", h2e)\n",
"\n",
" H = ecore * identity(2 * ncas)\n",
" # one-body term\n",
" for p in range(ncas):\n",
" for r in range(p, ncas):\n",
" H += t1e[p, r] * Exc[p][r - p]\n",
" # two-body term\n",
" for g in range(ng):\n",
" Lg = 0 * identity(2 * ncas)\n",
" for p in range(ncas):\n",
" for r in range(p, ncas):\n",
" Lg += Lop[p, r, g] * Exc[p][r - p]\n",
" H += 0.5 * Lg @ Lg\n",
"\n",
" return H.chop().simplify()"
]
},
{
"cell_type": "markdown",
"id": "ffcb2569-f832-4fa9-a8bb-d7626c2a233d",
"metadata": {},
"source": [
"Cargaremos los paquetes restantes para ejecutar el propio VQE, como el ansatz efficient\\_su2, y los minimizadores SciPy :\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "5f8dd8e0-6dfb-4be7-af15-1591f569201c",
"metadata": {},
"outputs": [],
"source": [
"# General imports\n",
"\n",
"# Pre-defined ansatz circuit and operator class for Hamiltonian\n",
"from qiskit.circuit.library import efficient_su2\n",
"from qiskit.quantum_info import SparsePauliOp\n",
"\n",
"# SciPy minimizer routine\n",
"from scipy.optimize import minimize\n",
"\n",
"# Plotting functions\n",
"\n",
"# Qiskit Runtime tools\n",
"from qiskit_ibm_runtime import QiskitRuntimeService\n",
"\n",
"service = QiskitRuntimeService()"
]
},
{
"cell_type": "markdown",
"id": "f94c1680-9c2a-4585-a63e-a49da4eb02f3",
"metadata": {},
"source": [
"Definiremos de nuevo la función de coste, pero ésta siempre toma como argumento un Hamiltoniano completamente construido y mapeado, por lo que nada cambia en esta función.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 34,
"id": "0a80ad8c-a9cb-4cab-835d-f65131b99c87",
"metadata": {},
"outputs": [],
"source": [
"def cost_func(params, ansatz, H, estimator):\n",
" pub = (ansatz, [H], [params])\n",
" result = estimator.run(pubs=[pub]).result()\n",
" energy = result[0].data.evs[0]\n",
" return energy\n",
"\n",
"\n",
"# def cost_func_sim(params, ansatz, H, estimator):\n",
"# energy = estimator.run(ansatz, H, parameter_values=params).result().values[0]\n",
"# return energy"
]
},
{
"cell_type": "markdown",
"id": "a068e645-8dfb-4e9c-b1b1-7dd936808188",
"metadata": {},
"source": [
"\n",
"\n",
"## Paso 2: Optimizar el problema para la ejecución cuántica\n",
"\n",
"Como el Hamiltoniano cambiará con cada nueva geometría, el transpiling del operador cambiará en cada paso. No obstante, podemos definir un gestor de pases general que se aplique en cada paso, específico para el hardware que queramos utilizar.\n",
"\n",
"Aquí utilizaremos el backend menos ocupado disponible. Utilizaremos ese backend como modelo para nuestro AerSimulator, permitiendo a nuestro simulador imitar, por ejemplo, el comportamiento del ruido del backend real. Estos modelos de ruido no son perfectos, pero pueden ayudarle a saber qué esperar del hardware real.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "fa07518f-a2c0-4f7a-b344-8d9576427478",
"metadata": {},
"outputs": [],
"source": [
"# Here, we select the least busy backend available:\n",
"backend = service.least_busy(operational=True, simulator=False)\n",
"print(backend)\n",
"# Or to select a specific real backend use the line below, and substitute 'ibm_strasbourg'\n",
"# for your chosen device. backend = service.get_backend('ibm_strasbourg')"
]
},
{
"cell_type": "code",
"execution_count": 36,
"id": "e67fc84a-431b-4efd-937f-49c8b7ac3abb",
"metadata": {},
"outputs": [],
"source": [
"# To run on a simulator:\n",
"# -----------\n",
"from qiskit_aer import AerSimulator\n",
"\n",
"backend_sim = AerSimulator.from_backend(backend)"
]
},
{
"cell_type": "markdown",
"id": "0636c55f-f46f-46ca-acaf-72bcc2f5f663",
"metadata": {},
"source": [
"Importamos el gestor de pases y los paquetes relacionados para ayudarnos a optimizar nuestro circuito. Este paso, y el anterior, son independientes del Hamiltoniano, por lo que no cambian respecto a la lección anterior.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 37,
"id": "8202332e-69ca-4049-af27-e98a77f15a5d",
"metadata": {},
"outputs": [],
"source": [
"from qiskit.transpiler import PassManager\n",
"from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager\n",
"from qiskit.transpiler.passes import (\n",
" ALAPScheduleAnalysis,\n",
" PadDynamicalDecoupling,\n",
" ConstrainedReschedule,\n",
")\n",
"from qiskit.circuit.library import XGate\n",
"\n",
"target = backend.target\n",
"pm = generate_preset_pass_manager(target=target, optimization_level=3)\n",
"pm.scheduling = PassManager(\n",
" [\n",
" ALAPScheduleAnalysis(target=target),\n",
" ConstrainedReschedule(\n",
" acquire_alignment=target.acquire_alignment,\n",
" pulse_alignment=target.pulse_alignment,\n",
" target=target,\n",
" ),\n",
" PadDynamicalDecoupling(\n",
" target=target,\n",
" dd_sequence=[XGate(), XGate()],\n",
" pulse_alignment=target.pulse_alignment,\n",
" ),\n",
" ]\n",
")"
]
},
{
"cell_type": "markdown",
"id": "4932ef5e-3c81-46a1-92cb-3082399734a0",
"metadata": {},
"source": [
"\n",
"\n",
"## Paso 3: Ejecutar utilizando Qiskit primitives.\n",
"\n",
"En el bloque de código que aparece a continuación, creamos una matriz para almacenar los resultados de cada paso de nuestro algoritmo de cálculo de distancias interatómicas $x$. Hemos elegido el intervalo de $x$ basándonos en nuestro conocimiento del valor experimental de la longitud de enlace en equilibrio: 0.74 angstroms. Primero lo ejecutaremos en un simulador, para lo cual importaremos nuestro Estimator ( BackendEstimator ) desde `qiskit.primitives`... Para cada paso geométrico, construimos el hamiltoniano y permitimos un número determinado de pasos de optimización (en este caso, 500) utilizando el optimizador «cobyla». En cada paso geométrico, almacenamos tanto la energía total como la energía electrónica. Debido al gran número de pasos del optimizador, esto puede tardar una hora o más. Quizás le interese modificar los datos que figuran a continuación para reducir el tiempo necesario.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 38,
"id": "4c41a221-8a02-4932-882b-5afcc98d1d8d",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"accuracy of Cholesky decomposition 1.1102230246251565e-15\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/home/porter284/.pyenv/versions/3.11.12/lib/python3.11/site-packages/scipy/_lib/pyprima/common/preproc.py:68: UserWarning: COBYLA: Invalid MAXFUN; it should be at least num_vars + 2; it is set to 34\n",
" warn(f'{solver}: Invalid MAXFUN; it should be at least {min_maxfun_str}; it is set to {maxfun}')\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = 1.316011435623847\n",
"The corresponding X is:\n",
"[2.32948769 5.39918229 3.03787975 4.11789904 4.97130735 2.68662232\n",
" 1.76573151 2.48982571 5.40431972 3.65780829 1.33792786 5.48472494\n",
" 6.18738702 1.78741883 0.78195251 2.96658955 1.35827677 5.599321\n",
" 4.54850148 1.0939048 4.26158726 0.52100721 0.82318 4.76796961\n",
" 3.75795507 3.8526447 5.51100375 5.91023075 2.61494836 1.79908918\n",
" 2.65937756 5.53964148]\n",
"\n",
"-0.44791260077615314\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: 1.316011435623847\n",
" x: [ 2.329e+00 5.399e+00 ... 2.659e+00 5.540e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 5.551115123125783e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = 0.7235003672327549\n",
"The corresponding X is:\n",
"[2.56282915 5.63369524 5.58059887 4.049643 4.2021266 3.06866011\n",
" 6.01619635 1.52520776 4.35403161 0.33673958 0.32623161 1.2179545\n",
" 2.84001371 3.98956684 4.89632562 1.38303588 1.96194695 2.13182089\n",
" 0.29739166 1.77895165 3.29151585 3.54355374 4.49626674 0.95756626\n",
" 0.87103927 4.53068385 1.31051302 0.37103108 1.02961355 3.13342311\n",
" 5.65815319 2.24770604]\n",
"\n",
"-0.5994426600672451\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: 0.7235003672327549\n",
" x: [ 2.563e+00 5.634e+00 ... 5.658e+00 2.248e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 5.551115123125783e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = 0.34960914928810116\n",
"The corresponding X is:\n",
"[5.44143165 6.75955835 1.56836472 3.09522093 4.67873235 1.67071481\n",
" 0.3056494 0.65998337 1.02197668 5.21162959 0.43690354 3.56522934\n",
" 4.56033119 1.90736037 0.40863891 2.87007312 3.2516952 5.90360196\n",
" 1.99057799 5.20726456 0.74710237 6.03179202 3.80685028 0.03844391\n",
" 5.88580196 3.62233258 3.98723567 2.50591888 5.44020267 2.2792993\n",
" 5.57102303 4.46548617]\n",
"\n",
"-0.7087452725518989\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: 0.34960914928810116\n",
" x: [ 5.441e+00 6.760e+00 ... 5.571e+00 4.465e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 2.220446049250313e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = 0.10594558882184543\n",
"The corresponding X is:\n",
"[5.35675483 2.26629567 1.45430546 5.56758296 5.76309509 0.73239338\n",
" 5.1216998 3.03258872 4.33624828 1.93197674 0.5292902 3.32274987\n",
" 3.43247633 0.81490741 0.48060245 1.9944799 5.67519646 5.12534057\n",
" 0.06510627 2.52989834 6.1699519 0.94828957 5.91634548 1.5994961\n",
" 4.27902164 2.3129213 1.82353095 2.10634209 1.43740426 4.06988733\n",
" 0.59624074 4.93925418]\n",
"\n",
"-0.7760164293781545\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: 0.10594558882184543\n",
" x: [ 5.357e+00 2.266e+00 ... 5.962e-01 4.939e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 1.1102230246251565e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = -0.06473600797229297\n",
"The corresponding X is:\n",
"[6.07735568 0.18019501 0.20743128 4.15445985 3.59388894 5.10047555\n",
" 6.09938474 6.54707528 3.36251167 2.05475223 3.67078456 5.96010605\n",
" 2.58589996 5.2723619 3.26352977 2.47432334 3.50289983 2.06620525\n",
" 6.0946056 1.22751903 0.97320057 2.19564095 5.73174941 2.05127682\n",
" 5.73805165 3.84046105 1.84816963 2.1247504 3.11106736 2.44136052\n",
" 3.39002685 0.81596991]\n",
"\n",
"-0.8207034521437214\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: -0.06473600797229297\n",
" x: [ 6.077e+00 1.802e-01 ... 3.390e+00 8.160e-01]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 5.551115123125783e-17\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = -0.19562982094782935\n",
"The corresponding X is:\n",
"[-0.02184462 3.67041038 7.25918653 5.89799546 0.63583624 1.84214506\n",
" 2.84059837 5.31485182 1.6053784 0.04556618 0.32018993 -0.03884066\n",
" 0.69131496 0.24203727 1.97397262 3.59723495 0.43355775 2.30131056\n",
" 4.63482292 3.9857415 4.32320753 4.55388437 2.18753433 5.99034987\n",
" 2.50489913 0.90650534 4.82518088 2.32954849 2.29901832 5.33658863\n",
" 5.91246716 3.2405013 ]\n",
"\n",
"-0.8571013345978292\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: -0.19562982094782935\n",
" x: [-2.184e-02 3.670e+00 ... 5.912e+00 3.241e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 1.1102230246251565e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = -0.2833766309947055\n",
"The corresponding X is:\n",
"[ 3.1700088 5.05055456 1.2545611 4.28751811 0.6255103 1.67526577\n",
" 5.48201473 4.83820497 7.34880059 5.99705431 4.2502643 0.32066274\n",
" 0.41001404 0.27271241 4.15682546 4.22393693 4.35148115 0.64538137\n",
" 5.26288622 5.03810489 4.62426621 4.74997689 1.09603919 0.34752466\n",
" 1.8116275 0.7474807 5.31754143 4.11181763 1.58797998 5.6299796\n",
" 3.0109383 -0.19062772]\n",
"\n",
"-0.8713513097947054\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: -0.2833766309947055\n",
" x: [ 3.170e+00 5.051e+00 ... 3.011e+00 -1.906e-01]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 1.1102230246251565e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = -0.3527503628484244\n",
"The corresponding X is:\n",
"[3.90513622 4.61398739 5.92552705 1.99953405 4.82157369 1.35702441\n",
" 2.77701782 5.73612247 4.22710527 1.83463189 0.45796297 4.62509318\n",
" 0.98998668 0.11666217 3.0234641 4.54298546 0.14034033 4.15635797\n",
" 1.41257357 4.48719602 2.39365535 0.19672041 5.0763044 1.86357581\n",
" 3.657757 4.60298344 2.49769577 1.88086199 3.00108725 1.84475841\n",
" 5.24047385 4.91142914]\n",
"\n",
"-0.8819275737684243\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: -0.3527503628484244\n",
" x: [ 3.905e+00 4.614e+00 ... 5.240e+00 4.911e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 2.7755575615628914e-17\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = -0.4022181851996095\n",
"The corresponding X is:\n",
"[6.09453981 3.5109422 3.37216019 4.94732621 1.25662002 5.89645164\n",
" 5.06403334 2.68073141 4.40385083 1.13638366 1.73347762 6.82932871\n",
" 1.15265014 2.07145964 4.36520459 1.14960341 1.62288871 4.32315915\n",
" 5.45622821 0.93554005 3.17418483 0.47230243 1.31535502 5.77698726\n",
" 2.04927925 2.50663538 5.9706002 5.4984681 2.9421232 1.56636313\n",
" 1.09394523 4.62582 ]\n",
"\n",
"-0.8832883769450639\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: -0.4022181851996095\n",
" x: [ 6.095e+00 3.511e+00 ... 1.094e+00 4.626e+00]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"accuracy of Cholesky decomposition 1.1102230246251565e-16\n",
"Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
"Number of function values = 34 Least value of F = -0.44423031870708934\n",
"The corresponding X is:\n",
"[4.05765050e+00 3.99144950e+00 3.13287593e+00 3.28855137e+00\n",
" 4.32613515e+00 4.91104512e+00 1.86521867e+00 2.18822879e+00\n",
" 6.01336171e+00 1.82501276e+00 2.64830637e+00 5.53045823e+00\n",
" 2.36110093e+00 3.98821703e+00 4.69013438e-01 4.38996815e+00\n",
" 7.78103801e-04 1.72994378e+00 2.24970934e+00 1.11978200e+00\n",
" 2.24846445e+00 4.90745512e+00 5.38474921e+00 5.03587994e+00\n",
" 3.54297277e+00 4.78147533e+00 1.25990218e+00 1.99168068e+00\n",
" 5.89203503e+00 1.77673987e+00 5.37848357e+00 5.60245198e-01]\n",
"\n",
"-0.8852113278070892\n",
" message: Return from COBYLA because the objective function has been evaluated MAXFUN times.\n",
" success: False\n",
" status: 3\n",
" fun: -0.44423031870708934\n",
" x: [ 4.058e+00 3.991e+00 ... 5.378e+00 5.602e-01]\n",
" nfev: 34\n",
" maxcv: 0.0\n",
"All energies have been calculated\n"
]
}
],
"source": [
"from qiskit.primitives import BackendEstimatorV2\n",
"\n",
"estimator = BackendEstimatorV2(backend=backend_sim)\n",
"\n",
"distances_sim = np.arange(0.3, 1.3, 0.1)\n",
"vqe_energies_sim = []\n",
"vqe_elec_energies_sim = []\n",
"\n",
"for dist in distances_sim:\n",
" xx = dist\n",
"\n",
" # Random initial state and efficient_su2 ansatz\n",
" H = build_hamiltonian(xx)\n",
" ansatz = efficient_su2(H.num_qubits)\n",
" ansatz_isa = pm.run(ansatz)\n",
" x0 = 2 * np.pi * np.random.random(ansatz_isa.num_parameters)\n",
" H_isa = H.apply_layout(ansatz_isa.layout)\n",
" nuclear_repulsion = ham_terms(xx)[0]\n",
"\n",
" res = minimize(\n",
" cost_func,\n",
" x0,\n",
" args=(ansatz_isa, H_isa, estimator),\n",
" method=\"cobyla\",\n",
" options={\"maxiter\": 20, \"disp\": True},\n",
" )\n",
"\n",
" # Note this returns the total energy, and we are often interested in the electronic energy\n",
" tot_energy = getattr(res, \"fun\")\n",
" electron_energy = getattr(res, \"fun\") - nuclear_repulsion\n",
" print(electron_energy)\n",
" vqe_energies_sim.append(tot_energy)\n",
" vqe_elec_energies_sim.append(electron_energy)\n",
"\n",
" # Print all results\n",
" print(res)\n",
"\n",
"print(\"All energies have been calculated\")"
]
},
{
"cell_type": "code",
"execution_count": 39,
"id": "bc6164ca-6909-4780-8009-6dc274c66268",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"np.float64(1.2000000000000004)"
]
},
"execution_count": 39,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"xx"
]
},
{
"cell_type": "markdown",
"id": "633a2a89-950f-4d13-a56e-6adc079245ea",
"metadata": {},
"source": [
"Los resultados de esta salida se discuten más adelante en la sección de post-procesamiento; por ahora, basta con señalar que la simulación se ha realizado correctamente. Ahora está listo para funcionar en hardware real. Estableceremos la resiliencia en `1`, indicando que se utilizará la mitigación de errores TREX. Ahora que estamos trabajando con hardware real, utilizaremos Qiskit Runtime, y las primitivas Runtime. Observe que tanto el bucle for relacionado con la geometría como los múltiples ensayos variacionales se encuentran dentro de la sesión.\n",
"\n",
"Dado que hay costes y límites de tiempo asociados a las ejecuciones en hardware real, hemos reducido el número de pasos de geometría y pasos del optimizador a continuación. Asegúrese de adaptar estos pasos en función de sus objetivos de precisión y sus límites de tiempo.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "b09744f5-87a9-4744-b495-cf993e5ffcb3",
"metadata": {},
"outputs": [],
"source": [
"# To continue running on real hardware use\n",
"from qiskit_ibm_runtime import Session\n",
"from qiskit_ibm_runtime import EstimatorV2 as Estimator\n",
"from qiskit_ibm_runtime import EstimatorOptions\n",
"\n",
"estimator_options = EstimatorOptions(resilience_level=1, default_shots=2000)\n",
"\n",
"distances = np.arange(0.5, 0.9, 0.1)\n",
"vqe_energies = []\n",
"vqe_elec_energies = []\n",
"\n",
"with Session(backend=backend) as session:\n",
" estimator = Estimator(mode=session, options=estimator_options)\n",
"\n",
" for dist in distances:\n",
" xx = dist\n",
"\n",
" # Random initial state and efficient_su2 ansatz\n",
"\n",
" H = build_hamiltonian(xx)\n",
" ansatz = efficient_su2(H.num_qubits)\n",
" ansatz_isa = pm.run(ansatz)\n",
" H_isa = H.apply_layout(ansatz_isa.layout)\n",
" nuclear_repulsion = ham_terms(xx)[0]\n",
" x0 = 2 * np.pi * np.random.random(ansatz_isa.num_parameters)\n",
"\n",
" res = minimize(\n",
" cost_func,\n",
" x0,\n",
" args=(ansatz_isa, H_isa, estimator),\n",
" method=\"cobyla\",\n",
" options={\"maxiter\": 50, \"disp\": True},\n",
" )\n",
"\n",
" # Note this returns the total energy, and we are often interested in the electronic energy\n",
" tot_energy = getattr(res, \"fun\")\n",
" electron_energy = getattr(res, \"fun\") - nuclear_repulsion\n",
" print(electron_energy)\n",
" vqe_energies.append(tot_energy)\n",
" vqe_elec_energies.append(electron_energy)\n",
"\n",
" # Print all results\n",
" print(res)\n",
"\n",
"print(\"All energies have been calculated\")"
]
},
{
"cell_type": "markdown",
"id": "5b6f1491-88bd-4fb1-a1fe-2e29ffa33f17",
"metadata": {},
"source": [
"\n",
"\n",
"## Paso 4: Postprocesamiento\n",
"\n",
"Tanto para el simulador como para el hardware real, podemos representar gráficamente las energías del estado básico calculadas para cada distancia interatómica y ver dónde se alcanza la energía más baja. Esa debería ser la distancia interatómica encontrada en la naturaleza, y de hecho se aproxima. Se podría obtener una curva más suave probando otros ansaetze, optimizadores y ejecutando el cálculo varias veces en cada paso de geometría y promediando sobre varias condiciones iniciales aleatorias.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "e81f3ead-27ac-415e-a9f2-64a51d4b7aa3",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# Here we can plot the results from this simulation.\n",
"plt.plot(distances_sim, vqe_energies_sim, label=\"VQE Energy\")\n",
"plt.xlabel(\"Atomic distance (Angstrom)\")\n",
"plt.ylabel(\"Energy\")\n",
"plt.legend()\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "1efe0d76-3c6c-4369-9ace-fc890084e676",
"metadata": {},
"source": [
"Obsérvese que no es probable que el simple aumento del número de pasos de optimización mejore los resultados del simulador, ya que todas las optimizaciones convergen realmente a la tolerancia requerida en menos del número máximo de iteraciones.\n",
"\n",
"Los resultados del hardware real son comparables, aparte de un rango ligeramente diferente de valores muestreados.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "de8f53cf-9547-4578-b6cb-20d2b5602ee0",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.plot(distances, vqe_energies, label=\"VQE Energy\")\n",
"plt.xlabel(\"Atomic distance (Angstrom)\")\n",
"plt.ylabel(\"Energy\")\n",
"plt.legend()\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "7e732366-425a-40e4-b6f2-d26f11f7d86b",
"metadata": {},
"source": [
"Además de esperar una longitud de enlace H2 de 0.74 Angstrom, la energía total debería ser -1.17 Hartrees. Vemos que los resultados del hardware real se acercaron más a estos valores que los del simulador. Esto se debe probablemente a que el ruido estaba presente (o simulado) en ambos casos, pero sólo en el caso del hardware real se empleó la mitigación de errores.\n",
"\n",
"\n",
"\n",
"### Cerrando\n",
"\n",
"Con esto concluye nuestro curso sobre VQE para química cuántica. Si está interesado en comprender parte de la teoría de la información subyacente utilizada en la computación cuántica, consulte el curso de John Watrous sobre los [Fundamentos de la Información Cuántica](/learning/courses/basics-of-quantum-information). Para ver un ejemplo breve adicional de un flujo de trabajo VQE, consulte nuestro [tutorial Estimación de la energía en estado base de la cadena de Heisenberg con VQE](/docs/tutorials/spin-chain-vqe). O navegue por los [tutoriales](/docs/tutorials) y [cursos](/learning) para encontrar más material educativo sobre la última tecnología en computación cuántica.\n",
"\n",
"No olvides hacer el examen de este curso. Una puntuación del 80% o superior le otorgará una insignia Credly, que se le enviará automáticamente por correo electrónico. Gracias por formar parte de la red IBM Quantum®\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "0ac6d249-d210-479d-a792-c8b4e94b8b88",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"1.3.2\n",
"0.35.0\n"
]
}
],
"source": [
"import qiskit\n",
"import qiskit_ibm_runtime\n",
"\n",
"print(qiskit.version.get_version_info())\n",
"print(qiskit_ibm_runtime.version.get_version_info())"
]
},
{
"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"
}
},
"nbformat": 4,
"nbformat_minor": 2
}