Quantum algorithms for equational reasoning

TL;DR

This paper introduces a quantum computational framework for equational reasoning that encodes all equivalent expressions in a quantum superposition, enabling the solution of problems previously infeasible for classical methods.

Contribution

It presents the quantum normal form reduction method, a novel approach to address exponential growth in equivalent expressions using quantum states and tensor networks.

Findings

Successfully solved instances with up to 10^28 equivalent expressions

Demonstrated quantum-inspired algorithms surpassing classical capabilities

Framework applicable to diverse fields like circuit design and data compression

Abstract

As a cornerstone of automated reasoning, equational reasoning finds equivalences between symbolic expressions and fuels advances across scientific disciplines. Yet, its potential remains limited by the exponential growth of equivalent expressions with increasing problem size. We introduce quantum normal form reduction, a quantum computational framework designed to address this challenge. We construct an efficiently implementable quantum Hamiltonian whose ground state encodes all equivalent expressions in a quantum superposition. By preparing and manipulating these states, we tackle fundamental problems in equational reasoning, including verifying and counting equivalent expressions and identifying structural properties of equivalence classes. We demonstrate a quantum-inspired version of the algorithm, using tensor networks to solve instances involving up to 10^28 equivalent expressions, far beyond the reach of classical graph exploration. This framework opens the path for quantum symbolic computation in areas from circuit design to data compression, computational group theory, linguistics, and macromolecular modeling, unlocking previously inaccessible problems.

Equations

$$S=\left[A|R\right],$$

$$A=\{{\alpha }_i,i\in [1,\dots ,d]\}$$

$$R=\{r_l|r_l:A^{\ast }\to A^{\ast },l\in [1,\dots ,n_r]\}$$

$$r_l=({\alpha }_{[j_1]}{\beta }_{[j_2]}\cdots \approx {\alpha }_{[j_1]}^{\prime }{\beta }_{[j_2]}^{\prime }\dots )$$

$$S=\left[\{a,b\}\mid \{a_1{b}_2\approx b_1{a}_2,a_2{b}_3\approx b_2{a}_3\}\right],$$

$$S^{\prime }=\left[\{a,b\}\mid \{a_1{a}_2\approx b_1{b}_2,,a_2{a}_3\approx b_2{b}_3\}\right],$$

$$|\psi ⟩=\sum _k{\psi }_k|{\omega }_k⟩,$$

$${\hat{r}}_l=|{\alpha }^{\prime }⟩{⟨\alpha |}_{j_1}\otimes |{\beta }^{\prime }⟩{⟨\beta |}_{j_2}\otimes \cdots +|\alpha ⟩{⟨{\alpha }^{\prime }|}_{j_1}\otimes |\beta ⟩{⟨{\beta }^{\prime }|}_{j_2}\otimes \cdots ,$$

$$|X_{S,\tilde{\omega }}⟩=\sum _{\omega \in X_{S,\tilde{\omega }}}\frac{1}{\sqrt{|X_{S,\tilde{\omega }}|}}|\omega ⟩.$$

$$E_{\text{D}}$$

$$=\sum _{({\omega }_k,{\omega }_{k^{\prime }})\in E}(|{\psi }_k{|}^2-{\psi }_{k^{\prime }}^{\ast }{\psi }_k-{\psi }_{k^{\prime }}{\psi }_k^{\ast }+|{\psi }_{k^{\prime }}{|}^2).$$

$$E_{\text{D}}=⟨\psi |{\hat{\mathcal{L}}}_S|\psi ⟩.$$

$${\hat{\mathcal{L}}}_S=\sum _{{\hat{r}}_l:r_l\in R}\left({\widehat{r_l}}^2-{\hat{r}}_l\right),$$

$${\hat{H}}_{S,\tilde{\omega }}(h)=(1-h){\hat{\mathcal{L}}}_S-h|\tilde{\omega }⟩⟨\tilde{\omega }|.$$

$$\underset{h\to 0}{\lim }|\varphi (h)⟩=|X_{S,\tilde{\omega }}⟩.$$

$$\hat{H}(t)={\hat{H}}_{S,\tilde{\omega }}\left(\frac{\tau -t}{\tau }\right)=\frac{t}{\tau }{\hat{\mathcal{L}}}_S-\left(\frac{\tau -t}{\tau }\right)|\tilde{\omega }⟩⟨\tilde{\omega }|$$

$$|\psi (t)⟩=|\varphi \left(\frac{\tau -t}{\tau }\right)⟩,$$

$$\hat{I}(t)=i{\hat{H}}_{S,\tilde{\omega }}\left(\frac{\tau -t}{\tau }\right),$$

$$F(X_{S_1,{\omega }_1},X_{S_2,{\omega }_2})=\frac{|X_{S_1,{\omega }_1}\cap X_{S_2,{\omega }_2}{|}^2}{|X_{S_1,{\omega }_1}|\cdot |X_{S_2,{\omega }_2}|},$$

$$|⟨X_{S_1,{\omega }_1}|X_{S_2,{\omega }_2}⟩{|}^2$$

$$=\frac{1}{|X_{S_1,{\omega }_1}|\cdot |X_{S_2,{\omega }_2}|}{\left|\sum _{\begin{array}{c}\omega \in X_{S_1,{\omega }_1}\\ {\omega }^{\prime }\in X_{S_2,{\omega }_2}\end{array}}⟨\omega |{\omega }^{\prime }⟩\right|}^2$$

$$=\frac{|X_{S_1,{\omega }_1}\cap X_{S_2,{\omega }_2}{|}^2}{|X_{S_1,{\omega }_1}|\cdot |X_{S_2,{\omega }_2}|}$$

$$=F,$$

