# The Complexity of Being Entangled

Stefano Baiguera, Shira Chapman, Giuseppe Policastro, Tal Schwartzman

► To cite this version:

Stefano Baiguera, Shira Chapman, Giuseppe Policastro, Tal Schwartzman. The Complexity of Being Entangled. Quantum, 2023, 8, pp.1472. 10.22331/q-2024-09-12-1472. hal-04276593

HAL Id: hal-04276593

https://hal.science/hal-04276593v1

Submitted on 2 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The Complexity of Being Entangled

Stefano Baiguera¹, Shira Chapman¹, Giuseppe Policastro², and Tal Schwartzman¹

¹Department of Physics, Ben-Gurion University of the Negev,
David Ben Gurion Boulevard 1, Beer Sheva 84105, Israel

²Laboratoire de Physique de l'École Normale Supérieure,
CNRS, PSL Research University and Sorbonne Universités,
24 rue Lhomond, 75005 Paris, France

Nielsen's approach to quantum state complexity relates the minimal num-
ber of quantum gates required to prepare a state to the length of geodesics
computed with a certain norm on the manifold of unitary transformations. For
a bipartite system, we investigate binding complexity, which corresponds to
norms in which gates acting on a single subsystem are free of cost. We reduce
the problem to the study of geodesics on the manifold of Schmidt coefficients,
equipped with an appropriate metric. Binding complexity is closely related to
other quantities such as distributed computing and quantum communication
complexity, and has a proposed holographic dual in the context of AdS/CFT.
For finite dimensional systems with a Riemannian norm, we find an exact re-
lation between binding complexity and the minimal Rényi entropy. We also
find analytic results for the most commonly used non-Riemannian norm (the
so-called $F_1$ norm) and provide lower bounds for the associated notion of state
complexity ubiquitous in quantum computation and holography. We argue that
our results are valid for a large class of penalty factors assigned to generators
acting across the subsystems. We demonstrate that our results can be bor-
rowed to study the usual complexity (not-binding) for a single spin for the case
of the $F_1$ norm which was previously lacking from the literature. Finally, we de-
rive bounds for multi-partite binding complexities and the related (continuous)
circuit complexity where the circuit contains at most 2-local interactions.

Stefano Baiguera: baiguera@post.bgu.ac.il

Shira Chapman: schapman@bgu.ac.il

Giuseppe Policastro: giuseppe.policastro@ens.fr

Tal Schwartzman: taljios@gmail.com
# 1 Introduction

Computational circuit complexity has been established as a useful tool of quantum computing theory for some time, see e.g., [1, 2]. Heuristically, complexity estimates the number of unitary gates (or time) required by an optimal circuit to produce a certain target state (typically entangled) starting from a reference state (typically factorized). As such, it is natural to ask whether there is a relation between the degree of entanglement and the complexity in producing a certain state [3]. In this work, we will show that under certain circumstances there is a relation between these two measures.

**Circuit complexity.** Circuit complexity plays a crucial role in classifying computational problems and in establishing the advantages of quantum over classical computation. From a practical point of view, finding the optimal circuit for a given task is an important goal for the compiler of any quantum computer.

In concrete circuits built in laboratories, quantum operations are usually performed using discrete elementary unitary operations (called *gates*) which act on one or two qubits at the same time. While it is not possible to cover all the unitary group with a finite set of gates, any unitary operator can be built with arbitrary precision by selecting a universal set of gates [4]. The implementation of those gates in a quantum circuit can vary in terms of difficulty, time, or errors. A proper definition of circuit complexity should take these variations into account. For generic circuits, finding the complexity is a hard problem.

Nielsen showed that lower and upper bounds can be obtained by looking at a continuous version of the problem [5–7]. The idea is to consider the unitary group, generate a continuous trajectory by a time-dependent Hamiltonian, and define a norm *F* which assigns different costs (called *penalty factors*) to the various directions along the tangent space of the group manifold. The shortest circuit is given by a geodesic in the space of unitaries, and the complexity will be the length of this geodesic, measured using the notion of distance induced by the norm with penalty factors. A similar definition of complexity can be introduced in the space of states by performing a quotient over the stabilizer at each point along the trajectory. A closely related notion to Nielsen’s approach is the framework of quantum optimal control, which studies the optimal manipulation of quantum dynamics under certain constraints [8–10].

**Complexity in holography.** Quantum complexity has gained much attention in holography after the observation that entanglement entropy in a conformal field theory (CFT) is insufficient to probe certain regions behind the horizons of dual black holes in Anti-de Sitter (AdS) space [11–14]. While entanglement entropy is unable to capture the long-time linear growth of the interior volume of black holes, it was proposed that quantum computational complexity could be the quantity dual to the volume of the Einstein-Rosen bridge (this is the so-called complexity=volume conjecture) [15, 16].<fnref n="1" /> Later on, the volume conjecture was supplemented by other holographic proposals for complexity [19–27], and their properties were compared by investigating their behavior in various backgrounds, see e.g., [28–53]. In parallel, efforts were made to put computational complexity on firmer ground from a quantum mechanical perspective and develop tools for calculating it in various setups. This led to several investigations in quantum mechanics (QM), quantum field theory (QFT), and CFT [26, 54–76]. For a comprehensive review, we refer the reader to [77].

<footnote n="1">¹See, however, [17, 18] for alternative boundary interpretations of the bulk volume of the Einstein-Rosen bridge.</footnote>
**Binding Complexity.** In this work we will investigate *Nielsen's binding complexity* (referred to, from now on, as binding complexity) [3], a specific version of Nielsen's geometric complexity. Given a partition of a system into multiple parts, binding complexity assigns non-vanishing cost only to operations that are non-local, *i.e.*, that act on multiple subsystems.<fnref n="2" /> Contrarily, it is free to act with operations that are local (*i.e.*, act within a single subsystem). While these requirements may seem unrealistic, we will find that the geodesics minimizing binding complexity do not present abrupt fluctuations along the free local directions during the time evolution. As a consequence, these geodesics would be a good approximation to the geodesics of an implementable circuit where the local operations have a small but non-zero cost.

There are several motivations to consider binding complexity. First of all, some of the main challenges for large-scale quantum computers (with many qubits) are to maintain coherence, *i.e.*, minimize the noise throughout the computation, and implement efficiently the necessary gates. Therefore, a prominent paradigm for future quantum computing is *distributed computing* [78–80], where the task of computation is distributed between multiple small quantum computers or nodes, and in each such computer, the processing is coherent to a good degree. In such a model, it is necessary to minimize the amount of gates between the different computers, which are often assumed to be less reliable or slower. This minimization is precisely the notion of binding complexity. Distributed computing is also related to *quantum communication complexity*, which measures the minimal information that distant parties have to exchange to perform a computational task [81]. It was argued in [3] that quantum communication complexity is upper bounded by the binding complexity. Secondly, we will show that binding complexity for a bipartite system leads to a geometric distance on the space of states which is fully determined by the entanglement spectrum. This has an operational interpretation in terms of the minimal length of a circuit required to turn a non-separable state into a separable one (see *e.g.*, [82], which is related to the case of the Riemannian norm discussed in section 4.3). Finally, binding complexity has been proposed to be dual to the interior volume of multi-boundary wormholes [3]. Thus one may hope to apply the methods developed in this work to unveil new features of the wormhole interior in multi-boundary wormhole geometries [83].

In [3], binding complexity was explored for Gaussian unitaries and states of a bosonic theory. Rough bounds for the binding complexity were also obtained in [84] using a discrete set of gates. Bounds on the minimal number of non-local gates required to optimize a quantum circuit were related in [85] to a sum (or maximum) of the entanglement entropy over all partitions. We will discuss the relation of [85] to our results in section 6.

In the present work, we compute binding complexity for bipartite finite dimensional systems and unveil its precise relation to the entanglement spectrum. One of the main novelties of our approach is that we will study the geometry of binding complexity using the so-called $F_1$ norm, which heuristically counts the number of gates required to perform a certain unitary operation, and find precise analytic results in this case. This quantity has received less attention in the literature because it is not a Riemannian norm, therefore certain analytic tools are not available in this case. In addition, we show that the binding complexity computed with the Riemannian $F_2$ norm is related to the minimal Rényi entropy of the reduced density matrix of the state. We argue that the $F_1$ binding complexity is, to some extent, independent of the penalty factors associated with the non-local generators acting between the two parties, up to a proportionality factor coming from the smallest of such penalty factors. We incidentally show that the methods developed to study the binding complexity of two qutrits can also be used to compute the $F_1$ (not

<footnote n="2">²Nonetheless, we will always require the gates to act at most on two subsystems at the same time.</footnote>
binding) complexity of a single qubit, which was missing from the previous Nielsen's complexity literature. Finally, we give lower bounds for multipartite binding complexity and demonstrate further inequalities for some more traditional notions of complexity using at most 2-local gates between different parts. The lower bounds are in terms of the bipartite binding complexity of various bipartitions of the state. A detailed summary of our results can be found at the beginning of section 7.

**Outline.** The paper is organized as follows. We begin in section 2 with a review of Nielsen's geometric approach to quantum computational complexity, and the definition of binding complexity. In section 3, we show that the computation of binding state complexity reduces to a problem defined on the manifold of Schmidt coefficients. In section 4, using this general technique, we compute the norm and the binding complexity for bipartite quantum-mechanical systems with generic dimensions. We apply the general method to study the binding complexity for two qutrits in section 5 and show that our results can be borrowed to study the usual notion of state complexity for a single qubit. In section 6 we derive bounds on some measures of complexity (and binding complexity) in the case of multipartite systems. We summarize the main results and discuss future directions in section 7. The appendices contain technical results. We derive a chain of inequalities for the cost functions in appendix A, which are required for the derivation of section 3. In appendix B, we consider the limit of a large number of Schmidt coefficients. In appendix C, we derive one of the bounds on the $F_1$ complexity in the Pauli basis mentioned in section 4.2. Appendix D contains an alternative derivation of certain results presented in section 4.3 by using the Euler-Arnold formalism.

# 2 Complexity geometry

We briefly review Nielsen's geometric approach to complexity [5-7] and set up the relevant notation. The idea of this approach is to replace circuits composed of discrete gates with continuous paths in a Lie group manifold. Complexity is then defined as the length of the shortest trajectory in a geometry defined over this manifold, once a notion of distance is introduced. While most of the discussions could apply to a generic Lie group, we will focus in the following on the special unitary group.

## 2.1 Unitary and state complexities

Given the special unitary group SU(N), we will be interested in finding the optimal way to construct a target unitary $U_T$ starting from the identity by means of a time-dependent Hamiltonian, which can be expanded in terms of a basis of Hermitian generators {$T_I$} as follows

$$H(t) = \sum_{I} Y^{I}(t)T_{I}, \qquad (1)$$

where $Y^I$ are called velocities and describe the tangent vector to a trajectory in the group manifold. Different ways to construct the target unitary are parameterized as curves in the space of unitary transformations with a path parameter $t \in [0, 1]$ such that $U(t=0) = \mathbf{1}$ and $U(t=1) = U_T$. The unitary at a point on the curve is expressed as

$$U(t) = \tilde{P} \exp \left( -i \int_{0}^{t} dt' H(t') \right), \qquad (2)$$

where $\tilde{P}$ denotes the path ordering such that the circuit is built from right to left. The time evolution of a generic curve is determined by the Schrödinger equation $\dot{U}(t) = -iH(t)U(t)$.
We define *unitary complexity* as the minimization of a cost function $F[\vec{Y}]$ introduced on the tangent space:

$$ C_F[U] = \min_{\{Y^I: U(0)=\mathbf{1}, U(1)=U_T\}} \int_0^1 dt F[Y^I(t)]. \quad (3) $$

The path in eq. (3) is subject to the boundary conditions imposing that it connects the identity with the target unitary.<fnref n="3" /> The cost function should satisfy certain properties so that it defines a notion of distance in the space of unitaries. A large class of cost functions that we will study in this work is given by

$$ F_{p,\vec{q}}[\vec{Y}] = \left( \sum_I q_I |Y^I|^p \right)^{\frac{1}{p}}, \quad (4) $$

where $q_I \ge 0$ are called *penalty factors*. Roughly speaking, penalties parametrize the difficulty of moving along the various directions of the group manifold and are supposed to model systems in which certain gates are more difficult to realize than others.

When $p=2$, the cost function (4) reduces to the distance induced by the Riemannian metric

$$ ds^2 = M_{IJ} \rho^I \rho^J, \quad \rho^I = \frac{i}{N} \mathrm{tr} (dU U^\dagger T^I), \quad M_{IJ} = \mathrm{diag}(q_1, \dots, q_{N^2-1}), \quad (5) $$

where $\rho^I$ are right-invariant forms, $M_{IJ}$ is the matrix containing the penalties and we have assumed that the generators are normalized according to $\mathrm{tr}(T^I T^J) = N \delta^{IJ}$ (that is $\rho^I = Y^I dt$). When $M_{IJ} = \delta_{IJ}$, the expression (5) corresponds to the Cartan-Killing metric on the unitary group, which is bi-invariant and has positive curvature. In the presence of non-trivial penalty factors, the metric is still right-invariant and therefore describes a homogeneous space, but now the curvature can also be negative [62]. Indeed, negative curvature is required to reproduce geodesic deviation [86] and the switchback effect [16, 87], which are the standard features of complexity. When $p=1$, the cost function (4) has the physical interpretation of counting the number of infinitesimal gates being used. For these reasons, the norms with $p=1, 2$ received special attention in the literature, see, e.g., [26, 55, 56, 62, 63, 66, 70–76], and we shall also focus on them here. Another property of the case $p=1$ is that the cost function is not smooth, and so we are not able to use the Euler Lagrange equations or calculus of variations to find geodesic trajectories. This makes its analysis more complicated, see e.g., [5] for a discussion on this point. In this paper, we develop efficient techniques to overcome this difficulty.

Let us now assume that the special unitary group $SU(N)$ acts on an $N$-dimensional Hilbert space $\mathcal{H}$, where $|\psi_R\rangle$ is a given reference state and $|\psi_T\rangle$ is a target state which we wish to construct. The *state complexity* is defined as the minimum of the unitary complexities computed for all the special unitary operators that connect the reference state to the target state:

