Le solveur d'équations quantiques variationnelles (VQE)
Cette leçon présente le résolveur quantique variationnel, explique son importance en tant qu'algorithme fondamental de l'informatique quantique et explore ses forces et ses faiblesses. La VQE en elle-même, sans méthodes complémentaires, n'est probablement pas suffisante pour les calculs quantiques modernes à l'échelle de l'utilité. Elle est néanmoins importante en tant qu'archétype de méthode hybride classique-quantique et constitue une base importante sur laquelle sont construits de nombreux algorithmes plus avancés.
Cette vidéo donne un aperçu de la VQE et des facteurs qui influencent son efficacité. Le texte ci-dessous donne plus de détails et met en œuvre le VQE à l'aide de Qiskit.
1. Qu'est-ce que le VQE?
Le solveur d'équations aux valeurs propres quantique variationnel est un algorithme qui combine l'informatique classique et l'informatique quantique pour accomplir une tâche. Le calcul de la VQE comporte quatre éléments principaux :
- Un opérateur : Souvent un hamiltonien, que nous appellerons , qui décrit une propriété de votre système que vous souhaitez optimiser. Une autre façon de le dire est de chercher le vecteur propre de cet opérateur qui correspond à la valeur propre minimale. Nous appelons souvent ce vecteur propre l'"état fondamental".
- Un "ansatz" (mot allemand signifiant "approche") : il s'agit d'un circuit quantique qui prépare un état quantique approchant le vecteur propre recherché. En réalité, l'ansatz est une famille de circuits quantiques, car certaines des portes de l'ansatz sont paramétrées, c'est-à-dire qu'elles sont alimentées par un paramètre que nous pouvons faire varier. Cette famille de circuits quantiques peut préparer une famille d'états quantiques proches de l'état fondamental.
- Un estimateur : un moyen d'estimer la valeur attendue de l'opérateur sur l'état quantique variationnel actuel. Parfois, ce qui nous importe vraiment, c'est simplement cette valeur attendue, que l'on appelle une fonction de coût. Parfois, nous nous intéressons à une fonction plus complexe qui peut néanmoins être écrite à partir d'une ou plusieurs valeurs attendues.
- Un optimiseur classique : un algorithme qui fait varier les paramètres pour tenter de minimiser la fonction de coût.
Examinons chacun de ces éléments plus en détail.
1.1 L'opérateur (hamiltonien)
Au cœur d'un problème VQE se trouve un opérateur qui décrit un système d'intérêt. Nous supposerons ici que la valeur propre la plus basse et le vecteur propre correspondant de cet opérateur sont utiles à des fins scientifiques ou commerciales. Il peut s'agir par exemple d'un hamiltonien chimique décrivant une molécule, de sorte que la valeur propre la plus basse de l'opérateur correspond à l'énergie de l'état fondamental de la molécule, et que l'état propre correspondant décrit la géométrie ou la configuration électronique de la molécule. L'opérateur peut également décrire le coût d'un certain processus à optimiser, et les états propres peuvent correspondre à des itinéraires ou à des pratiques. Dans certains domaines, comme la physique, un "hamiltonien" fait presque toujours référence à un opérateur décrivant l'énergie d'un système physique. Mais en informatique quantique, il est courant de voir les opérateurs quantiques qui décrivent un problème commercial ou logistique également appelés "hamiltoniens". Nous adopterons cette convention ici.
La mise en correspondance d'un problème physique ou d'optimisation avec des qubits est généralement une tâche non triviale, mais ces détails ne font pas l'objet de ce cours. Une discussion générale sur la cartographie d'un problème sur un opérateur quantique peut être trouvée dans L'informatique quantique en pratique. Un examen plus détaillé de la conversion des problèmes de chimie en opérateurs quantiques est présenté dans Quantum Chemistry with VQE (Chimie quantique avec VQE ).
Pour les besoins de ce cours, nous supposerons que la forme de l'hamiltonien est connue. Par exemple, un hamiltonien pour une molécule d'hydrogène simple (sous certaines hypothèses d'espace actif, et en utilisant le cartographe de Jordan-Wigner) est :
from qiskit.quantum_info import SparsePauliOp
hamiltonian = SparsePauliOp(
[
"IIII",
"IIIZ",
"IZII",
"IIZI",
"ZIII",
"IZIZ",
"IIZZ",
"ZIIZ",
"IZZI",
"ZZII",
"ZIZI",
"YYYY",
"XXYY",
"YYXX",
"XXXX",
],
coeffs=[
-0.09820182 + 0.0j,
-0.1740751 + 0.0j,
-0.1740751 + 0.0j,
0.2242933 + 0.0j,
0.2242933 + 0.0j,
0.16891402 + 0.0j,
0.1210099 + 0.0j,
0.16631441 + 0.0j,
0.16631441 + 0.0j,
0.1210099 + 0.0j,
0.17504456 + 0.0j,
0.04530451 + 0.0j,
0.04530451 + 0.0j,
0.04530451 + 0.0j,
0.04530451 + 0.0j,
],
)Notez que dans l'hamiltonien ci-dessus, il y a des termes comme ZZII et YYYY qui ne sont pas commutés entre eux. En d'autres termes, pour évaluer ZZII, nous devrions mesurer l'opérateur Z de Pauli sur le qubit 3 (entre autres mesures). Mais pour évaluer YYYY, nous devons mesurer l'opérateur Y de Pauli sur ce même qubit, le qubit 3. Il existe une relation d'incertitude entre les opérateurs Y et Z sur le même qubit; nous ne pouvons pas mesurer ces deux opérateurs en même temps. Nous reviendrons sur ce point plus loin, et en fait tout au long du cours.
L'hamiltonien ci-dessus est un opérateur matriciel . Il n'est pas difficile de diagonaliser l'opérateur pour trouver sa valeur propre de plus faible énergie.
import numpy as np
A = np.array(hamiltonian)
eigenvalues, eigenvectors = np.linalg.eigh(A)
print("The ground state energy is ", min(eigenvalues), "hartrees")Output:
The ground state energy is -1.1459778447627311 hartrees
Les résolveurs d'équations classiques de force brute ne peuvent pas décrire les énergies ou les géométries de très grands systèmes d'atomes, comme les médicaments ou les protéines. VQE est l'une des premières tentatives d'exploitation de l'informatique quantique pour résoudre ce problème.
Dans cette leçon, nous rencontrerons des hamiltoniens beaucoup plus grands que ceux mentionnés ci-dessus. Mais il serait inutile de repousser les limites de ce que VQE peut faire, avant de présenter certains des outils les plus avancés qui peuvent compléter ou remplacer VQE, plus tard dans ce cours.
1.2 Ansatz
Le mot "ansatz" signifie "approche" en allemand. Le pluriel correct en allemand est "ansätze", bien que l'on trouve souvent "ansatzes" ou "ansatze". Dans le contexte de l'EQV, un ansatz est le circuit quantique que vous utilisez pour créer une fonction d'onde multi-qubit qui se rapproche le plus de l'état fondamental du système que vous étudiez, et qui produit donc la valeur d'espérance la plus faible de votre opérateur. Ce circuit quantique contiendra des paramètres variationnels (souvent rassemblés dans le vecteur de variables ).
Un ensemble initial de valeurs des paramètres variationnels est choisi. Nous appellerons l'opération unitaire de l'ansatz sur le circuit . Par défaut, tous les qubits des ordinateurs quantiques IBM® sont initialisés à l'état . Lorsque le circuit fonctionne, l'état des qubits est le suivant
Si tout ce dont nous avions besoin était l'énergie la plus basse (en utilisant le langage des systèmes physiques), nous pourrions l'estimer en mesurant simplement l'énergie plusieurs fois et en prenant la plus basse. Mais en général, nous voulons aussi la configuration qui produit l'énergie ou la valeur propre la plus basse. L'étape suivante consiste donc à estimer la valeur d'espérance de l'hamiltonien, ce qui est réalisé grâce à des mesures quantiques. Il y a beaucoup de choses qui entrent en ligne de compte. Mais nous pouvons comprendre ce processus qualitativement en notant que la probabilité de mesurer une énergie (toujours en utilisant le langage des systèmes physiques) est liée à la valeur d'espérance par :
La probabilité est également liée au chevauchement entre l'état propre et l'état actuel du système :
Ainsi, en effectuant de nombreuses mesures des opérateurs de Pauli qui composent notre hamiltonien, nous pouvons estimer la valeur espérée de l'hamiltonien dans l'état actuel du système . L'étape suivante consiste à faire varier les paramètres et à essayer de s'approcher davantage de l'état d'énergie le plus bas (état fondamental) du système. En raison des paramètres variationnels dans l'ansatz, on l'appelle parfois forme variationnelle.
Avant de passer à ce processus variationnel, notez qu'il est souvent utile de commencer votre état à partir d'un état de "bonne supposition". Il se peut que vous en sachiez suffisamment sur votre système pour faire une meilleure estimation initiale que . Par exemple, il est courant d'initialiser les qubits à l'état Hartree-Fock dans les applications chimiques. Cette hypothèse de départ, qui ne contient aucun paramètre variationnel, est appelée état de référence. Appelons le circuit quantique utilisé pour créer l'état de référence . Chaque fois qu'il devient important de distinguer l'état de référence du reste de l'ansatz, utilisons : De manière équivalente
1.3 estimateur
Nous avons besoin d'un moyen d'estimer la valeur d'espérance de notre hamiltonien dans un état variationnel particulier . Si nous pouvions mesurer directement l'opérateur entier , il suffirait d'effectuer de nombreuses mesures (disons ) et de calculer la moyenne des valeurs mesurées :
Ici, le symbole nous rappelle que cette valeur attendue ne serait précisément correcte que dans la limite de . Mais avec des milliers de mesures effectuées sur un circuit, l'erreur d'échantillonnage de la valeur attendue est assez faible. D'autres considérations, telles que le bruit, posent problème pour les calculs très précis.
Cependant, il n’est généralement pas possible de mesurer l’ e en une seule fois. L’ e peut contenir plusieurs opérateurs de Pauli X, Y et Z non commutatifs. L'hamiltonien doit donc être décomposé en groupes d'opérateurs pouvant être mesurés simultanément; chaque groupe doit être évalué séparément, puis les résultats combinés pour obtenir une valeur d'espérance. Nous y reviendrons plus en détail lors de la prochaine leçon, lorsque nous aborderons la question de la mise à l'échelle des approches classiques et quantiques. Cette complexité des mesures est l'une des raisons pour lesquelles nous avons besoin d'un code très performant pour réaliser ce type d'estimation. Dans cette leçon et par la suite, nous utiliserons l'estimateur primitif « IBM Quantum » à cette fin.
1.4 Optimiseurs classiques
Un optimiseur classique est un algorithme classique conçu pour trouver les extrema d'une fonction cible (généralement un minimum). Ils parcourent l'espace des paramètres possibles à la recherche d'un ensemble qui minimise une fonction donnée. Elles peuvent être classées en deux grandes catégories : les méthodes basées sur le gradient, qui utilisent les informations relatives au gradient, et les méthodes sans gradient, qui fonctionnent comme des optimiseurs à boîte noire. Le choix d'un optimiseur classique peut avoir un impact significatif sur les performances d'un algorithme, en particulier en présence de bruit dans le matériel quantique. Les optimiseurs les plus populaires dans ce domaine sont Adam, AMSGrad et SPSA, qui ont donné des résultats prometteurs dans les environnements bruyants. Les optimiseurs plus traditionnels comprennent COBYLA et SLSQP.
Un processus courant (démontré dans la section 3.3 ) consiste à utiliser l'un de ces algorithmes comme méthode à l'intérieur d'un minimiseur tel que la fonction minimize de scipy. Elle prend comme arguments
- Une fonction à minimiser. Il s'agit souvent de la valeur énergétique attendue. Mais ces fonctions sont généralement appelées "fonctions de coût".
- Un ensemble de paramètres à partir desquels la recherche peut commencer. Souvent appelé ou .
- Arguments, y compris les arguments de la fonction de coût. Dans le cadre de l'informatique quantique avec Qiskit, ces arguments comprendront l'ansatz, l'hamiltonien et la primitive Estimator, qui sera abordée plus en détail dans la sous-section suivante.
- Une "méthode" de minimisation. Il s'agit de l'algorithme spécifique utilisé pour rechercher l'espace des paramètres. C'est ici que l'on spécifie, par exemple, COBYLA ou SLSQP.
- Options. Les options disponibles peuvent varier selon la méthode. Mais un exemple que pratiquement toutes les méthodes incluraient est le nombre maximum d'itérations de l'optimiseur avant de terminer la recherche : 'maxiter'.
À chaque étape itérative, la valeur espérée de l'hamiltonien est estimée en effectuant de nombreuses mesures. Cette énergie estimée est renvoyée par la fonction de coût, et le minimiseur met à jour les informations dont il dispose sur le paysage énergétique. Ce que fait exactement l'optimiseur pour choisir l'étape suivante varie d'une méthode à l'autre. Certains utilisent des gradients et sélectionnent la direction de la descente la plus raide. D'autres peuvent prendre en compte le bruit et exiger que le coût diminue dans une large mesure avant d'accepter que l'énergie réelle diminue dans cette direction.
# Example syntax for minimization
# from scipy.optimize import minimize
# res = minimize(cost_func, x0, args=(ansatz, hamiltonian, estimator), method="cobyla",
# options={'maxiter': 200})1.5 Le principe variationnel
Dans ce contexte, le principe de variation est très important; il stipule qu'aucune fonction d'onde variationnelle ne peut produire une valeur d'espérance d'énergie (ou de coût) inférieure à celle produite par la fonction d'onde de l'état fondamental. Mathématiquement,
Ceci est facile à vérifier si nous notons que l'ensemble de tous les états propres de forment une base complète pour l'espace de Hilbert. En d'autres termes, tout état, et en particulier , peut être écrit comme une somme pondérée (normalisée) de ces états propres de :
où sont des constantes à déterminer et . Nous laissons cet exercice au lecteur. Mais notez l'implication : l'état variationnel qui produit la valeur d'espérance de l'énergie la plus basse est la meilleure estimation de l'état fondamental réel.
Vérifiez votre compréhension
Vérifiez mathématiquement que pour tout état variationnel .
En utilisant l'expansion donnée de l'état variationnel en termes d'états propres énergétiques,
nous pouvons écrire la valeur espérée de l'énergie variationnelle comme suit
Pour tous les coefficients . On peut donc écrire
2. Comparaison avec le flux de travail classique
Supposons que nous nous intéressons à une matrice de N lignes et N colonnes. Supposons que votre matrice soit si grande qu'une diagonalisation exacte ne soit pas envisageable. Supposons en outre que vous en sachiez suffisamment sur votre problème pour pouvoir émettre quelques hypothèses sur la structure globale de l'état propre cible, et que vous souhaitiez sonder les états similaires à votre hypothèse initiale pour voir si votre coût/énergie peut être encore réduit. Il s'agit d'une approche variationnelle, qui est utilisée lorsque la diagonalisation exacte n'est pas possible.
2.1 Flux de travail classique
En utilisant un ordinateur classique, cela fonctionnerait comme suit :
- Faites un état supposé, avec quelques paramètres que vous ferez varier : . Bien que cette supposition initiale puisse être aléatoire, cela n'est pas conseillé. Nous voulons utiliser la connaissance du problème en question pour adapter notre supposition autant que possible.
- Calculer la valeur d'espérance de l'opérateur avec le système dans cet état :
- Modifier les paramètres variationnels et répéter : .
- Utilisez les informations accumulées sur le paysage des états possibles dans votre sous-espace variationnel pour faire des suppositions de plus en plus précises et vous rapprocher de l'état cible. Le principe variationnel garantit que notre état variationnel ne peut pas produire une valeur propre inférieure à celle de l'état fondamental cible. Ainsi, plus la valeur d'espérance est faible, meilleure est notre approximation de l'état fondamental :
Examinons la difficulté de chaque étape de cette approche. Le réglage ou la mise à jour des paramètres est facile à calculer; la difficulté réside dans la sélection de paramètres initiaux utiles et physiquement motivés. L'utilisation des informations accumulées lors des itérations précédentes pour mettre à jour les paramètres de manière à s'approcher de l'état fondamental est une tâche non triviale. Mais il existe des algorithmes d'optimisation classiques qui permettent de le faire de manière assez efficace. Cette optimisation classique n'est coûteuse que parce qu'elle peut nécessiter de nombreuses itérations; dans le pire des cas, le nombre d'itérations peut augmenter de façon exponentielle avec N. L'étape la plus coûteuse en termes de calcul est très certainement le calcul de la valeur d'espérance de votre matrice à partir d'un état donné :
La matrice doit agir sur le vecteur -élément, ce qui correspond à : opérations de multiplication dans le pire des cas. Cette opération doit être effectuée à chaque itération des paramètres. Pour les matrices de très grande taille, le coût de calcul est élevé.
2.2 Flux de travail quantique et groupes de Pauli de commutation
Imaginez maintenant que cette partie du calcul soit confiée à un ordinateur quantique. Au lieu de calculer cette valeur d'espérance, vous l'estimez en préparant l'état sur l'ordinateur quantique à l'aide de votre ansatz variationnel, puis en effectuant des mesures.
Cela peut sembler plus facile que cela ne l'est. n'est généralement pas facile à mesurer. Par exemple, il pourrait être composé de nombreux opérateurs X, Y et Z de Pauli non commutatifs. Mais peut être écrit comme une combinaison linéaire de termes, , dont chacun est facilement mesurable (par exemple, des opérateurs de Pauli ou des groupes d'opérateurs de Pauli commutés dans le sens du qubit). La valeur espérée de sur un certain état est la somme pondérée des valeurs espérées des termes constitutifs . Cette expression est valable pour n'importe quel état , mais nous l'utiliserons spécifiquement avec nos états variationnels .
où est une chaîne de Pauli comme IZZX…XIYX, ou plusieurs chaînes de ce type qui commutent entre elles. Une description de la valeur d'espérance plus proche des réalités de la mesure sur les ordinateurs quantiques est donc la suivante
Et dans le contexte de notre fonction d'onde variationnelle :
Chacun des termes peut être mesuré fois, ce qui donne des échantillons de mesure avec et renvoie une valeur attendue et un écart type . Nous pouvons additionner ces termes et propager les erreurs à travers la somme pour obtenir une valeur d'espérance globale et un écart type .
Cela ne nécessite pas de multiplication à grande échelle, ni de processus qui s'échelonne nécessairement comme . Au lieu de cela, il faut effectuer de multiples mesures sur l'ordinateur quantique. Si vous n'avez pas besoin d'un grand nombre d'entre eux, cette approche peut s'avérer efficace. C'est la partie quantique de VQE.
Mais examinons les raisons pour lesquelles cela pourrait ne pas être efficace. L'une des raisons d'effectuer de nombreuses mesures est de réduire l'incertitude statistique de vos estimations, pour des calculs de très haute précision. Une autre raison est le nombre de cordes de Pauli nécessaires pour couvrir l'ensemble de la matrice. Comme les matrices de Pauli (plus l'identité : X, Y, Z et I) couvrent l'espace de tous les opérateurs d'une dimension donnée, nous avons la garantie de pouvoir écrire la matrice qui nous intéresse comme une somme pondérée d'opérateurs de Pauli, comme nous l'avons fait précédemment.
où est une chaîne de Pauli agissant sur tous les qubits décrivant votre système comme IZZX…XIYX, ou plusieurs chaînes de ce type qui commutent entre elles. Rappelons que Qiskit utilise la notation little endian, dans laquelle l'opérateur de Pauli de la droite agit sur le qubit . Nous pouvons donc mesurer notre opérateur en mesurant une série d'opérateurs de Pauli.
Mais nous ne pouvons pas mesurer tous ces opérateurs de Pauli en même temps. Les opérateurs de Pauli (à l'exception de I) ne commutent pas entre eux s'ils sont associés au même qubit. Par exemple, nous pouvons mesurer IZIZ et ZZXZ simultanément, car nous pouvons mesurer I et Z simultanément pour le troisième qubit, et nous pouvons connaître I et X simultanément pour le premier qubit. Mais nous ne pouvons pas mesurer ZZZZ et ZZZX simultanément, car Z et X ne commutent pas, et agissent tous deux sur le qubit de rang 0. Les lecteurs avertis se souviendront peut-être que deux ensembles d'opérateurs de Pauli peuvent commuter en tant qu'ensemble, même si les mesures effectuées sur chaque qubit pris individuellement ne commutent pas. L'estimateur suppose des mesures de Pauli sous forme de produit tensoriel (via des rotations de base), correspondant à des opérateurs de regroupement qui commutent au niveau des qubits. Ainsi, pour estimer simultanément deux chaînes (A et B) d'opérateurs de Pauli à l'aide d'Estimator, les opérateurs de Pauli de chaque qubit dans A et B doivent commuter. Cela signifie que nous ne pouvons pas non plus mesurer ZZZZ et ZZXX simultanément.
Nous décomposons donc notre matrice en une somme de Paulis agissant sur différents qubits. Certains éléments de cette somme peuvent être mesurés en une seule fois; nous appelons cela un groupe de Paulis commutatifs. En fonction du nombre de mandats non pendulaires, nous pourrions avoir besoin de plusieurs groupes de ce type. Appelons le nombre de ces groupes de cordes de Pauli commutatives . Si est petit, cela pourrait bien fonctionner. Si a des millions de groupes, cela ne sera pas utile.
Les processus nécessaires au calcul de la valeur attendue sont regroupés dans la primitive « Estimator » de l' IBM Quantum. Pour en savoir plus sur Estimator, consultez la documentation de référence de l'API dans la documentation d' IBM Quantum®. On peut bien sûr utiliser directement la fonction « Estimator », mais celle-ci renvoie bien plus que la simple valeur propre d'énergie la plus faible. Par exemple, elle renvoie également des informations sur l'erreur-type d'ensemble. Ainsi, dans le cadre des problèmes de minimisation, on retrouve souvent un « estimateur » au sein d'une fonction de coût. Pour en savoir plus sur les données d'entrée et de sortie d'Estimator, consultez ce guide dans la documentation d' IBM Quantum.
Vous enregistrez la valeur d'espérance (ou la fonction de coût) pour l'ensemble des paramètres utilisés dans votre État, puis vous mettez à jour les paramètres. Au fil du temps, vous pouvez utiliser les valeurs d'espérance ou les valeurs de la fonction de coût que vous avez estimées pour obtenir une approximation du gradient de votre fonction de coût dans le sous-espace des états échantillonnés par votre ansatz. Il existe des optimiseurs classiques basés sur le gradient et des optimiseurs classiques sans gradient. Tous deux souffrent de problèmes potentiels de formation, tels que de multiples minima locaux et de vastes régions de l'espace des paramètres avec un gradient proche de zéro, appelées plateaux stériles.
2.3 Facteurs qui déterminent le coût de calcul
VQE ne résoudra pas tous vos problèmes de chimie quantique les plus difficiles. Non. Mais il ne s'agit pas d'être meilleur dans tous les calculs. Nous avons déplacé ce qui détermine le coût de calcul.
Nous sommes passés d'un processus dont la complexité dépend uniquement de la dimension de la matrice à un processus qui dépend de la précision requise et du nombre d'opérateurs de Pauli non commutatifs qui composent la matrice. Ce dernier point n'a pas d'analogue dans l'informatique classique.
Sur la base de ces dépendances, ce processus peut s'avérer utile pour les matrices peu denses ou les matrices comportant peu de chaînes de Pauli non commutatives. C'est le cas des systèmes de spins en interaction, par exemple. Pour les matrices denses, il peut être moins utile. Nous savons par exemple que les systèmes chimiques ont souvent des hamiltoniens qui impliquent des centaines, des milliers, voire des millions de chaînes de Pauli. Des travaux intéressants ont été menés pour réduire ce nombre de termes. Mais les systèmes chimiques peuvent être mieux adaptés à certains des autres algorithmes que nous aborderons dans ce cours.
Vérifiez votre compréhension
Considérons un hamiltonien sur quatre qubits qui contient les termes :
IIXX, IIXZ, IIZZ, IZXZ, IXXZ, ZZXZ, XZXZ, ZIXZ, ZZZZ, XXXX
Vous souhaitez classer ces termes en groupes de manière à ce que tous les termes d'un groupe puissent être mesurés simultanément. Quel est le plus petit nombre de groupes de ce type que l'on puisse constituer de manière à ce que tous les termes soient pris en compte?
On peut le faire en 4 groupes. Il convient de noter que ces solutions ne sont généralement pas uniques.
IIXX,XXXX,IIZZ,ZZZZIIXZ,IZXZ,ZIXZ,ZZXZIXXZXZXZ
Selon vous, qu'est-ce qui rend généralement difficile la chimie quantique avec la méthode VQE : le nombre de termes dans l'hamiltonien, ou la recherche d'une bonne approximation?
Il s'avère qu'il existe des réponses très optimisées pour les contextes chimiques. Le nombre de termes dans l'hamiltonien, et donc le nombre de mesures nécessaires, posent généralement plus de problèmes.
3. Exemple d'hamiltonien
Mettons cet algorithme en pratique en utilisant une petite matrice hamiltonienne afin de voir ce qui se passe à chaque étape. Nous utiliserons le cadre des modèles Qiskit :
- Étape 1 : Cartographier le problème en circuits et opérateurs quantiques - Étape 2 : Optimisation pour le matériel cible - Étape 3 : Exécution sur le matériel cible - Étape 4 : Post-traitement des résultats
3.1 Étape 1 : Transposer le problème en circuits quantiques et opérateurs
Nous utiliserons celle définie ci-dessus dans le contexte de la chimie. Nous commençons par quelques importations générales.
# General imports
import numpy as np
# SciPy minimizer routine
from scipy.optimize import minimize
# Plotting functions
import matplotlib.pyplot as pltLà encore, nous supposons que l'hamiltonien en question est connu. Nous utiliserons ici un hamiltonien extrêmement petit, car les autres méthodes abordées dans ce cours seront plus efficaces pour résoudre des problèmes plus importants.
from qiskit.quantum_info import SparsePauliOp
import numpy as np
hamiltonian = SparsePauliOp.from_list(
[("YZ", 0.3980), ("ZI", -0.3980), ("ZZ", -0.0113), ("XX", 0.1810)]
)
A = np.array(hamiltonian)
eigenvalues, eigenvectors = np.linalg.eigh(A)
print("The ground state energy is ", min(eigenvalues))Output:
The ground state energy is -0.702930394459531
Il existe de nombreux choix d'ansatz préfabriqués à Qiskit. Nous utiliserons efficient_su2.
# Pre-defined ansatz circuit and operator class for Hamiltonian
from qiskit.circuit.library import efficient_su2
# Note that it is more common to place initial 'h' gates outside the ansatz.
# Here we specifically wanted this layer structure.
ansatz = efficient_su2(
hamiltonian.num_qubits, su2_gates=["h", "rz", "y"], entanglement="circular", reps=1
)
num_params = ansatz.num_parameters
print("This circuit has ", num_params, "parameters")
ansatz.decompose().draw("mpl", style="iqp")Output:
This circuit has 4 parameters
Des réponses différentes auront des structures d'enchevêtrement différentes et des portes de rotation différentes. Celle présentée ici utilise des portes CNOT pour l'enchevêtrement, et des portes Y et des portes RZ paramétrées pour les rotations. Notez la taille de cet espace de paramètres; cela signifie que nous devons minimiser la fonction de coût sur 4 variables (les paramètres des portes RZ). Ce système peut être étendu, mais pas indéfiniment. L'exécution d'un problème similaire sur 4 qubits, en utilisant les 3 répétitions par défaut pour efficient_su2 , donne 16 paramètres variationnels.
3.2 Étape 2 : Optimisation pour le matériel cible
L'ansatz a été écrit en utilisant des portes familières, mais notre circuit doit être transposé pour utiliser les portes de base qui peuvent être mises en œuvre sur chaque ordinateur quantique. Nous sélectionnons le backend le moins sollicité.
# runtime imports
from qiskit_ibm_runtime import QiskitRuntimeService, Session
from qiskit_ibm_runtime import EstimatorV2 as Estimator
# To run on hardware, select the backend with the fewest number of jobs in the queue
service = QiskitRuntimeService()
backend = service.least_busy(operational=True, simulator=False)
print(backend)Output:
<IBMBackend('ibm_torino')>
Nous pouvons maintenant transpiler notre circuit pour ce matériel et visualiser notre ansatz transpilé.
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
target = backend.target
pm = generate_preset_pass_manager(target=target, optimization_level=3)
ansatz_isa = pm.run(ansatz)
ansatz_isa.draw(output="mpl", idle_wires=False, style="iqp")Output:
Notez que les portes utilisées ont changé et que les qubits de notre circuit abstrait ont été mis en correspondance avec des qubits numérotés différemment sur l'ordinateur quantique. Nous devons cartographier notre hamiltonien de manière identique pour que nos résultats soient significatifs.
hamiltonian_isa = hamiltonian.apply_layout(layout=ansatz_isa.layout)3.3 Étape 3 : Exécuter sur le matériel cible
3.3.1 Communication des valeurs
Nous définissons ici une fonction de coût qui prend comme arguments les structures que nous avons construites lors des étapes précédentes : les paramètres, l'ansatz et l'hamiltonien. Il utilise également Estimator, que nous n'avons pas encore défini. Nous intégrons du code permettant de suivre l'évolution de notre fonction de coût, afin de pouvoir vérifier son comportement en matière de convergence.
def cost_func(params, ansatz, hamiltonian, estimator):
"""Return estimate of energy from Estimator
Parameters:
params (ndarray): Array of ansatz parameters
ansatz (QuantumCircuit): Parameterized ansatz circuit
hamiltonian (SparsePauliOp): Operator representation of Hamiltonian
estimator (EstimatorV2): Estimator primitive instance
cost_history_dict: Dictionary for storing intermediate results
Returns:
float: Energy estimate
"""
pub = (ansatz, [hamiltonian], [params])
result = estimator.run(pubs=[pub]).result()
energy = result[0].data.evs[0]
cost_history_dict["iters"] += 1
cost_history_dict["prev_vector"] = params
cost_history_dict["cost_history"].append(energy)
print(f"Iters. done: {cost_history_dict['iters']} [Current cost: {energy}]")
return energy
cost_history_dict = {
"prev_vector": None,
"iters": 0,
"cost_history": [],
}Il est très avantageux de pouvoir choisir les valeurs initiales des paramètres en fonction de la connaissance du problème et des caractéristiques de l'état cible. Nous ne ferons aucune hypothèse sur ces connaissances et utiliserons des valeurs initiales aléatoires.
x0 = 2 * np.pi * np.random.random(num_params)# This required 13 min, 20 s QPU time on an Eagle processor, 28 min total time.
with Session(backend=backend) as session:
estimator = Estimator(mode=session)
estimator.options.default_shots = 10000
res = minimize(
cost_func,
x0,
args=(ansatz_isa, hamiltonian_isa, estimator),
method="cobyla",
options={"maxiter": 50},
)Output:
Iters. done: 1 [Current cost: 0.010575798722044727]
Iters. done: 2 [Current cost: 0.004040015974440895]
Iters. done: 3 [Current cost: 0.0020213258785942503]
Iters. done: 4 [Current cost: 0.18723082446726014]
Iters. done: 5 [Current cost: -0.2746792152068885]
Iters. done: 6 [Current cost: -0.3094547651648519]
Iters. done: 7 [Current cost: -0.05281985428356641]
Iters. done: 8 [Current cost: 0.00808560303514377]
Iters. done: 9 [Current cost: -0.0014821685303514388]
Iters. done: 10 [Current cost: -0.004759824281150161]
Iters. done: 11 [Current cost: 0.09942328705995292]
Iters. done: 12 [Current cost: 0.01092366214057508]
Iters. done: 13 [Current cost: 0.05017497496069776]
Iters. done: 14 [Current cost: 0.13028868414310696]
Iters. done: 15 [Current cost: 0.013747803514376994]
Iters. done: 16 [Current cost: 0.2583072432944498]
Iters. done: 17 [Current cost: -0.14422125655131562]
Iters. done: 18 [Current cost: -0.0004950150347678081]
Iters. done: 19 [Current cost: 0.00681082268370607]
Iters. done: 20 [Current cost: -0.0023377795527156544]
Iters. done: 21 [Current cost: 0.6027665591169237]
Iters. done: 22 [Current cost: 0.00596641373801917]
Iters. done: 23 [Current cost: -0.008318769968051117]
Iters. done: 24 [Current cost: -0.00026683306709265246]
Iters. done: 25 [Current cost: -0.007648222843450479]
Iters. done: 26 [Current cost: 0.004121086261980831]
Iters. done: 27 [Current cost: -0.004075019968051117]
Iters. done: 28 [Current cost: -0.004419369009584665]
Iters. done: 29 [Current cost: 0.213185460054037]
Iters. done: 30 [Current cost: -0.06505919572162797]
Iters. done: 31 [Current cost: -0.5334241316590271]
Iters. done: 32 [Current cost: 0.00218370607028754]
Iters. done: 33 [Current cost: 0.09579352143666908]
Iters. done: 34 [Current cost: -0.009274800319488819]
Iters. done: 35 [Current cost: -0.44395141360688106]
Iters. done: 36 [Current cost: 0.011747104632587858]
Iters. done: 37 [Current cost: -0.003344149361022364]
Iters. done: 38 [Current cost: 0.19138183916486304]
Iters. done: 39 [Current cost: 0.013513931813145209]
Nous pouvons examiner les résultats bruts.
resOutput:
message: Return from COBYLA because the trust region radius reaches its lower bound.
success: True
status: 0
fun: -0.5334241316590271
x: [ 1.024e+00 6.459e+00 3.625e+00 4.007e+00]
nfev: 39
maxcv: 0.0
3.4 Étape 4 : Post-traitement des résultats
Si la procédure se termine correctement, les valeurs de notre dictionnaire doivent être égales au vecteur de solution et au nombre total d'évaluations de la fonction, respectivement. Ceci est facile à vérifier :
cost_history_dictOutput:
{'prev_vector': array([1.02397956, 6.45886604, 3.62479262, 4.00744128]),
'iters': 39,
'cost_history': [np.float64(0.010575798722044727),
np.float64(0.004040015974440895),
np.float64(0.0020213258785942503),
np.float64(0.18723082446726014),
np.float64(-0.2746792152068885),
np.float64(-0.3094547651648519),
np.float64(-0.05281985428356641),
np.float64(0.00808560303514377),
np.float64(-0.0014821685303514388),
np.float64(-0.004759824281150161),
np.float64(0.09942328705995292),
np.float64(0.01092366214057508),
np.float64(0.05017497496069776),
np.float64(0.13028868414310696),
np.float64(0.013747803514376994),
np.float64(0.2583072432944498),
np.float64(-0.14422125655131562),
np.float64(-0.0004950150347678081),
np.float64(0.00681082268370607),
np.float64(-0.0023377795527156544),
np.float64(0.6027665591169237),
np.float64(0.00596641373801917),
np.float64(-0.008318769968051117),
np.float64(-0.00026683306709265246),
np.float64(-0.007648222843450479),
np.float64(0.004121086261980831),
np.float64(-0.004075019968051117),
np.float64(-0.004419369009584665),
np.float64(0.213185460054037),
np.float64(-0.06505919572162797),
np.float64(-0.5334241316590271),
np.float64(0.00218370607028754),
np.float64(0.09579352143666908),
np.float64(-0.009274800319488819),
np.float64(-0.44395141360688106),
np.float64(0.011747104632587858),
np.float64(-0.003344149361022364),
np.float64(0.19138183916486304),
np.float64(0.013513931813145209)]}
fig, ax = plt.subplots()
x = np.linspace(0, 10, 50)
# Define the constant function
constant = -0.7029
y_constant = np.full_like(x, constant)
ax.plot(
range(cost_history_dict["iters"]), cost_history_dict["cost_history"], label="VQE"
)
ax.set_xlabel("Iterations")
ax.set_ylabel("Cost")
ax.plot(y_constant, label="Target")
plt.legend()
plt.draw()Output:
IBM Quantum propose d'autres offres de perfectionnement liées à la VQE. Si vous êtes prêt à mettre VQE en pratique, consultez notre tutoriel : Estimation de l'énergie de l'état fondamental de la chaîne de Heisenberg avec VQE. Si vous souhaitez plus d'informations sur la création d'hamiltoniens moléculaires, consultez cette leçon dans notre cours sur la chimie quantique avec VQE. Si vous souhaitez mieux comprendre le fonctionnement des algorithmes variationnels tels que VQE, nous vous recommandons le cours Variational Algorithm Design.
Vérifiez votre compréhension
Dans cette section, nous avons calculé l'énergie de l'état fondamental à partir d'un hamiltonien. Si nous voulions appliquer ce principe à la détermination de la géométrie d'une molécule, par exemple, comment pourrions-nous l'étendre?
Nous devrions introduire des variables pour l'espacement inter-atomique et les angles entre les liaisons. Nous devrions les modifier. Pour chaque variation de ces éléments, nous produirons un nouvel hamiltonien (puisque les opérateurs décrivant l'énergie dépendent certainement de la géométrie). Pour chaque hamiltonien de ce type produit et mis en correspondance avec des qubits, nous devrions procéder à une optimisation comme celle effectuée ci-dessus. Parmi tous ces nombreux problèmes d'optimisation convergents, la géométrie qui produit l'énergie la plus faible serait celle adoptée par la nature. Cette opération est un peu plus complexe que celle présentée ci-dessus. Un tel calcul est effectué pour la molécule la plus simple, , ici.
4. Relation entre la VQE et d'autres méthodes
Dans cette section, nous passerons en revue les avantages et les inconvénients de l'approche originale de l'EQV et soulignerons ses liens avec d'autres algorithmes plus récents.
4.1 Les forces et les faiblesses du VQE
Certains points forts ont déjà été soulignés. Ils comprennent :
- Adaptation au matériel moderne : Certains algorithmes quantiques requièrent des taux d'erreur beaucoup plus faibles, ce qui les rapproche d'une tolérance aux pannes à grande échelle. Ce n'est pas le cas de VQE, qui peut être mis en œuvre sur les ordinateurs quantiques actuels.
- Circuits peu profonds : La VQE utilise souvent des circuits quantiques relativement peu profonds. La VQE est donc moins sensible aux erreurs de grille accumulées et se prête à de nombreuses techniques d'atténuation des erreurs. Bien entendu, les circuits ne sont pas toujours peu profonds; cela dépend de l'ansatz utilisé.
- Polyvalence : L'EQV peut (en principe) être appliquée à tout problème qui peut être considéré comme un problème de valeurs propres ou de vecteurs propres. Il existe de nombreuses mises en garde qui font que l'EQV n'est pas pratique ou désavantageuse pour certains problèmes. Certaines d'entre elles sont récapitulées ci-dessous.
Certaines faiblesses de l'EQV et les problèmes pour lesquels elle n'est pas pratique ont également été décrits ci-dessus. Exemples :
- Nature heuristique : La VQE ne garantit pas la convergence vers l'énergie correcte de l'état fondamental, car ses performances dépendent du choix de l'ansatz et des méthodes d'optimisation [1-2]. Si l'on choisit un ansatz médiocre qui n'offre pas l'intrication nécessaire à l'état fondamental souhaité, aucun optimiseur classique ne peut atteindre cet état fondamental.
- Nombreux paramètres potentiels : Un ansatz très expressif peut avoir tellement de paramètres que les itérations de minimisation prennent beaucoup de temps.
- Coût de calcul élevé : dans le cadre de la VQE, on utilise un estimateur pour estimer la valeur attendue de chaque terme de l'hamiltonien. La plupart des hamiltoniens d'intérêt comportent des termes qui ne peuvent pas être estimés simultanément. Cela peut rendre la méthode VQE très gourmande en ressources pour les grands systèmes dotés d'hamiltoniens complexes [1].
- Effets du bruit : Lorsque l'optimiseur classique recherche un minimum, des calculs bruyants peuvent l'embrouiller et l'éloigner du véritable minimum ou retarder sa convergence. Une solution possible consiste à tirer parti des techniques de pointe en matière d'atténuation et de suppression des erreurs [2-3] sur le site IBM.
- Plateaux stériles : Ces régions où les gradients s'évanouissent [2-3] existent même en l'absence de bruit, mais le bruit les rend plus problématiques, car la variation des valeurs attendues due au bruit peut être plus importante que la variation due à la mise à jour des paramètres dans ces régions stériles.
4.2 Relation avec d'autres approches
Adapt-VQE
L'algorithme ADAPT-VQE (Adaptive Derivative-Assembled Pseudo-Trotter Variational Quantum Eigensolver) est une amélioration de l'algorithme VQE original, conçu pour améliorer l'efficacité, la précision et l'évolutivité des simulations quantiques, en particulier dans le domaine de la chimie quantique.
L'algorithme VQE original décrit tout au long de cette leçon utilise un ansatz fixe prédéfini pour approximer l'état fondamental du système. Dans notre cas, nous avons utilisé efficient_su2, avec une seule répétition, en utilisant les portes de rotation Y et RZ. Bien que les paramètres des portes RZ aient changé, la structure de cet ansatz et les portes utilisées n'ont pas changé.
ADAPT-VQE répond aux limitations de VQE par la construction d'un ansatz adaptatif. Au lieu de commencer avec un ansatz fixe, ADAPT-VQE construit dynamiquement l'ansatz de manière itérative. À chaque étape, il sélectionne l'opérateur d'un ensemble prédéfini (comme les opérateurs d'excitation fermionique) qui présente le plus grand gradient par rapport à l'énergie. Cela garantit que seuls les opérateurs les plus importants sont ajoutés, ce qui permet d'obtenir un ansatz compact et efficace [4-6]. Cette approche peut avoir plusieurs effets bénéfiques :
- Réduction de la profondeur du circuit : En développant l'ansatz de manière incrémentielle et en se concentrant uniquement sur les opérateurs nécessaires, ADAPT-VQE minimise les opérations de porte par rapport aux approches VQE traditionnelles [5,7].
- Précision améliorée : La nature adaptative permet à ADAPT-VQE de récupérer davantage d'énergie de corrélation à chaque étape, ce qui la rend particulièrement efficace pour les systèmes fortement corrélés où la VQE traditionnelle a du mal à s'imposer [8,9].
- Évolutivité et résistance au bruit : L'ansatz compact réduit l'accumulation d'erreurs de porte, réduit la charge de calcul et limite le nombre de paramètres variationnels qui doivent être minimisés.
ADAPT-VQE n'est pas encore parfait. Dans certains cas, il peut être piégé ou ralenti par des minima locaux, et il peut souffrir d'une paramétrisation excessive. Elle peut également être assez gourmande en ressources, car elle nécessite le calcul de gradients et l'optimisation de paramètres avec de nombreuses structures de portes.
Estimation de phase quantique (QPE)
L'objectif du QPE est similaire à celui du VQE, mais sa mise en œuvre est très différente. L'EPQ nécessite des ordinateurs quantiques tolérants aux pannes en raison des circuits quantiques généralement profonds et du niveau élevé de cohérence qu'il requiert. Une fois que le QPE pourra être mis en œuvre, il sera plus précis que le VQE. Une façon de décrire la différence est de considérer la précision comme une fonction de la profondeur du circuit. Le QPE permet d'atteindre la précision avec des profondeurs de circuit s'échelonnant jusqu'à [10]. VQE nécessite échantillons pour atteindre la même précision [10,11].
Krylov, SQD, QSCI et autres dans ce cours
Le VQE a permis d'établir des algorithmes quantiques qui dépendent encore d'ordinateurs classiques, non seulement pour faire fonctionner l'ordinateur quantique, mais aussi pour des parties substantielles de l'algorithme. Plusieurs de ces algorithmes font l'objet du reste de ce cours. Nous donnons ici une explication succincte de quelques-uns d'entre eux, simplement pour les comparer et les opposer à VQE. Ils seront expliqués plus en détail dans les leçons suivantes.
Diagonalisation quantique de Krylov (KQD)
Les méthodes de sous-espace de Krylov permettent de projeter une matrice sur un sous-espace afin de réduire sa dimension et de la rendre plus facile à gérer, tout en conservant les caractéristiques les plus importantes. L'une des astuces de cette méthode consiste à générer un sous-espace qui conserve ces caractéristiques. Il s'avère que la génération de ce sous-espace est étroitement liée à une méthode bien établie sur les ordinateurs quantiques, appelée Trotterisation.
Il existe quelques variantes des méthodes de Krylov quantique, mais l'approche est généralement la suivante :
- Utiliser l'ordinateur quantique pour générer un sous-espace (le sous-espace de Krylov) par trotterisation
- Projeter la matrice d'intérêt sur ce sous-espace de Krylov
- Diagonaliser le nouvel hamiltonien projeté à l'aide d'un ordinateur classique
Diagonalisation quantique basée sur l'échantillonnage (SQD)
La diagonalisation quantique basée sur l'échantillonnage (SQD) est liée à la méthode de Krylov en ce sens qu'elle tente également de réduire la dimension d'une matrice à diagonaliser tout en préservant des caractéristiques essentielles. La SQD procède de la manière suivante :
- Commencez par une estimation initiale de votre état fondamental et préparez le système dans cet état fondamental.
- Utilisez Sampler pour échantillonner les chaînes de bits qui composent cet état.
- Utilisez la collection d'états de base de calcul de l'échantillonneur comme sous-espace sur lequel vous projetez votre matrice d'intérêt.
- Diagonaliser la plus petite matrice projetée à l'aide d'un ordinateur classique.
Elle est liée à la VQE en ce sens qu'elle exploite l'informatique classique et quantique pour des composants algorithmiques importants. Elles ont toutes deux en commun l'obligation de préparer une bonne estimation initiale ou un bon ansatz. Mais la répartition du travail entre l'ordinateur classique et l'ordinateur quantique dans la SQD ressemble davantage à celle de la méthode de Krylov.
En fait, la méthode de Krylov et la SQD ont récemment été combinées dans la méthode de diagonalisation quantique de Krylov basée sur l'échantillonnage (SKQD) [12].
Interaction de la configuration du sous-espace quantique
L' interaction de configuration sélectionnée quantique (QSCI )[13] est un algorithme qui produit un état fondamental approximatif d'un hamiltonien en échantillonnant une fonction d'onde d'essai afin d'identifier les états de base de calcul significatifs pour générer un sous-espace pour une diagonalisation classique. La SQD et la QSCI utilisent toutes deux un ordinateur quantique pour construire un sous-espace réduit. La force supplémentaire de QSCI réside dans la préparation des États, en particulier dans le contexte des problèmes de chimie. Il tire parti de diverses stratégies telles que l'utilisation d'états évoluant dans le temps [14] et d'un ensemble de réponses inspirées de la chimie. En se concentrant sur la préparation efficace des états, QSCI réduit les coûts de calcul quantique pour les hamiltoniens chimiques tout en maintenant une grande fidélité et en tirant parti de la robustesse au bruit des techniques d'échantillonnage d'états quantiques [15]. QSCI propose également une technique de construction adaptative qui fournit plus de réponses pour un meilleur résultat.
Le flux de travail par défaut de QSCI pour les problèmes de chimie est le suivant :
- Construisez l'hamiltonien moléculaire à l'aide du logiciel de votre choix (tel que SciPy ).
- Préparez un algorithme QSCI en sélectionnant un état initial approprié et un ansatz inspiré de la chimie avec un ensemble de paramètres présélectionnés.
- Échantillonnez des états de base significatifs et diagonalisez l'hamiltonien à l'aide d'un ordinateur classique pour obtenir l'énergie de l'état fondamental.
- On a souvent recours à la récupération de configuration [16] et à la post-sélection de symétrie [15] comme techniques de post-traitement.
- En option, le flux de travail de la QSCI adaptative comporte une boucle d'optimisation supplémentaire de step2 à step3, en utilisant davantage de réponses avec des états initiaux aléatoires.
Vérifiez votre compréhension
Qu'est-ce que la VQE a en commun avec toutes les autres méthodes énumérées ci-dessus (à l'exception de la QPE, qui n'est pas décrite en détail)?
Tous impliquent un état d'essai ou une fonction d'onde d'une sorte ou d'une autre. Tout fonctionne mieux lorsque l'estimation initiale de cet état d'essai est excellente.
Une autre réponse correcte est qu'elles sont toutes plus faciles à mettre en œuvre lorsque l'hamiltonien est facile à mesurer (peut être trié en un nombre relativement restreint de groupes d'opérateurs de Pauli commutés).
Qu'est-ce que la VQE a en commun avec aucune des autres méthodes énumérées ci-dessus?
Optimiseurs classiques. Aucun des autres n'utilise d'algorithmes d'optimisation classiques pour sélectionner les paramètres variationnels.
Références
[2] https://en.wikipedia.org/wiki/Variational_quantum_eigensolver
[3] https://journals.aps.org/prapplied/abstract/10.1103/PhysRevApplied.19.024047
[4] https://arxiv.org/abs/2111.05176
[6] https://inquanto.quantinuum.com/tutorials/InQ_tut_fe4n2_2.html
[7] https://www.nature.com/articles/s41467-019-10988-2
[8] https://arxiv.org/abs/2210.15438
[9] https://journals.aps.org/prresearch/abstract/10.1103/PhysRevResearch.6.013254
[10] https://arxiv.org/html/2403.09624v1
[11] https://www.nature.com/articles/s42005-023-01312-y
[13] https://arxiv.org/abs/1802.00171
[14] https://arxiv.org/abs/2103.08505
[15] https://arxiv.org/html/2501.09702v1
[16] https://quri-sdk.qunasys.com/docs/examples/quri-algo-vm/qsci/
[17] https://arxiv.org/abs/2412.13839