$$F=\{\begin{array}{ll}1& \text{if }{\omega }_1\text{ and }{\omega }_2\text{ are connected,}\\ 0& \text{otherwise,}\end{array}$$

$$S_A=\left[A,R=\{{\alpha }_i\approx {\alpha }_i^{\prime }|\alpha ,{\alpha }^{\prime }\in A,i\in [1,\dots ,L]\}\right],$$

$$|All⟩=|X_{S_A,{\omega }_1}⟩=\frac{1}{\sqrt{d^L}}\sum _{\omega }|\omega ⟩,$$

$$F=\frac{|X_{S_2,{\omega }_2}|}{d^L},$$

$$X_{S,\tilde{\omega },g}=\{\omega \in X_{S,\tilde{\omega }}\mid g(\omega )=1\},$$

$$\sum _{\omega \in X_{S,\tilde{\omega }}}\frac{1}{\sqrt{|X_{S,\tilde{\omega }}|}}|\omega ⟩\otimes |0⟩⟶\sum _{\omega \in X_{S,\tilde{\omega }}}\frac{1}{\sqrt{|X_{S,\tilde{\omega }}|}}|\omega ⟩\otimes |g(\omega )⟩.$$

$$\sum _{\omega \in X_{S,\tilde{\omega },g}}\frac{1}{\sqrt{|X_{S,\tilde{\omega },g}|}}|\omega ⟩.$$

$$p_1=\frac{|X_{S,\tilde{\omega },g}|}{|X_{S,\tilde{\omega }}|},$$

$$\hat{O}={\left(\widehat{\text{mathds}1}+{\hat{\sigma }}_X\right)}^{\otimes L}$$

Full text

Introduction

Equivalence relations group individual objects into categories based on shared structure or behavior, enabling reasoning at the level of entire classes rather than isolated instances. This abstraction is central across disciplines: biologists study species rather than individual organisms, mathematicians analyze functions instead of specific representations, physicists examine macrostates rather than microstates, and linguists consider languages beyond single sentences. In symbolic computation, equivalence relations can be encoded and manipulated algorithmically, allowing computers to automate such reasoning. A key approach is provided by term rewriting systems  (?, ?, ?): given one expression and a set of rewriting rules, all other equivalent expressions are generated by applying a sequence of rule-based substitutions.

Across scientific domains, it is well known that a few simple rewriting rules can encapsulate complex semantics and generate rich behaviors. In mathematics, term rewriting provides a unifying framework for encoding algebraic identities, logical equivalences, and inference rules. This underpins symbolic algorithms for solving equations  (?), performing computations in algebraic groups and monoids  (?) (Figure 1Aa) , and proving theorems  (?). In classical and quantum information processing, rewriting rules enable local substitutions of equivalent subcircuits (Figure 1Ad)  (?), which is crucial for verifying functional correctness and optimizing resource usage through circuit compilation  (?, ?, ?, ?). Rewriting also extends to formal grammars in the field of linguistics, which define the generative structure of natural language, and in computer science provide the syntactic foundation of programming languages, supporting tasks such as equivalence checking and compiler optimization  (?, ?, ?). In data compression, formal grammar-based encoding allows lossless compression of large datasets  (?). In automata theory, rewriting systems enable the systematic exploration of configuration spaces, being able to simulate any Turing machine  (?, ?). In biology and chemistry, rewriting frameworks capture both the structural and informational aspects of macromolecules such as DNA (Figure 1Ac), RNA, proteins, and polymers(Figure 1Ab)  (?, ?, ?).

The expressive power of rewriting systems is accompanied by intrinsic computational challenges. Iteratively applying rewriting rules to a single object generates an equivalence class that can grow combinatorially. Many foundational and applied scientific investigations in all previous domains require exploring such equivalence classes in full or in part. For example, verifying whether two circuits implement the same function reduces to determining whether they can be transformed into one another using a prescribed set of functionality-preserving rewrites  (?). In computational group theory, this corresponds to the word problem for finitely generated groups  (?): determining whether two sequences of generators produce the same group element. Related decision problems arise across disciplines. The grammar equivalence problem  (?, ?) — determining whether two formal grammars generate the same language — is central to both computational linguistics and programming language theory. In automata theory and statistical physics, a crucial question is to determine how many distinct states a given rewriting system or automaton can explore, or how many space-filling curves are admissible on a given lattice  (?, ?). This counting problem connects to entropy, complexity, and information content, with relevance in fields such as polymer thermodynamics  (?), automata-based models of physical systems  (?), and genetic diversity analysis  (?).

Due to the ability of rewriting systems to simulate arbitrary computations, the word problem is undecidable in general  (?). Only a narrow subset of the aforementioned questions can currently be answered using existing algorithms with feasible resource requirements, also for decidable problems, due to the exponential size of the equivalence classes. But a different form of computation, such as quantum computing, can offer a new theoretical perspective, or in some cases a computational advantage, in the vast field of equational reasoning. Quantum computation is a fundamentally new computational paradigm, enabling polynomial and exponential speedups for problems that are intractable for classical algorithms  (?, ?, ?).

Here, we extend quantum computation to automated symbolic reasoning by introducing quantum normal form reduction, a quantum computational framework tailored to address a variety of the aforementioned challenges. The key insight is that quantum mechanics enables the representation and manipulation of entire equivalence classes encoded as coherent quantum superpositions of the class members (Figure 1C), bypassing the need to sequentially represent individual elements (Figure 1B). We refer to these class-representing quantum states as orbit states. We demonstrate that these states can be prepared on a quantum computer, and that key questions can be efficiently answered by performing quantum operations on orbit states. For example, the word problem is solved by measuring overlap between two orbit states, and the number of elements of an equivalence class is inferred by measuring an appropriate observable.

As a quantum analogue to normal form reduction (see Supplementary Text for details), which is usually performed by exploiting the Knuth-Bendix algorithm  (?), the preparation of orbit states via quantum normal form reduction allows one to associate a unique quantum state to an entire class of equivalent expressions. We demonstrate that orbit states can be prepared as ground states of an appropriate sparse Hamiltonian, specifically the discrete Laplacian of the configuration graph generated by the action of the rewriting system  (?). This Hamiltonian is constructed as a linear combination of tensor products of local operators, each encoding a rule of the rewriting system. As a consequence, we show that the computational cost of simulating its action on a quantum device scales polynomially with the number of rules in the rewriting system itself and with the size of the rules, i.e., the number of characters on which each rule acts.

Orbit state preparation can be accomplished using current quantum optimization techniques, such as quantum annealing  (?, ?, ?, ?), optimal control  (?, ?, ?), quantum approximate optimization algorithms (QAOA)  (?, ?, ?) and imaginary time evolution  (?, ?, ?). Despite the exponential size of the equivalence classes, the amount of quantum memory required to represent the orbit state is polynomial, while the amount of time needed to prepare the state depends on the specific instance of the equivalence problem  (?). Once the state has been prepared, many global properties of the equivalence class can be measured efficiently. Solving the word problem and counting elements inside an equivalence class reduces to measuring the quantum fidelity of orbit states. This operation can be performed in polynomial time, for instance, through a swap test  (?), laying the groundwork for tackling previously inaccessible problems.

We employ Tensor Network (TN) methods  (?, ?, ?, ?, ?) to emulate the execution of the proposed quantum algorithm and to demonstrate its effectiveness. Tensor network methods enable efficient classical simulations of quantum processes, under specific conditions regarding the structure of quantum correlations and entanglement  (?). The TN implementation effectively defines a quantum-inspired classical algorithm that can be executed on existing classical hardware, already enabling efficient classical analysis, in some cases, beyond the current state of the art.

We focus on a toy term rewriting system in which an initial string of size $L$ is updated by substituting substrings of equal length. The string length $L$ corresponds to the quantum memory requirements of our device. We address both the word problem and the counting problem for input strings of different sizes. Despite the exponential size of the equivalence classes, we efficiently represent and manipulate orbit states generated from strings of length up to $L=100$, containing up to $10^{28}$ connected strings encoded in $1$ gigabyte of memory. The memory required to encode the same data as a list of strings would be approximately $10^{17}$ terabytes. This demonstrates the potential of our algorithm as a powerful compression method for sets of data that are related by rewriting relations.

Results

Term rewriting systems

A term rewriting system consists of a set of rewrite rules that operate on syntactic terms, capturing how complex structures evolve through rule-based replacements. Term rewriting can manipulate a variety of symbolic data structures, including strings, lattice configurations, and graphs. Since any data structure can be serialized and ultimately represented as a bit string, we focus on rewriting strings. As a consequence, rewriting rules that act locally on the original data structure might assume a non-local representation when applied to the serialized data structure. In this setting, we denote as a string rewriting system a term rewriting system in which all rewriting rules act locally on strings.

Thus, without loss of generality, we define an invertible term rewriting system $S$ as

[TABLE]

where

[TABLE]

is an alphabet of $d$ characters used to compose strings, $A^{*}$ is the set of all the possible strings over the alphabet $A$, and

[TABLE]

is a set of $n_{r}$ rewriting rules, i.e., functions mapping a string being in the set $A^{*}$ to another string in $A^{*}$ via replacement operations on a set of characters. In particular, a rule

[TABLE]

simultaneously replaces the character $\alpha$ at position $j_{1}$ with $\alpha^{\prime}$, the character $\beta$ at position $j_{2}$ with $\beta^{\prime}$, and so on. The equivalence symbol $\approx$ indicates that we restrict to invertible rules, meaning that the reverse transformation can also be applied. Here, we consider rules that preserve the length of strings. More general rules can be made length-preserving by the introduction of a blank character.

Invertible rules establish an equivalence relation among connected strings, then partitioning the space of strings into disjoint equivalence classes  (?). Each equivalence class consists of all strings that are mutually reachable via arbitrary sequences of rewriting rules from the set $S$. To label these classes, we choose a representative element $\tilde{\omega}$ within a given class and denote by $X_{S,\tilde{\omega}}$ the set of all strings that can be obtained by applying sequences of rules from $S$ to the initial string $\tilde{\omega}$. As an illustrative example, consider the rewriting system

[TABLE]

acting on strings of fixed length $L=3$. The space of all length-$3$ strings is partitioned by the rules of $S$ into four equivalence classes: $X_{S,aaa}=\{aaa\}$, $X_{S,aab}=\{aab,aba,baa\}$, $X_{S,abb}=\{abb,bab,bba\}$, and $X_{S,bbb}=\{bbb\}$, where we labeled each equivalence class by choosing the lexicographically smallest element as its representative. For this rewriting system, the word problem, namely, determining whether two given words are equivalent, has a positive answer for the pair $aab$ and $baa$, and a negative answer for the pair $aab$ and $aaa$. The counting problem, which asks for the number of words equivalent to a given input, yields $1$ for $aaa$ and $bbb$, and $3$ for $aab$ and $bab$.

More elaborate decision problems arise when comparing multiple rewriting systems. For instance, consider a second system

[TABLE]

and ask whether $S$ and $S^{\prime}$ generate the same equivalence class when applied to the same input word $aaa$. This question, which captures a simple instance of the grammar equivalence problem  (?, ?), has a negative answer in this case, since $X_{S^{\prime},aaa}=\{aaa,abb,bba\}\neq X_{S,aaa}$.

The rewriting systems in these examples exhibit particularly simple word and counting problems because their dynamics preserve an easily identifiable invariant. In the first system, this invariant is the total number of $a$ characters, while in the second it is the parity of that number. In contrast, equivalence relations relevant to more complex settings, such as equivalence problems for circuits, cannot generally be reduced to the conservation of a simple quantity. This increased complexity is largely due to the presence of context-dependent rules, that is, rules that apply only when a specific local pattern is matched. Within our formalism, the rule $aaa\approx aba$ constitutes an example of a context-dependent rewriting rule, since the central character $a$ is replaced only when surrounded by other characters $a$.

Quantum states and orbit states

Quantum computing is the use of controllable quantum many-body systems to process information and solve computational problems. We consider many-body systems composed of $L$ local subsystems, each associated with a finite set of $d$ discrete classical configurations, such as the two orientations of a spin in a magnetic field. When local measurements of all subsystems are performed on such a system, the readout yields a single classical configuration $\omega_{k}$ with $k\in[1,\dots,d^{L}]$, assigning to each local subsystem one of its classical configurations. Repeating the same quantum computation multiple times produces a statistical ensemble of outcomes, effectively sampling configurations $\omega_{k}$ with associated probabilities $p_{k}$. The goal of quantum computation is to manipulate the quantum system so that the solution to a problem of interest can be efficiently extracted from the measurement statistics.

What distinguishes quantum computation from classical probabilistic computation is the nature of the system’s evolution prior to measurement. While stochastic processes evolve probability distributions with non-negative real weights, quantum systems are described by complex probability amplitudes that may interfere constructively or destructively. In particular, the state of the system is represented by a normalized vector $\ket{\psi}$ in a Hilbert space constructed by associating a computational basis vector $\ket{\omega_{k}}$ to each classical configuration $\omega_{k}$,

[TABLE]

where the $\psi_{k}$ are complex amplitudes. The probability $p_{k}$ of observing the configuration $\omega_{k}$ upon measurement is given by the Born rule, $p_{k}=|\psi_{k}|^{2}$ and the measurement irreversibly projects the state onto the observed configuration, eliminating all other components of the superposition.

For implementing symbolic quantum computation, any string of $L$ characters from the alphabet $A$ is mapped to a configuration of a many-body quantum system. To this aim, we define a quantum system composed of $L$ local subsystems with $d$ internal states, or qudits, where $d=|A|$ is the size of the alphabet. Each internal state of the qudit is labeled by a character in $A$. Upon measuring such a quantum system in the computational basis, we sample a classical configuration $(\alpha_{k_{1}},\dots,\alpha_{k_{L}})$ encoding a string $\omega_{k}=\alpha_{k_{1}}\dots\alpha_{k_{L}}$, and we assert that the system is in the state $\ket{\omega_{k}}=\ket{\alpha_{k_{1}},\dots,\alpha_{k_{L}}}$. As we will show hereafter, before measurement, the state of the system is in general in a quantum superposition of states, specifically $\ket{\psi}=\sum_{k}\psi_{k}\ket{\omega_{k}}$. The Hilbert space has dimension $d^{L}$, corresponding to the maximum number of different strings that can be sampled.

Each rewriting rule $r_{l}\in R$ naturally defines a linear operator $\hat{r}_{l}$ that acts on a quantum superposition of strings by applying the rule $r$ to each computational basis state. For a rule $r_{l}=(\alpha_{[j_{1}]}\beta_{[j_{2}]}\dots\approx\alpha^{\prime}_{[j_{1}]}\beta^{\prime}_{[j_{2}]}\dots)$, the corresponding operator can be written in the bra-ket notation as:

[TABLE]

where $\ket{\alpha}\bra{\alpha^{\prime}}_{j}$ is the local operator that replace the character $\alpha^{\prime}$ with the character $\alpha$ in any computational basis state containing the character $\alpha^{\prime}$ at position $j$, while multiplies by [math] other computational basis states. The two terms in the latter expression correspond to the two directions in which the rule can be applied, and together ensure that the operator is Hermitian, as required for a valid quantum observable.

As shown in the previous section, the rewriting rules in the term rewriting system $S$ partition the space of all strings of a given size into disjoint equivalence classes. We encode each equivalence class $X_{S,\tilde{\omega}}$ as an equally weighted quantum superposition of all words $\omega\in X_{S,\tilde{\omega}}$, which we designate as the orbit state:

[TABLE]

As we will show, once the orbit state has been prepared, one can efficiently manipulate the entire equivalence class via quantum operations.

A parent Hamiltonian for orbit states

The dynamics of a closed quantum system is governed by a Hermitian operator acting on its state, known as the Hamiltonian. The Hamiltonian encodes the energetic structure of the system: its eigenvalues correspond to the possible energy levels, and its eigenvectors to the quantum states with definite energy. One important paradigm of quantum computation is to encode the solution of a computational problem into the lowest-energy eigenstate, or ground state, of a suitably designed Hamiltonian. The desired quantum state can then be prepared using physical processes such as adiabatic evolution.

We construct a parent Hamiltonian $\mathcal{L}_{S}$ whose degenerate ground states are precisely the orbit states defined in the previous section. While finding a parent Hamiltonian is generally a challenging task  (?, ?, ?), in this case, we systematically build $\mathcal{L}_{S}$ directly from the rewriting rules in $S$, exploiting the notion of discrete Laplacian on a graph  (?) (see Figure 2).

We define the graph $G=(V,E)$, where the vertices $V=\{\omega_{k}\}$ represent the possible strings of size $L$. The edges $E=\{(\omega_{k},\omega_{k^{\prime}})\}$ connect strings that are related by a rewriting rule $r$. By construction, the graph $G$ is a disconnected graph with the equivalence classes forming connected subgraphs. Note that the size of the graph $G$ increases exponentially with the string length $L$. Since each vertex represents a string and each string corresponds to a state of the computational basis, a quantum state $\ket{\psi}=\sum_{k}\psi_{k}\ket{\omega_{k}}$ can be represented as complex function $\psi:V\rightarrow\mathbb{C}$ that associates to each vertex $\omega_{k}$ of the graph the complex value $\psi_{k}$.

Orbit states correspond to the flattest complex functions on the graph $G$, as they remain constant in both modulus and phase on each connected component. For a function defined on a graph, flatness can be quantified as the minimization of differences between values assigned to adjacent vertices. The variation of $\psi$ along edges is measured by its discrete gradient $\vec{\nabla}\psi$ (in red in Figure 2), that is the function that associates to each edge $(\omega_{k},\omega_{k^{\prime}})\in E$ the difference $\psi_{k^{\prime}}-\psi_{k}$. The sum of $|\psi_{k^{\prime}}-\psi_{k}|^{2}$ on all edges measures the total variation of $\psi$, namely its Dirichlet energy  (?). It reads:

[TABLE]

In the Supplementary Text, we show that the Dirichlet energy associated with a function $\psi$ is given by the expectation value of a positive semi-definite graph Laplacian operator $\hat{\mathcal{L}}_{S}$ on the state $\ket{\psi}$:

[TABLE]

In this framework, orbit states are ground states of the Laplacian and have zero energy. The operator $\hat{\mathcal{L}}_{S}$ is sparse and defined by the rewriting rules $r_{l}\in R$ of the system $S$:

[TABLE]

where $\hat{r}_{l}$ is the operator encoding the rule $r_{l}$, as defined in Eq. (8).

Although both the size of the graph $G$ and the dimension of the Hilbert space grow exponentially with the string length $L$, the action of the Laplacian operator $\hat{\mathcal{L}}_{S}$ can nevertheless be simulated efficiently on a quantum computer. More precisely, for a term rewriting system specified by a number of rules that scales at most polynomially with $L$, the time evolution generated by $\hat{\mathcal{L}}_{S}$ can be approximated to fixed accuracy by a quantum circuit of depth scaling polynomially in $L$. This efficiency follows from the fact that the Laplacian decomposes into a sum of rule operators, each of which acts as a product of single-qudit operators (see Supplementary Text for details).

Preparing the orbit state

While the ground state space of $\hat{\mathcal{L}}_{S}$ is spanned by arbitrary superpositions of orbit states, our goal is to prepare a specific orbit state $\ket{X_{S,\tilde{\omega}}}$. A natural strategy is to drive the system toward the ground state of the Laplacian while constraining the dynamics to the subspace of the Hilbert space spanned solely by computational-basis states corresponding to strings reachable from $\tilde{\omega}$. To achieve this, we add to the Laplacian the projector $P_{\tilde{\omega}}=-\ket{\tilde{\omega}}\bra{\tilde{\omega}}$. This addition leads to the family of Hamiltonians:

[TABLE]

The off-diagonal component of this family of Hamiltonians is composed of the rule operators of the rewriting system. As a consequence, the dynamics generated by $\hat{H}_{S,\tilde{\omega}}(h)$ acting on the initial state $\ket{\tilde{\omega}}$ is confined to the subspace spanned by the equivalent strings in $X_{S,\tilde{\omega}}$. Within this subspace, the unique ground state is the orbit state $\ket{X_{S,\tilde{\omega}}}$.

To prepare the orbit state on a quantum computer or simulator, one possible strategy is to design an evolution that follows the instantaneous ground states $\ket{\phi(h)}$ of the Hamiltonian $\hat{H}_{S,\tilde{\omega}}(h)$ as $h$ transitions from 1 to 0. If the evolution is adiabatic, the final state at the end of this process is then

[TABLE]

As well as the Laplacian operator, the Hamiltonian $\hat{H}_{S,\tilde{\omega}}$ can be simulated with a quantum circuit involving a number of gates that scales polynomially with the number of rules in $S$ and their size (see Supplementary Text). Thus, the evolution can be implemented using state-of-the-art quantum algorithms such as Quantum Annealing (QA)  (?, ?, ?, ?), Optimal Control  (?, ?, ?), Quantum Approximate Optimization Algorithms  (?, ?, ?) or Imaginary Time Evolution (ITE)  (?, ?, ?) and shortcuts to adiabaticity  (?, ?).

In a QA-based approach, the initial state $\ket{\psi(t=0)}=\ket{\tilde{\omega}}$ evolves under the time-dependent Hamiltonian

[TABLE]

for $t:0\rightarrow\tau$. If $\tau$ is sufficiently large, the evolution adiabatically follows the ground state path

[TABLE]

whose final state is the target orbit state.

Here, to benchmark the proposed approach, we perform a TN simulation of the algorithm, focusing on a mixture of QA and ITE, which we refer to as Imaginary Quantum Annealing (IQA)  (?). In an IQA-based approach, the evolution is governed by the imaginary Hamiltonian

[TABLE]

and the system follows the adiabatic path by gradually suppressing the amplitude of the excited components generated by non-perfect adiabatic evolution. A key advantage of IQA is that, at any given time, the imaginary Hamiltonian dampens previously generated excitations, effectively mitigating errors and eventually outperforming standard QA  (?). In a tensor network (TN) simulation, damping can also suppress excitations that arise from approximating intermediate highly entangled states with finite bond dimension.

Quantum algorithms for equational reasoning

Many global properties of equivalence classes can be efficiently extracted via quantum operations on orbit states.

We begin by focusing on the reconstruction of the overlap function

[TABLE]

where $X_{S_{1},\omega_{1}}$ is the equivalence class generated by the action of a term rewriting system $S_{1}$ on the input word $\omega_{1}$, $X_{S_{2},\omega_{2}}$ is the equivalence class generated by the action of a possibly different term rewriting system $S_{2}$ on the input word $\omega_{2}$, and $|X|$ is the number of elements in $X$. The function $F(X_{S_{1},\omega_{1}},X_{S_{2},\omega_{2}})$ quantifies the squared size of the intersection between the two equivalence classes $X_{S_{1},\omega_{1}}$ and $X_{S_{2},\omega_{2}}$, normalized by the product of their sizes. This quantity measures the similarity between equivalence classes: it approaches zero when the overlap is negligible and reaches one when the two sets coincide. As we show here, the knowledge of $F$ enables the resolution of several important problems, such as the word problem, the counting problem, and the grammar-equivalence problem.

Measuring equivalence classes overlap via fidelity

Once two orbit states $\ket{X_{S_{1},\omega_{1}}}$ and $\ket{X_{S_{2},\omega_{2}}}$ are prepared in two memory registers of a digital quantum computer as shown in the previous section, their similarity can be quantified by measuring the squared magnitude of their overlap, i.e., $\left|\langle{X_{S_{1},\omega_{1}}}\ket{X_{S_{2},\omega_{2}}}\right|^{2}$. This function is also called quantum fidelity  (?). Using the definition in Eq. (9), the fidelity between orbit states $\ket{X_{S_{1},\omega_{1}}}$ and $\ket{X_{S_{2},\omega_{2}}}$ is reduced to the overlap function $F(X_{S_{1},\omega_{1}},X_{S_{2},\omega_{2}})$ as follows:

[TABLE]

since only the non-zero terms in the summation correspond to the strings that simultaneously belong to both equivalence classes. The ratio $F(X_{S_{1},\omega_{1}},X_{S_{2},\omega_{2}})$ can be efficiently measured using the swap test algorithm  (?), which requires only a polynomial number of gates relative to the length of the input words. Furthermore, as discussed in Supplementary Text, fidelity between orbit states can also be efficiently estimated on analog quantum simulators such as quantum annealers.

Word problem

When $S_{1}=S_{2}$, the strings $\omega_{1}$ and $\omega_{2}$ belong to the same equivalence class if and only if they are connected under the rewriting system. In this case,

[TABLE]

thus providing a solution to the word problem.

Counting problem

Let us consider the rewriting system

[TABLE]

which replaces an arbitrary character at any position with any other character. $S_{A}$ is thus capable of generating all the possible strings over the alphabet. Consequently, we can construct a uniform superposition over all such strings in the Hilbert space as an orbit state for $S_{A}$, i.e.,

[TABLE]

where $\omega_{1}$ is an arbitrary input string. In this special case, the function $F=|\langle X_{S_{A},\omega_{1}}\ket{X_{S_{2},\omega_{2}}}|^{2}$ reduces to

[TABLE]

which allows us to estimate the number of strings connected to $\omega_{2}$ under the action of $S_{2}$, thus solving the counting problem. Importantly, the uniform superposition state $\ket{All}$ can be efficiently implemented on digital quantum computers using Hadamard gates as commonly done in most quantum algorithms  (?).

Filtering

Another relevant operation on equivalence classes is the extraction of a specific subset of elements that satisfy a given condition. We name this operation filtering. For instance, one might be interested in generating all space-filling curves on a square lattice that exhibit inversion symmetry. The ability to filter elements is also essential for applications in formal grammars, where it may be necessary to exclude strings containing nonterminal symbols.

We define the target subset as

[TABLE]

where $g$ is a Boolean function that returns 1 if and only if the string $\omega$ satisfies the desired property. Once the orbit state has been prepared, subset extraction can be implemented by introducing an ancilla qubit initialized in the state $\ket{0}$, which is flipped conditionally based on the value of $g(\omega)$. This results in the following transformation:

[TABLE]

By measuring the ancilla qubit and post-selecting the outcome $\ket{1}$, the remaining system collapses into a quantum state that encodes a uniform superposition over the filtered subset:

[TABLE]

The probability of successfully measuring the ancilla in state $\ket{1}$ is given by

[TABLE]

which implies that, on average, the procedure must be repeated $\frac{|X_{S,\tilde{\omega}}|}{|X_{S,\tilde{\omega},g}|}$ times to obtain a successful outcome. When the size of the selected subset is not exponentially smaller than the full equivalence class, filtering can be carried out efficiently on a quantum computer.

Grammar equivalence problem

A formal grammar can be viewed as a special case of a string rewriting system in which attention is restricted to the equivalence class of strings reachable from a designated starting string $\tilde{\omega}$. The language generated by the grammar is obtained by filtering this equivalence class to retain only those strings that do not contain a designated set of symbols, called non-terminals. This construction reflects the interpretation of grammars as generative mechanisms for syntactically well-formed sentences  (?). The grammar equivalence problem then consists of determining whether two grammars generate the same language, that is, whether their filtered equivalence class coincide. In our language, if $g$ denotes a Boolean function that selects strings containing no non-terminal symbols, the grammar equivalence problem reduces to estimating the overlap $F$ between the filtered sets $X_{S_{1},\tilde{\omega},g}$ and $X_{S_{2},\tilde{\omega},g}$, generated by two rewriting systems $S_{1}$ and $S_{2}$ acting on the same initial string $\tilde{\omega}$.

Estimating classical expectation values

Finally, a wide range of statistical information can be extracted by measuring the system in the state $\ket{X_{S,\tilde{\omega}}}$, which results in uniformly sampling strings from the equivalence class $X_{S,\tilde{\omega}}$. This sampling process enables the estimation of expectation values of classical functions over strings. For example, we can estimate the probability of finding a particular character or substring at a specified position by counting its occurrences in the sampled strings.

Estimating quantum expectation values

We can also estimate expectation values that are not associated with any classical function. Indeed, the expectation value of any Hermitian operator $\langle\hat{O}\rangle$ that can be expressed as a sum of a polynomial number of tensor product operators can be estimated on a quantum computer via single qubit rotations and sampling. For example, in the case of a binary alphabet ($|A|=2$), the number of connected strings $|X_{S,\tilde{\omega}}|$ can be estimated by measuring the expectation value of the observable:

[TABLE]

where $\hat{\sigma}_{X}$ is the Pauli-X operator and $L$ is the string length. The operator $\hat{O}$ is the projector onto the equal superposition of all computational basis states, scaled by a factor of $2^{L}$. Its expectation value yields the fidelity between the orbit state and the uniform superposition state (in Eq. (22) ), multiplied by the size of the orbit state. This directly provides the total size of the equivalence class, i.e., a solution for the counting problem.

Estimating non-linear quantities

Finally, by preparing multiple copies of the orbit state, one can also estimate non-linear quantities, such as the $2$-Rényi entropy of subsystems  (?). This provides insights into the correlations between substrings located in different regions of the string. In particular, entropy serves as an indirect measure of computational complexity: when the Rényi entropy of a region tends toward zero, the structure of the equivalence class simplifies, becoming close to a Cartesian product of independent components. This factorization allows for a compact classical representation of the equivalence class and may lead to significantly faster classical algorithms for many tasks.

Tensor network implementations

In certain instances, the quantum algorithm introduced in this work admits efficient classical simulation via tensor network (TN) methods  (?, ?, ?, ?, ?). TN are numerical techniques that express quantum states as networks of tensors with contractions between indices. The size of each tensor index is determined by the local Hilbert space dimension $d$ or by an integer parameter $\chi$, referred to as the bond dimension. While a faithful representation of a generic quantum state typically demands a bond dimension that scales exponentially with system size $L$, this requirement can be significantly relaxed when entanglement between subsystems is limited. In such cases, the bond dimension may scale polynomially, enabling efficient classical representations of quantum states and allowing for the simulation of specific processes and measurements. This representation is especially efficient for states obeying an area-law for entanglement, such as ground states of local gapped Hamiltonians  (?, ?).

A Matrix Product State (MPS) is a one-dimensional instance of a tensor network, which represents a one-dimensional lattice of $L$ qudits using $L$ rank-3 tensors of shape $(\chi,d,\chi)$. The memory cost to store an MPS scales as $\mathcal{O}(L\cdot\chi^{2}\cdot d)$. We simulate the time evolution within the MPS formalism using a Time-Dependent Variational Principle (TDVP)  (?). A variety of relevant quantities — including fidelities, expectation values of observables, sampling, and subsystem entropies — can be efficiently computed, provided that $\chi$ grows sub-exponentially with system size  (?).

When the entanglement of the orbit state admits an efficient tensor network representation, the performance of the quantum algorithm can be assessed through classical simulation. In particular, the preparation of orbit states via imaginary quantum annealing is especially well-suited to tensor network methods, as it naturally suppresses excitations introduced by numerical errors and by approximations associated with a finite bond dimension.

Beyond serving as a testbed for the quantum normal form reduction algorithm, the tensor network representation naturally gives rise to a quantum-inspired algorithm. In the current computational landscape — where quantum hardware remains in its early stages while classical computing is highly mature — such an approach offers a novel and potentially advantageous framework for addressing the challenges inherent in term rewriting. Additionally, tensor networks enable the efficient computation of properties of quantum states that are otherwise difficult, or even exponentially hard, to extract on quantum hardware. Notable examples include the evaluation of von Neumann entropies and the estimation of exponentially small probability amplitudes  (?).

The primary limitation of the tensor network approach proposed here lies in the complexity of local correlations among strings within the equivalence class $X_{S,\omega}$. These correlations can become particularly intricate in certain relevant contexts, for example, in natural languages, where meaning emerges from the nuanced interplay among components of a sequence  (?, ?, ?). In this case, more sophisticated TN structures, such as Tree Tensor Networks  (?, ?, ?), could allow for an accurate simulation. An even more challenging scenario arises in rewriting systems that go beyond string rewriting. In one-dimensional string rewriting systems, the locality of the rules translates into a one-dimensional local structure for the Laplacian, which can yield a relatively simple entanglement pattern in the ground states and the possibility to approximate these states as MPS with bond dimension that scales polynomially with the system size.  (?). One-dimensional string rewriting appears in practice in DNA mutation modeling, formal languages, regex-based text processing, compiler peephole optimization, and grammar-based data compression. In contrast, when rewriting rules are non-local, such as large tree-level transformations used in compilers and symbolic algebra, the notion of one-dimensional distance between subsystems breaks down, making an efficient tensor network representation more challenging or, in some cases, inefficient. Nevertheless, recent advances in tensor network algorithms have demonstrated the ability to represent ground states of systems beyond one dimension through Projected Entangled Pair States  (?), Multi-scale Entanglement Renormalization Ansatz  (?), and Augmented Tree Tensor Networks  (?).

Numerical results

We simulate the quantum algorithm using tensor network methods and observe a favorable scaling of the required computational resources with the size of the input strings.

We consider the string rewriting system $S=\Big[A|R\Big]$, where $A=\{a,b\}$ and

[TABLE]

Despite its simplicity, the application of the aforementioned rewriting system to a string of length $L$ generates an exponential number of connected strings. We validate the quantum framework introduced in this work by solving a collection of instances of the word problem and the counting problem for the SRS defined in Eq. (Numerical results). Each problem instance is specified by two disjoint sets of four strings of length $L$:

[TABLE]

Within each set, all pairs of strings are mutually connected through a sequence of rewriting operations, i.e.,

[TABLE]

and

[TABLE]

No string from one set is connected to any string in the other, which guarantees that the corresponding equivalence classes are distinct, for example $X_{S,\omega_{1}}\neq X_{S,\omega_{5}}$.

We consider problem instances with string lengths $L$ ranging from $10$ to $100$. For each string $\tilde{\omega}$ belonging to $\Omega_{1}(L)$ or to $\Omega_{2}(L)$, we construct the corresponding orbit state by simulating the time evolution governed by the imaginary-time Hamiltonian:

[TABLE]

for a long enough annealing time $\tau$. Here, $\hat{\mathcal{L}}_{S}$ denotes the Laplacian operator associated with the rewriting system $S$.

To assess the quality of the prepared ground states, we evaluate both their flatness, quantified by the Dirichlet energy $E_{\text{D}}=\langle\psi|\hat{\mathcal{L}}_{S}|\psi\rangle$, and their non-connected probability $p_{\mathrm{NC}}$, defined as the probability of sampling a string which is not connected to the input string when performing measurements in the computational basis, that is, an error in our classification.

In Figure 3, we present the behavior of the Dirichlet energy $E_{D}$ and the non-connected probability $p_{\mathrm{NC}}$ as a function of the system size. We consider different input strings $\omega$, as well as varying bond dimensions $\chi$ and annealing times $\tau$. As expected, increasing the annealing time systematically lowers the Dirichlet energy of the final state. When the bond dimension is sufficiently large, the non-connected probability also decreases toward zero, indicating that the evolution converges to the orbit state. Conversely, if the bond dimension is too small, the projection onto the MPS manifold can steer the evolution toward alternative ground states of the Laplacian, which are typically linear combinations of orbit states and exhibit a high non-connected probability $p_{\mathrm{NC}}$. An extreme example of this behavior is the equal-amplitude superposition of all computational basis states, which is a ground state of the Laplacian and can be exactly represented with bond dimension $\chi=1$.

For the remainder of this analysis, we consider only MPS with both low non-connected probability ($p_{\mathrm{NC}}\leq\varepsilon_{NC}=0.05$) and low Dirichlet energy ($E_{D}\leq\varepsilon_{D}=0.0002$). We select the MPS generated in the shortest time among those that exceed the quality tolerance. This selection allows us to estimate the minimal computational time and memory required to classically encode orbit states via TN. As the cost of computing fidelities and expectation values is negligible compared to that of state preparation, the total runtime can be interpreted as the time-to-solution for both the word and counting problems.

The computational resources required to prepare orbit states via TN techniques are summarized in Figure 4. In the figure, we observe that increasing the annealing time from $\tau=1000$ to $\tau=3000$ suffices to scale the computation from systems of length $L=10$ to $L=100$. This trend suggests that the annealing time $\tau$ grows slowly with system size, supporting the possibility of efficient scaling and motivating further exploration toward a physical implementation of the algorithm. In the same regime, the bond dimension required to approximate orbit states above the quality threshold increases from $\chi=32$ to $\chi=256$.

Although the computational costs shown in Figure 4 are consistent with polynomial scaling, our analysis is currently limited to systems of size up to $L=100$, constrained by the computational effort required to simulate larger instances. As a result, Figure 4 alone may not fully capture the asymptotic behavior of the algorithm. To strengthen the case for a polynomial scaling for the string rewriting system in the exam, we turn to the entanglement structure of the orbit states, which directly impacts the bond dimension $\chi$ required for an accurate MPS approximation. As shown in Figure 5, the maximum entanglement entropy across all bipartitions increases logarithmically with both the system size (Panel A) and the size of the largest bipartition of the system (Panel B). Moreover, for each bipartition, the Schmidt singular values exhibit an exponential decay (Panel C). Together, these observations imply that the number of singular values required to achieve a high-fidelity approximation grows only polynomially with the bipartition size—and hence with the overall system size. This favorable entanglement structure supports the efficient representability of orbit states using MPS with a bond dimension that scales polynomially.

Since the Laplacian inherits locality from the rewriting rules, the observed polynomial scaling of annealing time, bond dimension, and computational cost can be deduced by the scaling of the spectral gap of the Laplacian. Indeed, as shown in the Section Computational complexity and final energy gap of the Supplementary Text, the number of time steps required to approximate the orbit state associated with an input word $\tilde{\omega}$ scales as $\mathcal{O}(\Delta_{\tilde{\omega}}^{-2})$, where $\Delta_{\tilde{\omega}}$ is the smallest nonzero eigenvalue of the Laplacian restricted to the corresponding equivalence class $X_{S,\tilde{\omega}}$. Each time step is simulated via a time-dependent variational principle algorithm that has computational complexity $\mathcal{O}(\chi^{3})$ and memory complexity $\mathcal{O}(\chi^{2})$. We also show that the errors accumulated during the imaginary-time annealing process are exponentially suppressed, implying that the bond dimension $\chi$ needed to reach good accuracy is mainly determined by the final stage of the evolution, i.e., in the vicinity of the orbit state. Since the Laplacian is a local operator, the bond dimension required to approximate the orbit states scales polynomially in the inverse gap, more precisely as $\mathcal{O}(\Delta_{\tilde{\omega}}^{-1/3})$  (?). In this way, the overall computational complexity of the proposed algorithm is controlled by the final spectral gap $\Delta_{\tilde{\omega}}$. The value of this gap depends on the structure of the underlying rewriting system. For the rewriting system considered in our example, Figure 6 reports the gaps of the Laplacian restricted to each equivalence class, for system sizes up to $20$. These gaps govern the average-case complexity when the input word is sampled uniformly at random. The smallest gap scales polynomially with system size, implying polynomial worst-case computational complexity.

Once orbit states are available, measuring fidelities and observables enables efficient solutions to both the word problem and the counting problem. The results of these tasks are presented in Figure 7.

Word problem

In Figure 7 Panel A), we illustrate the solution of the word problem via fidelity measurements between orbit states. Pairs of words that belong to the same equivalence class correspond to orbit states with high fidelity, whereas disconnected word pairs exhibit low fidelity. These results demonstrate that our algorithm correctly solves the word problem. Greedy strategies based on explicit graph exploration, such as breadth-first search, require exponential time and cannot reach comparable system sizes. In contrast, the Knuth-Bendix completion algorithm solves the word problem in polynomial time and remains faster within the system sizes explored here (see Supplementary Text). However, a definitive performance comparison would require extending the analysis to larger system sizes and a broader range of rewriting systems.