$$ C_F^{\mathrm{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{\{U \in G: |\psi_T\rangle = U |\psi_R\rangle\}} C_F[U]. \quad (6) $$

Recall that the space of quantum states is defined as the space of rays in the Hilbert state, i.e., vectors differing by a phase are identified. Elements of the special unitary group that leave a state invariant define an equivalence class

$$ V |\psi\rangle = e^{i\phi} |\psi\rangle \quad \Rightarrow \quad U' = UV \sim U. \quad (7) $$

<footnote n="3">³In general the trajectory in the space of unitaries will run along the time interval $t \in [0, t_f]$. If the cost function $F[\vec{Y}]$ is positive homogeneous of degree 1, we can always redefine the parameter $t$ in such a way that $t_f = 1$.</footnote>
We call $V$ the elements of the stabilizer of the state $|\psi\rangle$, which corresponds to the maximal subgroup $SU(N-1) \times U(1)$ of $SU(N)$.<fnref n="4" /> The equivalence class in eq. (7) defines a map from the unitary space to the quotient

$$ \pi: SU(N) \to \mathbb{CP}^{N-1} \equiv \frac{SU(N)}{SU(N-1) \times U(1)}. \qquad (8) $$

The procedure (6) induces a norm for the state complexity that depends only on the coordinates of the projective space. While the concept of complexity on a state requires the notion of unitary complexity, the opposite is not true: we can define the complexity of a unitary without reference to a particular representation.

For a given trajectory in the Hilbert space, the minimization over the stabilizer group can be done locally, and thus (4) induces a norm $F^{\text{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle]$ on the tangent space of $\mathcal{H}$. The Schrödinger equation $|\dot{\psi}(t)\rangle = -iH |\psi(t)\rangle$ constrains $2(N-1)$ real degrees of freedom out of the $N^2-1$ degrees of freedom of $H$. Those are the degrees of freedom that control the evolution in the Hilbert space. Minimizing over the remaining $(N-1)^2$ degrees of freedom, which belong to the tangent space of the stabilizer group with respect to $|\psi(t)\rangle$, gives the induced norm

$$ F^{\text{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle] \equiv \min_{\text{stab } |\psi(t)\rangle} F[\vec{Y}]. \qquad (9) $$

With this definition, we can re-express the distance between states as

$$ C_F^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{\substack{|\psi(t)\rangle \\ |\psi(0)\rangle = |\psi_R\rangle \\ |\psi(1)\rangle = |\psi_T\rangle}} \int_0^1 dt \, F^{\text{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle], \qquad (10) $$

where now we can study optimal trajectories directly in terms of a geometry on the space of states. When state and unitary complexity are defined using the cost function (4) with $p=2$, the projection $\pi$ in eq. (8) is a Riemannian submersion [66]. The induced metric obtained on $\mathbb{CP}^{N-1}$ is only left-invariant for trivial penalty factors [62], in which case it is proportional to the Fubini-Study metric.

## 2.2 Binding complexity

Let us now assume that the physical system is split into $n$ subsystems, i.e., the Hilbert space decomposes as

$$ \mathcal{H} = \mathcal{H}_{A_1} \otimes \cdots \otimes \mathcal{H}_{A_n}. \qquad (11) $$

This setting is relevant to compare complexity with multiparty entanglement. Heuristically, in a discrete model we define binding complexity as the minimal number of gates, acting on more than one subsystem at the same time, that are needed to build the target unitary $U_T$ starting from the identity [3]. At the same time, we will require gates to be bi-local, meaning that they can act at most on two subsystems at the same time. This is a special case of the unitary complexity defined above with specific choices of penalties. Let us focus on the case where the Hilbert space factorizes into two parts (we return to the general case in section 6)

$$ \mathcal{H} = \mathcal{H}_A \otimes \mathcal{H}_B, \qquad N_A \equiv \dim(\mathcal{H}_A) \le \dim(\mathcal{H}_B) \equiv N_B. \qquad (12) $$

<footnote n="4">⁴The stabilizer is defined as the subgroup of a group $G$ that leaves invariant an element $x \in X$, where $X$ is a set on which the group acts.</footnote>
Each subspace admits a basis of Hermitian traceless generators {$T_a^A, T_i^B$} for the Lie subalgebras $su(N_A)$ and $su(N_B)$, respectively.<fnref n="5" /> We refer to those as the *local* generators since they only act on one of the two sides of the system, and we will call local unitary (LU) the transformations they span (in the sense of acting locally within each part). We choose a set that is orthogonal with respect to the trace in the fundamental representation

$$ \mathrm{Tr}(T_a^A T_b^A) = \mathcal{N}_A \delta_{ab}, \quad \mathrm{Tr}(T_i^B T_j^B) = \mathcal{N}_B \delta_{ij}, \qquad (13) $$

where $\mathcal{N}_{A(B)}$ are normalization constants. We do not assume that the generators are normalized, since we will later want to use the Pauli and Gell-Mann bases for $SU(N_{A(B)})$ which oftentimes appear in the literature unnormalized.

The unitary operators that act on the entire system belong to the group $SU(N_A N_B)$, which also entails operations that can entangle the two parts. The latter are spanned by generators of the form $T_a^A \otimes T_i^B$, which we call *non-local*.<fnref n="6" /> The full Hamiltonian of the system can be decomposed as

$$ H(t) = Y_A^a(t) T_a^A \otimes \mathbf{1}^B + Y_B^i(t) \mathbf{1}^A \otimes T_i^B + Y^{ai}(t) T_a^A \otimes T_i^B. \qquad (14) $$

The splitting of the velocities into the subsets $Y^I \equiv \{Y_A^a, Y_B^i, Y^{ai}\}$ implies a similar decomposition of the penalty factors into the corresponding subsets $q_I \equiv \{q_a^A, q_i^B, q_{ai}\}$, respectively. Explicitly, the cost function (4) becomes

$$ F_{p,\vec{q}}[\vec{Y}] = \left( \sum_a q_a^A |Y_A^a|^p + \sum_i q_i^B |Y_B^i|^p + \sum_{a,i} q_{ai} |Y^{ai}|^p \right)^{\frac{1}{p}}. \qquad (15) $$

We consider the case where it is free to act on each subsystem separately, and therefore only the non-local generators contribute to the unitary complexity. The corresponding penalty factors are<fnref n="7" />

$$ q_a^A = q_i^B = \varepsilon \to 0, \quad q_{ai} \ge 1. \qquad (16) $$

We have chosen $q_{ai} \ge 1$ rather than any other finite value without loss of generality. This choice of penalty factors describes a framework where the implementation of gates that entangle two subsystems is much harder than performing operations on each of them separately, as we have discussed in the introduction.

Nielsen's binding complexity ($\mathcal{BC}$) is defined as the unitary complexity (3) with penalty factors (16), i.e.,<fnref n="8" />

$$ \mathcal{BC}_F[U] = \lim_{q_a^A=q_i^B \to 0} C_F[U]. \qquad (17) $$

We will focus on two norms of this kind:

<footnote n="5">⁵We will usually denote the generators of the special unitary group in subsystem A with indices $a, b \in \{1, \dots, N_A^2 - 1\}$ from the beginning of the alphabet, and the generators in B with indices $i, j \in \{1, \dots, N_B^2 - 1\}$ from the middle part of the alphabet. Capital indices $I, J \in \{1, \dots, (N_A N_B)^2 - 1\}$ will collectively refer to the special unitary group $SU(N_A N_B)$ of the full system. We will also denote with a sub(super)script A, B the subsystem where the local generators act.</footnote>

<footnote n="6">⁶In the nomenclature of [3, 84] these are the *relevant* or *straddling* gates.</footnote>

<footnote n="7">⁷Choosing vanishing penalty factors along certain directions in the group manifold has other interesting applications. For instance, in [88] such a choice was used to relate Nielsen complexity to Krylov complexity.</footnote>

<footnote n="8">⁸In reference [3], Nielsen's binding complexity corresponds to what we call the homogeneous norm defined by the cost function (20a). In this work, we extend the definition of binding complexity to include any choice of penalties for the non-local generators.</footnote>
• The homogeneous case assigns unit penalty factors to all the non-local generators

$$F_{p,\text{hom}}[\vec{Y}] = \left( \sum_{a,i} |Y^{ai}|^p \right)^{\frac{1}{p}}. \qquad (18)$$

• The universal-set case is defined by the cost function

$$F_{p,\text{set}}[\vec{Y}] = |Y^{11}|, \qquad (19)$$

where the set $\vec{Y}$ does not include non-local gates other than $Y^{11}$. This corresponds to the scenario where all the non-local generators are forbidden, except for a single one (with velocity $Y^{11}$) which, together with the local generators, can be shown to form a universal set of gates (the subscript set stands for universal set of gates). In this way, it is still possible to reach all the unitary operators in the group SU($N_A N_B$). We immediately see that $F_{p,\text{set}}$ is independent of $p$; therefore, it is sufficient to study one special case to get information about all the other values of $p$.

It is possible to recover these norms as formal limits of the cost function (15) as follows:

$$F_{p,\text{hom}}[\vec{Y}] = \lim_{\epsilon \to 0} \left( \sum_a \epsilon |Y_A^a|^p + \sum_i \epsilon |Y_B^i|^p + \sum_{a,i} |Y^{ai}|^p \right)^{\frac{1}{p}}, \qquad (20a)$$

$$F_{p>1,\text{set}}[\vec{Y}] = \lim_{\substack{\epsilon \to 0 \\ q_{ai} \to \infty}} \left( \sum_a \epsilon |Y_A^a|^p + \sum_i \epsilon |Y_B^i|^p + \sum_{(a,i) \neq (1,1)} q_{ai} |Y^{ai}|^p + |Y^{11}|^p \right)^{\frac{1}{p}}. \qquad (20b)$$

In the latter case, we relied on the results in [89], where it was shown that the formal $q_{ai} \to \infty$ limit forces the geodesics in the group manifold to live in the subspace of the tangent bundle where the directions associated to the penalties $q_{ai}$ are excluded.<fnref n="9" />

By applying the general definition of state complexity (6) to the case with penalty factors (16), we get the binding complexity in the space of states

$$\mathcal{BC}_F^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{\{U \in \text{SU}(N): \ |\psi_T\rangle = U |\psi_R\rangle\}} \mathcal{BC}_F[U]. \qquad (21)$$

Although the space of states has fewer coordinates than the unitary group manifold, it is not homogeneous and the minimization procedure quickly becomes difficult to perform with the increase of the dimensionality [67]. However, a special feature of the current setting is the decomposition of the Hilbert space into a product according to eq. (12). This fact, combined with the assignment of vanishing cost to local operators, will reduce the computation of binding state complexity to the evaluation of geodesics on the group manifold SO($N_A$), as we will show in section 3.

<footnote n="9">⁹In reference [89], the authors performed the limit $q_{ai} \to \infty$, while at the same time keeping $q_{ai} Y^{ai}$ finite (equivalently, they kept the Hamiltonian bounded in operator norm) in order to recover the quantum brachistochrone equation by starting from the Euler-Arnold equations [90]. The application of this prescription gives $q_{ai} |Y^{ai}|^p \to 0$ in eq. (20b) for all $p > 1$. In this work, we will not make explicit use of this limit, instead we will directly assume that certain directions in the tangent space of the group manifold are forbidden.</footnote>
## 3 The orthogonal group picture

Our main task is to compute the binding state complexity (21) for a bipartite system (12) and eventually relate this quantity to the entanglement spectrum. We will show that the relevant group manifold for this analysis is the orthogonal group $SO(N_A)$, defined with respect to the maximal number of Schmidt coefficients $N_A$ characterizing the size of the smaller of the two subsystems. This section is the main core of the work since it contains the general strategy that we will apply in order to study all the specific cases presented later.

As we are interested in the state complexity, the norms expressed before should be minimized over the stabilizer group of the state at any given point in the trajectory in the manifold of unitaries and eventually expressed in terms of the degrees of freedom of the state. Given the decomposition of the Hilbert space (12), a natural and convenient parametrization for any state along a trajectory in the quotient space is given by the Schmidt decomposition

$$|\psi(t)\rangle = \sum_{m=1}^{N_A} \lambda_m(t) |\phi_m(t)\chi_m(t)\rangle, \qquad (22)$$

where $\lambda_m(t)$ (called *Schmidt coefficients*) are real and positive numbers, and $|\phi_m(t)\rangle, |\chi_m(t)\rangle$ form orthonormal linearly-independent sets for the Hilbert spaces $\mathcal{H}_A$ and $\mathcal{H}_B$, respectively. With some abuse of notation, we will later refer to $(\phi_m, \chi_m)$ as *basis parameters*.<fnref n="10" /> Note that for the state to be normalized, the Schmidt coefficients should lie on a sphere, i.e., $\sum_m \lambda_m^2 = 1$, and thus only $N_A - 1$ of them are independent. Any entanglement monotone on a bipartite system, when evaluated on pure states, can be written as a function of the Schmidt coefficients [91].

In the context of binding complexity, local operations (which act only on one subsystem) play a special role. It is then convenient to introduce the following terminology: two states $|\psi\rangle, |\psi'\rangle$ belonging to the Hilbert space (12) are called local unitary equivalent (LUE), denoted as $|\psi\rangle \sim_{\text{LUE}} |\psi'\rangle$, if there exists a local unitary transformation mapping one into the other:

$$|\psi\rangle \underset{\text{LUE}}{\sim} |\psi'\rangle \Leftrightarrow |\psi'\rangle = U_A \otimes U_B |\psi\rangle. \qquad (23)$$

This relation defines an equivalence class whose representative is denoted with $[|\psi\rangle]$.

Starting from the unitary group manifold $SU(N_A N_B)$, we now argue that the evaluation of binding state complexity reduces to a problem defined in the special orthogonal space $SO(N_A)$. We refer to this setting as the $SO(N_A)$ picture. We will show in sections 3.1 and 3.2 that the problem eventually depends only on the Schmidt coefficients. Then we summarize the main steps of the procedure in section 3.3 and concretely apply them to a system composed of two qubits in section 3.4.

### 3.1 Reduction to a minimal set of parameters

We consider reference and target states defined by the Schmidt decompositions

$$|\psi_R\rangle = \sum_{m=1}^{N_A} \lambda_m |\phi_m \chi_m\rangle, \qquad |\psi_T\rangle = \sum_{m=1}^{N_A} \bar{\lambda}_m |\bar{\phi}_m \bar{\chi}_m\rangle, \qquad (24)$$

and we study the binding state complexity $BC_F^{|\psi_T\rangle, |\psi_R\rangle}$ induced by any $F$-norm of the form (15) with penalty factors (16). A trivial consequence of the zero cost assigned to

<footnote n="10">¹⁰Since $N_A \le N_B$, generally the set $|\chi_m\rangle$ does not form a complete orthonormal basis for $\mathcal{H}_B$, but rather an incomplete one.</footnote>
local generators is that the complexity of moving between LUE states vanishes

$$|\psi_R\rangle \underset{\text{LUE}}{\sim} |\psi_T\rangle \Rightarrow \mathcal{B}C_F^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = 0, \qquad (25)$$

because we can connect the states using a path generated by local operators only, which have zero cost. We depict a generic trajectory in the space of states given by eq. (22) in fig. 1(a). In the figure, we represent classes of LUE states as horizontal red lines. Those are null directions, along which the distance vanishes. It is important to remark that this statement is not only valid infinitesimally, but even when the parameters ($\phi_m, \chi_m$) of the basis undergo a finite change, see fig. 1(b).

<figure><img src="image_2.png" /><figcaption>Figure 1: (a) Continuous path in the space of states $\mathcal{B}$ given by eq. (22). The vertical line represents different values of the Schmidt coefficients $\lambda_m$ which parametrize the equivalence classes of LUE states, while the horizontal axis corresponds to variations of the parameters ($\phi_m, \chi_m$) of the basis. (b) The red lines represent the set of LUE states defined in eq. (23). Those are regions where the distance between two states vanishes, either when the trajectory is infinitesimal (path 1) or when there is a finite variation of the parameters of the basis (path 2).</figcaption></figure>

Consider now the induced norm $F^{\text{state}}$ on the Hilbert space, introduced in eq. (9). Along a generic path in the space of states defined by the Schmidt decomposition (22), we schematically denote the dependence of this norm on the parameters of the basis as<fnref n="11" />

$$F^{\text{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle] = g(\lambda_m, \phi_m, \chi_m, \dot{\lambda}_m, \dot{\phi}_m, \dot{\chi}_m). \qquad (26)$$

Let us focus on any $F$-norm with penalty factors (16). Since we assign vanishing penalty factors to local generators, which control the change in ($\phi_m, \chi_m$), the norm (26) cannot depend on the direction of the tangent vector in that subspace:

$$F^{\text{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle] = g(\lambda_m, \phi_m, \chi_m, \dot{\lambda}_m). \qquad (27)$$

However, notice that the line element depends on the parameters ($\phi_m, \chi_m$), implying that there is no translation invariance along these directions. Pictorially, this corresponds to the fact that the paths 2A and 2B in fig. 2 have different lengths, while the infinitesimally close paths 1A, 1B have the same length.

Binding state complexity, as defined in eq. (21), requires minimizing the distance among all the possible trajectories connecting the reference and target states. In particular, the identity (27) implies that the parameters ($\phi_m, \chi_m$) in the Schmidt decomposition (22) become non-dynamical degrees of freedom in the geodesic minimization. For this reason, the norm we obtain over the space of states, as per (10), can be minimized locally in terms

<footnote n="11">¹¹In writing the expression in this form (where it only depends on the parameters of the states) we have already performed the minimization over the stabilizer in eq. (9).</footnote>
<figure><img src="image_3.png" /><figcaption>Figure 2: Label 1 refers to infinitesimal paths, while label 2 to trajectories with finite variations of the parameters ($\phi_m, \chi_m$) of the basis. Consider reference and target states with infinitesimally close Schmidt coefficients $\lambda_m$ and $\lambda_m + d\lambda_m$. The length of the infinitesimally close trajectories 1A, 1B is the same. Instead, the length of the paths 2A, 2B is in general different. Therefore the metric is not invariant under translations along ($\phi_m, \chi_m$).</figcaption></figure>

of the parameters of the basis $\phi_m(t)$ and $\chi_m(t)$ to give a norm just in terms of the Schmidt coefficients $\lambda_m(t)$ and their differentials.<fnref n="12" />

We shall denote this norm as $\mathcal{BF}$. More formally, for any $F$-norm with vanishing cost associated with local generators, we define

$$ \mathcal{BF}[\vec{\lambda}, \dot{\vec{\lambda}}] = \min_{|\phi_n, \chi_n\rangle} F^{\text{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle] = \min_{\substack{|\phi_n, \chi_n\rangle \\ \text{stab } |\psi(t)\rangle}} F[\vec{Y}], \quad (28) $$

where $\min_{|\phi_n, \chi_n\rangle}$ is a minimization over the non-dynamical degrees of freedom in the possible bases in the Schmidt decomposition of the state $|\psi(t)\rangle$. The combined minimization on the right-hand side amounts, essentially, to minimizing over all control parameters $\vec{Y}$ and all basis parameters with the constraint of a prescribed change $\dot{\vec{\lambda}}(t)$ of the Schmidt coefficients. The last term in the definition (28) describes the procedure we are going to follow in the next sections.

Using the definition (28), the complexity (21) (see also eq. (10)) can be expressed as

$$ \mathcal{BC}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{\substack{\vec{\lambda}(t) \\ \vec{\lambda}(0) = \vec{\lambda}_R \\ \vec{\lambda}(1) = \vec{\lambda}_T}} \int_0^1 dt \mathcal{BF}[\vec{\lambda}, \dot{\vec{\lambda}}], \quad (29) $$

where $\vec{\lambda}_R$ and $\vec{\lambda}_T$ are the vectors of Schmidt coefficients of the reference and target states. We see that binding state complexity can be reduced to finding geodesics in the space of Schmidt coefficients. The boundary conditions for the reference and target states are given in terms of the equivalence classes introduced in eq. (23). In order to describe the trajectory in the full Hilbert space that takes us from the exact reference state to the exact target state, we will perform a finite LU transformation (of zero cost) at the beginning and at the end of the trajectory selected by the parameters of the basis that minimize the norm above. We discuss this further below.

## 3.2 Relation to the orthogonal group

Given the definition of $\mathcal{BF}$ in eq. (28), we study the dependence of the evolution of the Schmidt coefficients on the parameters of the basis and on the Hamiltonian. First, we find

<footnote n="12">¹²This is similar to partially solving the equations of motion for certain non-dynamical fields in QFT, whose momenta do not appear in the action, to obtain an on-shell action for the remaining fields.</footnote>
from the conservation of the norm under unitary time evolution that

$$ \begin{aligned} & \langle \phi_m(t) \chi_m(t) | \partial_t (|\phi_n(t) \chi_n(t)\rangle) + \partial_t (\langle \phi_n(t) \chi_n(t) |) |\phi_m(t) \chi_m(t)\rangle = \\ & = \delta_{mn} (\langle \phi_m(t) | \partial_t \phi_n(t) \rangle + \langle \partial_t \phi_n(t) | \phi_m(t) \rangle + \langle \chi_m(t) | \partial_t \chi_n(t) \rangle + \langle \partial_t \chi_n(t) | \chi_m(t) \rangle) = (30) \\ & = \partial_t (\langle \phi_n(t) | \phi_n(t) \rangle + \langle \chi_n(t) | \chi_n(t) \rangle) = 0. \end{aligned} $$

In terms of the state defined in eq. (22), the previous identity implies

$$ 2\dot{\lambda}_n(t) = \langle \phi_n(t) \chi_n(t) | \partial_t \psi(t) \rangle + \langle \partial_t \psi(t) | \phi_n(t) \chi_n(t) \rangle. \quad (31) $$

Using the Schrödinger equation and the general decomposition (14) of the Hamiltonian, we find

$$ \begin{aligned} \dot{\lambda}_n(t) &= \operatorname{Im}(\langle \phi_n(t) \chi_n(t) | H(t) | \psi(t) \rangle) = \sum_m \lambda_m(t) \operatorname{Im}(\langle \phi_n(t) \chi_n(t) | H(t) | \phi_m(t) \chi_m(t) \rangle) = \\ &= \sum_m \lambda_m(t) \operatorname{Im}(\langle \phi_n(t) \chi_n(t) | Y^{ai}(t) T_a^A \otimes T_i^B | \phi_m(t) \chi_m(t) \rangle). \end{aligned} \quad (32) $$

Note that local generators do not contribute to the change of the Schmidt coefficients, i.e., it is necessary to use operations involving both sides of the system (that create non-trivial entanglement) to modify them.<fnref n="13" /> In summary, BF is computed by minimizing $F[Y]$ over all the parameters of the basis and over all the degrees of freedom of $H$ that respect eq. (32).

We now proceed to show that the problem reduces from the group manifold $SU(N_A N_B)$ to $SO(N_A)$. The previous manipulations naturally related the Schmidt coefficients to a matrix $R$ such that

$$ \dot{\lambda}_m = \sum_m R_{mn} \lambda_n, \quad (33a) $$

$$ R_{mn}(t) = \operatorname{Im}(\langle \phi_m(t) \chi_m(t) | Y^{ai}(t) T_a^A \otimes T_i^B | \phi_n(t) \chi_n(t) \rangle). \quad (33b) $$

Since $R_{mn}(t)$ is real and antisymmetric, it is a generator of rotations of the special orthogonal group $SO(N_A)$ on the sphere of Schmidt coefficients. In other words, the constraint (33b) identifies a set of $\dim(SO(N_A)) = \frac{1}{2}N_A(N_A - 1)$ degrees of freedom, which are relevant for the evolution of the Schmidt coefficients in eq. (33a). Note however that this is less than the number of independent degrees of freedom in $Y_{ab}(t)$, which is $(N_A^2 - 1)(N_B^2 - 1)$, as long as $N_A \ge 2$. This means that multiple different values of the velocity vectors are associated with the same matrix $R_{mn}$.

As we mentioned earlier, around eq. (28), our final goal is to obtain a norm in terms of the Schmidt parameters. We can simplify the minimizations in eq. (28) by taking an intermediate step in which we minimize the norm $F$ over the velocities $\vec{Y}(t)$ and basis parameters $\phi_m(t), \chi_m(t)$ subject to the constraint that these parameters yield a fixed rotation generator $R_{mn}(t)$ at any point along the trajectory according to eq. (33b). In this way, we obtain a norm $BF(R)$ in terms of the parameters of the matrix $R$. By having to deal with $BF(R)$, instead of a norm on the Hamiltonian of $SU(N_A N_B)$, we say that we have reduced the problem to $SO(N_A)$.<fnref n="14" />

<footnote n="13">¹³This statement does not rely on the choice of penalty factors in the cost function, but it is simply a consequence of the time evolution of the Schmidt decomposition.</footnote>

<footnote n="14">¹⁴By referring to the reduction of the problem to $SO(N_A)$, we mean that the matrix $R_{mn}$ defined in eq. (32) and all the relevant degrees of freedom belong to are anti-symmetric and therefore belong to the algebra of $SO(N_A)$. However, it is worth emphasizing that $R_{mn}$ is defined through the velocities and the generators of the full $SU(N_A N_B)$ group. The precise problem that we are solving will be summarized in section 3.3 and the reduction to $SO(N_A)$ that we refer to is completed after step 3.</footnote>
We then perform an additional minimization over the set of matrices R which generate a given change $\vec{\lambda}$ in the Schmidt coefficients $\vec{\lambda}$ according to eq. (33a). This step amounts to minimizing over the stabilizer group of $\vec{\lambda}$ within SO($N_A$).

To perform the minimization of the basis parameters as in (28), we can re-express $R_{mn}$ in terms of a specific Hilbert space basis $|m\rangle$, $|\tilde{m}\rangle$ of our choice using $|\phi_m\chi_m\rangle = U_A \otimes U_B |m\tilde{m}\rangle$, where $U_A, U_B$ are unitary operators acting in the respective subsystems. This implies

$$R_{mn} = \mathrm{Im}(\langle m\tilde{m} | Y^{ai}(t)(\mathrm{Ad}_{UA})_a^b(\mathrm{Ad}_{UB})_i^j T_b^A \otimes T_j^B | n\tilde{n} \rangle), \quad (\mathrm{Ad}_U)_J^I T_I = U^\dagger T_J U. \quad (34)$$

The adjoint operators implement the LU transformations and can be reabsorbed into the velocities. In this way, the generator of rotations (33b) becomes

$$\tilde{Y}^{bj} = Y^{ai}(\mathrm{Ad}_{UA})_a^b(\mathrm{Ad}_{UB})_i^j, \quad R_{mn} = \mathrm{Im}(\langle m\tilde{m} | \tilde{Y}^{ai}(t)T_a^A \otimes T_i^B | n\tilde{n} \rangle). \quad (35)$$

In terms of the parameters in eq. (35), the cost function (18) becomes

$$F_{p,\mathrm{hom}}[\tilde{Y}] = \left( \sum_{b,j} |\tilde{Y}^{ai}(\mathrm{Ad}_{UA})_a^b(\mathrm{Ad}_{UB})_i^j|^p \right)^{\frac{1}{p}}. \quad (36)$$

Minimizing the above equation over all adjoint operators and all the $\tilde{Y}^{ai}$ that respect the constraint (35) will give a norm in terms of the parameters of R, that is:

$$\mathcal{BF}_{p,\mathrm{hom}}(R_{mn}) = \min_{\tilde{Y}, \mathrm{Ad}_U} \sum_{b,j} \left( |\tilde{Y}^{ai}(\mathrm{Ad}_{UA})_a^b(\mathrm{Ad}_{UB})_i^j|^p \right)^{\frac{1}{p}}. \quad (37)$$

The final result $\mathcal{BF}_{p,\mathrm{hom}}(\vec{\lambda}, \tilde{\lambda})$ is then obtained by minimizing (37) over all matrices R satisfying (33a).

We conclude this section by stating an important matrix inequality which will be useful for the evaluation of binding state complexity for the $F_{1,\mathrm{hom}}$ norm:<fnref n="15" />

$$\min_{\mathrm{Ad}_U} \sum_{b,j} |\tilde{Y}^{ai}(\mathrm{Ad}_{UA})_a^b(\mathrm{Ad}_{UB})_i^j| \ge \sum_a \sigma_a(\tilde{Y}) \ge \sum_a |\tilde{Y}^{aa}|, \quad (38)$$

where $\sum_a \sigma_a(\tilde{Y})$ is the nuclear norm, which sums the singular values of the matrix $\tilde{Y}$. The technical proof of this result is reported in Appendix A. When further minimized over all $\tilde{Y}$ with the appropriate boundary conditions (35), the first term on the left-hand side becomes the $\mathcal{BF}_{1,\mathrm{hom}}(R_{mn})$ norm in terms of the matrix R and the inequality (38) provides the following lower bound

$$\mathcal{BF}_{1,\mathrm{hom}}(R_{mn}) \ge \min_{\tilde{Y}: R_{mn} \text{ is fixed}} \sum_a |\tilde{Y}^{aa}|. \quad (39)$$

We will find in section 4.1 that this inequality is actually saturated, at least when considering the Gell-Mann basis for the generators of the unitary group, see eq. (55). In this way, we will be able to derive exact results for the binding state complexity with the $F_{1,\mathrm{hom}}$ norm.

<footnote n="15">¹⁵The index a goes from 1 to $N_A^2 - 1$, even if $\tilde{Y}$ is rectangular, see appendix A. Since the nuclear norm is invariant under rotations, $\tilde{Y}$ can be replaced in the rightmost-hand side by any $O_A \tilde{Y} O_B$, where O are orthogonal matrices, yielding other equally valid inequalities.</footnote>
## 3.3 Summary of the procedure

Based on the previous analysis, we outline the steps that will be followed in the remainder of the paper to compute binding state complexity for a bipartite system with Hilbert space (12):

1. Specify a basis for the generators of SU($N_A$), SU($N_B$) and define the reference and target states in the Hilbert space.

2. Impose that the rotation matrix $R_{mn}$ is given by eq. (35), and solve for as many velocities $\tilde{Y}^{ai}$ as the degrees of freedom encoded by this orthogonal matrix.<fnref n="16" />

3. Plug the solution inside the desired cost function of the form (15) with penalties (16) and minimize over the remaining velocities $\tilde{Y}^{ai}$ and the adjoint of the local operators $U^A, U^B$ to obtain $\mathcal{BF}(R_{mn})$.

4. Solve the evolution equation for the Schmidt coefficients (33a) and minimize over the remaining degrees of freedom inside the stabilizer of SO($N_A$) to obtain $\mathcal{BF}(\vec{\lambda}, \vec{\lambda})$.

5. Minimize $\int_0^1 \mathcal{BF}(\vec{\lambda}, \vec{\lambda}) dt$ with boundary conditions $\vec{\lambda}(0) = \vec{\lambda}_R, \vec{\lambda}(1) = \vec{\lambda}_T$ to find the complexity.

## 3.4 Example: two qubits

Before diving into the general analysis, we exemplify our method in a simple setting: a system composed of two qubits, split into two subsystems of one qubit each. In the following, we specify the step numbers according to the list outlined in section 3.3 above. We will treat in turn the cases of the $F_{1,\text{hom}}$ and $F_{2,\text{hom}}$ norms. Steps 1 and 2 are actually independent of the norm used, except for the fact that the norm $F_{2,\text{hom}}$ does not depend on the choice of basis of generators of the unitary group.

**Step 1.** The basis of generators for the SU(2) group acting on each qubit separately is given by the Pauli matrices, therefore the Hamiltonian of the full system reads

$$H(t) = Y_A^a(t) \sigma_a^A \otimes \mathbf{1}^B + Y_B^i(t) \sigma_i^B \otimes \mathbf{1}^A + Y^{ai}(t) \sigma_a^A \otimes \sigma_i^B. \quad (40)$$

From now until the end of the section, we will omit the superscripts referring to the subsystems {$A, B$}, since the generators for the two copies are the same. We denote the eigenstates of the Pauli matrix $\sigma_z$ by $\{|1\rangle \equiv |\uparrow\rangle, |2\rangle \equiv |\downarrow\rangle\}^{17}$ and define $\{|\bar{1}\rangle, |\bar{2}\rangle\} = \{|1\rangle, i|2\rangle\}$ which will serve as a convenient basis for the Hilbert space of the second qubit. We fix the reference state and allow for a generic Schmidt decomposition of the target state parametrized by a real angular coordinate $\theta$

$$|\psi_R\rangle = |11\rangle, \quad |\psi_T\rangle = \cos\theta |\bar{\phi}_1\bar{\chi}_1\rangle + \sin\theta |\bar{\phi}_2\bar{\chi}_2\rangle, \quad (41)$$

where $\bar{\phi}_m, \bar{\chi}_m$ (with $m = 1, 2$) denote a fixed basis determined by the Schmidt decomposition, which is related to the computational basis $|m\tilde{m}\rangle$ by LU transformations. These

<footnote n="16">¹⁶There always exists a solution to this equation if we allow a large enough set of non-local generators in our control Hamiltonian. However, the procedure may fail when the set of generators is not large enough. This will be the case for the universal-set of gates, which we discuss in detail in section 4.4. There, we will use a different strategy than the one described here.</footnote>

<footnote n="17">¹⁷The eigenstates of the matrix $\sigma_z$ are more commonly denoted by $\{|0\rangle, |1\rangle\}$, but here we have used indices starting from 1 to maintain consistency with the range of indices in eq. (24).</footnote>
states fix the boundary conditions for any trajectory. At each instant in the time evolution, the state can be expressed as

$$|\psi(t)\rangle = \cos\theta(t) |\phi_1(t)\chi_1(t)\rangle + \sin\theta(t) |\phi_2(t)\chi_2(t)\rangle, \quad (42)$$

identifying the Schmidt coefficients $\lambda_1(t) = \cos\theta(t)$ and $\lambda_2(t) = \sin\theta(t)$ in terms of a single parameter.

**Step 2.** As explained in section 3.2, we identify the relation between $R$, the generator of rotations of the Schmidt coefficients, and the Hamiltonian. This determines the degrees of freedom in the Hamiltonian which are relevant for the problem. The SO(2) rotation generator has only one degree of freedom that we parametrize as $R = v(t)i\sigma_y$ (here $\sigma_y$ is understood to act on the vector of Schmidt coefficients). The constraint on the control parameters of the Hamiltonian in terms of $R$ (35) reads

$$ \begin{aligned} R_{mn} = v(t)i(\sigma_y)_{mn} &= \operatorname{Im} (\langle \phi_m(t)\chi_m(t) | Y^{ai}(t)\sigma_a \otimes \sigma_i | \phi_n(t)\chi_n(t) \rangle) = \\ &= \operatorname{Im} (\langle m\tilde{m} | \tilde{Y}^{ai}(t)\sigma_a \otimes \sigma_i | n\tilde{n} \rangle) = \\ &= (\tilde{Y}_{11}(t) - \tilde{Y}_{22}(t))(\delta_{m1}\delta_{n2} - \delta_{m2}\delta_{n1}) = \\ &= (\tilde{Y}_{11}(t) - \tilde{Y}_{22}(t))i(\sigma_y)_{mn}. \end{aligned} \quad (43) $$

In the first line, we used the definition in eq. (33b). In the second line, we redefined the matrix of velocities using an LU transformation on the basis elements according to eq. (35), and we computed in the third line the non-vanishing elements. Finally, in the last line, we recognized the appearance of the Pauli generator $\sigma_y$.

Binding state complexity requires to minimize the $F$-norm under the constraint

$$\tilde{Y}_{11}(t) - \tilde{Y}_{22}(t) = v(t). \quad (44)$$

### $\mathcal{B}\mathcal{F}_{1,\text{hom}}$ norm

**Step 3.** Using the inequality (39), we find

$$ \mathcal{B}\mathcal{F}_{1,\text{hom}}(R_{mn}) \ge \min_{\tilde{Y}} (|\tilde{Y}_{11}| + |\tilde{Y}_{22}| + |\tilde{Y}_{33}|) = \min_{\tilde{Y}_{11}, \tilde{Y}_{33}} (|\tilde{Y}_{11}| + |-v(t)| + \tilde{Y}_{11}| + |\tilde{Y}_{33}|) = |v(t)|. \quad (45) $$

The lower bound, together with the constraint (43), can be attained, for instance, by the following Hamiltonian

$$H(t) = -v(t) \sigma_x \otimes \sigma_x, \quad (46)$$

and by choosing $U^A$ and $U^B$ to be the identity in eq. (35). We will, of course, require an additional LU transformation at $t = 1$ to reach the exact target state. This shows that the minimization problem reduced from the group manifold SU(2) to SO(2), with $\mathcal{B}\mathcal{F}_{1,\text{hom}}(R_{mn}) = |v(t)|$.

**Steps 4 and 5.** Finally, we should minimize the norm over the directions of the stabilizer of the Schmidt vector $\vec{\lambda}(t)$ under the constraint (33a). Since in this case $R$ and $\vec{\lambda}(t)$ depend on a single parameter, the minimization is trivial and we find $|v(t)| = |\dot{\theta}(t)|$, where $\theta(t)$ was introduced in eq. (42). We then conclude that $\mathcal{B}\mathcal{F}_{1,\text{hom}} = |\dot{\theta}(t)|$. One can argue that the geodesic for this norm does not flip the sign of $\dot{\theta}(t)$, since backtracking in $\theta$ will only elongate the integrated trajectory. Therefore a direct integration gives

$$\mathcal{B}\mathcal{C}_{1,\text{hom}}[|\psi_T\rangle, |\psi_R\rangle] = |\theta|, \quad \text{with } |\theta| \le \frac{\pi}{4}, \quad (47)$$
where the restriction on the angle comes from the fact that we can always use a free local unitary transformation to rotate $\theta \to \theta + \pi/2$. We notice that the maximally entangled state, precisely achieved when $\theta = \pi/4$, is also maximally complex. Finally, we remark that related results were obtained in [92] by evaluating the minimal time to entangle two qubits under certain constraints.

**Trajectory in the Hilbert space.** Since the binding complexity (6) associated with the norm (20a) is invariant under time reparametrizations (see footnote 3), we can actually choose $\theta(t) = \theta t$, such that the velocity is simply a constant along the trajectory $v(t) = \theta$. The trajectory of the state under the Hamiltonian evolution will first reach an intermediate state in the Hilbert space

$$|\tilde{\psi}_T\rangle = \vec{P}e^{-i\int_0^1 dt'H(t')} |\psi_R\rangle = e^{-iHt} |\psi_R\rangle \Big|_{t=1}, \quad H = -\theta \sigma_x \otimes \sigma_x, \quad (48)$$

where $|\tilde{\psi}_T\rangle$ is related to the exact target state $|\psi_T\rangle$ in eq. (41) by a LU transformation. Explicitly, the trajectory of the state generated by this Hamiltonian is given by

$$|\psi(t)\rangle = \cos(\theta t) |11\rangle + i \sin(\theta t) |22\rangle, \quad |\psi(0)\rangle = |\psi_R\rangle, \quad |\psi(1)\rangle = |\tilde{\psi}_T\rangle. \quad (49)$$

To reach the exact target state, a final costless (instantaneous) LU transformation is needed at $t=1$ such that $U \otimes V |m\tilde{m}\rangle = |\bar{\phi}_m \bar{\chi}_m\rangle$.

We observe that the Hamiltonian which generates the trajectories with minimal length is not unique. For example, another choice would be $H = a \sigma_x \otimes \sigma_y + b \sigma_y \otimes \sigma_x$, where $a$ and $b$ are positive and $a+b=\theta$.18 The trajectories obtained with this class of Hamiltonians will differ from the path in eq. (49) by the basis used for the Schmidt decomposition and will require a different LU transformation at $t=1$ to reach the exact target state.

### BF$_{2,hom}$ norm

**Step 3.** To obtain the complexity according to the $F_{2,hom}$ norm, see eq. (20a), we use the following chain of identities:

$$\mathcal{BF}_{2,hom} = \min_{\tilde{Y}} \sqrt{\tilde{Y}_{ai}\tilde{Y}_{ai}} = \min_{\tilde{Y}_{11}} \sqrt{\tilde{Y}_{11}^2 + (\tilde{Y}_{11} - v(t))^2} = \frac{1}{\sqrt{2}} |v(t)|. \quad (50)$$

The first step exploits the invariance of the norm under orthogonal transformations $Y \to O^A \tilde{Y} O^B$. The second step uses the constraint (44) and involves a minimization over all the velocities except $\tilde{Y}_{11}$, which is minimized in the last step. The minimum value is obtained when $\tilde{Y}_{11} = -\tilde{Y}_{22} = v(t)/2$. Choosing the reference and target states in eq. (41), we find that a unitary connecting them with cost (50) is given by

$$|\tilde{\psi}_T\rangle = e^{-iHt} |\psi_R\rangle \Big|_{t=1}, \quad H = -\frac{\theta}{2} (\sigma_x \otimes \sigma_x - \sigma_y \otimes \sigma_y). \quad (51)$$

where again a final LU transformation is needed to go from $|\tilde{\psi}_T\rangle$ to $|\psi_T\rangle$.

**Steps 4 and 5.** Steps 4 and 5 are trivial like in the previous case, because the problem has only one degree of freedom. We can again identify $v(t) = -\theta$ and compute the binding complexity

$$\mathcal{BC}_{2,hom}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \int_0^1 dt \frac{1}{\sqrt{2}} |\theta| = \frac{1}{\sqrt{2}} |\theta|. \quad (52)$$

This result is smaller by a factor of $\sqrt{2}$ with respect to the complexity associated with the $F_{1,hom}$ norm.

<footnote n="18">18 At the price of being pedantic, let us comment that this Hamiltonian defines the values of $Y^{ai}(t)$, not those of $\tilde{Y}^{ai}(t)$ in eq. (43).</footnote>
# 4 Results for arbitrary qudits

We apply the strategy outlined in section 3.3 to study the general case of two subsystems described by the unitary groups SU($N_A$) and SU($N_B$). One can interpret this setting as composed of two qudits.<fnref n="19" /> Alternatively, when $N_{A(B)} = 2^{n_{A(B)}}$, we can interpret the system as a spin chain split into two parts. We will be following the general discussion in section 3.1 and specifying the step numbers according to section 3.3. The results for the $F_{1,hom}$ norm depend on the choice of basis of generators used in decomposing the Hamiltonian (1). Therefore we conduct our analysis in the context of specific examples of bases. We first study the Gell-Mann basis of SU($N$) where we have better analytic control in section 4.1. We discuss the Pauli basis in section 4.2 where, in several cases, we can present bounds for the $F_{1,hom}$ complexity. Our results for the $F_{2,hom}$ norm, discussed in section 4.3, are (up to an overall constant) independent of the basis of generators and are therefore valid for any basis. In section 4.4, we demonstrate that our results for the binding complexity using the $F_1$ norm are, to some extent, insensitive to the assignment of penalty factors to the non-local generators.

**Step 1.** The minimal set of parameters required to study binding state complexity is given by the Schmidt coefficients and their differentials. We can then specify the reference and target states by the set of Schmidt coefficients

$$ \vec{\lambda}_{\text{ref}} = (1, 0, \dots, 0), \quad \vec{\lambda}_{\text{tar}} = (\bar{\lambda}_1, \bar{\lambda}_2, \dots, \bar{\lambda}_{N_A}), \quad \bar{\lambda}_1 \ge \bar{\lambda}_2 \ge \dots \ge \bar{\lambda}_{N_A}. \quad (53) $$

We will always take the reference state to be factorized, $|\psi_{\text{ref}}\rangle = |\psi_A\rangle \otimes |\psi_B\rangle$, so that it has only one non-vanishing Schmidt coefficient equal to 1. For the target state, we can always assume (modulo a local unitary operation) that the coefficients are positive real numbers and ordered from largest to smallest as indicated in (53).

**Choice of basis for the generators.** Consider the bases $\{T_a^A\}$ and $\{T_i^B\}$ for the algebras $\mathfrak{su}(N_A)$ and $\mathfrak{su}(N_B)$, containing $N_A^2 - 1$ and $N_B^2 - 1$ generators, respectively. The basis for the unitary algebra $\mathfrak{su}(N_A N_B)$ of the full system is built by taking the set of $(N_A N_B)^2 - 1$ generators

$$ T_a^A \otimes \mathbf{1}^B, \quad \mathbf{1}^A \otimes T_i^B, \quad T_a^A \otimes T_i^B. \quad (54) $$

In the following, we will consider two possibilities for the basis $T_{a(i)}^{A(B)}$ of each subsystem: (1) the basis of generalized Gell-Mann matrices or (2) tensor products of Pauli matrices. These two possibilities coincide for the case of SU(2) discussed above when $N_A = N_B = 2$. Let us start with the generalized Gell-Mann basis.

## 4.1 $\mathcal{BF}_{1,hom}$ norm with the Gell-Mann basis

The so-called Gell-Mann basis of generators of the special unitary group SU($N$) is a natural generalization of the Pauli matrices to describe qudits, the information units of higher dimensional quantum computing [93]. The basis is defined by

$$ T_{c(k,p)} \equiv \begin{cases} (T^x)_{k,p} \equiv |k\rangle \langle p| + |p\rangle \langle k| & \text{if } 1 \le k < p \le N \\ (T^y)_{k,p} \equiv i (|k\rangle \langle p| - |p\rangle \langle k|) & \text{if } 1 \le p < k \le N \\ (T^z)_{k,p} \equiv \sqrt{\frac{2}{k(k+1)}} (-k|k+1\rangle \langle k+1| + \sum_{i=1}^k |i\rangle \langle i|) & \text{if } 1 \le k = p < N, \end{cases} \quad (55) $$

<footnote n="19">¹⁹Qudits are invariant under the SU(d) group and model atomic systems with d excited states. Qutrits correspond to the case d = 3.</footnote>
where $|p\rangle$ denotes a fixed basis in the Hilbert space, and the three cases in the previous list generalize the Pauli matrices $\sigma_x, \sigma_y, \sigma_z$, respectively. The $N^2 - 1$ generators are labelled using the indices $k, p \in \{1, 2, \dots, N\}$, with the only constraint that we cannot have $k = p = N$. The generalized Gell-Mann matrices form an orthogonal basis normalized as $\mathrm{Tr}(T_{c(k,p)}T_{c(m,n)}) = 2\delta_{km}\delta_{pn}$. It is sometimes convenient to label the basis using a single index, defined by

$$c(k, p) = N(k - 1) + p. \qquad (56)$$

We build the basis for the algebra $su(N_A N_B)$ generating unitary transformations on the full system by combining the generators according to eq. (54). The bases of the Hilbert spaces of the systems A and B (12) will be denoted with $|p\rangle, |\tilde{p}\rangle$, and the associated Gell-Mann matrices will take the form (55) with or without tilde, respectively. Of course, as implied by this notation, the bases $|p\rangle, |\tilde{p}\rangle$, used here to define the Gell-Mann matrices will be the same fixed bases used in eq. (34)-(35) to define $\tilde{Y}$.

**Step 2.** We now set up the constraint equations in the $SO(N_A)$ picture (see section 3.2) starting from the rotation matrix $R_{mn}$ defined in eq. (35). Since $R_{mn}$ is antisymmetric, there is not any contribution of the diagonal generators (those with $m = n$), either in the subsystem A or B or both. Furthermore, since the expression (35) involves just the imaginary part, only non-local generators of the form $T^x \otimes T^y$ and $T^y \otimes T^x$ give rise to non-vanishing contributions.<fnref n="20" /> Assuming $m > n$, we get<fnref n="21" />

$$R_{mn} = \mathrm{Im} \left( \langle m \tilde{m} | \tilde{Y}_{c(a,b),c(i,j)} T_{c(a,b)}^A \otimes T_{c(i,j)}^B | n \tilde{n} \rangle \right) = \tilde{Y}_{c(m,n),c(n,m)} + \tilde{Y}_{c(n,m),c(m,n)}. \quad (57)$$

In this notation, the velocities $\tilde{Y}$ form a matrix with collective indices $c(a, b)$ and $c(i, j)$. Notice that $(T_{a,b}^y)^A$ form a subset of $\frac{1}{2}N_A(N_A - 1)$ Gell-Mann matrices that generate the orthogonal subgroup $SO(N_A)$ (they are angular momentum operators in Cartesian coordinates). Therefore, we can write

$$R = -i \sum_{m>n=1}^{N_A} R_{mn} (T_{m,n}^y)^A. \qquad (58)$$

We now perform for convenience an orthogonal transformation $O^B \in SO(N_B^2 - 1)$ to define a new matrix of velocities $Y' = \tilde{Y}O^B$, such that the constraints (57) are written only in terms of diagonal elements of $Y'$ as<fnref n="22" />

$$Y'_{c(m,n),c(m,n)} + Y'_{c(n,m),c(n,m)} = R_{mn}, \quad (m > n). \qquad (59)$$

The constraints (59) are the starting point to study binding state complexity in a system composed of two qudits (or equivalently a spin chain when $N_A$ is a power of two) using the Gell-Mann basis.

<footnote n="20">²⁰Note that in the two qubits example in section 3.4, we had, unlike here, a relative phase between the bases of the two systems, and therefore the generators that came into play were $\sigma_x \otimes \sigma_x$ and $\sigma_y \otimes \sigma_y$.</footnote>

<footnote n="21">²¹Assuming instead $n > m$ will give the same result up to an overall minus sign, as it should from the fact that $R_{mn}$ is anti-symmetric</footnote>

<footnote n="22">²²Explicitly, the relevant matrix is $O_{c(i,j),c(k,l)}^B = \delta_{i,l}\delta_{j,k}$, where here the indices $i, j, k, l$ run over the range $1, \dots, N_B$.</footnote>
**Step 3.** Using the inequalities (38)-(39) and the observation that the nuclear norm is invariant under an orthogonal transformation (applied separately on the left or right), i.e., $\sum_a \sigma_a(\tilde{Y}) = \sum_a \sigma_a(Y')$, we find the following bound

$$ \mathcal{BF}_{1,\text{hom}}(R_{mn}) \geq \min_{Y': R_{mn} \text{ is fixed}} \sum_a \sigma_a(Y') \geq \min_{Y': R_{mn} \text{ is fixed}} \sum_{m,n} |Y'_{c(m,n),c(m,n)}| \\ \geq \min_{Y': R_{mn} \text{ is fixed}} \sum_{m>n} |Y'_{c(m,n),c(m,n)}| + |Y'_{c(n,m),c(n,m)}| = \sum_{m>n} |R_{mn}|, \qquad (60) $$

where the velocities $Y'$ were defined below eq. (58) and the minimization is subject to the constraint (59). In going from the first to the second line, we used a trivial bound based on excluding $m=n$ terms from the summation. In the last step, we performed the explicit minimization over $Y'$ using the constraints (59). One can check that the lower bound, together with the constraint (57), can be achieved by using, for instance, the optimal Hamiltonian$^{23}$

$$ H = \sum_{m>n} R_{mn} T_{c(m,n)}^{A} \otimes T_{c(n,m)}^{B}, \qquad (61) $$

and by choosing $U^A$ and $U^B$ to be the identity in eq. (35). We will, of course, require an additional LU transformation at the end of the trajectory to reach the exact target state. Therefore the minimization yields

$$ \mathcal{BF}_{1,\text{hom}}(R_{mn}) = \sum_{m>n} |R_{mn}|. \qquad (62) $$

We note that $R_{mn}$ rotates a point along the plane spanned by the $m$ and $n$ Cartesian axes. This means that the problem now has the interpretation of the minimal sum of (absolute values of) rotation angles needed to reach a given point starting from the pole on a unit sphere, where the rotations are restricted to planes spanned by all pairs of Cartesian axes. To illustrate this point further, for a unit sphere embedded in three dimensions, the equivalent quantity is the minimal sum of absolute values of rotation angles, where only rotations around the $x$, $y$ and $z$ axes are allowed.

**Step 4.** To perform the minimization over the stabilizer group and find the state complexity, we need to solve eq. (33a), which in this setting reads

$$ \sum_{m<n} R_{nm}\lambda_m - \sum_{m>n} R_{mn}\lambda_m = \dot{\lambda}_n. \qquad (63) $$

There are $N_A - 1$ independent equations (recall that the Schmidt coefficients live on the sphere) and a total of $\frac{1}{2}N_A(N_A - 1)$ independent velocities $R_{mn}$. Let us enumerate the different velocity components in terms of a single index $R_\alpha$ with $\alpha \in \{1, \dots, \frac{1}{2}N_A(N_A - 1)\}$. We can, without loss of generality, solve the equations above for the first $N_A - 1$ velocities and express $R_\alpha$ for $\alpha \in \{1, \dots, N_A - 1\}$ in terms of the velocities $R_\beta$ for $\beta \ge N_A$ (if the equations cannot be solved for the first $N_A - 1$ velocities, just enumerate them differently). Inserting the solution back into the $F_{1,\text{hom}}$ norm from eq. (62) gives

$$ \mathcal{BF}_{1,\text{hom}}(\vec{\lambda}, \dot{\vec{\lambda}}) = \min_{\text{stab}} \sum_{m>n} |R_{mn}| = \min_{R_\beta} \left( \sum_{\alpha=1}^{N_A-1} |R_\alpha(\lambda_m, \dot{\lambda}_m, R_\beta)| + \sum_{\beta=N_A}^{N_A(N_A-1)/2} |R_\beta| \right), \qquad (64) $$

<footnote>$^{23}$One can check that in the two-qubit case presented in section 3.4, the Hamiltonian (61) degenerates to $H = \theta\sigma_x \otimes \sigma_y$, which connects the reference state $|\psi_R\rangle = |11\rangle$ to the generic target state $|\tilde{\psi}_T\rangle = \cos\theta|11\rangle + \sin\theta|22\rangle$, where we have chosen the basis $|1\rangle = |1\rangle$ and $|2\rangle = |2\rangle$. In particular, this $H$ is optimal because it gives the same binding state complexity as evaluated in eq. (47). This is a Hamiltonian in the class described below eq. (49).</footnote>
where the minimization “stab” refers to step 4 in section 3.3, and allows us to express the norm in terms of the Schmidt coefficients. To minimize over a single $R_β$, we notice that the right-hand side of eq. (64) is a piece-wise linear positive function, so the minimum must be at the points where one of the arguments of the absolute values vanish, i.e.,

$$R_{\beta} = 0 \lor R_{\alpha}(\lambda_{m}, \dot{\lambda}_{m}, R_{\beta}) = 0. \qquad (65)$$

Reiterating this process for all the $R_β$-s amounts to setting any selection of $(N_A - 1)(N_A - 2)/2$ of the above absolute values to zero. Of course, we then have to take the absolute minimum among all these choices. In other words, we should solve eq. (63) with only $N_A - 1$ non-vanishing velocities, when the rest of the velocities are set to zero. Then we substitute the solution, in terms of $\vec{\lambda}$ and $\dot{\vec{\lambda}}$, into the sum over the absolute values of all velocities and we finally take the minimum over all choices of those $N_A - 1$ non-vanishing velocities out of the $N_A(N_A - 1)/2$ velocities.

Let us now show that the minimum of (62) is obtained by setting to zero all the velocities other than $R_{m1}$ in eq. (58) subject to the constraint (63) with given $\vec{\lambda}$ and $\dot{\vec{\lambda}}$. In other words, the following matrix minimizes the norm, at least for points on the sphere that are not too far from the reference state $\vec{\lambda}_{ref}$ in eq. (53):

$$R_{\text{opt}} = -i \sum_{m>1}^{N_A} w_m (T_{m,1}^y)^A. \qquad (66)$$

The requirement to solve (63) fixes

$$\begin{align}
w_m &\equiv R_{m1} = \frac{\dot{\lambda}_m}{\lambda_1}, \quad (m > 1), \tag{67} \\
R_{mn} &= 0, \quad (m > n > 1). \nonumber
\end{align}$$

Let us first motivate this ansatz. At the start of the trajectory, all $(T^y)_{m,n}$ with $m, n > 1$ are members of the stabilizer subgroup of rotations of the vector of Schmidt coefficients. Therefore, the use of these generators in $R$ will not influence the constraint (63), but they will increase the norm and so there is no point in using those generators in $R$. When we are slightly away from the start of the trajectory, these generators are still close to the stabilizer, and so are less efficient than $(T^y)_{m,1}$ in modifying the Schmidt coefficients. Intuitively, this suggests that the same velocities will also be set to zero for a region on the Schmidt sphere that is close enough to the reference state. We will prove that this is the optimal solution also away from the pole, as long as $\lambda_2 + \lambda_3 \le \lambda_1$ (assuming that the Schmidt coefficients are ordered in decreasing order).

Let us now prove that eq. (66) gives the optimal matrix $R$ for the minimization of the norm. We compare the norm computed using eq. (66) to any other choice of $R$ subject to the constraint (63). A generic option for $R$ can be expressed as

$$R_{\text{gen}} = -i \sum_{m>n=1}^{N_A} \nu_{m,n} (T_{m,n}^y)^A, \qquad (68)$$

where $\nu_{m,n}$ form another set of velocities. By comparing $R_{\text{opt}}\vec{\lambda} = R_{\text{gen}}\vec{\lambda}$ defined in eqs. (66) and (68), we can express the velocities $w_m$ in terms of $\nu$ as

$$w_m = \frac{1}{\lambda_1} \left( \sum_{n<m} \nu_{m,n} \lambda_n - \sum_{n>m} \nu_{n,m} \lambda_n \right) \quad m > 1. \qquad (69)$$
The difference between the $F_{1,\text{hom}}$ norm using the two different rotation matrices reads

$$ \sum_{m>n=1}^{N_A} |\nu_{m,n}| - \sum_{m>1}^{N_A} |w_m| \geq \sum_{m>n>1}^{N_A} |\nu_{m,n}| \left( 1 - \frac{\lambda_n + \lambda_m}{\lambda_1} \right), \quad (70) $$

where we used the identity (69) and applied the triangle inequality to the second sum on the left-hand side following this replacement. We see that for a region on the sphere in which $\lambda_n + \lambda_m \le \lambda_1$, for any $m > n > 1$, this expression is positive (since we take the Schmidt coefficients in a decreasing order, this region amounts to $\lambda_2 + \lambda_3 \le \lambda_1$). Therefore the matrix $R_{\text{opt}}$ proposed in eq. (66) indeed minimizes the norm for this region on the sphere. The $\mathcal{BF}_{1,\text{hom}}$ norm then takes the following form in terms of the Schmidt coefficients $\lambda_m$ and $\dot{\lambda}_m$:

$$ \mathcal{BF}_{1,\text{hom}} dt = \sum_{m>1}^{N_A} |w_m| dt = \sum_{m>1}^{N_A} \frac{|d\lambda_m|}{\lambda_1} = \sum_{m>1}^{N_A} \frac{|d\lambda_m|}{\sqrt{1 - \sum_{n>1}^{N_A} \lambda_n^2}}, \quad (71) $$

where $\lambda_1$ is the maximal Schmidt coefficient, and in the third equality we used the solution for the velocities $w_m$ given in eq. (67). This cost function defines an infinitesimal distance in the space of Schmidt coefficients, and the length of its geodesics computes the binding state complexity.

**Step 5.** Next, we compute the geodesics connecting the reference and target states (53) in this geometry. Let us first show that for the shortest path, $d\lambda_m \ge 0$ for $m > 1$ at all points along the path. We will prove this by contradiction. Let us assume that, for a trajectory parameterized by $\bar{\lambda}(t)$ with $\bar{\lambda}(0) = \bar{\lambda}_{\text{ref}}$ defined in eq. (53), there exists some range of time $(t_2, t_4)$ during which $\dot{\lambda}_m(t) < 0$, while $\dot{\lambda}_m > 0$ for $t < t_2$. Then there must exist two times $t_1$ and $t_3$ such that $t_1 < t_2 < t_3$, for which $\lambda_m(t_1) = \lambda_m(t_3) \le \lambda_m(t)$ for $t \in (t_1, t_3)$. The factor of $(1 - \sum_{n>1}^{N_A} \lambda_n^2)^{-1/2}$ entering the norm (71) is a monotonically increasing function of the Schmidt coefficients, thus there is a shorter trajectory such that $\lambda_m(t) = \lambda_m(t_1) = \lambda_m(t_3)$ for $t \in [t_1, t_3]$, contradicting the assumption.

Since the Schmidt coefficients are arranged in decreasing order, the possible candidates for the geodesics are trajectories in a region defined by $\lambda_1(t) \ge \lambda_2(t) \ge \dots \ge \lambda_{N_A}(t)$, and have, as shown above, $\dot{\lambda}_m \ge 0$ for all $m > 1$. Consider the following piecewise smooth trajectory $\lambda_m^g(t)$:

$$ \begin{aligned} (1, 0, \dots, 0) &\rightarrow \left( \sqrt{1 - (N_A - 1)\bar{\lambda}_{N_A}^2}, \bar{\lambda}_{N_A}, \dots, \bar{\lambda}_{N_A} \right) \\ &\rightarrow \left( \sqrt{1 - \bar{\lambda}_{N_A}^2 - (N_A - 2)\bar{\lambda}_{N_A - 1}^2}, \bar{\lambda}_{N_A - 1}, \dots, \bar{\lambda}_{N_A - 1}, \bar{\lambda}_{N_A} \right) \quad (72) \\ &\rightarrow \dots \rightarrow (\bar{\lambda}_1, \bar{\lambda}_2, \dots, \bar{\lambda}_{N_A}), \end{aligned} $$

where in the $p$-th step, the differentials satisfy $d\lambda_m = d\lambda_n$ for $m, n \in \{2, 3, \dots, N_A - p + 1\}$ and $d\lambda_m = 0$ for $m > N_A - p + 1$. $^{24,25}$ We will show that $\lambda_m^g(t)$ is the geodesic.

First of all, let us explain how the trajectory works. The first Schmidt coefficient is always fixed by the normalization condition $\sum_m \lambda_m^2 = 1$. In the first step, all the other non-vanishing coefficients are increased at the same pace until the lowest coefficient reaches the target value $\bar{\lambda}_{N_A}$. After that, the last Schmidt coefficient stops evolving and we continue

<footnote>$^{24}$Note that for states where certain Schmidt coefficients do not appear in the final state, these Schmidt coefficients will also be turned off along all the trajectory.</footnote>

<footnote>$^{25}$The specification of the relation between the different differentials, rather than the precise time dependence, is enough to find the complexity due to the reparametrization invariance of the integrated norm.</footnote>
to increase the other ones at an equal rate until they reach the next-to-minimal value $\bar{\lambda}_{N_A-1}$. Then also the second-to-last coefficient stops evolving, and we repeatedly apply this method until the target state is reached.

We now prove that the trajectory (72) is the geodesic by showing that the difference between the length of any other trajectory and $\bar{\lambda}^g(t)$ is positive. Let us consider another trajectory $\bar{\lambda}(t)$ which we parametrize using the time $t = 1 - \lambda_1^g = 1 - \lambda_1$.<fnref n="26" />

With this parametrization, $t \in [0, t_f \equiv 1 - \bar{\lambda}_1]$ and $\lambda_{m>2}(0) = 0$. Let us compute the difference in the length of the two trajectories, i.e., the cost of $\lambda_m(t)$ minus the cost of $\lambda_m^g(t)$. We use integration by parts to find

$$ 
\begin{aligned}
\int_0^{t_f} dt \left[ \mathcal{B}\mathcal{F}_{1,\text{hom}}(\vec{\lambda}, \dot{\vec{\lambda}}) - \mathcal{B}\mathcal{F}_{1,\text{hom}}(\bar{\lambda}^g, \dot{\bar{\lambda}}^g) \right] &= \int_0^{t_f} dt \left( \sum_{m=2}^{N_A} \frac{\dot{\lambda}_m}{1-t} - \sum_{m=2}^{N_A} \frac{\dot{\lambda}_m^g}{1-t} \right) \\
&= \int_0^{t_f} dt \frac{1}{(1-t)^2} \sum_{m=2}^{N_A} (\lambda_m^g(t) - \lambda_m(t)).
\end{aligned}
\quad (73) $$

Because of the parametrization that we chose, the trajectories satisfy

• the normalization condition

$$ \sum_{m=2} \lambda_m^g(t) \lambda_m^g(t) = \sum_{m=2} \lambda_m(t) \lambda_m(t) = 1 - (1-t)^2. \quad (74) $$

• the inequality

$$ 0 \le \lambda_m^g(t), \lambda_m(t) \le \bar{\lambda}_m, \quad (\text{for } m > 1), \quad (75) $$

where $\bar{\lambda}_m$ is the endpoint of the trajectory.

• decreasing order

$$ \lambda_m^g(t) \le \lambda_{m-1}^g(t), \quad \lambda_m(t) \le \lambda_{m-1}(t), \quad (76) $$

chosen for convenience and without loss of generality.

Under these constraints, $\bar{\lambda}^g(t)$ maximizes the sum $\sum_{m=2}^{N_A} \lambda_m^g(t)$ and so the difference between the lengths in eq. (73) is always positive. To see why, consider the sum, with $\lambda_2$ replaced by the sphere constraint in eq. (74),

$$ \sum_{m=2}^{N_A} \lambda_m = \sum_{m=3}^{N_A} \lambda_m + \sqrt{1 - (1-t)^2 - \left(\sum_{m=3}^{N_A} \lambda_m^2\right)}. \quad (77) $$

Taking the gradient with respect to $\lambda_{m \ge 3}$ gives

$$ \partial_{\lambda_m} \left( \sum_{m=2}^{N_A} \lambda_m \right) = 1 - \frac{\lambda_m}{\sqrt{1 - (1-t)^2 - \left( \sum_{m=3}^{N_A} \lambda_m^2 \right)}} = 1 - \frac{\lambda_m}{\lambda_2}. \quad (78) $$

Since $\lambda_m \le \lambda_2$, the gradient is non negative, and vanishes only when $\lambda_m = \lambda_2$. This shows that the maximum of the sum is when the Schmidt coefficients for $m \ge 3$ are at their extreme values $\lambda_m(t) = \min(\lambda_{m-1}(t), \bar{\lambda}_m)$, which is exactly how the geodesic $\lambda_m^g(t)$ behaves.

<footnote n="26">²⁶Note that since the start and endpoints of the trajectory are the same and $\lambda_1$ is monotonic, see argument above (72), we have the freedom to choose a parametrization such that $\lambda_1 = \lambda_1^g$.</footnote>
Computing the length of the geodesic using the norm (71), we obtain the binding state complexity

$$ 
\begin{aligned}
\mathcal{BC}_{1,\text{hom}}^{\text{state}} &= \int_0^{\bar{\lambda}_{N_A}} \frac{(N_A - 1)dx}{\sqrt{1 - (N_A - 1)x^2}} + \int_{\bar{\lambda}_{N_A}}^{\bar{\lambda}_{N_A-1}} \frac{(N_A - 2)dx}{\sqrt{1 - \bar{\lambda}_{N_A}^2 - (N_A - 2)x^2}} \\
&\quad + \int_{\bar{\lambda}_{N_A-1}}^{\bar{\lambda}_{N_A-2}} \frac{(N_A - 3)dx}{\sqrt{1 - \bar{\lambda}_{N_A}^2 - \bar{\lambda}_{N_A-1}^2 - (N_A - 3)x^2}} + \cdots = \\
&= \sum_{m=1}^{N_A-1} \sqrt{m} \left[ \arcsin \left( \frac{\bar{\lambda}_{m+1}\sqrt{m}}{\sqrt{\sum_{n=1}^{m+1} \bar{\lambda}_n^2}} \right) - \arcsin \left( \frac{\bar{\lambda}_{m+2}\sqrt{m}}{\sqrt{\sum_{n=1}^{m+1} \bar{\lambda}_n^2}} \right) \right],
\end{aligned}
\quad (79) $$

where $\bar{\lambda}_{N_A+1} \equiv 0$. In the case of the maximally entangled state, the binding complexity reduces to $\sqrt{N_A - 1} \arccos(1/\sqrt{N_A})$. For a generic state of a spin chain consisting of $n_A$ spins, $\mathcal{BC}_{1,\text{hom}}^{\text{state}}$ scales exponentially as $\sqrt{N_A} = 2^{n_A/2}$, as we have checked numerically by considering random sequences of Schmidt coefficients sampled from a uniform distribution and then normalized. Eq. (79) is our main result in the case of the $F_{1,\text{hom}}$ norm in the Gell-Mann basis.

A Hamiltonian that builds a trajectory realizing the minimal cost can be constructed as follows

$$ H = \sum_{m>1}^{N_A} w_m(t) T_{c(m,1)}^A \otimes T_{c(1,m)}^B, \quad w_m \equiv \frac{\dot{\lambda}_m^g(t)}{\lambda_1^g(t)}, \quad (80) $$

where we used eqs. (61) and (67) with $\bar{\lambda}^g(t)$ the trajectory of the geodesic. This Hamiltonian preserves the basis for the Schmidt decomposition along the entire time evolution, and the state along the trajectory reads

$$ |\psi(t)\rangle = \sum_{m=1}^{N_A} \lambda_m^g(t) |mm\rangle. \quad (81) $$

Finally, to move from the state $|\psi(t=1)\rangle$ to the exact target state, an additional costless LU transformation is needed.

We notice that the trajectory built with the optimal Hamiltonian (80) does not have large fluctuations along the null directions (i.e., the red lines containing LUE states in fig. 1(a)), because a finite LU transformation is only performed as the final step at $t=1$ to move from an LUE target state to the exact one. Therefore the results obtained in the limit of vanishing cost for the local generators represent a good approximation for a more realistic setting, where one has a small but non-vanishing cost for any local gate.<fnref n="27" />

The above results are valid for systems with a finite-dimensional Hilbert space such as qubits and qudits. Recently, the study of complexity in continuous variable systems has raised significant interest. These studies served to find definitions of complexity in quantum field theory as a first attempt towards making the holographic complexity conjectures more precise, see e.g., [55, 56, 94]. The study of binding complexity in quantum field theories is outside the scope of this manuscript, but as a first step, we can consider the complexity for a very large number of Schmidt coefficients. In that case, the complexity formula (79) can be cast in the form of a continuous integral. We present a closed formula in this limit in appendix B.

<footnote n="27">²⁷See the discussion below eq. (110) regarding the $F_{1,\text{set}}$ norm, for which there are large fluctuations in the null directions and approximating the optimal circuits with non-vanishing cost for local gates is more involved.</footnote>
## 4.2 $\mathcal{BF}_{1,\text{hom}}$ norm with the Pauli basis

The tools developed in section 3 can also be applied to the case in which the generators of the unitary group are the generalized Pauli matrices. However, it turns out that the implementation of the constraints (35) and the minimization to perform are technically harder, making it difficult to determine the geodesics in the space of states for the $F_{1,\text{hom}}$ norm. For these reasons, we leave an exhaustive analysis of the $\mathcal{BC}_{1,\text{hom}}^{\text{state}}$ with the Pauli basis for future work, and instead derive bounds on the complexity and exact results for special cases.

**Generalized Pauli basis.** The generalized Pauli matrices give a basis of generators of the unitary group $U(N)$ which is natural in the case of a spin chain $N = 2^n$. They are defined as

$$T_{\vec{\mu}} = \sigma_{\mu_1}^1 \otimes \sigma_{\mu_2}^2 \otimes \dots \otimes \sigma_{\mu_n}^n, \qquad (82)$$

where each factor of the tensor product is a covariant two-dimensional Pauli matrix $\sigma_\mu = (\sigma_0, \sigma_i)$, with $\sigma_0 = 1$. The generalized Pauli matrices are orthogonal and normalized as $\text{Tr}(T_{\vec{\mu}}T_{\vec{\nu}}) = 2^n \delta_{\vec{\mu}\vec{\nu}}$. Assuming that the chain is split into two parts as in eq. (12), the basis for the full system is given by

$$T_I = T_{\vec{\mu}}^A \otimes T_{\vec{\nu}}^B = \sigma_{\mu_1}^{A,1} \otimes \dots \sigma_{\mu_n}^{A,n_A} \otimes \sigma_{\nu_1}^{B,1} \otimes \dots \otimes \sigma_{\nu_n}^{B,n_B}, \qquad (83)$$

The superscript $(A, \alpha)$ refers to the $\alpha$-th qubit in subsystem A, while $(B, \beta)$ to the $\beta$-th qubit in subsystem B. The index $I$ runs over the Pauli matrices of the combined system, cf. footnote 5. As we are interested in the generators of $SU(N_A N_B)$, we will exclude the case where all the matrices in the product are the two-dimensional identity matrix. The generators are normalized as $\text{Tr}(T_I T_J) = 2^{n_A+n_B} \delta_{IJ}$. In the case of the Pauli basis, the constraints (35) will take a more complicated form than those for the Gell-Mann basis in eq. (57) and they will involve, generically, a larger number of velocities. As a consequence, we won't be able to derive an exact result for $\mathcal{BF}_{1,\text{hom}}$ except for the special case where the smaller subsystem consists of a single qubit, $N_A = 2$; in the general case we can only prove bounds. Therefore we will not be specifying the step structure in this section.

**Upper bounds.** Consider the binding state complexity for a system composed of $n$ qubits, divided into two subregions with $n_A$ and $n_B$ qubits, respectively. We choose reference and target states described by the Schmidt coefficients in eq. (53) (i.e., the reference state is a product state). We take the path built with the Hamiltonian (80) (which was a geodesic for the geometry on the space of states defined by the choice of Gell-Mann matrices as generators), but now we compute its length according to the $\mathcal{BF}_{1,\text{hom}}^{\text{state}}$ norm for the Pauli basis. In appendix C, we show that the cost is the same as with the Gell-Mann basis norm, and as such, the result in eq. (79) gives an upper bound. However, as we show next, this upper bound can be rather weak in certain circumstances, and we can improve it in specific cases.

Let us now consider the special case in which the two systems are of equal size $n_A = n_B$ and focus on target states which are LUE to a product of entangled qubit pairs across the two subregions, i.e., 28

$$U_A \otimes U_B |\psi_T\rangle = \otimes_i (\cos(\theta_i) |\uparrow_{A_i} \uparrow_{B_i}\rangle + \sin(\theta_i) |\downarrow_{A_i} \downarrow_{B_i}\rangle). \qquad (84)$$

<footnote n="28">28 The discussion can be extended to $n_A < n_B$ by padding eq. (84) with, for example, tensor products of spin up $\otimes |\uparrow\rangle$ for the Hilbert spaces of the remaining unentangled qubits of the system B.</footnote>
The cost of a minimal circuit acting on two qubits (evaluated according to the $F_{1,hom}$ norm) was given in eq. (47). Extrapolating this result to a circuit $U_P(t)$ composed by the tensor product of the geodesics for each pair gives a total cost of $\sum_i |\theta_i|$. Since the generators acting on each pair of qubits form a subalgebra of $\mathfrak{su}(N_A N_B) = \mathfrak{su}(N_A^2)$, and since these subalgebras commute with each other, it is reasonable to expect that this cost is the binding state complexity

$$BC_{1,hom}^{\text{state}} = \sum_{i} |\theta_i|, \qquad (85)$$

which would then be additive for factorized states.<fnref n="29" /> Let us compare this bound to the previous one from the Gell-Mann basis expression for positive angles $0 < \theta_i < \pi/4$, where eq. (79) can be applied. In a small angle expansion the two results coincide, while for large angles they differ significantly. This can be seen by the scaling with $n_A$ of the complexity bounds for the maximally entangled state, which is linear for the bound (85), but exponential for the bound following from the Gell-Mann basis complexity (in this case of course the above bound from additivity is much stronger). The bounds as a function of the angle for the case of 6 entangled qubit pairs with equal angles are plotted in figure 3. We can see that they coincide for small angles, but deviate significantly for larger angles.

<figure><img src="image_4.png" /><figcaption>Figure 3: Comparison of the different bounds for the Pauli basis binding complexity of the state (84) in the case of equal angles and 6 entangled qubit pairs, as a function of the angle $\theta$. Red: the bound coming from the Gell-Mann optimal trajectory. Blue: the bound from eq. (85) coming from additivity. Orange: one-half the entanglement entropy of the state. The allowed region for the binding complexity is between the blue and orange curves, according to the lower bound of (90) when $c = 2$.</figcaption></figure>

**Exact Result.** Let us consider yet another special case, in which the subsystem A consists of a single qubit, and B of an arbitrary number. In this case, we can provide an exact result for the $BC_{1,hom}^{\text{state}}$ complexity with the Pauli basis. Here, the target state is always LUE to a state in which there is non-vanishing entanglement only between A and a single qubit in B. We can prove that all the other qubits (that we will call *ancillary qubits*) in subregion B do not give any advantage in minimizing a trajectory for the evaluation of binding state complexity according to the $BF_{1,hom}$ norm with the Pauli basis. Since there is only one independent Schmidt coefficient (two in total), the equation of motion (32) reduces to

$$\dot{\lambda}(t) = R_{12} \sqrt{1 - \lambda(t)^2}, \qquad (86a)$$

$$R_{12} = \mathrm{Im} \left( \langle \phi_1(t) \chi_1(t) | Y^{\mu_1 \nu_1 \dots \nu_{n_B}}(t) \sigma_{\mu_1} \otimes \sigma_{\nu_1} \otimes \dots \otimes \sigma_{\nu_{n_B}} | \phi_2(t) \chi_2(t) \rangle \right), \qquad (86b)$$

<footnote n="29">²⁹See [5], section IIIB, for a discussion on the additivity of geodesics lengths for tensor product target states in the context of Finsler metrics.</footnote>
where $|\phi_m(t)\rangle$ and $|\chi_m(t)\rangle$ are eigenstates for the subsystems A and B, respectively. Eq. (86b) provides a linear relation that can be solved for one of the Y components in terms of the others (and of $R_{12}$). The norm $\sum |Y|$ can then be minimized (for a fixed value of $R_{12}$) using a similar procedure to the one explained around eq. (64). That is, since the norm is a piecewise linear function, it will be minimized at points where the absolute values change signs. Therefore, the minimum of the norm will occur when only one of the $Y^{\mu_1\nu_1...\nu_{n_B}}(t)$ does not vanish. The norm will then be

$$ \mathcal{BF}_{1,\text{hom}} = \min_{\substack{\mu_1\nu_1\dots\nu_{n_B} \\ U_A \otimes U_B}} \left| \frac{\dot{\lambda}(t)}{\sqrt{1-\lambda(t)^2 \operatorname{Im}(\langle 1\tilde{1} | U_A \otimes U_B \sigma_{\mu_1} \otimes \sigma_{\nu_1} \otimes \dots \otimes \sigma_{\nu_{n_B}} U_A^\dagger \otimes U_B^\dagger | 2\tilde{2} \rangle)}} \right| = \frac{|\dot{\lambda}(t)|}{\sqrt{1-\lambda(t)^2}}, \quad (87) $$

where the minimization is done over all possible choices of the single non-vanishing $Y^{\mu_1\nu_1...\nu_{n_B}}$, and over all the possible LU transformations with respect to an arbitrary fixed orthonormal basis. We reach the final expression by observing that $\operatorname{Im}(\langle 1\tilde{1} | U_A \otimes U_B \cdot \sigma_{\mu_1} \otimes \sigma_{\nu_1} \otimes \dots \otimes \sigma_{\nu_{n_B}} \cdot U_A^\dagger \otimes U_B^\dagger | 2\tilde{2} \rangle)$ is upper bounded by the largest eigenvalue of $\sigma_{\mu_1} \otimes \sigma_{\nu_1} \otimes \dots \otimes \sigma_{\nu_{n_B}}$, which equals 1. This shows that the norm is the same as in the two-qubit case in section 3.4, and that the presence of ancillary qubits did not improve the complexity.

**Lower bound.** We can use the maximal entanglement rate studied in [95–97] to give lower bounds for the binding state complexity. Let us consider a trajectory in the space of states generated by a Hamiltonian $H(t) = H_A(t) + H_B(t) + H_{AB}(t)$, where $H_A$ and $H_B$ are composed by local generators in the subsystems A and B, while $H_{AB}$ is the interaction term between them. In such a case, an upper bound for the rate of change of the entanglement entropy $S_A$ between the two systems is given by [97]

$$ \frac{dS_A}{dt} \le c \log d ||H_{AB}(t)||. \quad (88) $$

Here $||*||$ is the operator norm (i.e., the maximal singular value), c is an $\mathcal{O}(1)$ constant (its optimal value has been numerically argued to be 2 in reference [97]), d is the dimension of the smaller subsystem out of A and B on which H has support.<fnref n="30" /> For instance, $d = N_A$ if H contains interactions between any of the qubits, and $d = 2$ if it only contains an interaction between a single pair of qubits. Since the operator norm of the generators, both in the Pauli and Gell-Mann bases, is equal to 1, we can use the triangle inequality to get $||H_{AB}|| = ||\sum_{a,i} Y^{ai} T_a^A \otimes T_i^B|| \le ||\sum_{a,i} |Y^{ai}| ||T_a^A \otimes T_i^B|| = \sum_{a,i} |Y^{ai}| = \mathcal{BF}_{1,\text{hom}}$. Integrating over time and remembering that the reference state is factorized, so it has zero entanglement entropy, gives

$$ \frac{S_A}{c \log d} \le \mathcal{BC}_1^{\text{state}}[|\psi_R\rangle, |\psi_T\rangle]. \quad (89) $$

We note that this lower bound is valid for both choices of the basis for the unitary group and for any kind of interactions we allow in $H(t)$, as long as the value of d is changed accordingly. The lower bound can be far from the binding complexity. For example, with the Gell-Mann basis, the binding complexity for the maximally entangled state grows as $\sqrt{N_A}$, while the entanglement entropy is $\log N_A$ and therefore for large $N_A$ the bound is

<footnote n="30">³⁰If the Hilbert space of the systems A and B can be factorised as $a \otimes \bar{a}$ and $b \otimes \bar{b}$, and $H_{AB}$ has the form of $1_{\bar{a}} \otimes H_{ab} \otimes 1_{\bar{b}}$, then $d = \min(\dim(a), \dim(b))$.</footnote>
far from being saturated. On the other hand with the Pauli basis, the binding complexity for this state is upper bounded by $\pi/4 \log(N_A)$, where $\log(N_A)$ is simply the number of spins, as can be seen from eq. (85). For the Pauli basis, the value of $d$ changes depending on whether we consider the $BC_{1,hom}^{state}$ or $BC_{1,sel}^{state}$ norms. For $BC_{1,hom}^{state}$, $d = N_A$, while for $BC_{1,sel}^{state}$, $d = 2$ if the allowed interaction has support only on one spin on each side. In section 4.4, we will actually show that $BC_{1,hom}^{state} = BC_{1,sel}^{state}$ and therefore for the Pauli basis,

$$ \frac{S_A}{c} \leq BC_{1,Pauli}^{state} [ |\psi_T\rangle, |\psi_R\rangle ] . \qquad (90) $$

We draw this lower bound in figure 3 and highlight in light blue the region between the most restrictive upper and lower bounds.

## 4.3 $\mathcal{BF}_{2,hom}$ norm in any basis

In this section, we compute the binding state complexity for a system of two qudits using the cost function $F_{2,hom}$ in eq. (18). The $F_{2,hom}$ norm can be expressed as

$$ F_{2,hom} = \sqrt{\sum_{a,i} |Y^{ai}|^2} = \sqrt{\frac{\mathrm{Tr}(Q(H)^2)}{\mathcal{N}_A \mathcal{N}_B}}, \qquad (91) $$

where $Q(H) = \sum_{a,i} Y^{ai} T_a^A \otimes T_i^B$ is the non-local part of the Hamiltonian, while $\mathcal{N}_A = \mathrm{Tr}(T_a^A T_b^A)\delta_{ab}$ and $\mathcal{N}_B = \mathrm{Tr}(T_i^B T_j^B)\delta_{ij}$ are the normalization factors for the subsystems A and B, respectively. Note that the only dependence of the $F_{2,hom}$ norm on the basis of generators comes from an overall normalization, as long as the generators are orthogonal. Therefore we can use the convenient Gell-Mann basis to compute the norm, and then rescale the result by an appropriate normalization in order to get a general expression valid for any basis. The normalization of the Gell-Mann basis for arbitrary rank of the unitary group is $\mathcal{N}_A = \mathcal{N}_B = 2$. Denoting the velocities of the Hamiltonian expanded in the Gell-Mann basis as $Y_G^{ai}$ and performing a convenient rescaling in terms of the normalization of the basis of generators, we find

$$ F_{2,hom} = \sqrt{\frac{4}{\mathcal{N}_A \mathcal{N}_B} \sum_{a,i} |Y_G^{ai}|^2}. \qquad (92) $$

The reference and target states are specified by the Schmidt coefficients in eq. (53). We consider, without loss of generality, the factorized reference state $|\psi_R\rangle = |1\bar{1}\rangle$ defined in the basis for the Hilbert space chosen in eq. (55). The constraint equations to get to the $SO(N_A)$ picture were previously studied in eqs. (57)-(59). This concludes **steps 1–2** of our procedure outlined in section 3.3.

**Step 3.** Since this norm is invariant under orthogonal transformations applied to the $Y$-s, we can use the velocities $Y'$ introduced below eq. (58) to find

$$ \begin{aligned} \mathcal{BF}_{2,hom}(R_{mn}) &= \frac{2}{\sqrt{\mathcal{N}_A \mathcal{N}_B}} \min_{Y'} \sqrt{\sum_{m,n} (Y'_{c(m,n),c(m,n)})^2} = \\ &= \frac{2}{\sqrt{\mathcal{N}_A \mathcal{N}_B}} \sqrt{\sum_{m>n} \frac{R_{mn}^2}{2}} = \frac{1}{\sqrt{\mathcal{N}_A \mathcal{N}_B}} \sqrt{\mathrm{Tr}(R^\mathrm{T} R)}, \end{aligned} \qquad (93) $$

where in the second equality we have explicitly performed the minimization over the redundant $Y'$ components.
**Step 4.** We define the orthonormal vectors

$$ \hat{r}_1 \equiv \vec{\lambda}(t), \qquad \hat{r}_2 \equiv \frac{\dot{\vec{\lambda}}(t)}{|\dot{\vec{\lambda}}(t)|}. \qquad (94) $$

and complete them to an orthonormal basis of $\mathbb{R}^{N_A}$. The frame associated with the new basis is related to the previous one by an orthogonal matrix $O$. In the rotated frame the constraint $R\vec{\lambda} = \vec{\lambda}$ fixes the following components of the matrix $\tilde{R} = ORO^T$

$$ \tilde{R}_{m,1} = -\tilde{R}_{1,m} = \delta_{m,2} \sqrt{\dot{\vec{\lambda}}(t)^2}. \qquad (95) $$

while all the other components of $\tilde{R}$ belong to the stabilizer group over which we have to minimize the norm to get the state complexity. The minimum is achieved when the latter components are set to zero, giving

$$ \mathcal{B}\mathcal{F}_{2,\text{hom}}(\vec{\lambda}, \dot{\vec{\lambda}}) = \frac{1}{\sqrt{\mathcal{N}_A \mathcal{N}_B}} \min_{\text{stab}} \sqrt{\text{Tr}(R^\text{T} R)} = \frac{1}{\sqrt{\mathcal{N}_A \mathcal{N}_B}} \min_{\text{stab}} \sqrt{\text{Tr}(\tilde{R}^\text{T} \tilde{R})} = \\ = \frac{1}{\sqrt{\mathcal{N}_A \mathcal{N}_B}} \sqrt{\tilde{R}_{1,2} \tilde{R}_{1,2} + \tilde{R}_{2,1} \tilde{R}_{2,1}} = \sqrt{\frac{2}{\mathcal{N}_A \mathcal{N}_B}} \sqrt{\dot{\vec{\lambda}}(t)^2}, \qquad (96) $$

where the “stab” minimization refers to step 4 in section 3.3, and the second step in the first line is a consequence of the rotational invariance of the $F_{2,\text{hom}}$ norm. It is implicit that the vector $\vec{\lambda}$ is normalized, and so this defines the Euclidean metric on a unit sphere.

**Step 5.** The geodesics are the arcs and the binding state complexity is the arc length from the reference to the target state:

$$ \mathcal{B}\mathcal{C}_{2,\text{hom}}^{\text{state}} = \sqrt{\frac{2}{\mathcal{N}_A \mathcal{N}_B}} \arccos(\bar{\lambda}_1) = \sqrt{\frac{2}{\mathcal{N}_A \mathcal{N}_B}} \arccos\left(e^{-\frac{1}{2}S_{\infty}(|\psi_T\rangle\rangle)}\right). \qquad (97) $$

Here $\bar{\lambda}_1$ is the maximal Schmidt coefficient of the target state defined in eq. (53), while $S_\infty(|\psi_T\rangle\rangle)$ is the minimal Rényi entropy of the reduced density matrix of the target state $|\psi_T\rangle\rangle$. Eq. (97) is our main result for the case of the $F_2$ norm. We note that up to a constant factor due to the normalization, this result is the Fubini-Study distance between the target state and the closest product state, which can be associated with the geometric measure of entanglement [98]. Here, in the context of complexity, this distance naturally receives an operational meaning (see [82] for a related discussion). The geometry generated by the $\mathcal{B}\mathcal{F}_{2,\text{hom}}$ norm and the corresponding geodesic are depicted in fig. 4.

For completeness, let us include the normalization constants for the different bases. For the Gell-Mann basis, the result reads

$$ \mathcal{B}\mathcal{C}_{2,\text{hom}}^{\text{state}} = \frac{1}{\sqrt{2}} \arccos(\bar{\lambda}_1) = \frac{1}{\sqrt{2}} \arccos\left(e^{-\frac{1}{2}S_{\infty}(|\psi_T\rangle\rangle)}\right), \qquad (98) $$

while for the Pauli basis of a system with $n_A$ and $n_B$ spins the result is

$$ \mathcal{B}\mathcal{C}_{2,\text{hom}}^{\text{state}} = \frac{1}{\sqrt{2^{n_A+n_B-1}}} \arccos(\bar{\lambda}_1) = \frac{1}{\sqrt{2^{n_A+n_B-1}}} \arccos\left(e^{-\frac{1}{2}S_{\infty}(|\psi_T\rangle\rangle)}\right). \qquad (99) $$

The two results coincide when $n_A = n_B = 1$, as the two bases become the same.
<figure><img src="image_5.png" /><figcaption>Figure 4: Unit sphere spanned by the vector of Schmidt coefficients, with the vertical axis representing the first direction. The red trajectory is the geodesic connecting the reference and target states in eq. (53).</figcaption></figure>

The Schmidt coefficients follow a trajectory along the great circle connecting the Schmidt vectors of the reference and target states. A particular Hamiltonian that realizes this trajectory and has the cost $BC_{2,hom}^{state}$ is

$$H = -i \frac{\arccos(\bar{\lambda}_1)}{\sqrt{1 - \bar{\lambda}_1^2}} \left( \sum_{m=1}^{N_A} \bar{\lambda}_m (|m\tilde{m}\rangle\langle 1\tilde{1}| - |1\tilde{1}\rangle\langle m\tilde{m}|) \right). \quad (100)$$

With this Hamiltonian, the basis of the Schmidt decomposition will remain in the computational basis $|ii\rangle$ until $t=1$, and the state will be related to the exact target state by a final costless LU transformation $|\psi_T\rangle = e^{-iS_0} \sum_{m=1}^{N_A} \bar{\lambda}_m |m\tilde{m}\rangle = e^{-iS_0} e^{-iH} |1\tilde{1}\rangle$, where $S_0$ is a time independent generator of LU transformations (i.e., take eq. (14) with $Y^{ai} = 0$).

Since the geometry associated with binding complexity contains null directions determined by LUE states (see the discussion in section 3.1), there is a degeneracy in the unitaries that take the reference state to the target state with minimal cost. Specifically for the case of the $F_2$ norm, it is possible to find a continuous trajectory that does not require a final LU transformation. We determine such geodesic in appendix D.1 by using the Euler-Arnold equations, and we show an alternative way to compute its cost (97) in appendix D.2.

A possible extension of the $F_{2,hom}$ norm to a multipartite system associates vanishing cost to each local transformation, and equal cost to any non-local generator. Related notions were considered in [82, 98]. We can speculate, extrapolating from the above bipartite case, that the state complexity of such generalization will be given in terms of the Fubini-Study distance to the closest multipartite pure factorized state. Notice that the above-mentioned extension will differ from the multipartite notion of complexity that will be discussed in section 6 because there, we will only allow non-local generators to act on a maximum of two subsystems at the same time.

## 4.4 Insensitivity of $BF_1$ binding complexity to penalty factors of non-local gates

We have obtained two results that give the binding state complexity for the $F_{1,hom}$ and $F_{2,hom}$ norms as two different functions of the Schmidt coefficients, see eqs. (79) and (97). One may be tempted to conjecture that by a suitable assignment of penalty factors, it
would be possible to obtain an arbitrary function of the Schmidt coefficients.<fnref n="31" /> We will now show that this is not the case, at least for the $\mathcal{BC}_1$ complexity. On the contrary, the binding complexity computed according to the $\mathcal{BF}_1$ norm is rather insensitive to the precise choice of penalty factors for the non-local gates.<fnref n="32" /> More precisely, we will argue that $\mathcal{BC}_1$ does not depend on the choice of penalties assigned to the non-local generators, up to an overall normalization stemming from the value of the minimal penalty, as long as the generators forming the optimal trajectory all lie on the same adjoint orbit of the $SU(N_A) \otimes SU(N_B)$ action as the cheapest generator. This condition is satisfied for the Gell-Mann and Pauli bases. The logic behind this is that if LU transformations are free, one could use them together with the cheapest non-local gate to create any other non-local gate.

To prove this statement, let us consider the norm

$$F_{1,\text{int}}[\vec{Y}] = \lim_{\substack{\varepsilon \to 0 \\ q_{ai} > 1}} \left( \sum_a \varepsilon |Y_A^a| + \sum_i \varepsilon |Y_B^i| + \sum_{(a,i) \neq (1, \tilde{1})} q_{ai} |Y^{ai}| + |Y^{11}| \right), \quad (101)$$

which represents an intermediate case between the homogeneous norm (20a) and the norm for a universal set of gates (20b). All three norms assign a vanishing cost to local operations, but the penalties $q_{ai}$ associated with non-local generators are 1 for $F_{1,\text{hom}}$, $q_{ai} > 1$ finite for $F_{1,\text{int}}$ and $q_{ai} = \infty$ for $F_{1,\text{set}}$. Since complexity increases with the penalties, we conclude that the binding state complexities associated with the three cost functions satisfy

$$\mathcal{BC}_{1,\text{hom}}^{\text{state}} \le \mathcal{BC}_{1,\text{int}}^{\text{state}} \le \mathcal{BC}_{1,\text{set}}^{\text{state}}. \quad (102)$$

We are going to show that these bounds are saturated, i.e., that all these complexities are the same as long as the conditions specified at the beginning of this subsection hold.

To make this argument formal, let us first notice that we can express a generic trajectory in the $SU(N_A N_B)$ group manifold as a circuit of infinitesimal steps

$$U = \tilde{P} \exp \left( -i \int_0^1 dt' H(t') \right) = \lim_{M \to \infty} \prod_{m=1}^M \exp \left( -\frac{i}{M} H \left( \frac{m}{M} \right) \right), \quad (103)$$

where $H(t) = Y_A^a(t)T_a^A \otimes \mathbf{1}^B + Y_B^i(t)\mathbf{1}^A \otimes T_i^B + Y^{ai}(t)T_a^A \otimes T_i^B$. Let us now consider a different circuit that creates the same unitary, in which we break the gate $\exp(-i/M H(m/M))$ into $(N_A N_B)^2 - 1$ consecutive gates

$$U = \lim_{M \to \infty} \prod_{m=1}^{M} \prod_{a} e^{-\frac{i}{M} Y_{A}^{a}(\frac{m}{M}) T_{a}^{A} \otimes \mathbf{1}^{B}} \prod_{i} e^{-\frac{i}{M} Y_{B}^{i}(\frac{m}{M}) \mathbf{1}^{A} \otimes T_{i}^{B}} \prod_{a,i} e^{-\frac{i}{M} Y^{ai}(\frac{m}{M}) T_{a}^{A} \otimes T_{i}^{B}}, \quad (104)$$

where the product over $a$ and $i$ is from 1 to $N_A^2 - 1$ and $N_B^2 - 1$ respectively, and we used the Baker-Campbell-Hausdorff formula in the limit $M \to \infty$. Notice that the cost of the circuits (103) and (104) is the same according to the $F_{1,\text{hom}}$ norm, but is a priori different using the $F_{p>1,\text{hom}}$ norm.<fnref n="33" /> This is because the cost of the circuit (104) is $\sum_{m,a,i} |Y^{ai}(\frac{m}{M})|$ regardless of $p$, as the Hamiltonian is composed out of a single generator at any given time step.

We now assume that the circuit in (103) is the optimal trajectory obtained using the $F_{1,\text{hom}}$ norm. Let us construct yet another circuit that creates the same unitary and has

<footnote n="31">³¹We thank Yaron Oz for raising this question.</footnote>

<footnote n="32">³²The insensitivity to penalty factors does not apply to other $\mathcal{BC}_{p>1}$ complexities, as we explain below.</footnote>

<footnote n="33">³³For this reason, the arguments presented in this subsection only apply to the $\mathcal{BF}_1$ norm.</footnote>
the same cost. This circuit will only use a single non-local generator $T_1^A \otimes T_1^B$, thus its cost will be the same whether the norm is $F_{1,hom}$ or $F_{1,\text{set}}$. To begin with, let us assume that for every non-local generator $T_a^A \otimes T_i^B$, there is an LU transformation such that

$$U_a^{A\dagger} T_1^A U_a^A \otimes U_b^{B\dagger} T_1^B U_b^B = T_a^A \otimes T_i^B. \quad (105)$$

Notice that this statement is about the individual systems, i.e., there exists an LU transformation that takes a $T_1^{A/B}$ to any other generator $T_I^{A/B}$. The construction is now straightforward: take the circuit in (104) and replace any occurrence of a non-local generator with the preferred non-local generator $T_1^A \otimes T_1^{B34}$

$$U = \lim_{M \to \infty} \prod_{m=1}^{M} \prod_{a} e^{-\frac{i}{M} Y_A^a (\frac{m}{M}) T_a^A} \prod_{i} e^{-\frac{i}{M} Y_B^i (\frac{m}{M}) T_i^B} \prod_{a,i} U_a^{A\dagger} U_i^{B\dagger} e^{-\frac{i}{M} Y^{ai} (\frac{m}{M}) T_1^A T_1^B} U_a^A U_i^B, \quad (106)$$

where we have omitted the tensor product notation for convenience. Since the cost of the LU transformations is zero, the cost of this circuit will be the same regardless of whether we used the $F_{1,hom}$ or $F_{1,\text{set}}$ norm. Due to the inequalities (102), we conclude that for $p=1$ all these binding complexities are the same, in particular this is true for $\mathcal{BC}^{\text{state}}_{1,\text{int}}$, which assigns arbitrary penalties $q_{ai} > 1$ to the non-local generators.

We are left to discuss the validity of our assumption (105). For the Pauli basis, take, without loss of generality, the single allowed generator to be $\sigma_x^{A,1}\sigma_x^{B,1}$. We show below that we can create any generator out of the set of generalized Pauli basis with $\sigma_x^{A,1}\sigma_x^{B,1}$ and an LU transformation. Let us focus on a system of $n$ spins and ask whether there is a $U \in SU(2^n)$ such that $U^\dagger\sigma_x^1 U = \sigma_{\mu_1}^1 \sigma_{\mu_2}^2 \dots \sigma_{\mu_n}^{n_a}$, for any sequence of $\mu_i$s except for all zeros, where the superscript marks the site, and the subscript marks which Pauli matrix is being considered (cf. the text below eq. (82)). Since we have not made any assumption about the number of spins, this argument will be valid for both systems A and B. Let us define $\bar{\sigma} = \sigma_{\mu_2}^2 \dots \sigma_{\mu_n}^n$. If $\mu_1 = y$ or $z$, then

$$e^{-i\frac{\pi}{4}\sigma_z^1\bar{\sigma}}\sigma_x^1 e^{i\frac{\pi}{4}\sigma_z^1\bar{\sigma}} = \sigma_y^1\bar{\sigma}, \quad e^{i\frac{\pi}{4}\sigma_y^1\bar{\sigma}}\sigma_x^1 e^{-i\frac{\pi}{4}\sigma_y^1\bar{\sigma}} = \sigma_z^1\bar{\sigma}. \quad (107)$$

To reach $\sigma_x^1\bar{\sigma}$, we can first perform the transformation to $\sigma_y^1\bar{\sigma}$ and then complement it with an additional LU transformation $e^{i\frac{\pi}{4}\sigma_z^1\bar{\sigma}}e^{-i\frac{\pi}{4}\sigma_z^1\bar{\sigma}} = \sigma_x^1\bar{\sigma}$. To reach $\mathbf{1}^1\bar{\sigma}$, we can use the swap $S = \exp(-i\frac{\pi}{4}(\sigma_x^1\sigma_x^2 + \sigma_y^1\sigma_y^2 + \sigma_z^1\sigma_z^2))$, which gives $S^\dagger\sigma_x^1 S = \sigma_x^2$, and then recursively use similar transformations on the subset of $n-1$ remaining spins to reach $\bar{\sigma}$ from $\sigma_x^2$. This concludes the argument for the equivalence of $\mathcal{BC}^{\text{state}}_{1,hom} = \mathcal{BC}^{\text{state}}_{1,\text{set}}$ in the Pauli basis.

For the Gell-Mann basis, one cannot hope to create from one generator and an LU transformation any other generator, as the diagonal and off-diagonal generators have in general different eigenvalues. If the single non-local generator is a tensor product of two off-diagonal generators $T_{c(n,m)}^A \otimes T_{c(k,p)}^B$ where $n \neq m, k \neq p$, then all the other off-diagonal generators can be reached with an LU transformation. To demonstrate that, we note the following relations for $T_{c(2,1)}$, where $p > k > 2$:

$$e^{-i\frac{\pi}{2}T_{c(p,2)}}T_{c(2,1)}e^{i\frac{\pi}{2}T_{c(p,2)}} = T_{c(p,1)}, \quad e^{-i\frac{\pi}{2}T_{c(2,p)}}T_{c(2,1)}e^{i\frac{\pi}{2}T_{c(2,p)}} = T_{c(1,p)}, \quad (108a)$$

$$e^{i\frac{\pi}{2}T_{c(p,1)}}T_{c(2,1)}e^{-i\frac{\pi}{2}T_{c(p,1)}} = T_{c(p,2)}, \quad e^{i\frac{\pi}{2}T_{c(1,p)}}T_{c(2,1)}e^{-i\frac{\pi}{2}T_{c(1,p)}} = T_{c(2,p)}, \quad (108b)$$

$$e^{i\frac{\pi}{2}T_{c(2,p)}}e^{i\frac{\pi}{2}T_{c(1,k)}}T_{c(2,1)}e^{-i\frac{\pi}{2}T_{c(1,k)}}e^{-i\frac{\pi}{2}T_{c(2,p)}} = T_{c(p,k)}, \quad (108b)$$

$$e^{i\frac{\pi}{2}T_{c(2,p)}}e^{i\frac{\pi}{2}T_{c(k,1)}}T_{c(2,1)}e^{-i\frac{\pi}{2}T_{c(k,1)}}e^{-i\frac{\pi}{2}T_{c(2,p)}} = T_{c(k,p)}.$$

<footnote n="34">³⁴It would be sufficient to show this only for the generators forming the optimal trajectory, but since we do not know the optimal trajectory in the Pauli basis, showing it for all generators will be the strategy there.</footnote>
These relations show that any off-diagonal generator can be reached from $T_{c(2,1)}$ by an LU transformation. This is not special to $T_{c(2,1)}$, since from the above identities we learn that all the off-diagonal generators are on the same adjoint orbit.

Fortunately, the geodesic circuit with the Gell-Mann basis, which is generated by the Hamiltonian in eq. (80), only uses off-diagonal generators. Thus, the previous argument is still valid if we limit the generality of the allowed single non-local generator such that it is off-diagonal. For the geodesic obtained using the Gell-Mann basis, the associated circuit in eq. (104) will take the form

$$U = \lim_{M \to \infty} \prod_{m=1}^{M} \prod_{p>1}^{N_A} \exp \left( -\frac{i}{M} R_{p1} \left( \frac{m}{M} \right) T_{c(p,1)}^{A} \otimes T_{c(1,p)}^{B} \right). \quad (109)$$

For instance, if the allowed single non-local generator is $T_{c(2,1)}^A \otimes T_{c(2,1)}^B$, we can construct the circuit as

$$U = \lim_{M \to \infty} \prod_{m=1}^{M} \prod_{p>1}^{N_A} e^{-i \frac{\pi}{2} \left[ T_{c(p,2)}^{A} \otimes \mathbf{1}^{B} + \mathbf{1}^{A} \otimes T_{c(2,p)}^{B} \right]} e^{-\frac{i}{M} R_{p1} \left( \frac{m}{M} \right) T_{c(2,1)}^{A} \otimes T_{c(2,1)}^{B}} e^{i \frac{\pi}{2} \left[ T_{c(p,2)}^{A} \otimes \mathbf{1}^{B} + \mathbf{1}^{A} \otimes T_{c(2,p)}^{B} \right]} \quad (110)$$

We still need to discuss an additional subtlety. The derivation given above applies when the local generators have exactly zero cost. However, in any realistic experimental set-up, the penalty assigned to local generators is $\epsilon \neq 0$, in which case the circuit (110) would have divergent cost in the limit $M \to \infty$. However, it is possible to remedy this issue by considering the unitary (110), but with finite $M$. If the additional cost and the error coming from the Trotterization of a Hamiltonian [99, 100] can be kept small, this will describe an efficient approximation for the target unitary. To estimate the error and cost of such a circuit, consider a simple case in which the Hamiltonian (103) is constant. This is the form of the geodesic in the Gell-Mann basis if all but the largest Schmidt coefficients are equal. The Trotterization error can be estimated as $O(BC^2/M)$ (compare to eq. (16) in [100]), and the extra cost to the binding complexity will be $O(MN_A\epsilon)$. Depending on the engineer's budget, states with small enough $BC$ could be reached to a good approximation. Furthermore, as long as $\epsilon$ is small enough, the cost of the above circuit will provide a good estimate of the states' complexity.

# 5 Binding complexity of two qutrits and complexity of a single qubit

The simplest system on which one can ask questions of Nielsen's complexity is a single qubit. However, much less is known for the single-qubit complexity with an $F_1$ norm. It turns out that the results we obtained earlier in this paper for the binding complexity can be applied to this problem, too. More specifically, we will show that the binding complexity for a system of two qutrits in the Gell-Mann basis is related to the $F_1$ Nielsen's complexity for a single qubit.

In the case of two qutrits, we expand the SO(3) generator $R$ introduced in eq. (33b) in terms of the rotations around each axis

$$R(t) = v_1(t)L_1 + v_2(t)L_2 + v_3(t)L_3, \quad (111)$$

where $(L_i)_{jk} = \epsilon_{ijk}$ with $\epsilon$ the Levi-Civita symbol. The relation to the Gell-Mann generators introduced in eq. (55) is $L_1 = i(T^y)_{32}$, $L_2 = -i(T^y)_{31}$, $L_3 = i(T^y)_{21}$. According to the result (62), in the Gell-Mann basis, the minimization yields the norm

$$\mathcal{BF}_{1,\text{hom}}(R_{mn}) = |v_1(t)| + |v_2(t)| + |v_3(t)|. \quad (112)$$
When $\lambda_1 \ge \lambda_2 + \lambda_3$, as we have shown in eq. (71), the metric can be written in terms of the Schmidt coefficients as

$$ \mathcal{B}\mathcal{F}_{1,\text{hom}} dt = \frac{|d\lambda_2| + |d\lambda_3|}{|\lambda_1|}, \qquad (113) $$

and the optimal trajectory from $\bar{\lambda}(0) = (1, 0, 0)$ to $\bar{\lambda}(1) = (\bar{\lambda}_1, \bar{\lambda}_2, \bar{\lambda}_3)$, depicted in fig. 5, has the following cost determined according to eq. (79):

$$ \mathcal{B}\mathcal{C}_{1,\text{hom}}^{\text{state}} = \sqrt{2} \arccos \left( \sqrt{1 - 2\bar{\lambda}_3^2} \right) + \arccos \left( \frac{\bar{\lambda}_3\bar{\lambda}_2 + \bar{\lambda}_1\sqrt{1 - 2\bar{\lambda}_3^2}}{1 - \bar{\lambda}_3^2} \right). \qquad (114) $$

For the maximally entangled state, the complexity is $\sqrt{2} \arccos(1/\sqrt{3})$.

<figure><img src="image_6.png" /><figcaption>Figure 5: Piecewise trajectory (72) in the space of Schmidt coefficients for the case of two qutrits. Binding complexity is given by the sum of the lengths of the red and blue curves in this geometry.</figcaption></figure>

Since the so(3) and su(2) algebras are isomorphic, and since the nature of the complexity problem is completely algebraic, eq. (112) can also be interpreted as the unitary (not binding) complexity norm for a single qubit using the basis of Pauli matrices. This is because a trajectory in SU(2) can be written as

$$ U(t) = \tilde{P} \exp\left(-i \int_{0}^{t} dt' H(t')\right), \quad H(t') = v_1(t')\sigma_1 + v_2(t')\sigma_2 + v_3(t')\sigma_3, \qquad (115) $$

for which the norm $F_1 = \sum_i |\text{Tr}(\dot{U}dU\sigma_i)|/2$ is equal to (112). The density matrix of a single spin in a pure state can be written as

$$ \rho = |\phi\rangle\langle\phi| = \frac{1}{2}(1 + x_1\sigma_1 + x_2\sigma_2 + x_3\sigma_3), \quad \text{where } |\vec{x}| = 1. \qquad (116) $$

The vector $\vec{x}$ is real, while the condition $|\vec{x}| = 1$ defines the Bloch sphere and ensures that the state is pure. The six points $x_i = \pm 1$ correspond to the eigenstates of the Pauli matrices $\sigma_i$. To some extent, as we discuss below, the complexity problems can be mapped to each other by identifying the vector of Schmidt coefficients $\bar{\lambda}$ with the position on the Bloch sphere $\vec{x}$. Starting from the equations of motion $\dot{\rho} = -i[H, \rho]$, one can find the state complexity norm by applying similar methods to those discussed in section 4.1 (see eq. (63) with the map $R_{32} = 2v_1$, $R_{31} = -2v_2$, $R_{21} = 2v_3$). The minimization over the stabilizer group with respect to $\rho$ implies that either $v_1$, $v_2$ or $v_3$ vanish, giving

$$ F_1 dt = \frac{1}{2} \min \left( \frac{|dx_2| + |dx_3|}{|x_1|}, \frac{|dx_3| + |dx_1|}{|x_2|}, \frac{|dx_1| + |dx_2|}{|x_3|} \right). \qquad (117) $$
Certain points on the sphere are related by LU transformations, which are free in the framework of binding complexity. Therefore, any reference and target states can be mapped to the region $\lambda_1 \ge \lambda_2 \ge \lambda_3 \ge 0$. For the complexity of a single qubit, we do not have this freedom. Nevertheless, the norm in (117) is invariant under a sign flip of each coordinate, and under permutations of the coordinates. Therefore, we can map any reference state which is an eigenstate of one of the Pauli matrices to $(x_1 = 1, x_2 = 0, x_3 = 0)$, and any target state to another state that has $\bar{x}_2 \ge \bar{x}_3 \ge 0$. It is then possible to map the single qubit complexity to the binding complexity of two qutrits if the inequality $\bar{x}_1 \ge \bar{x}_2$ is also satisfied. Since we could only compute the binding complexity for two qutrits when $\bar{\lambda}_1 \ge \bar{\lambda}_2 + \bar{\lambda}_3$, then the result in eq. (114) corresponds to twice the complexity of a single qubit in the regime $\bar{x}_1 \ge \bar{x}_2 + \bar{x}_3$, provided that we identify $\bar{\lambda}_i = \bar{x}_i$, that is

$$C_{1,qubit}^{\text{state}} = \frac{1}{\sqrt{2}} \arccos \left( \sqrt{1 - 2\bar{x}_3^2} \right) + \frac{1}{2} \arccos \left( \frac{\bar{x}_3\bar{x}_2 + \bar{x}_1\sqrt{1 - 2\bar{x}_3^2}}{1 - \bar{x}_3^2} \right). \quad (118)$$

This expression is directly applicable when the reference state is an eigenstate of the Pauli matrix $\sigma_i$. More generally, when the reference state is an eigenstate of $\sigma_i$ for some given $i$, the state complexity of the single qubit is given by half of eq. (114), with $\bar{\lambda}_1 = |\bar{x}_i|$, $\bar{\lambda}_2 = \max(|\bar{x}_j|, |\bar{x}_k|)$ and $\bar{\lambda}_3 = \min(|\bar{x}_j|, |\bar{x}_k|)$, where $j, k$ label the other two directions in the Bloch sphere. This solution is valid as long as the inequality $|\bar{x}_i| \ge |\bar{x}_j| + |\bar{x}_k|$, defining two separate connected regions on the sphere (around the positive and negative $i$-th axis), is satisfied, and as long as the reference and target states both belong to the same region out of the two.

It is interesting to notice that the minimization over the stabilizer group, which yielded the norm for the state complexity, amounts to setting one of the velocities $v_i$ to zero (see eq. (65)). In the context of the complexity of a single qubit, this gives a trajectory in SU(2) that does not use the $\sigma_i$ Pauli matrix associated with the velocity $v_i$. Thus, giving $\sigma_i$ a larger penalty will not change the results. It was argued in [62] that in certain experimental two-level systems the $\sigma_z$ gate is much harder to implement (see their discussion below eq. (2.4)). In order to estimate the complexity, the authors therefore proposed to use a Riemannian $F_2$ norm with a very large penalty assigned to the $\sigma_z$ direction. In our results we see that we can alternatively use the $F_1$ cost function with reference state $|\uparrow\rangle$ and the same effect of suppressing the use of the $\sigma_z$ gate will be achieved, even without the use of penalty factors. We expect that a similar effect will come into play for larger systems and intend to investigate this issue further in the future.

# 6 Geometrically local complexity and multipartite binding complexity

We now consider a multipartite system with a factorizable Hilbert space of the form (11). As we mentioned below eq. (11), oftentimes one of the constraints imposed on complexity (even outside the context of binding complexity) is that only generators between pairs of qubits are allowed. This type of locality constraint plays an important role in circuit complexity in realistic quantum computation scenarios, where many-qubit gates are typically harder to implement. Motivated by these considerations, in this section, we study Nielsen's state complexity according to the $F_1$ norm in eq. (4), using the following nomenclature adopted in [85]:

* Local complexity $C_{1,l}^{\text{state}}$ allows gates acting on one part or any pair of parts of the system.

* Geometrically-local complexity $C_{1,g}^{\text{state}}$ allows gates acting only on nearest-neighbour parts or within each part.
These notions of complexity are depicted in fig. 6. As a side comment, let us mention that the authors of [101] have noted that the notion of geometrically-local complexity seems incompatible with the holographic notions of complexity. On the other hand, many works have taken locality (as opposed to geometrical-locality) as a guiding principle in looking for a quantum mechanical dual of holographic complexity, see e.g., [3, 62, 63, 86].

In [85], the authors derived lower bounds on the notion of local and geometrically-local complexity in terms of the entanglement entropies of different bipartitions of the target state. We will now use the techniques developed in this work to derive new lower bounds, in terms of binding complexity, which improve the lower bounds of [85].

<figure><img src="image_7.png" /><figcaption>Figure 6: Pictorial representation of a spin chain where the vertical purple lines represent a separation between parts. (a) Circuit with gates acting on at most two parts of a system. (b) Circuit with gates acting only on nearest-neighbour parts or within each part. Both circuits can be used to study local complexity (120), but only case (b) is a circuit that can be used to study geometrically-local complexity (121). When the penalties associated with operations within a subsystem vanish (in this example, gate $G_3$ for both circuits), the definitions reduce to multipartite binding complexity, see eqs. (122a) and (122b), respectively.</figcaption></figure>

Let us define the set-up. The system is split into $n$ identical parts, each of them has an $N$-dimensional Hilbert space on which the special unitary space $SU(N)$ can act. Therefore, the possible circuits in the composite system belong to $SU(N^n)$.<fnref n="35" /> We denote the set of generators for each part as $T_i^A$, where the subscript $i$ runs over the list of generators $\{1, \dots, N^2 - 1\}$ within each subsystem, while the superscript $A$ refers to the site on which they act (i.e., a generator $T_i^A$ is the identity on all sites except for the $A$-th site).<fnref n="36" /> The full set of generators for the system is made by taking products of these $T_i^A$ for the different parts. Suppose that we construct our control Hamiltonian and unitary circuits from gates which act on at most two parties, that is

$$H_l(t) = \sum_{A=1}^{n} \sum_{i=1}^{N^2-1} Y_i^A T_i^A + \sum_{A=1}^{n-1} \sum_{B=A+1}^{n} \sum_{i,j=1}^{N^2-1} Y_{ij}^{AB} T_i^A T_j^B, \quad t \in [0, 1], \qquad (119)$$

$$H_g(t) = \sum_{A=1}^{n} \sum_{i=1}^{N^2-1} Y_i^A T_i^A + \sum_{A=1}^{n-1} \sum_{i,j=1}^{N^2-1} Y_{ij}^{AA+1} T_i^A T_j^{A+1}, \quad t \in [0, 1],$$

where the two Hamiltonians differ by the restriction of the non-local generators of $H_g(t)$ to nearest-neighbour-parts interactions. We will denote our reference and target states as $|\psi_R\rangle, |\psi_T\rangle$. The notions of complexity introduced in this section correspond to

* The local complexity

$$C_{1,l,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{H_l(t)} \int_0^1 dt \left( \sum_{A=1}^n \sum_{i=1}^{N^2-1} |Y_i^A| + \sum_{A=1}^{n-1} \sum_{B=A+1}^n \sum_{i,j=1}^{N^2-1} |Y_{ij}^{AB}| \right). \quad (120)$$

<footnote n="35">³⁵This applies, for example, to a subdivision of a spin chain to $n$ equal parties, where $N = 2^{\# \text{ spins}}$.</footnote>

<footnote n="36">³⁶The notation has changed from the previous sections in which $A$ denoted one of the two parties $A$ and $B$, and not a running index.</footnote>
• The geometrically-local complexity

$$C_{1,g,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{H_g(t)} \int_0^1 dt \left( \sum_{A=1}^n \sum_{i=1}^{N^2-1} |Y_i^A| + \sum_{A=1}^{n-1} \sum_{i,j=1}^{N^2-1} |Y_{ij}^{A,A+1}| \right). \quad (121)$$

In both cases, the minimization is performed over all the control Hamiltonians of the form (119) which connect the reference and target states. We also added the subscript “hom” referring to the homogeneous choice of penalties equal to unity for all the generators (including the local ones).

It is useful for the following to introduce also a multipartite notion of binding complexity $\mathcal{MBC}[|\psi_T\rangle, |\psi_R\rangle]$. This can be defined as a special case of the (geometrically-) local complexity, when the penalty factors assigned to local gates vanish:

$$\mathcal{MBC}_{1,l,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{H_l(t)} \int_0^1 dt \sum_{A=1}^{n-1} \sum_{B=A+1}^{n} \sum_{i,j=1}^{N^2-1} |Y_{ij}^{AB}|, \quad (122a)$$

$$\mathcal{MBC}_{1,g,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] = \min_{H_g(t)} \int_0^1 dt \sum_{A=1}^{n-1} \sum_{i,j=1}^{N^2-1} |Y_{ij}^{A,A+1}|. \quad (122b)$$

These notions of binding complexity are pictorially depicted in fig. 6. Since complexity decreases when we decrease the penalties, we immediately get

$$\mathcal{MBC}_{1,l,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] \le C_{1,l,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle], \quad \mathcal{MBC}_{1,g,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle] \le C_{1,g,\text{hom}}^{\text{state}}[|\psi_T\rangle, |\psi_R\rangle]. \quad (123)$$

From now until the end of this section, we will omit the explicit dependence of the complexity on the reference and target states. In order to derive bounds on the (geometrically-) local complexity and the multipartite binding complexity, we introduce some useful definitions:

• The bipartite binding complexity $\mathcal{BC}_{1,\text{hom}}^{\text{state}}[\bar{A}]$, computed using the norm (18) with the two subsystems given by the union of the subsystems $\{1, \dots, \bar{A}\}$ in one part, and $\{\bar{A} + 1, \dots, N\}$ in the other part. This is the kind of binding complexity that we studied, e.g., in section 4.1.

• The bipartite binding complexity $\mathcal{BC}_{1,g}^{\text{state}}[\bar{A}]$, where all the non-local gates which have support on both sides of the partitioning, are not allowed, except from nearest-neighbour interactions of the form $T_i^{\bar{A}} T_j^{\bar{A}+1}$. This is similar to $\mathcal{BC}_{1,\text{set}}^{\text{state}}$, except that we allow several generators $T_i^{\bar{A}} T_j^{\bar{A}+1}$ with support on parties $\bar{A}$ and $\bar{A} + 1$, instead of just one.

The previous notions of bipartite binding complexity are pictorially represented in fig. 7. We note that those quantities refer to the bipartite binding complexity and therefore, all generators that have support on just one side of the partitioning are free. Of course, this implies that

$$\mathcal{BC}_{1,g}^{\text{state}}[\bar{A}] \ge \mathcal{BC}_{1,\text{hom}}^{\text{state}}[\bar{A}], \quad (124)$$

because for the circuit that minimizes $\mathcal{BC}_{1,g}^{\text{state}}$ the cost function is the same, but on the right-hand side we can build circuits using a larger set of generators, which allows to find in principle a better geodesic compared to the left-hand side.
<figure><img src="image_8.png" /><figcaption>Figure 7: Pictorial representation of a spin chain where the vertical purple lines represent a separation between parts $P_1, \dots, P_4$. (a) Circuit with gates acting on at most two parts of a system. (b) Circuit with gates acting only on nearest-neighbour parts or within each part. Consider a bipartition of the system into $A = \{P_1, P_2, P_3\}$ and $B = \{P_4\}$. Both circuits can be used to study the bipartite binding complexity $BC_{1,\text{hom}}^{\text{state}}(P_3)$, but only case (b) is a circuit that can be used to study $BC_{1,g}^{\text{state}}(P_3)$ because $G_1$ in subfigure (a) does not act on nearest-neighbor-parts. The complexity is evaluated according to a norm only assigning non-vanishing cost to gates connecting the parts A and B among the allowed gates.</figcaption></figure>

Consider now the following term involving nearest-neighbour parts in the expression for the geometrically-local complexity (122b)

$$ \int_0^1 dt \sum_{i,j=1}^{N^2-1} |Y_{ij}^{\bar{A},\bar{A}+1}|, \qquad (125) $$

where the value of the control functions $Y$ in this expression are those determined by the minimization done in eq. (122b). Intuitively, when local generators (with respect to the bipartitioning) have no cost, the minimal circuit can be further optimized by using those gates. In other words, $BC_{1,g}^{\text{state}}[\bar{A}]$ minimizes the cost of producing the target state with the norm of (125). If we consider a sub-optimal trajectory given by the minimization done in eq. (122b), we will certainly obtain a larger value for this norm. This formally corresponds to the inequality

$$ \int_0^1 dt \sum_{i,j=1}^{N^2-1} |Y_{ij}^{\bar{A},\bar{A}+1}| \geq BC_{1,g}^{\text{state}}[\bar{A}]. \qquad (126) $$

Summing over all parties gives on the left-hand side the definition (122b), and on the right-hand side a bound involving the bipartite notion of binding complexity:

$$ \mathcal{MC}_{1,g}^{\text{state}} \geq \sum_{A=1}^{n-1} \mathcal{BC}_{1,g}^{\text{state}}[A]. \qquad (127) $$

If we now combine the inequalities (123), (124) and (127), we obtain a lower bound on geometric complexity:

$$ C_{1,g,\text{hom}}^{\text{state}} \geq \sum_{A=1}^{n-1} \mathcal{BC}_{1,\text{hom}}^{\text{state}}[A]. \qquad (128) $$

Similar arguments (applied to the subset of terms in the expression (122a) corresponding to the binding complexity $BC_{1,\text{hom}}^{\text{state}}[A]$), lead to a lower bound on local complexity:

$$ C_{1,l,\text{hom}}^{\text{state}} \geq \max_{A} \mathcal{BC}_{1,\text{hom}}^{\text{state}}[A]. \qquad (129) $$
The bounds above are an improvement of the results obtained by [85] using the notion of entanglement power; the results of [85] are recovered by also using the inequality (89).

The bounds we derived in terms of bipartite binding complexity are tighter than the ones of [85]. To illustrate this, consider a system consisting of a three-qubit chain, with reference and target states

$$|\psi_R\rangle = |000\rangle, \quad |\psi_T\rangle = \cos \alpha |000\rangle + \sin \alpha |111\rangle, \qquad (130)$$

where we assume $0 \le \alpha \le \pi/4$. To use (128), when $A = 1$ we partition the system between the first qubit and the rest, and when $A = 2$, the partition is between the third qubit and the rest. In both partitions, the Schmidt decomposition is straightforward, and the coefficients are $\vec{\lambda} = (\cos \alpha, \sin \alpha)$. In section 4.2 we argued that the binding complexity for such partitions is given by the binding complexity between two qubits having the same Schmidt coefficients. Using eq. (47), we find that the binding complexity for each partition is given by $\alpha$, therefore the bound in (128) becomes

$$C_{1,g,\mathrm{hom}}^{\mathrm{state}}[|\psi\rangle] \ge \sum_{A=1}^{n-1} \mathcal{B}C_{1,\mathrm{hom}}^{\mathrm{state}}[A] = 2\alpha \ge -\cos^2 \alpha \log(\cos^2 \alpha) - \sin^2 \alpha \log(\sin^2 \alpha). \quad (131)$$

The rightmost term in this inequality is the bound derived in reference [85], which we then proved to be less restrictive, in this case, than (128).

# 7 Discussion

## 7.1 Summary of results

In this paper, we studied binding state complexity, a notion of Nielsen's geometric complexity for a multipartite system (11) such that gates acting on each subsystem separately (we refer to such gates as local) have zero cost. Here we focused, mainly, on the bipartite case (12). The freedom to act with local generators implies that the norm $F^{\mathrm{state}}$ induced on the space of states can be reduced to a norm that only depends on the Schmidt coefficients $\vec{\lambda}$ of the state and on their derivatives $\vec{\lambda}$ along a trajectory connecting the reference and target states:

$$F_p^{\mathrm{state}}[|\psi(t)\rangle, |\dot{\psi}(t)\rangle] \to \mathcal{BF}_p[\vec{\lambda}, \dot{\vec{\lambda}}]. \qquad (132)$$

To find this norm, we took an intermediate step in which we showed how the problem of finding the geodesics in SU($N_A N_B$) reduces to the manifold SO($N_A$), which is the rotation group of the Schmidt coefficients. The norm on SO($N_A$) is in terms of the parameters of $R$, the generator of rotations, and the equations of motion for the Schmidt coefficients read $\vec{\lambda}(t) = R\vec{\lambda}(t)$. An additional minimization over the stabilizer of $\vec{\lambda}$ inside SO($N_A$) is required to find the norm in terms of the Schmidt coefficients. Finally, the optimal trajectory whose length computes the binding state complexity is found by minimizing the length according to this last norm.

The main results of the work are the **exact analytic expressions** for binding state complexity using the homogeneous norm defined in eq. (20a). We summarize them in table 1 and in the text below:

* The $\mathcal{BF}_1$ norm depends on the choice of a basis of generators of the unitary group. In table 1 we report the results obtained using the generalized Gell-Mann generators (55). The optimal trajectory corresponds to a set of rotations along planes spanned by pairs of Cartesian axes to get from $\vec{\lambda}(0) = (1, 0, \dots, 0)$ to $\vec{\lambda}(1)$, see eq. (72).
• The $\mathcal{BF}_2$ norm is independent of the basis of generators in the unitary group, and turns out to be the standard round metric of the unit sphere. The binding state complexity is the arc length from the reference to the target state.

<table><caption>Table 1: Results for the homogeneous norm (20a). First row: based on eqs. (71), (79). $N_A$ represents the dimension of the smaller Hilbert space in the decomposition (12), while $\lambda_1$ is the maximal Schmidt coefficient. Second row: based on eqs. (96), (97). $N_{A(B)}$ are overall factors related to the normalization of the generators, see eq. (13), and it is understood that $|\vec{\lambda}(t)|^2 = 1$. $S_\infty(|\psi_T\rangle)$ is the minimal Rényi entropy of the reduced density matrix of the target state $|\psi_T\rangle$.</caption><thead><tr><th></th><th><i>B</i><i>F</i><sub>p,hom</sub> norm</th><th>Binding state complexity</th></tr></thead><tbody><tr><td><i>p</i> = 1 (Gell-Mann basis)</td><td>&sum;<sub><i>m</i> &gt; 1</sub><sup><i>N</i><sub>A</sub></sup> |<i>&lambda;</i><sub><i>m</i></sub>| / &lambda;<sub>1</sub></td><td>Eq. (79)</td></tr><tr><td><i>p</i> = 2 (any basis)</td><td>(2 / (<i>N</i><sub>A</sub><i>N</i><sub>B</sub>) &ddot;<i>&lambda;</i>(<i>t</i>))<sup>1/2</sup> / &radic;(2 / (<i>N</i><sub>A</sub><i>N</i><sub>B</sub>) arccos(e<sup>-1/2 <i>S</i><sub>&infin;</sub>(|<i>&psi;</i><sub><i>T</i></sub>)</sup>))</td><td></td></tr></tbody></table>

For the $\mathcal{BF}_{1,hom}$ with the Pauli basis of generators (83), we gave an exact result for a system of any number of qubits which is bipartitioned to a single qubit and the rest, and we found bounds in the general case.

We argued that the binding complexity computed by any $F_1$ norm of the form (101) is the same, up to an overall normalization, as long as the generators forming the optimal trajectory all lie on the same adjoint orbit of the $SU(N_A) \otimes SU(N_B)$ action as the cheapest generator. The latter condition is satisfied for the Gell-Mann and Pauli bases. As a consequence, the binding state complexity is rather insensitive to the penalty factors chosen for the non-local gates. This provides a huge class of systems that are accessible analytically. The freedom to choose the penalty factors for non-local generators may help in implementing these circuits in an experimental setting, with the caveat that local gates can have a small but non-vanishing cost, as discussed at the end of section 4.4.

We have shown that the binding state complexity $\mathcal{BC}_{1,hom}^{\text{state}}$ applied to a system composed of two qutrits is equivalent to the $F_1$ state complexity (not binding) for a single qubit, for target states that belong to a certain region around the reference state. This extends previous studies in the literature that focused on the Riemannian $F_2$ norm for a single spin [62].

Finally, we inferred lower bounds for the multipartite binding complexity and the geometrically local complexity [85], in terms of the bipartite binding complexity for different bipartitions of a given system.

## 7.2 Future directions

The approach presented in this work allowed us to pursue various directions that were previously hindered by the fact that the $F_1$ norm is not a Riemannian metric. Several questions remain open, which we intend to study in the future:

1. **Generalized Pauli basis.** The generalized Pauli basis is natural when considering systems of spins on a lattice that entail a notion of spatial locality in which each spin resides at a different point, and therefore it is commonly used in the literature on Nielsen's complexity. We showed in section 4.3 that the binding complexity using the $\mathcal{BF}_2$ norm in the generalized Pauli-basis differs from the one obtained using the Gell-Mann basis just by a normalization factor. The $\mathcal{BF}_1$ norm case, however, is more involved. We derived several bounds in section 4.2 for the complexity associated
with the $\mathcal{BF}_1$ norm in the Pauli-basis, but it would be desirable to find the optimal trajectory and the exact complexity in the Pauli basis for an arbitrary number of spins. A discrete version of this problem was considered in [84].

2. **Extension of our results for other $\mathcal{BF}$ norms.** We were able to calculate $\mathcal{BC}_{1,\text{hom}}$ only in a region of the group manifold close to the reference state, see section 4.1. We would like to extend the analysis to the whole space. Furthermore, here we focused on the $\mathcal{BF}_1$ and $\mathcal{BF}_2$ norms, finding that the binding complexity for the $p=2$ case is a function of the Rényi min-entropy. All the other $\mathcal{BF}_p$ norms will result in binding complexities that depend on the Schmidt coefficients. It is reasonable to assume that one can obtain all the Schmidt coefficients from the knowledge of binding complexities for all the $\mathcal{BF}_p$ norms, as one can obtain them from knowing all the Rényi entropies. We plan to obtain explicit results for these relations.

3. **Curvature properties of binding complexity geometry.** The Gell-Mann basis for the unitary group $SU(N)$ has the advantage of clearly distinguishing between elements of the stabilizer of the reference state, and elements of the coset space.<fnref n="37" /> In the framework of Nielsen's complexity, it was argued that geodesic deviation and the switchback effect, which are essential features of complexity geometry, require regions with negative sectional curvature [86, 87]. In the Riemannian case, a negatively curved geometry can be obtained whenever the commutator structure of the algebra satisfies a relation of the kind $[\mathcal{e}\mathcal{a}\mathcal{y}, \mathcal{e}\mathcal{a}\mathcal{y}] = \mathcal{h}\mathcal{a}\mathcal{r}\mathcal{d}$, where $\mathcal{e}\mathcal{a}\mathcal{y}$ and $\mathcal{h}\mathcal{a}\mathcal{r}\mathcal{d}$ refer to a smaller or bigger penalty factor inside the cost function, respectively. This statement was argued to be a consequence of Pythagorean's theorem in curved space in reference [62], and it was later confirmed for the state complexity of a single qutrit in [66]. In the case of the $\mathcal{BF}_{p,\text{set}}$ norm, we assign penalties 0 to local generators, 1 to a single non-local generator $T_{11} \equiv T_1^A \otimes T_1^B$ and the formal cost $\infty$ to all the other non-local generators, which are not allowed. As a consequence, the algebra naturally admits a subset of commutation relations of the form $[\textit{local}, T_{11}] = \textit{non-local}$, which satisfies the above-mentioned pattern $[\mathcal{e}\mathcal{a}\mathcal{y}, \mathcal{e}\mathcal{a}\mathcal{y}] = \mathcal{h}\mathcal{a}\mathcal{r}\mathcal{d}$. We therefore expect that the sectional curvature is negative along the planes generated by any local generator and the single allowed gate $T_{11}$. Since binding state complexity involves non-trivial limits on the penalty factors, it would be interesting to check whether these statements on the negative sectional curvature are still valid. Furthermore, it is interesting to understand the implications of the above structure of commutators in the context of the switchback effect. Let us further remind the reader that in the context of a single qubit, we have shown in section 5 that the circuit (as opposed to binding) complexity associated with the $F_1$ norm provides a natural way of suppressing hard gates which was achieved in previous works by the assignment of penalty factors. It is interesting to study under which conditions the $F_1$ norm can also provide an alternative mechanism for geodesic deviation.

4. **Characterizing multipartite entanglement.** In reference [3], the authors suggested that the multipartite version of binding complexity can distinguish whether a multipartite state will be entangled after tracing over one of the parts. The authors use as an example the GHZ state which is separable after the tracing and has an estimated binding complexity of $\mathcal{O}(n)$, while the $W$ state which is not and has an estimated binding complexity of $\mathcal{O}(n^2)$, where $n$ is the number of spins. We plan to

<footnote n="37">³⁷The Gell-Mann matrices introduced in eq. (55) are proportional to the set of generators $E_{ij}$, $F_{ij}$ defined in eqs. (I.6.60)-(I.6.61) on page 208 of reference [102]. The latter form a basis that distinguishes between the generators of the maximal subgroup $SU(N-1) \times U(1)$ and the coset space.</footnote>
consider the multipartite version in spin systems and understand more generally the relation between multipartite binding complexity and the robustness of multipartite entanglement.

5. **Generalized complexity measures.** Binding state complexity can be interpreted as an example of a larger class of complexity measures in which one does not require to reach a specific target state, but only an equivalence class of states that satisfy certain properties. In this case, the required property is based on the entanglement structure encoded by the Schmidt coefficients. Other examples could require that the equivalence class of target states has a fixed vacuum expectation value of certain operators, i.e., $\langle\psi_T|O_i|\psi_T\rangle = Q_i$. This seems natural because usually we cannot distinguish states with absolute precision.

6. **Holographic interpretation.** It remains to be investigated if there is a holographic dual to the quantities we considered. In [3], it was conjectured that the binding complexity is dual to the volume of the interior of a multi-boundary wormhole (see e.g., [103]), minimized over all the boundary times (see their discussion in sections 5 and 7). While in this work we mostly considered the bipartite case, which is conjectured to be related to the two-boundary wormhole, we should still be able to compare the binding complexity to the expected results from the gravity side in that specific case. A holographic computation of the volume of a three-dimensional wormhole with three boundaries has been performed in [83]. It would be interesting to compare this case to the multipartite binding complexity, partially considered in section 6 and which we hope to investigate further in the near future. Another interesting candidate for the holographic dual may be deduced as follows. Given the relation between the $\mathcal{BF}_2$ binding complexity and the minimal Rényi entropy, the holographic dual in that case may be obtained along the lines of [104, 105]. It would be interesting to compare the proposal obtained this way to the previous one.

7. **Relation to kinematic space.** For the vacuum state in a two-dimensional CFT, there is a conjectured relation between volume complexity and entanglement via the Crofton formula [106, 107]. Furthermore, reference [108] claims that the minimal Rényi entropy of a quantum chain at the critical point is related to the entanglement entropy. Combining these observations, one may argue that the volume complexity could be expressed as an integral over the binding complexity of an appropriate trajectory in kinematic space. We intend to investigate this link more deeply.

8. **Subsystem complexity.** Binding complexity is naturally defined when there is a partition of the Hilbert space into subsystems. In this case, the state of each subsystem separately is generically mixed. Other interesting notions of complexity have also been defined in the context of mixed states, such as the purification or spectrum complexity proposed in reference [109].<fnref n="38" /> It would be interesting to compare these notions of complexity to the binding complexity of a state which is the purification of such mixed states (where the bipartition is between the system and its purification).

9. **Ancillary qubits.** We would like to understand the role of ancillary qubits in binding complexity. In the Gell-Mann basis the result (79) is unaffected by enlarging the dimension of $\mathcal{H}_B$. Therefore there is no influence of ancillary qubits on the results in the Gell-Mann basis. Furthermore, in section 4.2, we discussed an exact result for the complexity in the Pauli basis for two systems, one of which (system A) contains

<footnote n="38">³⁸These definitions were originally suggested in the context of holographic complexity for subregions [106, 107, 109–116] and were further developed for Gaussian states, see e.g., [114, 117, 118].</footnote>
a single spin and the other (system B) an arbitrary number. In this case, we were able to prove that when the single spin in system A was entangled with a single spin in system B, the extra spins in system B, which we referred to as ancillary, did not participate in the optimal trajectory and did not improve the result for the complexity. It would be interesting to understand if this statement is generally valid in the Pauli basis for the addition of spins which are unentangled both in the reference and the target state for arbitrary numbers of spins in the systems A and B, and if a similar statement can be made for binding complexity with more than two parties.

## Acknowledgements

We are happy to thank Vijay Balasubramanian, Ramy Brustein, Hugo Camargo, Nat Levine, Christian Northe, Yaron Oz, Eran Palti, Daniel Rohrlich, Michael Walter and Erez Zohar for valuable discussions. We are particularly grateful to Adam Chapman for providing us with the proof of the inequalities in appendix A. The work of SB, SC and TS is supported by the Israel Science Foundation (grant No. 1417/21), by the German Research Foundation through a German-Israeli Project Cooperation (DIP) grant “Holography and the Swampland”, by Carole and Marcus Weinstein through the BGU Presidential Faculty Recruitment Fund and by the ISF Center of Excellence for theoretical high energy physics. SB is supported by an Azrieli fellowship funded by the Azrieli foundation. SB gratefully acknowledges the Perimeter Institute for Theoretical Physics, Brandeis University and the Simons Center for Geometry and Physics at Stony Brook University, where part of the research in this paper was carried out. Research at Perimeter Institute is supported by the Government of Canada through the Department of Innovation, Science and Economic Development and by the Province of Ontario through the Ministry of Colleges and Universities. SB and GP thank the organizers and participants of the workshop YITP-T-23-01 “Quantum information, quantum matter and quantum gravity” in Kyoto for hospitality and for interesting discussions.

## A An inequality for the $F_1$ cost function

This appendix is devoted to the proof of the inequality (38). We derive the following chain of relations:

$$ \begin{aligned} & \min_{\mathrm{Ad}_U} \sum_{b,j} |\tilde{Y}^{ai}(\mathrm{Ad}_{UA})_a^b(\mathrm{Ad}_{UB})_i^j| \ge \min_O \sum_{b,j} |\tilde{Y}^{ai}(O^A)_a^b(O^B)_i^j| = \\ & = \sum_a \sigma_a(\tilde{Y}) \ge \sum_a |\tilde{Y}^{aa}|. \end{aligned} \quad (133) $$

The first minimum in the first line is over all possible adjoint matrices. $O^{A(B)}$ denote general orthogonal matrices and the second minimum in the first line is with respect to these matrices. $\sigma_a(\tilde{Y})$ are the singular values of the matrix $\tilde{Y}_{ai}$ and the first quantity in the second line is the so-called nuclear or Schatten 1-norm [119]. Let us go through the derivation of the previous steps:

* The inequality in the first line follows from the fact that the adjoints $\mathrm{Ad}_U$ are specific orthogonal matrices, thus the complete set of orthogonal matrices is larger than those in the adjoint. Note that here we are comparing the adjoint of $\mathrm{SU}(N_A)$ to $\mathrm{SO}(N_A^2-1)$ (not to $\mathrm{SO}(N_A)$).
• In going from the first to the second line, we used the inequality$^{39} \sum_a \sigma_a(T) \le \sum_{a,b} |T_{a,b}|$ with $T = (O^A)^T \tilde{Y} O^B$, together with the observation that the nuclear norm is invariant under orthogonal transformations. The inequality is saturated when performing the minimization over the orthogonal matrices, once we choose $O^A, O^B$ according to the singular value decomposition of $\tilde{Y}^{ai}$.

• The last inequality (see e.g., [121]) can be proven by applying von Neumann's trace inequality $|\text{Tr}(\tilde{Y}B)| \le \sum_i \sigma_i(\tilde{Y})\sigma_i(B)$. This is done by choosing $B$ as a diagonal matrix with elements $\pm 1$, where the sign is chosen such that $|\text{Tr}(\tilde{Y}B)| = \sum_a |\tilde{Y}^{aa}|$ recovers the last inequality in (133). The singular values of the matrix $B$ are of course $\sigma_i(B) = 1$. Notice that this argument can also be applied to rectangular $\tilde{Y}$ matrices. In fact, one can pad $\tilde{Y}$ with zeroes to define a square matrix $\hat{Y}$, and then perform the same steps outlined above to obtain the last inequality for $\hat{Y}$. The non-zero singular values and the non-zero diagonal elements of $\hat{Y}$ are equal to the ones of $\tilde{Y}$.

Once we recognize that the first term in (133) is the definition of the homogeneous $F$-norm in eq. (36) with $p = 1$, and upon an additional minimization over $\tilde{Y}$ subject to the constraint (35), the chain of relations (133) can be written as (39).

## B Complexity with a large number of Schmidt coefficients

In this appendix, we focus on a quantum-mechanical system with a large number of Schmidt coefficients. In this case, we may recast the binding complexity (79) in the form of an integral. The integral is generally easier to compute than the original discrete sum. Furthermore, as mentioned at the end of section 4.1, this can be viewed as a first step towards the study of binding complexity in QFT.

First, let us re-write eq. (79) in reverse order of summation,$^{40}$

$$ \mathcal{BC}^{\text{state}}_{1,\text{hom}} = - \sum_{k=1}^{N_A-1} \int_{\bar{\lambda}_{k+1}}^{\bar{\lambda}_{k+2}} \frac{kd\lambda}{\sqrt{1 - \sum_{j=k+2}^{N_A} \bar{\lambda}_j^2 - k\lambda^2}}. \quad (134) $$

The index $k$ counts the number of Schmidt coefficients larger than $\bar{\lambda}_{k+1}$. In the limit of an infinite number of Schmidt coefficients $N_A \to \infty$, we obtain a distribution $\bar{\lambda}(k)$. We will treat $\bar{\lambda}_1$ as a normalization that does not have to be a part of the continuous distribution. This means that the distribution represents all the Schmidt coefficients apart from $\bar{\lambda}_1$ and that $\bar{\lambda}_1^2 + \int_0^\infty \bar{\lambda}(k)^2 dk = \bar{\lambda}_1^2 + \int_0^\bar{\lambda}(0) \rho(\bar{\lambda}) \bar{\lambda}^2 d\bar{\lambda} = 1$, where $\rho(\bar{\lambda})$ is the density of Schmidt coefficients. This type of separation can be seen for the entanglement spectrum derived in [122] for a one-dimensional system near the scale-invariant regime.

In the continuum, the index $k$ becomes a continuous function of $\bar{\lambda}$, such that $k(\bar{\lambda}) = \int_{\bar{\lambda}}^{\bar{\lambda}(0)} \rho(\lambda) d\lambda$, and the sum $\sum_{j=k+2}^\infty \bar{\lambda}_j^2$ turns to $\int_0^{\bar{\lambda}(k)} \rho(\lambda) \lambda^2 d\lambda$. With these replacements the complexity reads

$$ \mathcal{BC}^{\text{state}}_{1,\text{hom}} = \int_0^{\bar{\lambda}(0)} d\lambda \frac{k(\lambda)}{\sqrt{\bar{\lambda}_1^2 + \int_{\bar{\lambda}}^{\bar{\lambda}(0)} \rho(x) x^2 dx - \lambda^2 k(\lambda)}} = \int_0^{\bar{\lambda}(0)} d\lambda \frac{k(\lambda)}{\sqrt{\bar{\lambda}_1^2 + 2 \int_{\bar{\lambda}}^{\bar{\lambda}(0)} xk(x) dx}}. \quad (135) $$

<footnote>$^{39}$In turn, this inequality can be proven using the following property of the singular values $\sum_a \sigma_a(M+N) \le \sum_a \sigma_a(M) + \sigma_a(N)$ [120] if we decompose the matrix $T$ into matrix units of the form $t_{a,b}\delta_{a,b}$ and apply this inequality to the sum. In this case $|t_{a,b}|$ is the only non-zero singular value of each matrix unit.</footnote>

<footnote>$^{40}$Note that we are using the form of eq. (79) before performing the integration. The reason is that we will be expanding under the assumption that successive Schmidt coefficients are close to each other and therefore the integrands themselves will appear naturally in our expressions.</footnote>
With an additional change of variables, the above formula can be rewritten as

$$ \mathcal{BC}_{1,\text{hom}}^{\text{state}} = -\int_{0}^{\infty} dk \frac{k\tilde{\lambda}'(k)}{\sqrt{1 - \int_{k}^{\infty} \tilde{\lambda}(x)^2 dx - k\tilde{\lambda}(k)^2}} = -\int_{0}^{\infty} dk \frac{k\tilde{\lambda}'(k)}{\sqrt{\bar{\lambda}_{1}^{2} + \int_{0}^{k} \tilde{\lambda}(x)^2 dx - k\tilde{\lambda}(k)^2}}. \quad (136) $$

We note that this result can be obtained by considering a continuous version of the norm of (71) of the form $\mathcal{BF}_{1,\text{hom}}dt = \int_{0}^{\infty} dk \frac{|d\tilde{\lambda}(k,t)|}{\bar{\lambda}_{1}(t)}$, and by then evaluating it on a trajectory that is analogous to the discrete one. A caveat should be made regarding the range of validity of (135) and (136). We remind the reader that the result in eq. (134) was derived in the region $\bar{\lambda}_{1} \ge \bar{\lambda}_{2} + \bar{\lambda}_{3}$, where we obtained the exact form of the norm $\mathcal{BF}_{1,\text{hom}}$. A natural guess for the continuum version of the previous constraint is $\lambda_{1} \ge 2\tilde{\lambda}(0)$. When the spectrum distribution does not have a separation between the largest Schmidt coefficient and the continuous part, i.e., $\int_{0}^{\infty} \tilde{\lambda}(k)^2 dk = 1$, the expected constraint reads $\tilde{\lambda}(0) \ge 2\tilde{\lambda}(0)$, and is clearly violated. Nevertheless, we suspect that the result (71) may in fact be valid even outside the range dictated by the inequality $\bar{\lambda}_{1} \ge \bar{\lambda}_{2} + \bar{\lambda}_{3}$, even for the discrete case, and in that case the continuum result could perhaps be used for distributions without a spectrum separation.

The expressions above, besides being suitable for a continuous distribution, can help in approximating the discrete $\mathcal{BC}_{1,\text{hom}}^{\text{state}}$ by a simpler integral expression. In this case, we have to explain how to turn our discrete distribution into a continuous function over which we can integrate. We do this as follows. Let $\tilde{\lambda}_{n}$, where $n$ is an integer that goes from 1 to $\infty$, be a set of ordered Schmidt coefficients. Define a smooth monotonic function $\tilde{\lambda}(x)$ that interpolates the set starting from $n \ge 2$, at points that differ by a unit interval, i.e.,

$$ \tilde{\lambda}(x_0 + m) = \bar{\lambda}_{m+2}, \quad m \ge 0. \quad (137) $$

In addition, set $x_0$ such that $\int_{0}^{\infty} \tilde{\lambda}(x)^2 dx = \sum_{n=2}^{\infty} \bar{\lambda}_n^2 = 1 - \bar{\lambda}_1^2$. Approximating the integrals of (136) (or (135)) by a sum using the values of the integrands at integer values gives an expression that is close to the one obtained by approximating the integrals of the exact result in (134).

The error in approximating (134) by (136) (or (135)) can be determined by standard numerical techniques, and will depend on the difference between adjacent Schmidt coefficients and the exact choice of the interpolating function $\tilde{\lambda}(x)$. We note that if $\bar{\lambda}_k - \bar{\lambda}_{k+1}$ is not sufficiently small for some specific value of $k$, one could combine the continuous approximation with the exact calculation. A particular case of interest would be if there is a finite number $N_A$ Schmidt coefficients, where $\bar{\lambda}_{N_A}$ is of order $\frac{1}{\sqrt{N_A}}$. In this case, the boundaries of integration in (136) could be taken from 0 to $N_A - 2$ and an additional term of value $\sqrt{N_A - 1} \arcsin[\bar{\lambda}_{N_A} \sqrt{N_A - 1}]$ coming from the $N_A - 1$ term of (134), should be added.

Let us exemplify the above-mentioned continuous approximation in the case of the thermofield double state of two harmonic oscillators. We start from a discrete system whose Hilbert space decomposes according to eq. (12), with Schmidt coefficients$^{41}$

$$ \bar{\lambda}_n = \sqrt{e^\beta - 1} e^{-\frac{\beta n}{2}} \quad (n \ge 1), \quad \text{such that} \quad \sum_{n=1}^{\infty} \bar{\lambda}_n^2 = 1. \quad (138) $$

<footnote>$^{41}$Notice that the Schmidt coefficients do not satisfy $\bar{\lambda}_{1} \ge \bar{\lambda}_{2} + \bar{\lambda}_{3}$ when $\beta$ is sufficiently small, and therefore we are uncertain whether the binding complexity should indeed be computed using (134). Nevertheless, this serves as a good example for the continuous approximation (136). Our methods can be applied to systems that have an entanglement gap in the thermodynamic limit, for instance, those discussed in [93].</footnote>
In particular, $\bar{\lambda}_1$ is the largest among the previous Schmidt coefficients. According to the discussion around eq. (137), in the continuum we define a new distribution related to the discrete one by

$$ \tilde{\lambda}(x) = \sqrt{e^{\beta} - 1} e^{-\frac{\beta(x-x_0)}{2}}, \quad \text{with } x_0 \text{ such that } \bar{\lambda}_1^2 + \int_0^\infty dx \tilde{\lambda}(x)^2 = 1. \quad (139) $$

We can now compare the exact result for the binding complexity in the discrete case (79) using the Schmidt coefficients (138) with the continuous approximation (136) using the distribution (139). In the regime $\beta \to 0$, the binding complexity is of order $1/\sqrt{\beta}$, while the relative error found by this procedure can be verified numerically to be smaller than $\beta\sqrt{\beta}$.

## C Gell-Mann trajectory length in the Pauli basis

In this appendix, we prove the claim of section 4.2 that if we use the trajectory optimizing the $\mathcal{B}\mathcal{F}_{1,\text{hom},\text{Gell-Mann}}$ norm and evaluate its length according to the $\mathcal{B}\mathcal{F}_{1,\text{hom},\text{Pauli}}$ norm, we recover the same length. More explicitly, we claim that the $\mathcal{B}\mathcal{F}_{1,\text{hom},\text{Pauli}}$ norm evaluated on the Hamiltonian (80)

$$ H = \sum_{m>1}^{N_A=2^{n_A}} w_m(t) T_{c(m,1)}^A T_{c(1,m)}^B, \quad w_m = \frac{\dot{\lambda}_m^g(t)}{\lambda_1^g(t)}, \quad (140) $$

reads

$$ \mathcal{B}\mathcal{F}_{1,\text{hom},\text{Pauli}} = \sum_{m>1}^{N_A} |w_m(t)|, \quad (141) $$

which is the same as the result for the Gell-Mann basis (71).

To demonstrate this, let us use the definition (55) for the Gell-Mann basis, in particular

$$ T_{c(m,1)} = i(|m\rangle\langle 1| - |1\rangle\langle m|), \quad T_{c(1,m)} = |m\rangle\langle 1| + |1\rangle\langle m|. \quad (142) $$

Here it is understood that the basis vectors used in this definition are

$$ |m\rangle^A = |\vec{m}\rangle^A = |m_1 m_2 \dots m_{n_A}\rangle, \quad |m\rangle^B = |\vec{m}\rangle^B = |m_1 m_2 \dots m_{n_A} \uparrow \dots \uparrow\rangle, \quad (143) $$

where $m_i \in \{|\uparrow\rangle, |\downarrow\rangle\}$ represent spins up or down, and the basis vectors for the system B have been padded with $|\uparrow\rangle$ to reach the correct length. The different generators in the Pauli basis can be enumerated by vectors $\vec{\mu}, \vec{\nu}$ as in the right-hand side of eq. (83)

$$ T_{\vec{\mu}\vec{\nu}} = \sigma_{\mu_1} \otimes \dots \otimes \sigma_{\mu_{n_A}} \otimes \sigma_{\nu_1} \otimes \dots \otimes \sigma_{\nu_{n_B}}. \quad (144) $$

To find the control functions, we should take the normalized trace of each generator against the Hamiltonian for our trajectory, that is

$$ Y_{\vec{\mu}\vec{\nu}} = \frac{1}{2^{n_A+n_B}} \text{tr}[H \cdot T_{\vec{\mu}\vec{\nu}}]. \quad (145) $$

Using the following properties of the trace and tensor products $\text{tr}(A \otimes B) = \text{tr}(A) \cdot \text{tr}(B)$ and $(a_1 \otimes b_1) \cdot (a_2 \otimes b_2) = (a_1 \cdot a_2) \otimes (b_1 \cdot b_2)$, as well as the cyclic property of the trace, we can rewrite the above control functions as

$$ Y_{\vec{\mu}\vec{\nu}} = \frac{i}{2^{n_A+n_B}} \sum_{\vec{m} \neq \vec{1}} \tilde{v}_{\vec{m},1} \langle \vec{1} | \sigma_{\vec{\mu}} | \vec{m} \rangle - \text{c.c.}_A \langle \vec{1} | \sigma_{\vec{\nu}} | \vec{m} \rangle + \text{c.c.}_B, \quad (146) $$
where we used $\vec{1}$ to denote the state with all spins up, “c.c.” stands for complex conjugation and $\sigma_{\vec{\mu}}$, $\sigma_{\vec{\nu}}$ collectively denote the generalized Pauli matrices in the two subsystems A, B separately. By splitting these generalized Pauli generators into tensor products of standard Pauli matrices, we compute the matrix elements as follows:

$$ \langle \vec{1} | \sigma_{\vec{\mu}} | \vec{m} \rangle = \prod_{k=1}^{n_A} (\delta_{m_k, \uparrow}(\delta_{\mu_k, 0} + \delta_{\mu_k, 3}) + \delta_{m_k, \downarrow}(\delta_{\mu_k, 1} - i\delta_{\mu_k, 2})). \quad (147) $$

We now use the above machinery to evaluate

$$ \mathcal{B}\mathcal{F}_{1,\text{hom},\text{Pauli}} = \sum_{\vec{\mu}\vec{\nu}} |Y_{\vec{\mu}\vec{\nu}}| = \frac{1}{2^{n_A+n_B}} \sum_{\vec{m} \neq \vec{1}} |w_m| \sum_{\vec{\mu}\vec{\nu}} |\langle \vec{1} | \sigma_{\vec{\mu}} | \vec{m} \rangle - \text{c.c.}| \times |\langle \vec{1} | \sigma_{\vec{\nu}} | \vec{m} \rangle + \text{c.c.}|. \quad (148) $$

We see that for any given $\vec{m}$, we are left with a combinatorial exercise where we have to sum products similar to those in eq. (147) over all values of $\vec{\mu}$ and $\vec{\nu}$. For the term in the first absolute value not to vanish when subtracting the complex conjugate, we need an odd number of $\mu_k = 2$ components, in this case the “c.c.” simply doubles the contribution. Similarly, for the second absolute value not to vanish, we need an even number of $\nu_k = 2$ components. In this case, the contribution of “c.c.” is again a factor of 2. For specific values of $\vec{\mu}, \vec{\nu}, \vec{m}$, we have at most one out of all the products of delta functions which survives. Therefore, every single contribution to this sum is equal to 4 and we just need to understand how many of those contributions do not vanish for a given value of $\vec{m}$.

What is the number of non-vanishing contributions? Assume to have a certain vector $\vec{m}$ with $m_{\downarrow}$ number of spins down. Notice that the number of spins down is the same for $|m\rangle^A$ and $|m\rangle^B$. For every spin up, we have two options: either $\delta_{\mu_k, 0}$ or $\delta_{\mu_k, 3}$ will contribute for this value of $k$, so in total this yields a factor of $2^{n_A-m_{\downarrow}} \times 2^{n_B-m_{\downarrow}}$ possible contributions. Then we have to choose the location of an odd (even) number of $\delta_{\mu_k, 2}$ in order for the first (second) term not to vanish after subtracting (adding) the complex conjugation. These are chosen out of the down arrows, resulting in a multiplicative factor $\sum_{k(\text{odd/even})}^{m_{\downarrow}} (m_{\downarrow}) = 2^{m_{\downarrow}-1}$ for each system. Multiplying everything together we obtain $4 \cdot 2^{n_A-m_{\downarrow}} \cdot 2^{n_B-m_{\downarrow}} \cdot 2^{m_{\downarrow}-1} \cdot 2^{m_{\downarrow}-1} = 2^{n_A+n_B}$, which precisely cancels the factor of $2^{n_A+n_B}$ in the denominator, and we get the same result (141) as advertised in the main text.

# D $\mathcal{B}\mathcal{F}_2$ homogeneous norm from the Euler-Arnold equations

In this appendix, we compute binding state complexity for a system of two qudits, whose Hilbert space dimensions are $N_A$ and $N_B$, see eq. (12). For convenience, we assume that $N_A \le N_B$. The complexity is computed according to the $F_{2,\text{hom}}$ norm defined in eq. (20a) using an arbitrary basis of the unitary algebra. The following calculation is based on the Euler-Arnold equations and provides a cross-check of the result (97).

## D.1 Geodesics from the Euler-Arnold equations

To compare more easily with the literature, we use the notation of [6]. Given the space of unitary operators $SU(N \equiv N_A N_B)$ acting on the qudits, we define the length $\ell$ of a trajectory $U(t)$ (2) generated by a Hamiltonian $H(t)$ as

$$ \ell \equiv \int_{0}^{1} dt \sqrt{\langle H(t), H(t) \rangle}, \quad (149) $$

where $\langle \cdot, \cdot \rangle$ is a suitable Riemannian metric independent of the group point $U$. Explicitly, we define the norm at the origin as

$$ \langle H, K \rangle \equiv \frac{1}{N} \mathrm{Tr}(H \mathcal{G}(K)). \quad (150) $$
Here $H$ and $K$ are traceless Hermitian Hamiltonians, $\mathcal{G}$ is a positive-definite linear operator acting on the algebra $\mathfrak{su}(N)$ as $\mathcal{G}(T_I) = q_I T_I$, where $T_I$ are the generators of the algebra and $q_I$ the corresponding penalty factors, and $\mathcal{N}$ is a normalization constant chosen such that $\mathrm{Tr}(T^I T^J) = \mathcal{N} \delta^{IJ}$. Eq. (149) is equivalent to the definition of the $F_2$ norm according to eq. (4). Unitary complexity is given by the length (149) of the geodesic $U(t)$ with boundary conditions $U(0) = \mathbf{1}$ and $U(t) = U_T$, the latter being a target unitary. This corresponds to the definition (3).

We now consider a bipartition of the two-qudit system and decompose the Hamiltonian as in (14). We define an $F_2$ norm corresponding to the following choice of $\mathcal{G}$

$$ \mathcal{G}_2(H) = \varepsilon S + Q, \qquad (151) $$

where $S$ and $Q$ represent respectively the local and non-local parts of $H$:

$$ S = Y_A^a T_a^A \otimes \mathbf{1}^B + Y_B^i T_a^A \otimes T_i^B, \quad Q = Y^{ai} T_a^A \otimes T_i^B, \quad H(t) = S(t) + Q(t). \qquad (152) $$

The geodesic equations in the space of unitaries with inner product (150) are the Euler-Arnold equations [90]$^{42}$

$$ \dot{H} = -i\mathcal{G}^{-1}([H, \mathcal{G}(H)]). \qquad (153) $$

When working with the operator $\mathcal{G}_2$, we can use the commutators $[S, S] = [Q, Q] = 0$ and the fact that $[S, Q]$ is non-local to show that the Euler-Arnold equations take the form

$$ \dot{S} = 0, \qquad \dot{Q} = -i(1 - \varepsilon)[S, Q]. \qquad (154) $$

The $F_{2,\mathrm{hom}}$ norm in eq. (20a) corresponds to (150) with the operator $\mathcal{G}_{2,\mathrm{hom}} \equiv \lim_{\varepsilon \to 0} \mathcal{G}_2$. The solution for the Hamiltonian in this case reads

$$ H(t) = \lim_{\varepsilon \to 0} (S_0 + e^{-i(1-\varepsilon)S_0t} Q_0 e^{i(1-\varepsilon)S_0t}) = S_0 + e^{-iS_0t} Q_0 e^{iS_0t}, \qquad (155) $$

where $S_0$ and $Q_0$ are the values of the local and non-local parts of the Hamiltonian at $t=0$, respectively. By solving the Schrödinger equation $\dot{U} = -iHU$ for the unitary along the curve, we get

$$ U(t) = \lim_{\varepsilon \to 0} e^{-i(1-\varepsilon)S_0t} e^{-i(1-(1-\varepsilon))S_0t - iQ_0t} = e^{-iS_0t} e^{-iQ_0t}. \qquad (156) $$

The length of this geodesic in the space of unitaries is

$$ \ell_{2,\mathrm{hom}} = \int_0^1 dt \lim_{\varepsilon \to 0} \sqrt{\frac{\varepsilon}{\mathcal{N}} \mathrm{Tr}(S(t)S(t)) + \frac{1}{\mathcal{N}} \mathrm{Tr}(Q(t)Q(t))} = \sqrt{\mathrm{Tr}\left(\frac{Q_0Q_0}{\mathcal{N}}\right)}, \qquad (157) $$

where we used the result (155) and the time-independence of $Q_0$ to directly perform the trivial integration. Binding unitary complexity $\mathcal{BC}_{2,\mathrm{hom}}$ is defined as the minimal length (157) of the global geodesic connecting the identity with the target unitary $U_T$. We immediately observe that this corresponds to minimizing $\mathrm{Tr}(Q_0Q_0)$, independently of the initial value $S_0$ of the local part of the Hamiltonian.

<footnote>$^{42}$Euler-Arnold equations have been recently used to study complexity in the space of unitaries, see e.g., [63, 68, 69, 71, 73].</footnote>
## D.2 Computation of binding state complexity

Let us now move to the computation of binding state complexity. First of all, we specify the reference and target states in the Hilbert space (12) of the two qudits

$$ egin{aligned} |\psi_R\rangle &= |1\tilde{1}\rangle, \\ |\psi_T\rangle &= \sum_{m=1}^{N_A} \bar{\lambda}_m |\bar{\phi}_m \bar{\chi}_m\rangle = U_A \otimes U_B \sum_{m=1}^{N_A} \bar{\lambda}_m |m\tilde{m}\rangle, \quad \bar{\lambda}_1 \ge \bar{\lambda}_2 \ge \dots \ge \bar{\lambda}_{N_A}, \end{aligned} \quad (158) $$

where $|m\tilde{m}\rangle$ is a fixed basis in the full Hilbert space. The state at a generic time along a geodesic trajectory is given by

$$ |\psi(t)\rangle = e^{-iS_0 t} e^{-iQ_0 t} |1\tilde{1}\rangle. \quad (159) $$

The globally minimal geodesic is obtained by minimizing eq. (157) over all the $Q_0$ that bring the reference state to a state with the same Schmidt coefficients $\bar{\lambda}_m$ as specified in eq. (158). Since $S_0$ can be chosen without affecting the length of the trajectory, we are free to act with LU operators (which modify the basis of the Schmidt decomposition) until we precisely reach the target state.

Let us further comment that there are Hamiltonians relating the reference and the target state continuously, i.e., without a final LU transformation. If the target state $|\psi_T\rangle$ and the final state of the geodesic trajectory $|\tilde{\psi}_T\rangle = e^{-iHt}|\psi_R\rangle|_{t=1}$ are related by $|\psi_T\rangle = e^{-iS_0}|\tilde{\psi}_T\rangle$ where $S_0$ is a time-independent generator of LU transformations, then the unitary $U(t)$ such that

$$ U(t) = e^{-iS_0 t} e^{-iQ_0 t}, \quad (160) $$

will have the same binding complexity as $e^{-iQ_0 t}$. In summary, we have shown by a different method than the one used in section 3.1 that binding state complexity is defined in terms of equivalence classes of LUE reference and target states, see eq. (23).

We are now going to show that the length of the geodesic in the space of states can be expressed in terms of the Fubini-Study distance on the projective space $\mathbb{CP}^{N-1}$. The Fubini-Study distance is proportional to the state complexity $C_{\text{IP}}^{\text{state}}$, defined according to eq. (6), and associated with the standard inner-product norm on the unitary group $G = U(N)$, which assigns penalties 1 to all the generators (either local or non-local), i.e.,

$$ egin{aligned} C_{\text{IP}}^{\text{state}}(|\psi_T\rangle, |\psi_R\rangle) &\equiv \min_{\{U \in \text{SU}(N): |\psi_T\rangle = U |\psi_R\rangle\}} C_{\text{IP}}[U], \\ F_{\text{IP}}[\vec{Y}] &= \sqrt{\frac{\text{Tr}(H^2)}{N}} = \sqrt{\sum_a |Y_A^a|^2 + \sum_i |Y_B^i|^2 + \sum_{a,i} |Y^{ai}|^2}, \end{aligned} \quad (161) $$

where IP stands for inner product. The norm on the projective space $\mathbb{CP}^{N-1}$ obtained after the minimization over the stabilizer as in eq. (9) is (proportional to) the Fubini-Study metric [123, 124]$^{43}$

$$ F_{\text{IP}}^{\text{state}} dt = \sqrt{\frac{2}{N}} ds_{\text{FS}}^2 \quad (162) $$

where

$$ ds_{\text{FS}}^2 = \langle d\psi|d\psi\rangle - \langle d\psi|\psi\rangle \langle \psi|d\psi\rangle, \quad (163) $$

<footnote>$^{43}$The minimization can be done by expressing the elements of $H$ in a Hilbert space basis for which the first and second basis vectors are $|\psi_1\rangle = |\psi(t)\rangle$ and $|\psi_2\rangle = [|\dot{\psi}(t)\rangle - \langle\psi(t)|\dot{\psi}(t)\rangle|\psi(t)\rangle] / \sqrt{ds_{\text{FS}}^2}$. In this basis, $H_{1n}$ and $H_{n1}$ for any $n$ are fixed by the Schrödinger equation, and the rest of the elements belong to the stabilizer over which we minimize. A further minimization over the time derivative of the global phase of the state $|\psi(t)\rangle$ will result in (162).</footnote>
and the associated distance reads

$$d_{\text{FS}}(|\psi_T\rangle, |\psi_R\rangle) = \arccos |\langle\psi_T|\psi_R\rangle|. \qquad (164)$$

From this, we conclude that

$$C_{\text{IP}}^{\text{state}} = \sqrt{\frac{2}{N}} d_{\text{FS}}. \qquad (165)$$

Now the main idea of this section is to define a metric over equivalence classes of states by minimizing the above metric over LUE reference and target states:

$$C_{2,LUE}^{\text{state}}(\{|\psi_T\rangle\}, \{|\psi_R\rangle\}) = \min_{\text{LUE}} C_{\text{IP}}^{\text{state}}(|\psi_T\rangle, |\psi_R\rangle) \qquad (166)$$

where the curly brackets denote equivalence classes of states. Let us now explain how all this relates to the binding complexity defined in the main text, where zero cost is associated with local generators. A trivial observation is that cost functions with bigger penalty factors imply a bigger complexity, leading to $\mathcal{BC}_{2,hom}^{\text{state}} \le C_{\text{IP}}^{\text{state}}$. Next, by minimizing this inequality over LUE states, we find that $\mathcal{BC}_{2,hom}^{\text{state}} \le C_{2,LUE}^{\text{state}}$, because the value of $\mathcal{BC}_{2,hom}^{\text{state}}$ does not change within the LUE class of states. On the other hand, taking the Hamiltonian of the optimal circuit for the $\mathcal{BC}_{2,hom}^{\text{state}}$ and setting $S_0 = 0$ will result in a different circuit connecting the reference state to a target state in the same LUE class. This circuit has the same cost according to either the inner-product metric, or the distance defined by eq. (157), therefore

$$\mathcal{BC}_{2,hom}^{\text{state}} = C_{2,LUE}^{\text{state}}. \qquad (167)$$

Let us now exploit this result to find the binding complexity by actually computing the geodesic length using the norm (166). The state complexity $C_{2,LUE}^{\text{state}}$ is nothing but the length of a global geodesic in the geometry (162), minimized over all LUE target states. We then find

$$C_{2,LUE}^{\text{state}} = \sqrt{\frac{2}{N}} \min_{\{|\psi_T\rangle\}, \{|\psi_R\rangle\}} d_{\text{FS}}(|\psi_T\rangle, |\psi_R\rangle), \qquad (168)$$

where $d_{\text{FS}}$ denotes the Fubini-Study distance between a given representative of the equivalence classes of the reference and target states. Recall that the Fubini-Study metric is left-invariant, and that we are minimizing modulo LUE states. For these reasons, we can perform arbitrary local unitary transformations on the reference state

$$U_A \otimes U_B |1\bar{1}\rangle = |\phi\chi\rangle, \qquad (169)$$

and equivalently express the complexity as

$$C_{2,LUE}^{\text{state}} = \sqrt{\frac{2}{N}} \min_{\phi, \chi} \arccos \left( \sum_{m=1}^{N_A} \langle \phi \chi | \bar{\lambda}_m | m \tilde{m} \rangle \right), \qquad (170)$$

where the definition (158) has been used. By performing a change of basis $|\chi\rangle = \sum_m a_m^* |m\rangle$, such that $\sum_m |a_m|^2 = 1$, we find the chain of relations

$$|\langle\phi\chi|\bar{\lambda}_m|m\tilde{m}\rangle| = |\langle\phi|a_m\bar{\lambda}_m|m\rangle| \le |\langle\phi|\phi\rangle| |\langle n|a_m a_n^*\bar{\lambda}_m\bar{\lambda}_n|m\rangle| = \sqrt{|a_m\bar{\lambda}_m|^2} \le \bar{\lambda}_1. \qquad (171)$$

In the second step we used the Cauchy-Schwarz inequality, while in the last step we used the fact that $\bar{\lambda}_1$ is the largest Schmidt coefficient. Since the arccos is a monotonically decreasing function of its argument, we conclude that

$$\arccos \left( \sum_{m=1}^{N_A} \langle \phi \chi | \bar{\lambda}_m | m \tilde{m} \rangle \right) \ge \arccos(\bar{\lambda}_1). \qquad (172)$$
The minimum is actually achieved by constructing a unitary (2) with a Hamiltonian which has the following non-local part at $t = 0$

$$Q_0 = -i \frac{\arccos(\bar{\lambda}_1)}{\sqrt{1 - \bar{\lambda}_1^2}} \left( \sum_{m=1}^{N_A} \bar{\lambda}_m (|m\tilde{m}\rangle \langle 1\tilde{1}| - |1\tilde{1}\rangle \langle m\tilde{m}|) \right). \quad (173)$$

Plugging the explicit form of $Q_0$ inside eq. (157), we obtain that the length of this geodesic is indeed proportional to the Fubini-Study distance, with the proportionality factor determined by eq. (165). Using $N = N_A N_B$,<fnref n="44" /> this implies that the binding state complexity reads

$$\mathcal{BC}_{2,\text{hom}}^{\text{state}} = C_{2,\text{IP}}^{\text{state}} = \sqrt{\frac{2}{N_A N_B}} \arccos(\bar{\lambda}_1), \quad (174)$$

which matches the result (97).

To conclude, we compare this procedure, which uses the Euler-Arnold equations, to the $SO(N_A)$ method presented in section 3.

**Advantages:**

* The use of Euler-Arnold equations may be more familiar to the reader.

* The trajectory in the unitary space is simple and smooth, see eq. (159).

**Disadvantages:**

* The reduction of the problem to the space of Schmidt coefficients is not manifest until the end of the computation.

* This method only applies to the $F_2$ norm.

## References

[1] John Watrous. Quantum Computational Complexity. 4 2008. DOI: 10.1007/978-0-387-30440-3_428.

[2] Scott Aaronson. The Complexity of Quantum States and Transformations: From Quantum Money to Black Holes. 7 2016. DOI: 10.48550/arXiv.1607.05256.

[3] Vijay Balasubramanian, Matthew DeCross, Arjun Kar, and Onkar Parrikar. Binding Complexity and Multiparty Entanglement. JHEP, 02:069, 2019. DOI: 10.1007/JHEP02(2019)069.

[4] Michael A. Nielsen and Isaac L. Chuang. Quantum Computation and Quantum Information: 10th Anniversary Edition. Cambridge University Press, 2010. DOI: 10.1017/CBO9780511976667.

[5] Michael A. Nielsen. A geometric approach to quantum circuit lower bounds. Quant. Inf. Comput., 6(3):213-262, 2006. DOI: 10.26421/QIC6.3-2.

[6] Mark R. Dowling and Michael A. Nielsen. The geometry of quantum computation. Quant. Inf. Comput., 8(10):0861-0899, 2008. DOI: 10.26421/QIC8.10-1.

[7] Michael A. Nielsen, Mark R. Dowling, Mile Gu, and Andrew C. Doherty. Quantum computation as geometry. Science, 311(5764):1133-1135, feb 2006. DOI: 10.1126/science.1121541. URL https://doi.org/10.1126/2Fscience.1121541.

<footnote n="44">⁴⁴Here, $N_A$ and $N_B$ are the normalization constants of the traces of pairs of generators in $SU(N_A)$ and $SU(N_B)$ describing each of the qudits, respectively, as defined in the text below eq. (150).</footnote>
[8] Alberto Carlini, Akio Hosoya, Tatsuhiko Koike, and Yosuke Okudaira. Time-optimal quantum evolution. *Phys. Rev. Lett.*, 96:060503, Feb 2006. DOI: 10.1103/PhysRevLett.96.060503.

[9] J Werschnik and EKU Gross. Quantum optimal control theory. *Journal of Physics B: Atomic, Molecular and Optical Physics*, 40(18):R175, 2007. DOI: 10.1088/0953-4075/40/18/R01.

[10] Christiane P. Koch et al. Quantum optimal control in quantum technologies. Strategic report on current status, visions and goals for research in Europe. *EPJ Quant. Technol.*, 9(1):19, 2022. DOI: 10.1140/epjqt/s40507-022-00138-x.

[11] Thomas Hartman and Juan Maldacena. Time Evolution of Entanglement Entropy from Black Hole Interiors. *JHEP*, 05:014, 2013. DOI: 10.1007/JHEP05(2013)014.

[12] Leonard Susskind. Entanglement is not enough. *Fortsch. Phys.*, 64:49-71, 2016. DOI: 10.1002/prop.201500095.

[13] Ben Freivogel, Ro Jefferson, Laurens Kabir, Benjamin Mosk, and I-Sheng Yang. Casting Shadows on Holographic Reconstruction. *Phys. Rev. D*, 91(8):086013, 2015. DOI: 10.1103/PhysRevD.91.086013.

[14] Vijay Balasubramanian, Borun D. Chowdhury, Bartlomiej Czech, and Jan de Boer. Entwinement and the emergence of spacetime. *JHEP*, 01:048, 2015. DOI: 10.1007/JHEP01(2015)048.

[15] Leonard Susskind. Computational Complexity and Black Hole Horizons. *Fortsch. Phys.*, 64:24-43, 2016. DOI: 10.1002/prop.201500092. [Addendum: Fortsch.Phys. 64, 44-48 (2016), DOI: 10.1002/prop.201500093]".

[16] Douglas Stanford and Leonard Susskind. Complexity and Shock Wave Geometries. *Phys. Rev. D*, 90(12):126007, 2014. DOI: 10.1103/PhysRevD.90.126007.

[17] Alexandre Belin, Aitor Lewkowycz, and Gábor Sárosi. Complexity and the bulk volume, a new York time story. *JHEP*, 03:044, 2019. DOI: 10.1007/JHEP03(2019)044.

[18] E. Rabinovici, A. Sánchez-Garrido, R. Shir, and J. Sonner. A bulk manifestation of Krylov complexity. *JHEP*, 08:213, 2023. DOI: 10.1007/JHEP08(2023)213.

[19] Adam R. Brown, Daniel A. Roberts, Leonard Susskind, Brian Swingle, and Ying Zhao. Complexity, action, and black holes. *Phys. Rev. D*, 93(8):086006, 2016. DOI: 10.1103/PhysRevD.93.086006.

[20] Adam R. Brown, Daniel A. Roberts, Leonard Susskind, Brian Swingle, and Ying Zhao. Holographic Complexity Equals Bulk Action? *Phys. Rev. Lett.*, 116(19): 191301, 2016. DOI: 10.1103/PhysRevLett.116.191301.

[21] Josiah Couch, Willy Fischler, and Phuc H. Nguyen. Noether charge, black hole volume, and complexity. *JHEP*, 03:119, 2017. DOI: 10.1007/JHEP03(2017)119.

[22] Alexandre Belin, Robert C. Myers, Shan-Ming Ruan, Gábor Sárosi, and Antony J. Speranza. Does Complexity Equal Anything? *Phys. Rev. Lett.*, 128(8):081602, 2022. DOI: 10.1103/PhysRevLett.128.081602.

[23] Alexandre Belin, Robert C. Myers, Shan-Ming Ruan, Gábor Sárosi, and Antony J. Speranza. Complexity equals anything II. *JHEP*, 01:154, 2023. DOI: 10.1007/JHEP01(2023)154.

[24] Eivind Jørstad, Robert C. Myers, and Shan-Ming Ruan. Complexity=anything: singularity probes. *JHEP*, 07:223, 2023. DOI: 10.1007/JHEP07(2023)223.

[25] Johanna Erdmenger, Mario Flory, Marius Gerbershagen, Michal P. Heller, and Anna-Lena Weigel. Exact Gravity Duals for Simple Quantum Circuits. *SciPost Phys.*, 13(3):061, 2022. DOI: 10.21468/SciPostPhys.13.3.061.

[26] Johanna Erdmenger, Anna-Lena Weigel, Marius Gerbershagen, and Michal P. Heller.
From complexity geometry to holographic spacetime. Phys. Rev. D, 108(10):106020, 2023. DOI: 10.1103/PhysRevD.108.106020.

[27] A. Ramesh Chandra, Jan de Boer, Mario Flory, Michal P. Heller, Sergio Hörtner, and Andrew Rolph. Cost of holographic path integrals. SciPost Phys., 14(4):061, 2023. DOI: 10.21468/SciPostPhys.14.4.061.

[28] Luis Lehner, Robert C. Myers, Eric Poisson, and Rafael D. Sorkin. Gravitational action with null boundaries. Phys. Rev. D, 94(8):084046, 2016. DOI: 10.1103/Phys-RevD.94.084046.

[29] Jose L. F. Barbon and Eliezer Rabinovici. Holographic complexity and spacetime singularities. JHEP, 01:084, 2016. DOI: 10.1007/JHEP01(2016)084.

[30] Shira Chapman, Hugo Marrochio, and Robert C. Myers. Complexity of Formation in Holography. JHEP, 01:062, 2017. DOI: 10.1007/JHEP01(2017)062.

[31] Rong-Gen Cai, Shan-Ming Ruan, Shao-Jiang Wang, Run-Qiu Yang, and Rong-Hui Peng. Action growth for AdS black holes. JHEP, 09:161, 2016. DOI: 10.1007/JHEP09(2016)161.

[32] Dean Carmi, Shira Chapman, Hugo Marrochio, Robert C. Myers, and Sotaro Sugishita. On the Time Dependence of Holographic Complexity. JHEP, 11:188, 2017. DOI: 10.1007/JHEP11(2017)188.

[33] Shira Chapman, Hugo Marrochio, and Robert C. Myers. Holographic complexity in Vaidya spacetimes. Part I. JHEP, 06:046, 2018. DOI: 10.1007/JHEP06(2018)046.

[34] Shira Chapman, Hugo Marrochio, and Robert C. Myers. Holographic complexity in Vaidya spacetimes. Part II. JHEP, 06:114, 2018. DOI: 10.1007/JHEP06(2018)114.

[35] Shira Chapman, Dongsheng Ge, and Giuseppe Policastro. Holographic Complexity for Defects Distinguishes Action from Volume. JHEP, 05:049, 2019. DOI: 10.1007/JHEP05(2019)049.

[36] Paolo Braccia, Aldo L. Cotrone, and Erik Tonni. Complexity in the presence of a boundary. JHEP, 02:051, 2020. DOI: 10.1007/JHEP02(2020)051.

[37] Yoshiki Sato and Kento Watanabe. Does Boundary Distinguish Complexities? JHEP, 11:132, 2019. DOI: 10.1007/JHEP11(2019)132.

[38] Alice Bernamonti, Federico Galli, Juan Hernandez, Robert C. Myers, Shan-Ming Ruan, and Joan Simón. First Law of Holographic Complexity. Phys. Rev. Lett., 123 (8):081601, 2019. DOI: 10.1103/PhysRevLett.123.081601.

[39] Alice Bernamonti, Federico Galli, Juan Hernandez, Robert C. Myers, Shan-Ming Ruan, and Joan Simón. Aspects of The First Law of Complexity. J. Phys. A, 53:29, 2020. DOI: 10.1088/1751-8121/ab8e66.

[40] Shira Chapman, Damián A. Galante, and Eric David Kramer. Holographic complexity and de Sitter space. JHEP, 02:198, 2022. DOI: 10.1007/JHEP02(2022)198.

[41] Roberto Auzzi, Stefano Baiguera, Sara Bonansea, Giuseppe Nardelli, and Kristian Toccacelo. Volume complexity for Janus AdS₃ geometries. JHEP, 08:045, 2021. DOI: 10.1007/JHEP08(2021)045.

[42] Stefano Baiguera, Sara Bonansea, and Kristian Toccacelo. Volume complexity for the nonsupersymmetric Janus AdS5 geometry. Phys. Rev. D, 104(8):086030, 2021. DOI: 10.1103/PhysRevD.104.086030.

[43] Roberto Auzzi, Stefano Baiguera, Sara Bonansea, and Giuseppe Nardelli. Action complexity in the presence of defects and boundaries. JHEP, 02:118, 2022. DOI: 10.1007/JHEP02(2022)118.

[44] Roberto Emparan, Antonia Micol Frassino, Martin Sasieta, and Marija Tomašević. Holographic complexity of quantum black holes. JHEP, 02:204, 2022. DOI: 10.1007/JHEP02(2022)204.
[45] Eivind Jørstad, Robert C. Myers, and Shan-Ming Ruan. Holographic complexity in $dS_{d+1}$. JHEP, 05:119, 2022. DOI: 10.1007/JHEP05(2022)119.

[46] Roberto Auzzi, Giuseppe Nardelli, Gabriel Pedde Ungureanu, and Nicolo Zenoni. Volume complexity of dS bubbles. Phys. Rev. D, 108(2):026006, 2023. DOI: 10.1103/PhysRevD.108.026006.

[47] Takanori Anegawa, Norihiro Iizuka, Sunil Kumar Sake, and Nicolò Zenoni. Is action complexity better for de Sitter space in Jackiw-Teitelboim gravity? JHEP, 06:213, 2023. DOI: 10.1007/JHEP06(2023)213.

[48] Takanori Anegawa and Norihiro Iizuka. Shock waves and delay of hyperfast growth in de Sitter complexity. JHEP, 08:115, 2023. DOI: 10.1007/JHEP08(2023)115.

[49] Stefano Baiguera, Rotem Berman, Shira Chapman, and Robert C. Myers. The cosmological switchback effect. JHEP, 07:162, 2023. DOI: 10.1007/JHEP07(2023)162.

[50] Sergio E. Aguilar-Gutierrez, Michal P. Heller, and Silke Van der Schueren. Complexity = Anything Can Grow Forever in de Sitter. 5 2023. DOI: 10.48550/arXiv.2305.11280.

[51] Sergio E. Aguilar-Gutierrez. C=Anything and the switchback effect in Schwarzschild-de Sitter space. JHEP, 03:062, 2024. DOI: 10.1007/JHEP03(2024)062.

[52] Amin Akhavan and Farzad Omidi. On the Role of Counterterms in Holographic Complexity. JHEP, 11:054, 2019. DOI: 10.1007/JHEP11(2019)054.

[53] Farzad Omidi. Regularizations of Action-Complexity for a Pure BTZ Black Hole Microstate. JHEP, 07:020, 2020. DOI: 10.1007/JHEP07(2020)020.

[54] Pawel Caputa, Nilay Kundu, Masamichi Miyaji, Tadashi Takayanagi, and Kento Watanabe. Liouville Action as Path-Integral Complexity: From Continuous Tensor Networks to AdS/CFT. JHEP, 11:097, 2017. DOI: 10.1007/JHEP11(2017)097.

[55] Ro Jefferson and Robert C. Myers. Circuit complexity in quantum field theory. JHEP, 10:107, 2017. DOI: 10.1007/JHEP10(2017)107.

[56] Shira Chapman, Michal P. Heller, Hugo Marrochio, and Fernando Pastawski. Toward a Definition of Complexity for Quantum Field Theory States. Phys. Rev. Lett., 120 (12):121602, 2018. DOI: 10.1103/PhysRevLett.120.121602.

[57] Rifath Khan, Chethan Krishnan, and Sanchita Sharma. Circuit Complexity in Fermionic Field Theory. Phys. Rev. D, 98(12):126001, 2018. DOI: 10.1103/Phys-RevD.98.126001.

[58] Arpan Bhattacharyya, Pawel Caputa, Sumit R. Das, Nilay Kundu, Masamichi Miyaji, and Tadashi Takayanagi. Path-Integral Complexity for Perturbed CFTs. JHEP, 07:086, 2018. DOI: 10.1007/JHEP07(2018)086.

[59] Shira Chapman, Jens Eisert, Lucas Hackl, Michal P. Heller, Ro Jefferson, Hugo Marrochio, and Robert C. Myers. Complexity and entanglement for thermofield double states. SciPost Phys., 6(3):034, 2019. DOI: 10.21468/SciPostPhys.6.3.034.

[60] Hugo A. Camargo, Pawel Caputa, Diptarka Das, Michal P. Heller, and Ro Jefferson. Complexity as a novel probe of quantum quenches: universal scalings and purifications. Phys. Rev. Lett., 122(8):081601, 2019. DOI: 10.1103/PhysRevLett.122.081601.

[61] Dongsheng Ge and Giuseppe Policastro. Circuit Complexity and 2D Bosonisation. JHEP, 10:276, 2019. DOI: 10.1007/JHEP10(2019)276.

[62] Adam R. Brown and Leonard Susskind. Complexity geometry of a single qubit. Phys. Rev. D, 100(4):046020, 2019. DOI: 10.1103/PhysRevD.100.046020.

[63] Vijay Balasubramanian, Matthew Decross, Arjun Kar, and Onkar Parrikar. Quantum Complexity of Time Evolution with Chaotic Hamiltonians. JHEP, 01:134, 2020. DOI: 10.1007/JHEP01(2020)134.

[64] Shira Chapman and Hong Zhe Chen. Charged Complexity and the Thermofield Double State. JHEP, 02:187, 2021. DOI: 10.1007/JHEP02(2021)187.
[65] Pawel Caputa and Javier M. Magan. Quantum Computation as Gravity. *Phys. Rev. Lett.*, 122(23):231302, 2019. DOI: 10.1103/PhysRevLett.122.231302.

[66] Roberto Auzzi, Stefano Baiguera, G. Bruno De Luca, Andrea Legramandi, Giuseppe Nardelli, and Nicolò Zenoni. Geometry of quantum complexity. *Phys. Rev. D*, 103 (10):106021, 2021. DOI: 10.1103/PhysRevD.103.106021.

[67] Reginald J. Caginalp and Samuel Leutheusser. Complexity in One- and Two-Qubit Systems. 10 2020. DOI: 10.48550/arXiv.2010.15099.

[68] Mario Flory and Michal P. Heller. Geometry of Complexity in Conformal Field Theory. *Phys. Rev. Res.*, 2(4):043438, 2020. DOI: 10.1103/PhysRevResearch.2.043438.

[69] Mario Flory and Michal P. Heller. Conformal field theory complexity from Euler-Arnold equations. *JHEP*, 12:091, 2020. DOI: 10.1007/JHEP12(2020)091.

[70] Nicolas Chagnet, Shira Chapman, Jan de Boer, and Claire Zukowski. Complexity for Conformal Field Theories in General Dimensions. *Phys. Rev. Lett.*, 128(5):051601, 2022. DOI: 10.1103/PhysRevLett.128.051601.

[71] Pablo Basteiro, Johanna Erdmenger, Pascal Fries, Florian Goth, Ioannis Matthaiakakis, and René Meyer. Quantum complexity as hydrodynamics. *Phys. Rev. D*, 106(6):065016, 2022. DOI: 10.1103/PhysRevD.106.065016.

[72] Adam R. Brown, Michael H. Freedman, Henry W. Lin, and Leonard Susskind. Universality in long-distance geometry and quantum complexity. *Nature*, 622(7981): 58–62, 2023. DOI: 10.1038/s41586-023-06460-3.

[73] Vijay Balasubramanian, Matthew DeCross, Arjun Kar, Yue (Cathy) Li, and Onkar Parrikar. Complexity growth in integrable and chaotic models. *JHEP*, 07:011, 2021. DOI: 10.1007/JHEP07(2021)011.

[74] Adam R. Brown. Polynomial Equivalence of Complexity Geometries. *Quantum*, 8: 1391, 2024. DOI: 10.22331/q-2024-07-02-1391.

[75] Ben Craps, Marine De Clerck, Oleg Evnin, and Philip Hacker. Integrability and complexity in quantum spin chains. *SciPost Phys.*, 16:041, 2024. DOI: 10.21468/SciPostPhys.16.2.041. URL https://scipost.org/10.21468/SciPostPhys.16.2.041.

[76] Ben Craps, Marine De Clerck, Oleg Evnin, Philip Hacker, and Maxim Pavlov. Bounds on quantum evolution complexity via lattice cryptography. *SciPost Phys.*, 13(4):090, 2022. DOI: 10.21468/SciPostPhys.13.4.090.

[77] Shira Chapman and Giuseppe Policastro. Quantum computational complexity from quantum information to black holes and back. *Eur. Phys. J. C*, 82(2):128, 2022. DOI: 10.1140/epjc/s10052-022-10037-1.

[78] Rodney Van Meter, W. J. Munro, Kae Nemoto, and Kohei M. Itoh. Arithmetic on a distributed-memory quantum multicomputer. *ACM Journal on Emerging Technologies in Computing Systems*, 3(4):1–23, jan 2008. DOI: 10.1145/1324177.1324179. URL https://doi.org/10.1145/2F1324177.1324179.

[79] Robert Beals, Stephen Brierley, Oliver Gray, Aram W. Harrow, Samuel Kutin, Noah Linden, Dan Shepherd, and Mark Stather. Efficient distributed quantum computing. *Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences*, 469(2153):20120686, may 2013. DOI: 10.1098/rspa.2012.0686. URL https://doi.org/10.1098/rspa.2012.0686.

[80] Marcello Caleffi, Michele Amoretti, Davide Ferrari, Daniele Cuomo, Jessica Illiano, Antonio Manzalini, and Angela Sara Cacciapuoti. Distributed Quantum Computing: a Survey. 12 2022. DOI: 10.48550/arXiv.2212.10609.

[81] Harry Buhrman, Richard Cleve, Serge Massar, and Ronald de Wolf. Nonlocality and communication complexity. *Reviews of Modern Physics*, 82(1):665–698, mar 2010. DOI: 10.1103/revmodphys.82.665. URL https://doi.org/10.1103/2Frevmodphys.82.665.
[82] Łukasz Rudnicki. Quantum speed limit and geometric measure of entanglement. *Physical Review A*, 104(3), sep 2021. DOI: 10.1103/physreva.104.032417. URL https://doi.org/10.1103/2Fphysreva.104.032417.

[83] Hamed Zolfi. Complexity and Multi-boundary Wormholes in 2 + 1 dimensions. *JHEP*, 04:076, 2023. DOI: 10.1007/JHEP04(2023)076.

[84] Yuxuan Zhang. Straddling-gates problem in multipartite quantum systems. *Phys. Rev. A*, 105(6):062430, 2022. DOI: 10.1103/PhysRevA.105.062430.

[85] J. Eisert. Entangling Power and Quantum Circuit Complexity. *Phys. Rev. Lett.*, 127 (2):020501, 2021. DOI: 10.1103/PhysRevLett.127.020501.

[86] Adam R. Brown and Leonard Susskind. Second law of quantum complexity. *Phys. Rev. D*, 97(8):086015, 2018. DOI: 10.1103/PhysRevD.97.086015.

[87] Adam R. Brown, Leonard Susskind, and Ying Zhao. Quantum Complexity and Negative Curvature. *Phys. Rev. D*, 95(4):045010, 2017. DOI: 10.1103/Phys-RevD.95.045010.

[88] Ben Craps, Oleg Evnin, and Gabriele Pascuzzi. A Relation between Krylov and Nielsen Complexity. *Phys. Rev. Lett.*, 132(16):160402, 2024. DOI: 10.1103/Phys-RevLett.132.160402.

[89] Xiaoting Wang, Michele Allegra, Kurt Jacobs, Seth Lloyd, Cosmo Lupo, and Masoud Mohseni. Quantum brachistochrone curves as geodesics: Obtaining accurate minimum-time protocols for the control of quantum systems. *Physical Review Letters*, 114(17), apr 2015. DOI: 10.1103/physrevlett.114.170501. URL https://doi.org/10.1103/2Fphysrevlett.114.170501.

[90] Vladimir Arnold. Sur la géométrie différentielle des groupes de Lie de dimension infinie et ses applications à l'hydrodynamique des fluides parfaits. *Annales de l'Institut Fourier*, 16(1):319-361, 1966. DOI: 10.5802/aif.233. URL http://www.numdam.org/articles/10.5802/aif.233/.

[91] Guifre Vidal. On the characterization of entanglement. *J. Mod. Opt.*, 47:355, 2000. DOI: 10.1080/09500340008244048.

[92] W. Dur, G. Vidal, J. I. Cirac, N. Linden, and S. Popescu. Entanglement capabilities of nonlocal Hamiltonians. *Phys. Rev. Lett.*, 87:137901, 2001. DOI: 10.1103/Phys-RevLett.87.137901.

[93] Yuchen Wang, Zixuan Hu, Barry C. Sanders, and Sabre Kais. Qudits and high-dimensional quantum computing. *Frontiers in Physics*, 8, nov 2020. DOI: 10.3389/fphy.2020.589504. URL https://doi.org/10.3389/2Fphy.2020.589504.

[94] Minyong Guo, Juan Hernandez, Robert C. Myers, and Shan-Ming Ruan. Circuit Complexity for Coherent States. *JHEP*, 10:011, 2018. DOI: 10.1007/JHEP10(2018)011.

[95] Sergey Bravyi. Upper bounds on entangling rates of bipartite hamiltonians. *Phys. Rev. A*, 76:052319, Nov 2007. DOI: 10.1103/PhysRevA.76.052319. URL https://link.aps.org/doi/10.1103/PhysRevA.76.052319.

[96] Karel Van Acoleyen, Michaël Mariën, and Frank Verstraete. Entanglement rates and area laws. *Physical Review Letters*, 111(17), oct 2013. DOI: 10.1103/phys-revlett.111.170501. URL https://doi.org/10.1103/2Fphysrevlett.111.170501.

[97] Michaël Mariën, Koenraad M. R. Audenaert, Karel Van Acoleyen, and Frank Verstraete. Entanglement rates and the stability of the area law for the entanglement entropy. *Communications in Mathematical Physics*, 346(1):35-73, 2016. DOI: 10.1007/s00220-016-2709-5.

[98] Tzu-Chieh Wei and Paul M. Goldbart. Geometric measure of entanglement and applications to bipartite and multipartite quantum states. *Physical Review A*, 68