Counting

In Figure 7 Panel B), we estimate the number of connected words for each input word. Words belonging to the same set are associated with the same number of connected words. We reconstruct connected sets containing up to $10^{28}$ words of length $100$. Storing this amount of information sequentially would require approximately $10^{17}$ terabytes of memory, underscoring the compression power of the quantum-inspired representation. For comparison, we also report the exact counts obtained via greedy enumeration for $L\leq 30$, which agree with our estimates. This baseline was computed using greedy graph exploration, i.e., breadth-first search with memoization  (?), which constitutes the state-of-the-art alternative for exact enumeration. Its computational cost grows polynomially with the connected set size—and hence exponentially with $L$—making it unfeasible for large-scale instances. In contrast, our quantum normal form reduction technique enables the enumeration of connected sets at scales inaccessible to graph exploration.

Discussion

We have introduced quantum normal form reduction, a general paradigm for automating equational reasoning on quantum computers. Our approach leverages the ability of quantum systems to encode and manipulate exponentially large sets of semantically equivalent symbolic expressions as a single quantum state, an orbit state. Orbit states are prepared as the ground states of suitable sparse Hamiltonians, which can be efficiently simulated on quantum devices when the rewriting system encoding the equivalence relations contains a polynomial number of rules in the string size $L$. Quantum optimization techniques are employed to prepare orbit states. As it is common in many optimization problems, the computational cost of these procedures depends on the specific problem instance.

We have simulated our algorithm using tensor network techniques, demonstrating its effectiveness in solving the word problem and the counting problem for a toy rewriting system. The results obtained from these simulations suggest the potential for a quantum advantage in equational reasoning, as well as for the development of novel quantum-inspired algorithms capable of outperforming classical approaches. While our results are promising, the present TN emulations are currently limited to strings up to $100$ characters. This limitation primarily stems from the growth in classical computational time associated with the increasing entanglement of orbit states, which appears to scale polynomially within the investigated size range. Further investigation is necessary to corroborate any general claim regarding quantum or quantum-inspired speedups. Interestingly, tensor network algorithms have been shown to outperform state-of-the-art solvers for a related counting problem, namely counting the solutions of SAT problems  (?, ?). Potential connections between these algorithms and the methods introduced in this work remain to be explored. The comparison with state-of-the-art classical algorithms (see Supplementary text) clarifies in which limits quantum normal form reduction can be regarded as a quantum extension of classical methods, and how insights from these approaches may be leveraged to improve the proposed quantum approach.

The tensor-network approach demonstrates the possibility of using computational tools from many-body physics in equational reasoning. For instance, one may employ density matrix renormalization group methods  (?) to approximate a random ground state of the Laplacian operator, i.e., a random superposition of orbit states. The fidelity between an input word and a random ground state obtained in this way would define a sound and complete equational hash, i.e., a function constant within each equivalence class but different across distinct classes. Once constructed, this function allows one to solve the word problem for any pair of words simply by comparing their hash values, without explicitly reconstructing the corresponding orbit states.

Among the many potential real-world applications outlined in the introduction, the design of quantum algorithms for formal language processing represents a particularly significant future direction. Since formal grammars capture the underlying structure of both human and programming languages, this could open new perspectives for the development of quantum algorithms in language processing  (?, ?) and software design. Moreover, grammar-based algorithms  (?) for lossless compression of classical data could be further improved through the quantum techniques introduced here.

Another promising future direction is the development of quantum algorithms for optimizing cost functions within an equivalence class. A preliminary example is the design of quantum algorithms for quantum circuit compilation  (?), where the set of circuits implementing the same unitary operator is explored through quantum dynamics generated by a set of Hermitian operators that encode rewriting rules.

Finally, our findings underscore the effectiveness of TN as a compressed representation for extensive datasets structured by equivalence rules, a property of increasing relevance in the era of massive data generation, and suggest promising directions toward real-world applications.

Beyond tensor network simulation, the proposed algorithm can also be implemented using quantum circuits (see Supplementary text). This circuit-based implementation relies on multi-controlled gates, which can be efficiently realized on universal quantum computers using a linear number of elementary gates and ancilla qubits  (?). Moreover, such gates may be natively supported on certain hardware platforms, including Rydberg-atom arrays and superconducting circuits  (?).

Materials and Methods

The graph Laplacian operator

Here we construct the graph Laplacian operator $\hat{\mathcal{L}}_{S}$ whose expectation value corresponds to the Dirichlet energy in Eq. (A parent Hamiltonian for orbit states).

First, we define the discrete Laplacian matrix

[TABLE]

where the degree $\deg(\omega_{k})$ is the number of edges attached to the vertex $\omega_{k}$, and the multiplicity $\text{mult}(\omega_{k},\omega_{k^{\prime}})$ is the number of edges connecting $\omega_{k}$ and $\omega_{k^{\prime}}$. The Dirichlet energy of $\psi$ can be expressed via the Laplacian matrix as

[TABLE]

By introducing the Laplacian operator

[TABLE]

we can finally express the Dirichlet energy as the expectation value of $\hat{\mathcal{L}}_{S}$ on the quantum state $\ket{\psi}$:

[TABLE]

By construction, the Dirichlet energy is always larger than or equal to zero, and is zero only for a function that is constant on each subgraph (see Eq. (A parent Hamiltonian for orbit states) ). Thus, the Laplacian operator is a positive semi-definite operator having the orbit states as degenerate ground states with zero energy.

As the graph $G$ is generated by the action of the rewriting rules, the Laplacian operator is a function of the rewriting rules in $R$. In particular, the entire Laplacian operator is the sum of Laplacian operators $\hat{\mathcal{L}}_{r}$ of the graphs induced by each single rule $\hat{r}$, since the degeneracy of a vertex is the sum of the degeneracies introduced by each rule, and the multiplicity of an edge is the sum of multiplicities. The Laplacian induced by the rule $r$ can be written as

[TABLE]

In this expression, the first two terms are diagonal in the computational basis and, for each basis state, count the number of configurations connected via rule $r$, thereby contributing to the vertex degrees in Eq. (34). The last two terms contribute to the off-diagonal structure of Eq. (34), assigning a weight of $1$ to pairs of basis states connected by rule $r$ and [math] otherwise. Finally, we write the Laplacian operator of the whole rewriting system as the sum of the Laplacians of each rule:

[TABLE]

The ground state of $\hat{\mathcal{L}}_{S}$ is also a ground state for each rule-associated operator $\hat{\mathcal{L}}_{r}$, so that $\hat{\mathcal{L}}_{S}$ is frustration-free.

Instances generation for the numerical experiment

In the numerical experiment, we consider different instances for the word problem and counting. Each instance consists of two disjoint sets of four strings of length $L$. Pairs of strings within the same set are connected through some sequence of rewriting operations, while no string from one set is connected to any string in the other. Each set of connected strings is generated by applying the Knuth-Bendix algorithm to reduce a randomly sampled initial string to distinct but connected strings, i.e., the normal forms produced by the Knuth-Bendix algorithms for different orderings (see Supplementary Text for details about normal forms and the Knuth-Bendix algorithms). We use the computer algebra system GAP  (?) to run the Knuth-Bendix algorithm. Additional random applications of the rewriting rules in $S$ are then performed to diversify the strings within each equivalence class. The construction of normal forms also allowed us to verify that no string in $\Omega_{1}$ is connected to any string in $\Omega_{2}$.

Simulation details

The imaginary quantum annealing evolution is simulated using the MPS formalism and a Time-Dependent Variational Principle (TDVP)  (?) via the tensor network emulator Quantum TEA Leaves  (?).

The Hamiltonian dynamics are discretized in time steps of size $\delta_{t}=0.5$. Note that the Trotter error introduced by $\delta_{t}$ decreases with the annealing time $\tau$. Indeed, increasing $\tau$ reduces the change in the Hamiltonian at each time step, thereby decreasing the magnitude of the commutator term arising from the Baker–Campbell–Hausdorff formula. For each input string, we perform simulations for annealing times $\tau\in\{1000,2000,3000\}$ and bond dimensions $\chi\in\{32,64,128,256\}$. Increasing $\tau$ and $\chi$ enhances the fidelity of the resulting orbit state with respect to the ideal target, but at the cost of increased computational resources.

All simulations were performed on a virtual machine equipped with 20 AMD EPYC 7413 CPUs and 128 GB of memory, using parallel execution across groups of eight input strings with identical length, bond dimension, and annealing time.

Acknowledgments

We acknowledge Nicola Assolini, Diego di Bernardo, Alessandra Di Pierro, Aleks Kissinger, Sergii Strelchuk and David Yu Yuan for useful discussions and valuable feedback.

Author Contributions:

Conceptualization: DR, MB, IS, DJ, SM. Methodology: DR, MB, DJ. Software: DR, DJ. Validation: DR. Formal analysis: DR. Investigation: DR. Resources: DJ. Data curation: DR. Writing - original draft: DR. Writing - review & editing: MB, IS, DJ, SM. Visualization: DR. Supervision: IS, DJ, SM. Project administration: IS, SM. Funding acquisition: SM.

Funding:

The research leading to these results has received funding from the following organizations: European Union via Italian Research Center on HPC, Big Data and Quantum Computing (NextGenerationEU Project No. CN00000013), project EuRyQa (Horizon 2020), project PASQuanS2 (Quantum Technologies Flagship); Italian Ministry of University and Research (MUR) via: Quantum Frontiers (the Departments of Excellence 2023-2027); the World Class Research Infrastructure - Quantum Computing and Simulation Center (QCSC) of Padova University; Istituto Nazionale di Fisica Nucleare (INFN): iniziativa specifica IS-QUANTUM; the German Federal Ministry of Education and Research (BMBF) via the project QRydDemo. We acknowledge computational resources from Cloud Veneto, as well as computation time on Cineca’s Leonardo machine.

Competing interests:

The authors declare that they have no competing interests.

Data and materials availability:

Data, software, simulation outputs, and all instructions needed to reproduce the results of this paper are available on Zenodo  (?).

Supplementary materials

Supplementary Text

References and Notes

Supplementary Materials for

Quantum algorithms for equational reasoning

Davide Rattacaso1,2∗, Daniel Jaschke1,2,3,4, Marco Ballarin1,2,5,

Ilaria Siloi1,2, Simone Montangero1,2

1Dipartimento di Fisica e Astronomia “G. Galilei” & Padua Quantum Technologies Research Center,

Università degli Studi di Padova, Italy I-35131, Padova, Italy.

2INFN, Sezione di Padova, via Marzolo 8, I-35131, Padova.

3Institute for Complex Quantum Systems, Ulm University, Albert-Einstein-Allee 11, 89069 Ulm, Germany.

4Current affiliation: PlanQC GmbH, Lichtenbergstr. 8, 85748 Garching, Germany.

5Current affiliation: Quantinuum, Partnership House, Carlisle Place, London SW1P 1BX, United Kingdom.

∗Corresponding author. Email: davide.rattacaso@unipd.it

Simulating the Hamiltonian

Here, we analyze the quantum resources required for simulating the Hamiltonian $\hat{H}_{S,\tilde{\omega}}$, specifically in terms of the number of quantum gates.

We do not consider native qudit platforms  (?), which might eventually prove advantageous for operating with rewriting systems whose alphabet contains more than two elements. Instead, we focus on standard qubit-based computation. To this end, we note that any rewriting system defined on an alphabet of size $d$ can be translated into a rewriting system over a binary alphabet, as happens when strings are manipulated in standard (Boolean) classical computers. In particular, each symbol of the original alphabet can be encoded as a binary string of length $\left\lceil\log_{2}{d}\right\rceil$, i.e., its binary logarithm rounded to the nearest bigger integer. This corresponds to the binary representation of the symbol index. For example, for the alphabet $A=\{a,b,c\}$, we can assign: $a\rightarrow(0,0)$, $b\rightarrow(0,1)$, and $c\rightarrow(1,0)$. Similarly, a string of $L$ characters over $A$ is mapped to a binary string of length $L\cdot\left\lceil\log_{2}{d}\right\rceil$. Each rewriting rule of length $l$ is likewise mapped to a rule of length $l\cdot\left\lceil\log_{2}{d}\right\rceil$, while the total number $n_{r}$ of rewriting rules remains unchanged.

With this binary encoding in place, we can restrict our analysis to simulating $\hat{H}_{S,\tilde{\omega}}$ over the binary alphabet $A=\{0,1\}$, which is compatible with standard qubit-based quantum devices.

The simulation typically relies on a discretization of the time evolution, allowing the approximation of the evolution operator as a product of exponentials of the individual Hamiltonian terms  (?). The Trotterized evolution, or any analogous optimization ansatz such as QAOA, can be implemented by the application of the operator

[TABLE]

where the $\alpha_{i}$ and $\beta_{i}$ are real coefficients scaling as the time step $\delta_{\tau}=\tau/N$, the indices $i$ run over the $N$ time steps, and $\tau$ is the total evolution time. The operator consists of one term involving the exponential of the projector $\ket{\tilde{\omega}}\bra{\tilde{\omega}}$, and a sequence of terms derived from the Laplacian. The number of Laplacian-related operators is twice the number $n_{r}$ of rewriting rules in the system, due to the squared and linear terms for each rule.

Thus, the entire evolution can be realized by implementing as a circuit consisting of $N\cdot(2n_{r}+1)$ operators

[TABLE]

where $\mathbf{b^{\prime}}$ and $\mathbf{b^{\prime\prime}}$ are binary substrings of size $l$ (for operators encoding rules) or $L$ (for the operator encoding the projection on the word $\tilde{\omega}$), and $N$ is the total number of time steps.

The operators in Eq. (S2) can be implemented in a quantum circuit (see Figure S1). First, we use a combination of X gates and a multi-controlled NOT gate to flip an ancilla qubit, conditioned on the local quantum state matching $\ket{\mathbf{b^{\prime}}}$. Next, we apply two-qubit gates with the control qubit being the ancilla to transform the state $\ket{\mathbf{b^{\prime}}}$ into $\mathrm{e}^{-i\theta}\ket{\mathbf{b^{\prime\prime}}}$. Finally, we reverse the ancilla operation by applying the same $X$ gates and multi-controlled NOT gate, effectively uncomputing the ancilla. Thus, when the operator $W$ to be implemented corresponds to a rewriting rule $r$, it requires $\mathcal{O}(w_{r})$ gates, where $w_{r}$ denotes the number of characters the rule acts upon. The implementation involves two multi-controlled gates. Similarly, if $W$ represents the evolution generated by the projection onto the input state, it requires $\mathcal{O}(L)$ gates, with $L$ being the system size. In this case as well, only two multi-controlled gates are needed. Multi-controlled gates can be implemented efficiently on universal quantum computers using a linear number of elementary gates and ancilla qubits  (?). Alternatively, they may be natively supported on certain hardware platforms, such as Rydberg-atom arrays and superconducting circuits  (?). In both scenarios, implementing the operator $W$ requires a number of long-range two-qubit gates that scales linearly with the length of the binary string $\mathbf{b^{\prime}}$.

Overall, simulating the Trotterized evolution over $N$ steps requires

[TABLE]

gates, where $w$ is the maximum number of characters affected by any rule in the rewriting system. These resources are linear in the number of bits needed to describe the problem instance classically—that is, the input string and the set of rewriting rules.

As a consequence of the Baker–Campbell–Hausdorff formula, the first-order Trotterized simulation in Eq. (S1) incurs, at each time step, a local error scaling as $\delta_{\tau}^{2}=\mathcal{O}(\tau^{2}/N^{2})$. Higher-order Trotter–Suzuki decompositions can be employed to systematically reduce this error. Restricting here to the first-order approximation, the cumulative error over the full time evolution scales as $\mathcal{O}(\tau^{2}/N)$. Therefore, in order to approximate the continuous-time dynamics within a target Trotter error $\epsilon$, one requires $N=\mathcal{O}(\tau^{2}/\epsilon)$ discrete time steps.

When the objective is ground-state preparation, and assuming a simple adiabatic evolution as the heuristic method of choice, the total evolution time required to reach the target ground state is expected to scale as $\Delta_{\min}^{-2}$  (?), where $\Delta_{\min}$ denotes the minimum energy gap of the Hamiltonian restricted to the dynamically accessible subspace. This subspace is spanned by computational basis states that are equivalent under the rewriting rules, to which the dynamics is constrained by construction. Combining the adiabatic time scaling with the first-order Trotterization error bound yields a total number of Trotter steps scaling as $N=\mathcal{O}(\Delta_{\min}^{-4})$. The resulting circuit depth can thus be bounded by

[TABLE]

This scaling can be easily improved through a variety of techniques, including higher-order Trotterization schemes, counterdiabatic driving and other shortcuts to adiabaticity, as well as heuristics better suited to near-term quantum hardware, such as the Quantum Approximate Optimization Algorithm.

Fidelity between orbit states on quantum annealers

Here, we show how the fidelity between orbit states can be measured on quantum annealers.

We prepare the orbit states $\ket{X_{S^{\prime},\omega^{\prime}}}$ and $\ket{X_{S^{\prime\prime},\omega^{\prime\prime}}}$ on a quantum annealer. To this aim, we initialize the annealer respectively in the states $\ket{\omega^{\prime}}$ and $\ket{\omega^{\prime\prime}}$. We evolve these initial states with the Hamiltonians $H_{S^{\prime},\omega^{\prime}}(t)$ and $H_{S^{\prime\prime},\omega^{\prime\prime}}(t)$ for large enough times $\tau^{\prime}$ and $\tau^{\prime\prime}$. In the adiabatic regime, we obtain

[TABLE]

and

[TABLE]

where $\mathcal{T}$ is the time-ordering operator.

Considering that $\mathcal{T}\left[\mathrm{e}^{-i\int_{0}^{\tau^{\prime}}dtH_{S^{\prime},\omega^{\prime}}(t)}\right]^{\dagger}=\mathcal{T}\left[\mathrm{e}^{-i\int_{0}^{\tau^{\prime}}dtH_{S^{\prime},\omega^{\prime}}(\tau^{\prime}-t)}\right]$, the fidelity between the orbit states is

[TABLE]

Now, we define the state

[TABLE]

This state can be prepared on the quantum annealer by evolving the initial state $\ket{\omega^{\prime\prime}}$ first with Hamiltonian $H_{S^{\prime\prime},\omega^{\prime\prime}}(t)$ for a time $\tau^{\prime\prime}$, and then with the Hamiltonian $H_{S^{\prime},\omega^{\prime}}(\tau^{\prime}-t)$ for a time $\tau^{\prime}$.

Once $\ket{\omega^{\prime\prime}_{S^{\prime},S^{\prime\prime}}}$ has been prepared, the fidelity in Eq. (Fidelity between orbit states on quantum annealers) is measured as the expectation value of the projector $P_{\omega^{\prime}}=\ket{\omega^{\prime}}\bra{\omega^{\prime}}$ onto the computational basis state $\ket{\omega^{\prime}}$, i.e.,

[TABLE]

This quantity is the probability of sampling $\omega^{\prime}$. It can be estimated by performing $N_{s}$ shots of the experiment and counting the frequency of outcomes corresponding to the configuration $\omega^{\prime}$. The error on this estimation scales as $\mathcal{O}(1/\sqrt{N_{s}})$.

Comparison to classical approaches

Here, we analyze connections and differences between quantum normal form reduction and the main classical approaches to the word problem and the counting problem. This comparison clarifies in which limits quantum normal form reduction can be regarded as a quantum extension of classical methods, and how insights from these approaches may be leveraged to improve quantum normal form reduction.

Graph exploration

The most direct classical approach is based on explicit graph-exploration algorithms, such as breadth-first search with memoization  (?). These algorithms sequentially explore the connected component of the rewriting graph containing the input string, allowing one to test equivalence and to enumerate or count all connected words. Since the size of the connected subgraph typically grows exponentially with the size of the rewriting system or the input, such approaches quickly become computationally infeasible.

Adiabatically preparing the ground state of the Laplacian associated with the rewriting system corresponds to exploring the graph of connected words, but in quantum superposition rather than via sequential traversal. From this perspective, quantum normal form reduction is more naturally compared to random walks on graphs, whose equilibrium distribution is uniform over connected components. In both the classical random-walk and the quantum settings, the computational complexity is governed by the bottlenecks of the graph, which determine the spectral gap of the Laplacian. Overall, the possibility of a quadratic quantum speedup relative to classical random walks has been extensively studied in the literature, although its realization depends sensitively on structural properties of the specific rewriting system  (?).

Canonical form reduction

More sophisticated approaches to the word problem are based on canonical form reduction  (?), that is, on constructing a procedure that maps all elements of the same equivalence class to a unique representative. This allows one to solve the word problem by comparing the normal forms associated with the input words. Quantum normal form reduction can be viewed as a quantum analogue of this paradigm, since it maps all equivalent strings to a single coherent quantum superposition. Unlike classical canonical forms, this superposition can address the counting problem, since it encodes information about the entire equivalence class.

Completion-based methods, most notably the Knuth–Bendix algorithm  (?), provide a systematic way to obtain such canonical forms by transforming the original rewriting system into one that is both terminating and confluent. Completion can fail or not terminate, limiting the applicability of this approach. In contrast to quantum normal form reduction, the Knuth–Bendix procedure requires the user to specify a suitable reduction ordering, and its success is highly sensitive to this choice. In Section Knuth-Bendix algorithm of this Supplementary Text, we compare the performance of the Knuth–Bendix algorithm with our approach for a specific case.

Automata-based methods

Automata-based methods consist of constructing an accepting automaton, i.e., a finite-state automaton that accepts exactly the strings connected from a given input, thereby solving the word problem  (?). Since each accepting path through the automaton corresponds to a distinct string, dynamic programming techniques can be used to count connected words.

The tensor-network representation of orbit states provides a bridge between quantum normal form reduction and automata-based constructions. Indeed, matrix product states are equivalent to weighted finite-state automata that compute functions on strings, where the bond dimension controls the number of internal automaton configurations. Thus, orbit states can be interpreted as accepting automata with exponentially large expressive power, since a generic quantum state corresponds to a matrix product state with a bond dimension exponential in the system size. Quantum normal form reduction defines a systematic procedure for implicitly constructing such automata.

Reduction to SAT

Boolean satisfiability (SAT)–based methods typically encode the predicate “$\omega_{1}$ rewrites to $\omega_{2}$ within $k$ rewrite steps” as a Boolean formula that is satisfiable if and only if such a derivation exists  (?). This encoding enables the use of highly optimized SAT solvers to address instances of the word problem. Furthermore, reductions to SAT allow one to enumerate satisfying assignments using blocking clauses, although counting an exponentially large set of connected strings in this way generally requires an exponential number of solver calls. Approximate counting can be achieved with a polynomial number of SAT solver invocations using hashing-based techniques. These methods estimate the number of satisfying assignments of a Boolean formula by randomly partitioning the solution space in buckets using XOR constraints, and then exactly counting the solutions in a randomly selected bucket whose expected size is bounded by a fixed threshold  (?).

Since the Laplacian can be written as a frustration-free sum of positive semi-definite operators $r^{2}-r$, preparing orbit states is a quantum analogue of SAT-based approaches, but with non-commuting clauses  (?). While the annealing-based scheme proposed here should therefore be regarded as a heuristic that may be effective in typical instances, exploring connections with modern SAT solver heuristics could improve this approach.

Unlike the solution of a SAT encoding, the ground state of the Laplacian is unique and naturally encodes a coherent superposition of all classical solutions. Estimating the number of connected strings nevertheless still requires repeated runs of the algorithm to sample the overlap between the orbit state and the uniform superposition over connected strings  (?). An analogue of approximate counting would be to restrict the Hilbert space to randomly chosen buckets, implemented by adding suitable penalty terms to the Laplacian operator.

Knuth-Bendix algorithm

The Knuth-Bendix algorithm  (?) is a state-of-the-art classical approach for solving the word problem. This procedure transforms the original rewriting system $S$ into a new and non-invertible system $S_{C}$, whose rules increase a specified total order on strings, such as lexicographic order. The new system is equivalent to the original one, meaning that two strings are connected by $S$ if and only if they are connected by $S_{C}$. Moreover, the transformed system $S_{C}$ is both terminating, meaning that no infinite sequence of rule applications is possible, and confluent, meaning that any sequence of valid rule applications yields the same final result, regardless of the order in which rules are applied.

For a rewriting system that is both confluent and terminating, two strings $\omega^{\prime}$ and $\omega^{\prime\prime}$ are equivalent if and only if repeated application of rewriting rules — regardless of the order in which they are applied — reduces both to the same string. This unique representative is called the normal form. Reduction to normal form thus enables efficient resolution of the word problem, provided that such a system $S_{C}$ can be constructed.

The Knuth-Bendix algorithm is not guaranteed to terminate: in general, the construction of $S_{C}$ may fail, reflecting the undecidability of the word problem for arbitrary rewriting systems. However, a fair comparison with the Knuth-Bendix algorithm for the purposes of this work must account for the finite size of the strings. This constraint limits the maximum length of the rewriting rules generated during the execution of the algorithm, ensuring termination within finite time and memory resources that depend on $L$.

As a benchmark for our quantum algorithm, we use the computer algebra system GAP  (?) to run the Knuth-Bendix algorithm for the string rewriting system in Eq. (Numerical results). The algorithm is executed for both the shortlex and recursive orderings, and for both possible permutations of the alphabet, $(a,b)$ and $(b,a)$ (see GAP documentation for further details). All executions were performed on a virtual machine equipped with 6 Intel(R) Core(TM) i5-8500 CPUs and 16 GB of memory. The total execution time until termination, as well as the memory footprint of the resulting confluent rewriting system $S_{C}$, are reported in Figure S2 for string lengths up to $L=400$. Across different choices of ordering, the computational time scales asymptotically as $\mathcal{O}(L^{\sim 6.4})$. The number of rules in $S_{C}$ grows as $\mathcal{O}(L^{\sim 2.1})$, while the total memory required to store the system — measured by the cumulative length of all rules — scales as $\mathcal{O}(L^{\sim 3.0})$.

Computational complexity of the imaginary quantum annealing

Here, we bound the computational complexity of imaginary quantum annealing (IQA) with respect to the energy gap and the fidelity susceptibility of the driving Hamiltonian ground state.

As in the main text, we consider IQA discretized in $N$ steps. We fix the time duration $\delta_{\tau}$ of each step, so that the total annealing time is $\tau=\delta_{\tau}N$.

The system’s Hamiltonian at the step $s$ is

[TABLE]

where $\delta=1/N$ is the variation of the Hamiltonian parameter per step.

We call $\ket{\psi_{s}}$ the state of the system at step $s\in[0,\dots,N]$, $\ket{0_{s}}$ the ground state of the Hamiltonian $\hat{H}_{s}$, and $E_{n,s}$ the $n$-th energy level of $\hat{H}_{s}$.

We measure the error in approximating the ground state at a step $s$ as the infidelity $1-F_{s}$ between $\ket{\psi_{s}}$ and $\ket{0_{s}}$, where the fidelity $F_{s}$ is defined as

[TABLE]

The state evolution at each step is

[TABLE]

so that the fidelity is

[TABLE]

The numerator of the last equation can be written as

[TABLE]

Let $P_{s}=\ket{0_{s}}\bra{0_{s}}$ be the projector on the ground state $\ket{0_{s}}$, and $P_{s}^{\perp}=1-\ket{0_{s}}\bra{0_{s}}$ its orthogonal complement. We have:

[TABLE]

so that the denominator becomes

[TABLE]

where at the second line we exploited the equation $P_{s+1}H_{s+1}P_{s+1}^{\perp}=0$, and, at the third line, we consider that $\bra{\psi_{s}}P_{s+1}\mathrm{e}^{-2H_{s+1}\delta_{\tau}}P_{s+1}\ket{\psi_{s}}=\bra{\psi_{s}}P_{s+1}\ket{\psi_{s}}\mathrm{e}^{-2E_{0,s+1}\delta_{\tau}}$ and $\bra{\psi_{s}}P_{s+1}^{\perp}\mathrm{e}^{-2H_{s+1}\delta_{\tau}}P_{s+1}^{\perp}\ket{\psi_{s}}\leq\bra{\psi_{s}}P_{s+1}^{\perp}\ket{\psi_{s}}||P_{s+1}^{\perp}\mathrm{e}^{-2H_{s+1}\delta_{\tau}}P_{s+1}^{\perp}||_{\text{op}}\leq\bra{\psi_{s}}P_{s+1}^{\perp}\ket{\psi_{s}}\mathrm{e}^{-2E_{1,s+1}\delta_{\tau}}$, where $||A||_{\text{op}}$ is the operator norm of $A$.

Substituting the corresponding terms in Eq. (S13) with the results of Eq. (S14) and Eq. (Computational complexity of the imaginary quantum annealing), we obtain

[TABLE]

where $\Delta_{s+1}$ is the first energy gap of $\hat{H}_{s+1}$.

Thus, the infidelity can be bounded as

[TABLE]

Now, we want to relate the infidelity at the step $s+1$ to the infidelity at the step $s$.

First, we decompose $\ket{\psi_{s}}$ on its parallel and orthogonal component with respect to $\ket{0_{s}}$, that is, $\ket{\psi_{s}}=\sqrt{F_{s}}\ket{0_{s}}+\mathrm{e}^{i\theta_{s}}\sqrt{1-F_{s}}\ket{e_{s}}$ for some complex phase $\theta_{s}$. Thus we have

[TABLE]

At this point, we introduce the fidelity susceptibility $f_{s}$ for the ground state path of $H_{s}$ as a measure of the infinitesimal variation of the ground state  (?), i.e.:

[TABLE]

which implies

[TABLE]

where $\ket{\varepsilon}$ is the orthogonal part of $\ket{0_{s+1}}$ with respect to $\ket{0_{s}}$ and $\theta^{\prime}_{s}$ is a complex phase.

Equation (S21) holds whenever the ground state path $\ket{0_{s}}$ is derivable. Being the Hamiltonian $H_{s}$ derivable, this happens whenever there is no level crossing between the first and the second energy level. In the case under exam here, where $H_{s}$ has negative non-diagonal entries, the Perron–Frobenius theorem guarantees the uniqueness of the ground state at each $s$ and therefore the absence of level crossing. Considering together Eq. (S19) and Eq. (S21) we obtain

[TABLE]

Exploiting the triangular inequality and considering that for small enough $\delta$, that is for large $N$, the modulus of the second addend is smaller than the modulus of the first, we obtain

[TABLE]

Now we substitute Eq. (S23) in Eq. (Computational complexity of the imaginary quantum annealing) to obtain the relationship between fidelity at successive steps for large $N$:

[TABLE]

Finally, we iteratively apply the last equation to obtain the evolution of the infidelity.

Considering that at the step $s=0$ the system is in the exact ground state of the Hamiltonian, so that $F_{0}=1$, we have

[TABLE]

and for $s=2$, considering small $\delta$, we have

[TABLE]

Iterating this process, we get the following upper-bound for the final fidelity:

[TABLE]

Equation (S27) can be interpreted as follows. At each time step, an error proportional to the fidelity susceptibility of the ground state path is accumulated. This error is exponentially damped down during all the remaining evolution at a rate proportional to the instantaneous first energy gap.

Now, we recall that by construction, the operator $\hat{H}_{s}$ can be written in block diagonal form, where each block acts on a connected subgraph, i.e., on a single equivalence class. The dynamic is thus restricted to the subspace spanned by the equivalence class of the input word $\tilde{\omega}$. Thus, given the input word $\tilde{\omega}$, the fidelity susceptibility refers to the path of ground states in the corresponding sector of the Hilbert space, as well as the first energy gap.

Computational complexity and final energy gap

Equation (S27) is a bound that depends on the specific path of ground states in the exam, and, ultimately, on the input state $\ket{\tilde{\omega}}$. A less strong bound can be derived that does not depend on $\tilde{\omega}$ but only on the first energy gap of the Laplacian in the block where $\tilde{\omega}$ belongs. In the following, we call this gap the final gap $\Delta_{\tilde{\omega}}$.

First of all, we observe that in the last part of the dynamics, the accumulated error is damped at a rate that only depends on the gap of the Hamiltonian for $s\approx N$. Since the Hamiltonian norm is bounded, this gap can not deviate too much from $\Delta_{\tilde{\omega}}$. This observation can be placed on a more rigorous mathematical footing using Weyl’s perturbation theorem  (?), which establishes that the rate of change of eigenvalues is bounded by the operator norm of the variation of the Hamiltonian as

[TABLE]

We can upper-bound the Hamiltonian norm as

[TABLE]

where we used the triangular inequality for the operator norm, $n_{r}$ is the number of operators $\hat{r}^{2}-\hat{r}$ in the Laplacian, and each operator $\hat{r}^{2}-\hat{r}$ has norm $2$. Thus, Eq. (S28) becomes

[TABLE]

which implies

[TABLE]

Let $s_{*}$ be the minimum value of $s$ for which the last expression is positive, that is

[TABLE]

For $s<s_{*}$ we have

[TABLE]

and for $s\geq s*$

[TABLE]

We substitute these two bounds in Eq. (S27) to obtain

[TABLE]

For large $N$, the second summation can be written as an integral. With a change of variable, we obtain

[TABLE]

For large $N$, we can use Laplace’s method for approximating the integral by expanding the exponent around the endpoints of the domain, thus obtaining:

[TABLE]

Now we bound the fidelity susceptibility at the end of the process with respect to the final energy gap  (?) as

[TABLE]

which, in our case, considering the bound on the Laplacian norm in Eq. (S29), becomes

[TABLE]

Replacing the bound on the final fidelity in Eq. (S37) we obtain

[TABLE]

Then, the number of time steps needed to get a final infidelity $1-F_{N}\leq\epsilon$ is

[TABLE]

Thus, when the number of rules in the rewriting system grows polynomially, the computational complexity scales polynomially with the inverse of the final minimum energy gap. This result holds even if the fidelity susceptibility diverges exponentially during the evolution, as a consequence of an exponentially closing gap. By contrast, this advantage is lost in real-time quantum annealing: due to the unitary nature of the dynamics, errors accumulated at intermediate times cannot be dissipated at the end of the process.