QuantumCore

Lesson 44 · 100 min · advanced

Hardware-Efficient Ansatz

Loading progress…

The hardware-efficient ansatz (HEA) builds the trainable unitary W(θ)W(\theta) of a variational circuit from the gates a device executes natively: single-qubit RYR_Y and RZR_Z rotations on every qubit, followed by a ladder of nearest-neighbour CNOTs, repeated LL times. This lesson counts its resources (2nL2nL parameters, (n−1)L(n-1)L CNOTs, a ladder of depth exactly n−1n-1), proves that the all-zero initialisation outputs ∣0n⟩\lvert0^n\rangle because the ladder permutes basis states by prefix parities, tracks basis states to show that a single non-zero angle in the last layer prepares the GHZ state while the same angle in the first layer only reaches fidelity 14\tfrac14, derives the Haar fidelity density (d−1)(1−F)d−2(d-1)(1-F)^{d-2} behind the expressibility measure of Sim, Johnson and Aspuru-Guzik, and relates expressibility to trainability and to the connectivity of real hardware.

Prerequisites: Variational Quantum Circuits: Architecture

Objectives and key results

  • Write the layer Wl(θ)=E∏qRZ(θl,q,1)RY(θl,q,0)W_l(\theta)=E\prod_qR_Z(\theta_{l,q,1})R_Y(\theta_{l,q,0}) of a hardware-efficient ansatz and draw it as a circuit with the CNOT ladder E=CNOT(n−2,n−1)⋯CNOT(0,1)E=\mathrm{CNOT}(n-2,n-1)\cdots\mathrm{CNOT}(0,1).

  • Count parameters (2nL2nL), CNOTs ((n−1)L(n-1)L) and depth, prove that the ladder has depth exactly n−1n-1 and compute the depth of a scheduled LL-layer circuit.

  • Prove that the ladder maps ∣x⟩\lvert x\rangle to the prefix parities of xx, hence E∣0n⟩=∣0n⟩E\lvert0^n\rangle=\lvert0^n\rangle and the zero initialisation outputs ∣0n⟩\lvert0^n\rangle.

  • Track basis states through a parameterised circuit to find angles that prepare a target state, and explain why the layer in which a rotation sits matters.

  • Define expressibility as DKL(P^PQC(F) ∥ PHaar(F))D_{KL}\big(\hat P_{\mathrm{PQC}}(F)\,\Vert\,P_{\mathrm{Haar}}(F)\big) and derive PHaar(F)=(d−1)(1−F)d−2P_{\mathrm{Haar}}(F)=(d-1)(1-F)^{d-2} from Gaussian vectors and the Gamma–Beta relation.

  • Explain the tension between expressibility and trainability, and estimate the SWAP overhead that non-native connectivity would add.

definition · Hardware-efficient layer
W(θ)=WL−1(θ)⋯W0(θ),Wl(θ)=E∏q=0n−1RZ(q)(θl,q,1) RY(q)(θl,q,0),E=CNOT(n−2,n−1)⋯CNOT(0,1)W(\theta)=W_{L-1}(\theta)\cdots W_0(\theta),\qquad W_l(\theta)=E\prod_{q=0}^{n-1}R_Z^{(q)}(\theta_{l,q,1})\,R_Y^{(q)}(\theta_{l,q,0}),\qquad E=\mathrm{CNOT}(n-2,n-1)\cdots\mathrm{CNOT}(0,1)
formula · Resource counts on a line
#θ=2nL,#CNOT=(n−1)L,depth(E)=n−1,depth(W)≤L(n+1)\#\theta=2nL,\qquad\#\mathrm{CNOT}=(n-1)L,\qquad\mathrm{depth}(E)=n-1,\qquad\mathrm{depth}(W)\le L(n+1)
identity · Ladder action and identity at zero
E∣x0x1⋯xn−1⟩=∣x0, x0⊕x1, …, x0⊕⋯⊕xn−1⟩,W(0)∣0n⟩=EL∣0n⟩=∣0n⟩E\lvert x_0x_1\cdots x_{n-1}\rangle=\lvert x_0,\ x_0\oplus x_1,\ \dots,\ x_0\oplus\cdots\oplus x_{n-1}\rangle,\qquad W(0)\lvert0^n\rangle=E^L\lvert0^n\rangle=\lvert0^n\rangle
theorem · Haar fidelity density
PHaar(F)=(d−1)(1−F)d−2,Pr⁡(F≤f)=1−(1−f)d−1,E[F]=1d,E[F2]=2d(d+1)P_{\mathrm{Haar}}(F)=(d-1)(1-F)^{d-2},\qquad\Pr(F\le f)=1-(1-f)^{d-1},\qquad\mathbb E[F]=\frac1d,\qquad\mathbb E[F^2]=\frac{2}{d(d+1)}
definition · Expressibility (Sim et al. 2019)
Expr=DKL(P^PQC ∥ PHaar)=∑j=0B−1P^jln⁡P^jQj,Qj=(1−jB)d−1−(1−j+1B)d−1\mathrm{Expr}=D_{KL}\big(\hat P_{\mathrm{PQC}}\,\Vert\,P_{\mathrm{Haar}}\big)=\sum_{j=0}^{B-1}\hat P_j\ln\frac{\hat P_j}{Q_j},\qquad Q_j=\Big(1-\frac jB\Big)^{d-1}-\Big(1-\frac{j+1}B\Big)^{d-1}
theorem · Cost concentration for Haar-random states
Eψ∼Haar[⟨ψ∣Z0∣ψ⟩]=0,Varψ∼Haar[⟨ψ∣Z0∣ψ⟩]=12n+1\mathbb E_{\psi\sim\mathrm{Haar}}\big[\langle\psi\vert Z_0\vert\psi\rangle\big]=0,\qquad\mathrm{Var}_{\psi\sim\mathrm{Haar}}\big[\langle\psi\vert Z_0\vert\psi\rangle\big]=\frac1{2^n+1}

Layers built from native gates#

The previous lesson split a variational quantum circuit into a feature map U(x)U(x), a trainable unitary W(θ)W(\theta) and an observable HH. It left open how W(θ)W(\theta) should be chosen. Two design goals pull in different directions. Physics-inspired ansätze, such as the unitary coupled-cluster ansatz of quantum chemistry, encode knowledge of the problem but compile into long sequences of gates that a noisy device cannot execute faithfully. Hardware-efficient ansätze, introduced for the variational eigensolver experiments of Kandala et al. (2017), take the opposite view: use only the gates the processor implements natively and only couplings between qubits that are physically connected, stack as many layers as the coherence time allows, and let the optimiser find the angles. This lesson studies the simplest member of that family, which is also the default template of most QML libraries: RYR_Y and RZR_Z rotations on every qubit followed by a ladder of CNOTs between neighbours.

Definition

Hardware-efficient ansatz (HEA)

Let n≥1n\ge1 qubits be arranged on a line, L≥1L\ge1 the number of layers, and θ=(θl,q,k)\theta=(\theta_{l,q,k}) with l∈{0,…,L−1}l\in\{0,\dots,L-1\}, q∈{0,…,n−1}q\in\{0,\dots,n-1\}, k∈{0,1}k\in\{0,1\} the parameter array. Write RP(q)(α)R^{(q)}_P(\alpha) for RP(α)=e−iαP/2R_P(\alpha)=e^{-i\alpha P/2} acting on qubit qq (tensored with identities). The CNOT ladder and the layer are E=CNOT(n−2,n−1)⋯CNOT(1,2) CNOT(0,1),Wl(θ)=E∏q=0n−1RZ(q)(θl,q,1) RY(q)(θl,q,0),E=\mathrm{CNOT}(n-2,n-1)\cdots\mathrm{CNOT}(1,2)\,\mathrm{CNOT}(0,1),\qquad W_l(\theta)=E\prod_{q=0}^{n-1}R_Z^{(q)}(\theta_{l,q,1})\,R_Y^{(q)}(\theta_{l,q,0}), with E=IE=I for n=1n=1, and the hardware-efficient ansatz with its output state is W(θ)=WL−1(θ)⋯W1(θ)W0(θ),∣ψ(θ)⟩=W(θ)∣0n⟩.W(\theta)=W_{L-1}(\theta)\cdots W_1(\theta)W_0(\theta),\qquad\lvert\psi(\theta)\rangle=W(\theta)\lvert0^n\rangle . In time order, every layer applies RYR_Y and then RZR_Z to each qubit, then CNOT(0,1)\mathrm{CNOT}(0,1), then CNOT(1,2)\mathrm{CNOT}(1,2), and so on down to CNOT(n−2,n−1)\mathrm{CNOT}(n-2,n-1). Rotations on different qubits act on different tensor factors and commute, so the order of the factors in ∏q\prod_q is irrelevant; on one qubit the order matters, because RZR_Z and RYR_Y do not commute.

Circuit · 3 qubits

Quantum circuit on 3 qubitsRY(0) on q0; RZ(0) on q0; RY(0) on q1; RZ(0) on q1; RY(0) on q2; RZ(0) on q2; X on q1 controlled by q0; X on q2 controlled by q1; RY(0) on q0; RZ(0) on q0; RY(0) on q1; RZ(0) on q1; RY(0) on q2; RZ(0) on q2; X on q1 controlled by q0; X on q2 controlled by q1q0q1q2RY(0)RZ(0)RY(0)RZ(0)RY(0)RZ(0)RY(0)RZ(0)RY(0)RZ(0)RY(0)RZ(0)

The template of the lab: n=3n=3 qubits, L=2L=2 layers, all 1212 angles set to 00. Each layer is a column of RYR_Y, a column of RZR_Z, then CNOT(0,1)\mathrm{CNOT}(0,1) and CNOT(1,2)\mathrm{CNOT}(1,2). In the lab every angle is a named trainable parameter θl,q,k\theta_{l,q,k}.

Theorem

Parameter count, CNOT count and depth

For n≥2n\ge2 qubits and LL layers the hardware-efficient ansatz has (a) 2nL2nL parameters; (b) (n−1)L(n-1)L CNOT gates, each between neighbouring qubits q,q+1q,q+1, so a linear coupling map needs no SWAP gates; (c) a CNOT ladder EE of depth exactly n−1n-1; (d) layers of depth exactly 2+(n−1)=n+12+(n-1)=n+1, so that the circuit executed layer after layer has depth L(n+1)L(n+1). Here a moment is a set of gates acting on pairwise disjoint qubits, and the depth of a gate sequence is the smallest number of consecutive moments into which it can be scheduled such that two gates sharing a qubit keep their order.

Proof (8 steps)
  1. Step 1:
    #θ=∑l=0L−1∑q=0n−1#{k∈{0,1}}=L⋅n⋅2=2nL\#\theta=\sum_{l=0}^{L-1}\sum_{q=0}^{n-1}\#\{k\in\{0,1\}\}=L\cdot n\cdot2=2nL
    Each layer has one RYR_Y and one RZR_Z angle per qubit, and every angle is an independent parameter.
  2. Step 2:
    #CNOT=∑l=0L−1#{(q,q+1):0≤q≤n−2}=(n−1)L,∣q−(q+1)∣=1\#\mathrm{CNOT}=\sum_{l=0}^{L-1}\#\{(q,q+1):0\le q\le n-2\}=(n-1)L,\qquad\lvert q-(q+1)\rvert=1
    The ladder contains one CNOT per neighbouring pair. On a line every such pair is physically coupled, so each CNOT is native and no routing is required.
  3. Step 3:
    m(CNOT(k,k+1))=k+1,k=0,…,n−2 ⇒ depth(E)≤n−1m\big(\mathrm{CNOT}(k,k+1)\big)=k+1,\quad k=0,\dots,n-2\ \Rightarrow\ \mathrm{depth}(E)\le n-1
    An explicit schedule: one CNOT per moment, in ladder order. It respects the order of all gates and uses n−1n-1 moments.
  4. Step 4:
    CNOT(k−1,k) and CNOT(k,k+1) share qubit k ⇒ m(CNOT(k,k+1))≥m(CNOT(k−1,k))+1\mathrm{CNOT}(k-1,k)\ \text{and}\ \mathrm{CNOT}(k,k+1)\ \text{share qubit }k\ \Rightarrow\ m\big(\mathrm{CNOT}(k,k+1)\big)\ge m\big(\mathrm{CNOT}(k-1,k)\big)+1
    Two gates on a common qubit cannot be in the same moment, and the earlier gate must stay earlier.
  5. Step 5:
    m(CNOT(0,1))≥1 ⇒ m(CNOT(k,k+1))≥k+1 ⇒ depth(E)≥m(CNOT(n−2,n−1))≥n−1m\big(\mathrm{CNOT}(0,1)\big)\ge1\ \Rightarrow\ m\big(\mathrm{CNOT}(k,k+1)\big)\ge k+1\ \Rightarrow\ \mathrm{depth}(E)\ge m\big(\mathrm{CNOT}(n-2,n-1)\big)\ge n-1
    Induction on kk with the previous inequality. Together with the explicit schedule, depth(E)=n−1\mathrm{depth}(E)=n-1: consecutive CNOTs of the ladder share a qubit, so no two of them can run in parallel.
  6. Step 6:
    moment12345q0RYRZ∙q1RYRZ⊕∙q2RYRZ⊕∙q3RYRZ⊕\begin{array}{c|ccccc}\text{moment}&1&2&3&4&5\\\hline q_0&R_Y&R_Z&\bullet&&\\q_1&R_Y&R_Z&\oplus&\bullet&\\q_2&R_Y&R_Z&&\oplus&\bullet\\q_3&R_Y&R_Z&&&\oplus\end{array}
    Moment diagram of one layer for n=4n=4 (∙\bullet control, ⊕\oplus target): all RYR_Y fit into moment 11 and all RZR_Z into moment 22 because they act on distinct qubits; the ladder is a staircase occupying moments 33 to n+1=5n+1=5.
  7. Step 7:
    RY(0)→RZ(0)→CNOT(0,1)→CNOT(1,2)→⋯→CNOT(n−2,n−1) ⇒ depth(Wl)=2+(n−1)=n+1R_Y^{(0)}\to R_Z^{(0)}\to\mathrm{CNOT}(0,1)\to\mathrm{CNOT}(1,2)\to\cdots\to\mathrm{CNOT}(n-2,n-1)\ \Rightarrow\ \mathrm{depth}(W_l)=2+(n-1)=n+1
    The diagram schedules a layer in n+1n+1 moments. Conversely, in this chain of n+1n+1 gates each gate shares a qubit with the next one, so the chain needs n+1n+1 distinct, ordered moments.
  8. Step 8:
    depthlayer by layer(W)=∑l=0L−1depth(Wl)=L(n+1)\mathrm{depth}_{\text{layer by layer}}(W)=\sum_{l=0}^{L-1}\mathrm{depth}(W_l)=L(n+1)
    If every layer is executed as a separate block (as a compiler does when layers are separated by barriers), the depths of the blocks add. The next proposition shows that a scheduler which starts each gate as soon as its qubits are free can overlap consecutive layers.

∎

Proposition

As-soon-as-possible scheduling overlaps consecutive layers

Schedule each gate of W(θ)W(\theta) in the first moment after the last moment used on any of its qubits (as-soon-as-possible, ASAP, scheduling, which is how the platform's depth counter packs gates). Then depthASAP(W)={L(n+1),n≤2,n+1+4(L−1),n≥3,\mathrm{depth}_{\mathrm{ASAP}}(W)=\begin{cases}L(n+1),&n\le2,\\ n+1+4(L-1),&n\ge3,\end{cases} and no valid schedule is shorter. For n=3n=3 both formulas give 4L4L (depth 88 for the lab template); for n≥4n\ge4 and L≥2L\ge2 the ASAP depth is strictly smaller than L(n+1)L(n+1), because the rotations of the next layer start on the top qubits while the ladder is still running down the chain.

Proof (9 steps)
  1. Step 1:
    m(g)=1+max⁡q∈supp(g)m(prevq(g)),m(∅)=0m(g)=1+\max_{q\in\mathrm{supp}(g)}m\big(\mathrm{prev}_q(g)\big),\qquad m(\varnothing)=0
    Definition of ASAP: prevq(g)\mathrm{prev}_q(g) is the previous gate on qubit qq. By induction over the gate sequence, every valid schedule places each gate at a moment ≥m(g)\ge m(g), so ASAP is optimal and its depth is the length of the longest chain of gates in which consecutive gates share a qubit.
  2. Step 2:
    n=1: m=1,2,…,2L;n=2: layer l occupies moments 3l+1, 3l+2, 3l+3n=1:\ m=1,2,\dots,2L;\qquad n=2:\ \text{layer }l\text{ occupies moments }3l+1,\,3l+2,\,3l+3
    For n=1n=1 all 2L2L rotations act on one qubit. For n=2n=2 each layer ends with CNOT(0,1)\mathrm{CNOT}(0,1), which touches both qubits, so the next layer cannot start before it: 3L=L(n+1)3L=L(n+1).
  3. Step 3:
    n≥3, l=0:m(RY)=1, m(RZ)=2,c0(k):=m(CNOT(k,k+1))=3+k,T0:=c0(n−2)=n+1n\ge3,\ l=0:\quad m(R_Y)=1,\ m(R_Z)=2,\quad c_0(k):=m\big(\mathrm{CNOT}(k,k+1)\big)=3+k,\quad T_0:=c_0(n-2)=n+1
    The first layer starts on idle qubits, so it is scheduled exactly as in the moment diagram. Note c0(k)=T0−(n−2−k)c_0(k)=T_0-(n-2-k).
  4. Step 4:
    hypothesis: cl(k)=Tl−(n−2−k) ⇒ fl(q)=Tl−(n−2−q) (q≤n−2),fl(n−1)=Tl\text{hypothesis: }c_l(k)=T_l-(n-2-k)\ \Rightarrow\ f_l(q)=T_l-(n-2-q)\ (q\le n-2),\qquad f_l(n-1)=T_l
    fl(q)f_l(q) is the last moment of layer ll on qubit qq: for q≤n−2q\le n-2 the last gate on qq is CNOT(q,q+1)\mathrm{CNOT}(q,q+1), which comes after CNOT(q−1,q)\mathrm{CNOT}(q-1,q); for q=n−1q=n-1 it is CNOT(n−2,n−1)\mathrm{CNOT}(n-2,n-1).
  5. Step 5:
    r(q):=m(RZ(q) in layer l+1)=fl(q)+2={Tl−n+4+q,q≤n−2,Tl+2,q=n−1r(q):=m\big(R_Z^{(q)}\text{ in layer }l+1\big)=f_l(q)+2=\begin{cases}T_l-n+4+q,&q\le n-2,\\T_l+2,&q=n-1\end{cases}
    On each qubit the next layer's RYR_Y follows immediately, at fl(q)+1f_l(q)+1, and its RZR_Z at fl(q)+2f_l(q)+2.
  6. Step 6:
    cl+1(0)=1+max⁡(r(0),r(1))=Tl−n+6,cl+1(k)=1+max⁡(cl+1(k−1), r(k+1))=Tl−n+6+k(1≤k≤n−3)c_{l+1}(0)=1+\max\big(r(0),r(1)\big)=T_l-n+6,\qquad c_{l+1}(k)=1+\max\big(c_{l+1}(k-1),\,r(k+1)\big)=T_l-n+6+k\quad(1\le k\le n-3)
    ASAP for the ladder of layer l+1l+1. Qubit 1≤n−21\le n-2 because n≥3n\ge3. In the inductive step both arguments of the maximum equal Tl−n+5+kT_l-n+5+k.
  7. Step 7:
    Tl+1=cl+1(n−2)=1+max⁡(cl+1(n−3), r(n−1))=1+max⁡(Tl+3, Tl+2)=Tl+4T_{l+1}=c_{l+1}(n-2)=1+\max\big(c_{l+1}(n-3),\,r(n-1)\big)=1+\max(T_l+3,\,T_l+2)=T_l+4
    The last CNOT of the ladder. Moreover cl+1(k)=Tl−n+6+k=Tl+1−(n−2−k)c_{l+1}(k)=T_l-n+6+k=T_{l+1}-(n-2-k) for all kk, which is the hypothesis for layer l+1l+1.
  8. Step 8:
    depthASAP(W)=TL−1=n+1+4(L−1),fL−1(q)≤TL−1 ∀q\mathrm{depth}_{\mathrm{ASAP}}(W)=T_{L-1}=n+1+4(L-1),\qquad f_{L-1}(q)\le T_{L-1}\ \forall q
    Unrolling Tl+1=Tl+4T_{l+1}=T_l+4 from T0=n+1T_0=n+1. No gate ends later than the last CNOT of the last layer. For n=3n=3 this is 4L=L(n+1)4L=L(n+1); for n≥4n\ge4 it is smaller than L(n+1)L(n+1) whenever L≥2L\ge2.
  9. Step 9:
    moment123456789q0Y0Z0∙Y1Z1∙q1Y0Z0⊕∙Y1Z1⊕∙q2Y0Z0⊕∙Y1Z1⊕∙q3Y0Z0⊕Y1Z1⊕\begin{array}{c|ccccccccc}\text{moment}&1&2&3&4&5&6&7&8&9\\\hline q_0&Y_0&Z_0&\bullet&Y_1&Z_1&&\bullet&&\\q_1&Y_0&Z_0&\oplus&\bullet&Y_1&Z_1&\oplus&\bullet&\\q_2&Y_0&Z_0&&\oplus&\bullet&Y_1&Z_1&\oplus&\bullet\\q_3&Y_0&Z_0&&&\oplus&Y_1&Z_1&&\oplus\end{array}
    ASAP moment diagram for n=4n=4, L=2L=2 (YlY_l, ZlZ_l: the rotations of layer ll): depth 9=n+1+49=n+1+4, one moment less than L(n+1)=10L(n+1)=10. Qubit 00 starts layer 11 in moment 44, while CNOT(1,2)\mathrm{CNOT}(1,2) of layer 00 is still running.

∎

Quick check

How many trainable parameters does a hardware-efficient ansatz with n=5n=5 qubits and L=4L=4 layers have?

tolerance ±0.5

Quick check

What is the depth of a hardware-efficient ansatz with n=6n=6 qubits and L=3L=3 layers when the layers are executed one after the other as separate blocks?

tolerance ±0.5

The ladder as a permutation and the zero initialisation#

Before training, all angles must be given values. The most common choice in practice is θ=0\theta=0. To understand what the circuit computes there, we first describe the CNOT ladder exactly: it never creates superpositions, it only relabels the computational basis states, and the relabelling is a simple parity rule.

Lemma

The CNOT ladder computes prefix parities

For n≥2n\ge2 and every bit string x=x0x1⋯xn−1x=x_0x_1\cdots x_{n-1}, E∣x0x1⋯xn−1⟩=∣y0y1⋯yn−1⟩,yk=x0⊕x1⊕⋯⊕xk.E\lvert x_0x_1\cdots x_{n-1}\rangle=\lvert y_0y_1\cdots y_{n-1}\rangle,\qquad y_k=x_0\oplus x_1\oplus\cdots\oplus x_k . Hence E∣0n⟩=∣0n⟩E\lvert0^n\rangle=\lvert0^n\rangle, E∣10n−1⟩=∣1n⟩E\lvert10^{n-1}\rangle=\lvert1^n\rangle, E∣1n⟩=∣1010⋯ ⟩E\lvert1^n\rangle=\lvert1010\cdots\rangle (yk=1y_k=1 exactly for even kk), and EE is a permutation matrix.

Proof (7 steps)
  1. Step 1:
    CNOT(c,t)∣⋯xc⋯xt⋯ ⟩=∣⋯xc⋯(xc⊕xt)⋯ ⟩\mathrm{CNOT}(c,t)\lvert\cdots x_c\cdots x_t\cdots\rangle=\lvert\cdots x_c\cdots(x_c\oplus x_t)\cdots\rangle
    Definition CNOT(c,t)=∣0⟩⟨0∣c⊗It+∣1⟩⟨1∣c⊗Xt\mathrm{CNOT}(c,t)=\lvert0\rangle\langle0\rvert_c\otimes I_t+\lvert1\rangle\langle1\rvert_c\otimes X_t: the target bit is flipped exactly when the control bit is 11; all other bits are unchanged.
  2. Step 2:
    x(j):=label after CNOT(0,1),…,CNOT(j−1,j);claim: xk(j)={x0⊕⋯⊕xk,k≤j,xk,k>jx^{(j)}:=\text{label after }\mathrm{CNOT}(0,1),\dots,\mathrm{CNOT}(j-1,j);\qquad\text{claim: }x^{(j)}_k=\begin{cases}x_0\oplus\cdots\oplus x_k,&k\le j,\\x_k,&k>j\end{cases}
    Each CNOT maps a basis state to a basis state, so it suffices to follow the label. We prove the claim by induction on the number jj of ladder gates applied.
  3. Step 3:
    j=0:x(0)=x,x0(0)=x0j=0:\quad x^{(0)}=x,\qquad x^{(0)}_0=x_0
    Base case: no gate applied; the prefix parity of length one is x0x_0 itself.
  4. Step 4:
    xj+1(j+1)=xj(j)⊕xj+1(j)=(x0⊕⋯⊕xj)⊕xj+1,xk(j+1)=xk(j) (k≠j+1)x^{(j+1)}_{j+1}=x^{(j)}_j\oplus x^{(j)}_{j+1}=(x_0\oplus\cdots\oplus x_j)\oplus x_{j+1},\qquad x^{(j+1)}_k=x^{(j)}_k\ (k\ne j+1)
    Inductive step: gate j+1j+1 is CNOT(j,j+1)\mathrm{CNOT}(j,j+1), with control bit xj(j)x^{(j)}_j (a prefix parity by hypothesis) and target bit xj+1(j)=xj+1x^{(j)}_{j+1}=x_{j+1} (still untouched).
  5. Step 5:
    E∣x⟩=∣x(n−1)⟩=∣y⟩,yk=x0⊕⋯⊕xkE\lvert x\rangle=\lvert x^{(n-1)}\rangle=\lvert y\rangle,\qquad y_k=x_0\oplus\cdots\oplus x_k
    After all n−1n-1 gates every position holds its prefix parity.
  6. Step 6:
    x=0n⇒y=0n;x=10n−1⇒y=1n;x=1n⇒yk=(k+1) mod 2x=0^n\Rightarrow y=0^n;\qquad x=10^{n-1}\Rightarrow y=1^n;\qquad x=1^n\Rightarrow y_k=(k+1)\bmod2
    Every prefix of 0n0^n has parity 00; every prefix of 10n−110^{n-1} contains exactly one 11; the prefix of 1n1^n of length k+1k+1 has parity (k+1) mod 2(k+1)\bmod2. For n=3n=3: E∣111⟩=∣101⟩E\lvert111\rangle=\lvert101\rangle.
  7. Step 7:
    x0=y0,xk=yk−1⊕yk (k≥1) ⇒ x↦y bijectivex_0=y_0,\quad x_k=y_{k-1}\oplus y_k\ (k\ge1)\ \Rightarrow\ x\mapsto y\ \text{bijective}
    The map can be inverted (consecutive differences), so EE sends the 2n2^n basis vectors to the 2n2^n basis vectors bijectively: it is a permutation matrix, in particular real and unitary.

∎

E=(1000000001000000000100000010000000000010000000010000010000001000),000↦000, 001↦001, 010↦011, 011↦010,100↦111, 101↦110, 110↦100, 111↦101E=\begin{pmatrix}1&0&0&0&0&0&0&0\\0&1&0&0&0&0&0&0\\0&0&0&1&0&0&0&0\\0&0&1&0&0&0&0&0\\0&0&0&0&0&0&1&0\\0&0&0&0&0&0&0&1\\0&0&0&0&0&1&0&0\\0&0&0&0&1&0&0&0\end{pmatrix},\qquad\begin{aligned}&000\mapsto000,\ 001\mapsto001,\ 010\mapsto011,\ 011\mapsto010,\\&100\mapsto111,\ 101\mapsto110,\ 110\mapsto100,\ 111\mapsto101\end{aligned}
E = CNOT(1,2) CNOT(0,1) on three qubits8×8

Unitary verified: U†U=IU^\dagger U = I (max deviation 0.0e+0)

Theorem

Identity at zero

For all nn and LL: Wl(0)=EW_l(0)=E for every layer, W(0)=ELW(0)=E^L, and W(0)∣0n⟩=∣0n⟩.W(0)\lvert0^n\rangle=\lvert0^n\rangle . So the all-zero initialisation outputs the input state ∣0n⟩\lvert0^n\rangle, and the fidelity of the initial circuit with a target ∣ϕ⟩\lvert\phi\rangle is ∣⟨ϕ∣0n⟩∣2\lvert\langle\phi\vert0^n\rangle\rvert^2; for ∣GHZn⟩=12(∣0n⟩+∣1n⟩)\lvert\mathrm{GHZ}_n\rangle=\tfrac1{\sqrt2}(\lvert0^n\rangle+\lvert1^n\rangle) it is 12\tfrac12. As an operator, however, W(0)=ELW(0)=E^L is not the identity in general: EE has order 2⌈log⁡2n⌉2^{\lceil\log_2n\rceil}, the smallest power of two that is ≥n\ge n, so W(0)=IW(0)=I exactly when LL is a multiple of that order.

Proof (8 steps)
  1. Step 1:
    RY(0)=cos⁡0 I−isin⁡0 Y=I,RZ(0)=cos⁡0 I−isin⁡0 Z=IR_Y(0)=\cos0\,I-i\sin0\,Y=I,\qquad R_Z(0)=\cos0\,I-i\sin0\,Z=I
    RP(α)=cos⁡α2 I−isin⁡α2 PR_P(\alpha)=\cos\frac\alpha2\,I-i\sin\frac\alpha2\,P at α=0\alpha=0.
  2. Step 2:
    Wl(0)=E∏q=0n−1I⋅I=E,W(0)=E⋯E⏟L=ELW_l(0)=E\prod_{q=0}^{n-1}I\cdot I=E,\qquad W(0)=\underbrace{E\cdots E}_{L}=E^L
    Every rotation factor of every layer is the identity.
  3. Step 3:
    E∣0n⟩=∣0n⟩ ⇒ EL∣0n⟩=EL−1(E∣0n⟩)=EL−1∣0n⟩=⋯=∣0n⟩E\lvert0^n\rangle=\lvert0^n\rangle\ \Rightarrow\ E^{L}\lvert0^n\rangle=E^{L-1}\big(E\lvert0^n\rangle\big)=E^{L-1}\lvert0^n\rangle=\cdots=\lvert0^n\rangle
    The lemma (y=0ny=0^n for x=0nx=0^n), or directly: every CNOT of the ladder sees a control bit 00 and acts trivially. Induction on LL.
  4. Step 4:
    ∣⟨GHZn∣0n⟩∣2=∣12⟨0n∣0n⟩+12⟨1n∣0n⟩∣2=12\big\lvert\langle\mathrm{GHZ}_n\vert0^n\rangle\big\rvert^2=\Big\lvert\tfrac1{\sqrt2}\langle0^n\vert0^n\rangle+\tfrac1{\sqrt2}\langle1^n\vert0^n\rangle\Big\rvert^2=\tfrac12
    ∣0n⟩\lvert0^n\rangle and ∣1n⟩\lvert1^n\rangle are orthonormal. In the lab the template starts at fidelity 12\tfrac12 before any angle is changed.
  5. Step 5:
    y=Ax over F2,A=(111⋮⋱11⋯1)=I+Ny=Ax\ \text{over}\ \mathbb F_2,\qquad A=\begin{pmatrix}1&&&\\1&1&&\\\vdots&&\ddots&\\1&1&\cdots&1\end{pmatrix}=I+N
    By the lemma, EE acts on labels as the lower-triangular all-ones matrix AA; NN is strictly lower triangular, so Nn=0N^n=0. ELE^L acts as ALA^L, hence EL=IE^L=I if and only if AL=IA^L=I.
  6. Step 6:
    A2m=(I+N)2m=I+N2mover F2A^{2^m}=(I+N)^{2^m}=I+N^{2^m}\quad\text{over }\mathbb F_2
    II and NN commute, and the binomial coefficients (2mi)\binom{2^m}{i} with 0<i<2m0<i<2^m are even (by induction, (I+N)2m=((I+N)2m−1)2=(I+N2m−1)2=I+2N2m−1+N2m=I+N2m(I+N)^{2^{m}}=\big((I+N)^{2^{m-1}}\big)^2=(I+N^{2^{m-1}})^2=I+2N^{2^{m-1}}+N^{2^m}=I+N^{2^m}).
  7. Step 7:
    (Nj)j,0=1  (1≤j≤n−1) ⇒ N2m=0  ⟺  2m≥n(N^j)_{j,0}=1\ \ (1\le j\le n-1)\ \Rightarrow\ N^{2^m}=0\iff2^m\ge n
    (Nj)a,b(N^j)_{a,b} counts mod 22 the strictly decreasing index paths a→⋯→ba\to\dots\to b of length jj; from jj to 00 there is exactly one, so Nj≠0N^j\ne0 for j<nj<n, while Nn=0N^n=0.
  8. Step 8:
    A2m∗=I,m∗=⌈log⁡2n⌉;ord(A)∣2m∗, A2m∗−1≠I ⇒ ord(E)=2⌈log⁡2n⌉A^{2^{m^\ast}}=I,\quad m^\ast=\lceil\log_2n\rceil;\qquad\mathrm{ord}(A)\mid2^{m^\ast},\ A^{2^{m^\ast-1}}\ne I\ \Rightarrow\ \mathrm{ord}(E)=2^{\lceil\log_2n\rceil}
    For n≥2n\ge2: the order divides 2m∗2^{m^\ast}, so it is a power of two, and it is not a smaller power since 2m∗−1<n2^{m^\ast-1}<n. Examples: n=2n=2: E=CNOTE=\mathrm{CNOT}, order 22; n=3,4n=3,4: order 44 (for the lab, W(0)=E2≠IW(0)=E^2\ne I); n=5,…,8n=5,\dots,8: order 88.

∎

Quick check

A hardware-efficient ansatz on n=4n=4 qubits with L=3L=3 layers is initialised with all angles equal to 00. Which statement is correct?

Select one answer

Reaching the GHZ state with one angle#

The lab asks for angles of the 33-qubit, 22-layer template that prepare ∣GHZ⟩=12(∣000⟩+∣111⟩)\lvert\mathrm{GHZ}\rangle=\tfrac1{\sqrt2}(\lvert000\rangle+\lvert111\rangle). The template starts at ∣000⟩\lvert000\rangle with fidelity 12\tfrac12. Because the ladder permutes basis states, the effect of a few non-zero angles can be computed by following a handful of basis states through the circuit, without multiplying 8×88\times8 matrices. We write θl,q,k\theta_{l,q,k} for the angle that the lab calls t_l_q_k\texttt{t\_}l\texttt{\_}q\texttt{\_}k.

Derivation

One rotation in the last layer prepares GHZ

  1. Step 1:
    θ1,0,0=π2,θl,q,k=0 otherwise\theta_{1,0,0}=\tfrac\pi2,\qquad\theta_{l,q,k}=0\ \text{otherwise}
    The candidate: a single RYR_Y angle on qubit 00 in layer 11, the last layer.
  2. Step 2:
    W0(0)∣000⟩=E∣000⟩=∣000⟩W_0(0)\lvert000\rangle=E\lvert000\rangle=\lvert000\rangle
    Layer 00 has only zero angles, so it acts as EE, which fixes ∣000⟩\lvert000\rangle (identity at zero).
  3. Step 3:
    RY(0)(π2)∣000⟩=(cos⁡π4∣0⟩+sin⁡π4∣1⟩)∣00⟩=12(∣000⟩+∣100⟩)R_Y^{(0)}\big(\tfrac\pi2\big)\lvert000\rangle=\Big(\cos\tfrac\pi4\lvert0\rangle+\sin\tfrac\pi4\lvert1\rangle\Big)\lvert00\rangle=\tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert100\rangle\big)
    The first column of RY(α)=(cos⁡α2−sin⁡α2sin⁡α2cos⁡α2)R_Y(\alpha)=\begin{pmatrix}\cos\frac\alpha2&-\sin\frac\alpha2\\\sin\frac\alpha2&\cos\frac\alpha2\end{pmatrix} at α=π2\alpha=\frac\pi2; all other rotations of layer 11 are II. This is ∣+00⟩\lvert{+}00\rangle.
  4. Step 4:
    CNOT(0,1): ∣000⟩↦∣000⟩,  ∣100⟩↦∣110⟩ ⇒ 12(∣000⟩+∣110⟩)\mathrm{CNOT}(0,1):\ \lvert000\rangle\mapsto\lvert000\rangle,\ \ \lvert100\rangle\mapsto\lvert110\rangle\ \Rightarrow\ \tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert110\rangle\big)
    Control bit x0x_0: the first term has x0=0x_0=0 and is unchanged; in the second term x1x_1 flips. By linearity the CNOT acts term by term.
  5. Step 5:
    CNOT(1,2): ∣000⟩↦∣000⟩,  ∣110⟩↦∣111⟩ ⇒ 12(∣000⟩+∣111⟩)=∣GHZ⟩\mathrm{CNOT}(1,2):\ \lvert000\rangle\mapsto\lvert000\rangle,\ \ \lvert110\rangle\mapsto\lvert111\rangle\ \Rightarrow\ \tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert111\rangle\big)=\lvert\mathrm{GHZ}\rangle
    Control bit x1x_1: 00 in the first term, 11 in the second. Equivalently, E∣100⟩=∣111⟩E\lvert100\rangle=\lvert111\rangle by the prefix-parity lemma.
  6. Step 6:
    F=∣⟨GHZ∣ψ(θ)⟩∣2=1F=\big\lvert\langle\mathrm{GHZ}\vert\psi(\theta)\rangle\big\rvert^2=1
    The output equals the target exactly. This is the lab's reference solution.
  7. Step 7:
    θ1,0,0=a:∣ψ⟩=cos⁡a2∣000⟩+sin⁡a2∣111⟩,F(a)=12(cos⁡a2+sin⁡a2)2=1+sin⁡a2\theta_{1,0,0}=a:\quad\lvert\psi\rangle=\cos\tfrac a2\lvert000\rangle+\sin\tfrac a2\lvert111\rangle,\qquad F(a)=\tfrac12\Big(\cos\tfrac a2+\sin\tfrac a2\Big)^2=\frac{1+\sin a}{2}
    The same computation for a general angle, with 2sin⁡a2cos⁡a2=sin⁡a2\sin\frac a2\cos\frac a2=\sin a. The grader's threshold F≥0.99F\ge0.99 means sin⁡a≥0.98\sin a\ge0.98, i.e. a∈[arcsin⁡0.98, π−arcsin⁡0.98]≈[1.3705, 1.7711]a\in[\arcsin0.98,\ \pi-\arcsin0.98]\approx[1.3705,\ 1.7711] (mod 2π2\pi). For a=−π2a=-\frac\pi2 the output is 12(∣000⟩−∣111⟩)\tfrac1{\sqrt2}(\lvert000\rangle-\lvert111\rangle) and F=0F=0.
  8. Step 8:
    θ1,0,1=φ:RZ(φ)∣0⟩+∣1⟩2=e−iφ/2∣0⟩+eiφ∣1⟩2 ⇒ F=∣1+eiφ2∣2=cos⁡2φ2\theta_{1,0,1}=\varphi:\quad R_Z(\varphi)\tfrac{\lvert0\rangle+\lvert1\rangle}{\sqrt2}=e^{-i\varphi/2}\tfrac{\lvert0\rangle+e^{i\varphi}\lvert1\rangle}{\sqrt2}\ \Rightarrow\ F=\Big\lvert\frac{1+e^{i\varphi}}{2}\Big\rvert^2=\cos^2\frac\varphi2
    An additional RZR_Z angle on qubit 00 in the last layer becomes the relative phase of the GHZ state; the global phase e−iφ/2e^{-i\varphi/2} drops out of the fidelity, and ∣1+eiφ∣2=2+2cos⁡φ=4cos⁡2φ2\lvert1+e^{i\varphi}\rvert^2=2+2\cos\varphi=4\cos^2\frac\varphi2.

Result

θ1,0,0=π2 ⇒ W(θ)∣000⟩=12(∣000⟩+∣111⟩)\boxed{\theta_{1,0,0}=\tfrac\pi2\ \Rightarrow\ W(\theta)\lvert000\rangle=\tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert111\rangle\big)}
Derivation

The same rotation in the first layer only reaches fidelity 1/4

  1. Step 1:
    θ0,0,0=π2,θl,q,k=0 otherwise\theta_{0,0,0}=\tfrac\pi2,\qquad\theta_{l,q,k}=0\ \text{otherwise}
    The tempting alternative: put the rotation in the first layer.
  2. Step 2:
    W0(θ)∣000⟩=12(∣000⟩+∣111⟩)W_0(\theta)\lvert000\rangle=\tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert111\rangle\big)
    Layer 00 now performs exactly the computation of the previous derivation: after the first layer the register already holds GHZ.
  3. Step 3:
    ∏qRZ(q)(0)RY(q)(0)=I\prod_{q}R_Z^{(q)}(0)R_Y^{(q)}(0)=I
    All rotations of layer 11 are zero, but its CNOT ladder is still applied.
  4. Step 4:
    CNOT(0,1): ∣000⟩↦∣000⟩,  ∣111⟩↦∣101⟩\mathrm{CNOT}(0,1):\ \lvert000\rangle\mapsto\lvert000\rangle,\ \ \lvert111\rangle\mapsto\lvert101\rangle
    In ∣111⟩\lvert111\rangle the control x0=1x_0=1 flips x1x_1 from 11 to 00.
  5. Step 5:
    CNOT(1,2): ∣000⟩↦∣000⟩,  ∣101⟩↦∣101⟩ ⇒ ∣ψ⟩=12(∣000⟩+∣101⟩)\mathrm{CNOT}(1,2):\ \lvert000\rangle\mapsto\lvert000\rangle,\ \ \lvert101\rangle\mapsto\lvert101\rangle\ \Rightarrow\ \lvert\psi\rangle=\tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert101\rangle\big)
    Now the control x1x_1 is 00 in both terms. In one line: E∣111⟩=∣1, 1⊕1, 1⊕1⊕1⟩=∣101⟩E\lvert111\rangle=\lvert1,\,1\oplus1,\,1\oplus1\oplus1\rangle=\lvert101\rangle.
  6. Step 6:
    ⟨GHZ∣ψ⟩=12(⟨000∣000⟩+⟨000∣101⟩+⟨111∣000⟩+⟨111∣101⟩)=12(1+0+0+0)=12\langle\mathrm{GHZ}\vert\psi\rangle=\tfrac12\big(\langle000\vert000\rangle+\langle000\vert101\rangle+\langle111\vert000\rangle+\langle111\vert101\rangle\big)=\tfrac12(1+0+0+0)=\tfrac12
    Expand both states in the orthonormal computational basis; only the ∣000⟩\lvert000\rangle components overlap.
  7. Step 7:
    F=∣12∣2=14F=\big\lvert\tfrac12\big\rvert^2=\tfrac14
    This is the value detected by the lab's diagnostic. The state 12(∣000⟩+∣101⟩)=12(∣00⟩+∣11⟩)02⊗∣0⟩1\tfrac1{\sqrt2}(\lvert000\rangle+\lvert101\rangle)=\tfrac1{\sqrt2}(\lvert00\rangle+\lvert11\rangle)_{02}\otimes\lvert0\rangle_1 is a Bell pair between qubits 00 and 22, with qubit 11 unentangled.

Result

θ0,0,0=π2 ⇒ W(θ)∣000⟩=12(∣000⟩+∣101⟩),F=14\boxed{\theta_{0,0,0}=\tfrac\pi2\ \Rightarrow\ W(\theta)\lvert000\rangle=\tfrac1{\sqrt2}\big(\lvert000\rangle+\lvert101\rangle\big),\qquad F=\tfrac14}

The general lesson is that every ladder after the last non-trivial rotation keeps acting. Earlier ladders relabel the basis states that the later rotations see, and later ladders relabel the result. A rotation in an early layer is therefore not "the same" as a rotation in the last layer. Keeping θ0,0,0=π2\theta_{0,0,0}=\tfrac\pi2, no choice of the layer-11 angles can repair the damage: those rotations act on the GHZ state as a local unitary V0⊗V1⊗V2V_0\otimes V_1\otimes V_2, and E(V0⊗V1⊗V2)∣GHZ⟩=∣GHZ⟩E(V_0\otimes V_1\otimes V_2)\lvert\mathrm{GHZ}\rangle=\lvert\mathrm{GHZ}\rangle would require (V0⊗V1⊗V2)∣GHZ⟩=E−1∣GHZ⟩=12(∣000⟩+∣100⟩)=∣+00⟩(V_0\otimes V_1\otimes V_2)\lvert\mathrm{GHZ}\rangle=E^{-1}\lvert\mathrm{GHZ}\rangle=\tfrac1{\sqrt2}(\lvert000\rangle+\lvert100\rangle)=\lvert{+}00\rangle (invert the prefix parities: xk=yk−1⊕ykx_k=y_{k-1}\oplus y_k maps 111111 to 100100). That is impossible, because a local unitary conjugates every one-qubit reduced state and so preserves its entropy, which is 11 bit for each qubit of GHZ and 00 for the product state ∣+00⟩\lvert{+}00\rangle. The meter below runs the reference solution: every single-qubit cut of the GHZ state has entanglement entropy 11 bit, the maximum for one qubit.

Quick check

In the lab template you set θ1,0,0=π2\theta_{1,0,0}=\tfrac\pi2 and, in addition, the RZR_Z angle θ1,0,1=π\theta_{1,0,1}=\pi on the same qubit and layer. What is the fidelity with ∣GHZ⟩\lvert\mathrm{GHZ}\rangle?

tolerance ±0.000001

Expressibility and the Haar fidelity distribution#

A layer of the HEA can be repeated as often as the hardware allows. How much does each additional layer add? Sim, Johnson and Aspuru-Guzik (2019) proposed to measure how uniformly an ansatz covers the state space: draw two independent parameter vectors, compute the fidelity of the two output states, and compare the distribution of this fidelity with the one obtained for Haar-random states, i.e. states drawn from the unique probability distribution on unit vectors that is invariant under every unitary. A circuit whose states cover the state space uniformly reproduces the Haar distribution; a circuit confined to a small region produces fidelities that pile up near 11.

Definition

Fidelity distribution and expressibility

Let ∣ψ(θ)⟩\lvert\psi(\theta)\rangle be the output of a parameterised circuit with PP parameters on nn qubits, d=2nd=2^n, and let θ,θ′\theta,\theta' be independent and uniform on [0,2π)P[0,2\pi)^P. The fidelity distribution PPQC(F)P_{\mathrm{PQC}}(F) is the law of F=∣⟨ψ(θ)∣ψ(θ′)⟩∣2∈[0,1]F=\lvert\langle\psi(\theta)\vert\psi(\theta')\rangle\rvert^2\in[0,1]. Split [0,1][0,1] into BB bins Ij=[jB,j+1B)I_j=\big[\frac jB,\frac{j+1}B\big), j=0,…,B−1j=0,\dots,B-1 (the last bin closed), sample NN pairs and let P^j\hat P_j be the fraction of fidelities in IjI_j. With the Haar bin probabilities Qj=Pr⁡Haar(F∈Ij)Q_j=\Pr_{\mathrm{Haar}}(F\in I_j), the expressibility is the Kullback–Leibler divergence Expr=DKL(P^ ∥ Q)=∑j=0B−1P^jln⁡P^jQj,0ln⁡0:=0.\mathrm{Expr}=D_{KL}\big(\hat P\,\Vert\,Q\big)=\sum_{j=0}^{B-1}\hat P_j\ln\frac{\hat P_j}{Q_j},\qquad0\ln0:=0 . Smaller values mean a fidelity distribution closer to Haar, i.e. a more expressive circuit.

Lemma

The Kullback–Leibler divergence is nonnegative

If Qj>0Q_j>0 for all jj and P^\hat P is a probability vector, then DKL(P^ ∥ Q)≥0D_{KL}(\hat P\,\Vert\,Q)\ge0, with equality if and only if P^=Q\hat P=Q.

Proof (4 steps)
  1. Step 1:
    g(t)=t−1−ln⁡t,g′(t)=1−1t,g(1)=0 ⇒ ln⁡t≤t−1  (t>0), equality iff t=1g(t)=t-1-\ln t,\quad g'(t)=1-\tfrac1t,\quad g(1)=0\ \Rightarrow\ \ln t\le t-1\ \ (t>0),\ \text{equality iff }t=1
    gg decreases on (0,1)(0,1) and increases on (1,∞)(1,\infty), so its minimum is g(1)=0g(1)=0.
  2. Step 2:
    −DKL=∑j:P^j>0P^jln⁡QjP^j≤∑j:P^j>0P^j(QjP^j−1)=∑j:P^j>0Qj−1≤0-D_{KL}=\sum_{j:\hat P_j>0}\hat P_j\ln\frac{Q_j}{\hat P_j}\le\sum_{j:\hat P_j>0}\hat P_j\Big(\frac{Q_j}{\hat P_j}-1\Big)=\sum_{j:\hat P_j>0}Q_j-1\le0
    Apply the inequality with t=Qj/P^jt=Q_j/\hat P_j; the bins with P^j=0\hat P_j=0 contribute nothing by the convention 0ln⁡0=00\ln0=0, and the QjQ_j over a subset of bins sum to at most 11.
  3. Step 3:
    DKL=0  ⟺  Qj=P^j ∀j:P^j>0 and ∑j:P^j>0Qj=1  ⟺  P^=QD_{KL}=0\iff Q_j=\hat P_j\ \forall j:\hat P_j>0\ \text{and}\ \sum_{j:\hat P_j>0}Q_j=1\iff\hat P=Q
    Both inequalities must be equalities; the second forces Qj=0Q_j=0 outside the support of P^\hat P, which together with the first gives P^=Q\hat P=Q.
  4. Step 4:
    Qj>0 ∀j ⇒ DKL<∞Q_j>0\ \forall j\ \Rightarrow\ D_{KL}<\infty
    Every term is finite; the theorem below shows that all Haar bin probabilities are positive.

∎

Lemma

Normalised complex Gaussian vectors are Haar-random states

Let g=(g1,…,gd)∈Cdg=(g_1,\dots,g_d)\in\mathbb C^d have independent real and imaginary parts, all distributed as N(0,12)\mathcal N(0,\tfrac12). Then ∣ψ⟩=g/∥g∥\lvert\psi\rangle=g/\lVert g\rVert is Haar-distributed: U∣ψ⟩U\lvert\psi\rangle has the same law as ∣ψ⟩\lvert\psi\rangle for every unitary UU. Consequently, for independent Haar-random ∣ψ⟩,∣ϕ⟩\lvert\psi\rangle,\lvert\phi\rangle, F=∣⟨ϕ∣ψ⟩∣2 =d ∣g1∣2∣g1∣2+⋯+∣gd∣2.F=\lvert\langle\phi\vert\psi\rangle\rvert^2\ \overset{d}{=}\ \frac{\lvert g_1\rvert^2}{\lvert g_1\rvert^2+\cdots+\lvert g_d\rvert^2}.

Proof (6 steps)
  1. Step 1:
    p(g)=∏k=1d1πe−(Re gk)21πe−(Im gk)2=π−de−∥g∥2p(g)=\prod_{k=1}^d\frac1{\sqrt\pi}e^{-(\mathrm{Re}\,g_k)^2}\frac1{\sqrt\pi}e^{-(\mathrm{Im}\,g_k)^2}=\pi^{-d}e^{-\lVert g\rVert^2}
    The density of N(0,12)\mathcal N(0,\frac12) is 1πe−x2\frac1{\sqrt\pi}e^{-x^2}; independence multiplies the 2d2d densities on Cd≅R2d\mathbb C^d\cong\mathbb R^{2d}.
  2. Step 2:
    ∥Ug∥=∥g∥ ⇒ g↦Ug is an orthogonal map of R2d,∣det⁡∣=1\lVert Ug\rVert=\lVert g\rVert\ \Rightarrow\ g\mapsto Ug\ \text{is an orthogonal map of }\mathbb R^{2d},\quad\lvert\det\rvert=1
    UU is real-linear on R2d\mathbb R^{2d} and preserves the Euclidean norm ∥g∥2=∑k(Re gk)2+(Im gk)2\lVert g\rVert^2=\sum_k(\mathrm{Re}\,g_k)^2+(\mathrm{Im}\,g_k)^2, so it is orthogonal and its Jacobian has modulus 11.
  3. Step 3:
    pUg(h)=p(U†h) ∣det⁡U†∣=π−de−∥U†h∥2=π−de−∥h∥2=p(h)p_{Ug}(h)=p(U^\dagger h)\,\lvert\det U^\dagger\rvert=\pi^{-d}e^{-\lVert U^\dagger h\rVert^2}=\pi^{-d}e^{-\lVert h\rVert^2}=p(h)
    Change of variables for the invertible linear map g↦Ugg\mapsto Ug: the vector UgUg has the same density as gg.
  4. Step 4:
    Pr⁡(g=0)=0,Ug∥g∥=Ug∥Ug∥ =d g∥g∥\Pr(g=0)=0,\qquad U\frac{g}{\lVert g\rVert}=\frac{Ug}{\lVert Ug\rVert}\ \overset{d}{=}\ \frac{g}{\lVert g\rVert}
    Normalisation is defined almost surely and commutes with UU; equal laws of UgUg and gg give equal laws of the normalised vectors.
  5. Step 5:
    μ unitarily invariant probability measure on S2d−1⊂Cd ⇒ μ=μHaar\mu\ \text{unitarily invariant probability measure on }S^{2d-1}\subset\mathbb C^d\ \Rightarrow\ \mu=\mu_{\mathrm{Haar}}
    The unitary group acts transitively on the unit sphere, and a transitive action of a compact group admits exactly one invariant probability measure; it is by definition the distribution of Haar-random states.
  6. Step 6:
    F=∣⟨ϕ∣ψ⟩∣2=∣⟨0n∣Vϕψ⟩∣2,Vϕ∣ϕ⟩=∣0n⟩,Vϕ∣ψ⟩=d∣ψ⟩F=\lvert\langle\phi\vert\psi\rangle\rvert^2=\lvert\langle0^n\vert V_\phi\psi\rangle\rvert^2,\qquad V_\phi\lvert\phi\rangle=\lvert0^n\rangle,\qquad V_\phi\lvert\psi\rangle\overset{d}{=}\lvert\psi\rangle
    Condition on ∣ϕ⟩\lvert\phi\rangle and choose a unitary VϕV_\phi that maps it to the first basis vector; by invariance and independence, Vϕ∣ψ⟩V_\phi\lvert\psi\rangle is again Haar-random. Hence FF has the law of ∣ψ1∣2=∣g1∣2/∥g∥2\lvert\psi_1\rvert^2=\lvert g_1\rvert^2/\lVert g\rVert^2.

∎

Lemma

Squared moduli of complex Gaussians are exponential

The variables Ek=∣gk∣2E_k=\lvert g_k\rvert^2, k=1,…,dk=1,\dots,d, are independent with Pr⁡(Ek>t)=e−t\Pr(E_k>t)=e^{-t} for t≥0t\ge0: each is Exp(1)=Gamma(1,1)\mathrm{Exp}(1)=\mathrm{Gamma}(1,1), with density e−te^{-t}.

Proof (4 steps)
  1. Step 1:
    Ek=(Re gk)2+(Im gk)2E_k=(\mathrm{Re}\,g_k)^2+(\mathrm{Im}\,g_k)^2
    Each EkE_k is a function of the pair (Re gk,Im gk)(\mathrm{Re}\,g_k,\mathrm{Im}\,g_k), and these pairs are independent for different kk, so the EkE_k are independent.
  2. Step 2:
    Pr⁡(Ek>t)=∬a2+b2>t1πe−(a2+b2) da db\Pr(E_k>t)=\iint_{a^2+b^2>t}\frac1\pi e^{-(a^2+b^2)}\,\mathrm da\,\mathrm db
    Integrate the joint density of the real and imaginary parts over the exterior of the disc of radius t\sqrt t.
  3. Step 3:
    =1π∫02π ⁣dα∫t∞e−r2 r dr=1π⋅2π⋅[−12e−r2]t∞=e−t=\frac1\pi\int_0^{2\pi}\!\mathrm d\alpha\int_{\sqrt t}^{\infty}e^{-r^2}\,r\,\mathrm dr=\frac1\pi\cdot2\pi\cdot\Big[-\tfrac12e^{-r^2}\Big]_{\sqrt t}^{\infty}=e^{-t}
    Polar coordinates a=rcos⁡αa=r\cos\alpha, b=rsin⁡αb=r\sin\alpha with da db=r dr dα\mathrm da\,\mathrm db=r\,\mathrm dr\,\mathrm d\alpha.
  4. Step 4:
    fEk(t)=−ddte−t=e−t=t1−1e−tΓ(1)f_{E_k}(t)=-\frac{\mathrm d}{\mathrm dt}e^{-t}=e^{-t}=\frac{t^{1-1}e^{-t}}{\Gamma(1)}
    Differentiate the survival function; Γ(1)=1\Gamma(1)=1, so this is the Gamma(1,1)\mathrm{Gamma}(1,1) density.

∎

Lemma

Gamma ratios are Beta-distributed (the Dirichlet aggregation)

Let X∼Gamma(a,1)X\sim\mathrm{Gamma}(a,1) and Y∼Gamma(b,1)Y\sim\mathrm{Gamma}(b,1) be independent, with densities xa−1e−xΓ(a)\frac{x^{a-1}e^{-x}}{\Gamma(a)} and yb−1e−yΓ(b)\frac{y^{b-1}e^{-y}}{\Gamma(b)} on (0,∞)(0,\infty), a,b>0a,b>0. Then T=X+YT=X+Y and U=XX+YU=\frac X{X+Y} are independent, T∼Gamma(a+b,1)T\sim\mathrm{Gamma}(a+b,1), and U∼Beta(a,b)U\sim\mathrm{Beta}(a,b) with density fU(u)=Γ(a+b)Γ(a)Γ(b) ua−1(1−u)b−1,0<u<1.f_U(u)=\frac{\Gamma(a+b)}{\Gamma(a)\Gamma(b)}\,u^{a-1}(1-u)^{b-1},\qquad0<u<1 . In particular, by induction, E2+⋯+Ed∼Gamma(d−1,1)E_2+\cdots+E_d\sim\mathrm{Gamma}(d-1,1) for independent Ek∼Gamma(1,1)E_k\sim\mathrm{Gamma}(1,1).

Proof (8 steps)
  1. Step 1:
    fX,Y(x,y)=xa−1yb−1e−(x+y)Γ(a)Γ(b),x,y>0f_{X,Y}(x,y)=\frac{x^{a-1}y^{b-1}e^{-(x+y)}}{\Gamma(a)\Gamma(b)},\qquad x,y>0
    Independence: the joint density is the product of the marginal densities.
  2. Step 2:
    (u,t)↦(x,y)=(ut, (1−u)t),(0,1)×(0,∞)→(0,∞)2,t=x+y, u=xx+y(u,t)\mapsto(x,y)=\big(ut,\,(1-u)t\big),\qquad(0,1)\times(0,\infty)\to(0,\infty)^2,\qquad t=x+y,\ u=\tfrac x{x+y}
    A smooth bijection with the stated inverse.
  3. Step 3:
    ∂(x,y)∂(u,t)=det⁡(tu−t1−u)=t(1−u)+tu=t\frac{\partial(x,y)}{\partial(u,t)}=\det\begin{pmatrix}t&u\\-t&1-u\end{pmatrix}=t(1-u)+tu=t
    Jacobian of the change of variables.
  4. Step 4:
    fU,T(u,t)=fX,Y(ut,(1−u)t)⋅t=ua−1(1−u)b−1⏟h(u)⋅ta+b−1e−tΓ(a)Γ(b)⏟k(t)f_{U,T}(u,t)=f_{X,Y}\big(ut,(1-u)t\big)\cdot t=\underbrace{u^{a-1}(1-u)^{b-1}}_{h(u)}\cdot\underbrace{\frac{t^{a+b-1}e^{-t}}{\Gamma(a)\Gamma(b)}}_{k(t)}
    (ut)a−1((1−u)t)b−1 t=ua−1(1−u)b−1ta+b−1(ut)^{a-1}\big((1-u)t\big)^{b-1}\,t=u^{a-1}(1-u)^{b-1}t^{a+b-1} and e−(ut+(1−u)t)=e−te^{-(ut+(1-u)t)}=e^{-t}.
  5. Step 5:
    fU,T=h(u) k(t) ⇒ U⊥T,fT(t)∝ta+b−1e−t ⇒ T∼Gamma(a+b,1)f_{U,T}=h(u)\,k(t)\ \Rightarrow\ U\perp T,\qquad f_T(t)\propto t^{a+b-1}e^{-t}\ \Rightarrow\ T\sim\mathrm{Gamma}(a+b,1)
    A joint density that factorises into a function of uu times a function of tt makes the variables independent with marginals proportional to the factors; the only normalised density proportional to ta+b−1e−tt^{a+b-1}e^{-t} is the Gamma(a+b,1)\mathrm{Gamma}(a+b,1) density.
  6. Step 6:
    1=∫01h(u) du∫0∞k(t) dt=∫01h(u) du⋅Γ(a+b)Γ(a)Γ(b) ⇒ ∫01ua−1(1−u)b−1du=Γ(a)Γ(b)Γ(a+b)1=\int_0^1h(u)\,\mathrm du\int_0^\infty k(t)\,\mathrm dt=\int_0^1h(u)\,\mathrm du\cdot\frac{\Gamma(a+b)}{\Gamma(a)\Gamma(b)}\ \Rightarrow\ \int_0^1u^{a-1}(1-u)^{b-1}\mathrm du=\frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)}
    The joint density integrates to 11, and ∫0∞ta+b−1e−tdt=Γ(a+b)\int_0^\infty t^{a+b-1}e^{-t}\mathrm dt=\Gamma(a+b) by definition of Γ\Gamma. This evaluates the Beta integral B(a,b)B(a,b).
  7. Step 7:
    fU(u)=h(u)∫01h=Γ(a+b)Γ(a)Γ(b)ua−1(1−u)b−1f_U(u)=\frac{h(u)}{\int_0^1h}=\frac{\Gamma(a+b)}{\Gamma(a)\Gamma(b)}u^{a-1}(1-u)^{b-1}
    Normalise the factor hh.
  8. Step 8:
    E2+⋯+Em+1∼Gamma(m,1) ⇒ (E2+⋯+Em+1)+Em+2∼Gamma(m+1,1)E_2+\cdots+E_{m+1}\sim\mathrm{Gamma}(m,1)\ \Rightarrow\ (E_2+\cdots+E_{m+1})+E_{m+2}\sim\mathrm{Gamma}(m+1,1)
    Induction on the number of summands with a=ma=m, b=1b=1, using that Em+2E_{m+2} is independent of the earlier ones. The vector (E1,…,Ed)/∑kEk(E_1,\dots,E_d)/\sum_kE_k is the Dirichlet(1,…,1)(1,\dots,1) distribution, and the lemma aggregates its coordinates.

∎

Theorem

Fidelity density of Haar-random states

For independent Haar-random pure states in Cd\mathbb C^d (d=2nd=2^n for nn qubits), F=∣⟨ϕ∣ψ⟩∣2∼Beta(1,d−1)F=\lvert\langle\phi\vert\psi\rangle\rvert^2\sim\mathrm{Beta}(1,d-1): PHaar(F)=(d−1)(1−F)d−2,Pr⁡(F≤f)=1−(1−f)d−1,0≤f≤1.P_{\mathrm{Haar}}(F)=(d-1)(1-F)^{d-2},\qquad\Pr(F\le f)=1-(1-f)^{d-1},\qquad0\le f\le1 . The Haar bin probabilities are Qj=(1−jB)d−1−(1−j+1B)d−1>0Q_j=\big(1-\frac jB\big)^{d-1}-\big(1-\frac{j+1}B\big)^{d-1}>0, and E[F]=1d\mathbb E[F]=\frac1d, E[F2]=2d(d+1)\mathbb E[F^2]=\frac2{d(d+1)}. For one qubit (d=2d=2) the density is uniform on [0,1][0,1].

Proof (8 steps)
  1. Step 1:
    F=dE1E1+S,S=E2+⋯+EdF\overset{d}{=}\frac{E_1}{E_1+S},\qquad S=E_2+\cdots+E_d
    The Gaussian lemma, with Ek=∣gk∣2E_k=\lvert g_k\rvert^2.
  2. Step 2:
    E1∼Gamma(1,1),S∼Gamma(d−1,1),E1⊥SE_1\sim\mathrm{Gamma}(1,1),\qquad S\sim\mathrm{Gamma}(d-1,1),\qquad E_1\perp S
    The exponential lemma and the induction of the Gamma–Beta lemma; SS depends only on g2,…,gdg_2,\dots,g_d.
  3. Step 3:
    PHaar(F)=Γ(d)Γ(1)Γ(d−1)F0(1−F)d−2=(d−1)!(d−2)!(1−F)d−2=(d−1)(1−F)d−2P_{\mathrm{Haar}}(F)=\frac{\Gamma(d)}{\Gamma(1)\Gamma(d-1)}F^{0}(1-F)^{d-2}=\frac{(d-1)!}{(d-2)!}(1-F)^{d-2}=(d-1)(1-F)^{d-2}
    The Gamma–Beta lemma with a=1a=1, b=d−1b=d-1, and Γ(m)=(m−1)!\Gamma(m)=(m-1)! for integers.
  4. Step 4:
    Pr⁡(F≤f)=∫0f(d−1)(1−u)d−2 du=[−(1−u)d−1]0f=1−(1−f)d−1\Pr(F\le f)=\int_0^f(d-1)(1-u)^{d-2}\,\mathrm du=\Big[-(1-u)^{d-1}\Big]_0^f=1-(1-f)^{d-1}
    Integrate the density; the antiderivative of (d−1)(1−u)d−2(d-1)(1-u)^{d-2} is −(1−u)d−1-(1-u)^{d-1}.
  5. Step 5:
    Qj=Pr⁡(F≤j+1B)−Pr⁡(F≤jB)=(1−jB)d−1−(1−j+1B)d−1>0Q_j=\Pr\Big(F\le\tfrac{j+1}B\Big)-\Pr\Big(F\le\tfrac jB\Big)=\Big(1-\frac jB\Big)^{d-1}-\Big(1-\frac{j+1}B\Big)^{d-1}>0
    Difference of the CDF at the bin edges (single points have probability 00). It is positive because 1−jB>1−j+1B≥01-\frac jB>1-\frac{j+1}B\ge0 and s↦sd−1s\mapsto s^{d-1} is strictly increasing on [0,1][0,1] for d≥2d\ge2.
  6. Step 6:
    E[F]=∫01Pr⁡(F>f) df=∫01(1−f)d−1 df=1d\mathbb E[F]=\int_0^1\Pr(F>f)\,\mathrm df=\int_0^1(1-f)^{d-1}\,\mathrm df=\frac1d
    For a variable with values in [0,1][0,1], E[F]=∫01Pr⁡(F>f) df\mathbb E[F]=\int_0^1\Pr(F>f)\,\mathrm df (write F=∫011f<F dfF=\int_0^1\mathbf 1_{f<F}\,\mathrm df and exchange expectation and integral). For nn qubits this is 2−n2^{-n}, the mean kernel value of the concentration analysis in the quantum-kernel lesson.
  7. Step 7:
    E[F2]=∫012fPr⁡(F>f) df=2∫01f(1−f)d−1 df=2∫01(1−s)sd−1 ds=2(1d−1d+1)=2d(d+1)\mathbb E[F^2]=\int_0^12f\Pr(F>f)\,\mathrm df=2\int_0^1f(1-f)^{d-1}\,\mathrm df=2\int_0^1(1-s)s^{d-1}\,\mathrm ds=2\Big(\frac1d-\frac1{d+1}\Big)=\frac2{d(d+1)}
    The same argument with F2=∫012f 1f<F dfF^2=\int_0^12f\,\mathbf 1_{f<F}\,\mathrm df, then the substitution s=1−fs=1-f.
  8. Step 8:
    d=2:PHaar(F)=1⋅(1−F)0=1,Qj=1B,E[F]=12,E[F2]=13d=2:\quad P_{\mathrm{Haar}}(F)=1\cdot(1-F)^0=1,\qquad Q_j=\frac1B,\qquad\mathbb E[F]=\frac12,\quad\mathbb E[F^2]=\frac13
    One qubit: every fidelity in [0,1][0,1] is equally likely, which the next derivation confirms geometrically.

∎

Derivation

Check for one qubit: Archimedes on the Bloch sphere

  1. Step 1:
    F=∣⟨ϕ∣ψ⟩∣2=12(1+r⃗ϕ⋅r⃗ψ)F=\lvert\langle\phi\vert\psi\rangle\rvert^2=\tfrac12\big(1+\vec r_\phi\cdot\vec r_\psi\big)
    For pure qubit states Tr(ρσ)=12(1+r⃗⋅s⃗)\mathrm{Tr}(\rho\sigma)=\frac12(1+\vec r\cdot\vec s) (quantum-kernel lesson, Bloch lemma).
  2. Step 2:
    r⃗ϕ=e⃗z ⇒ F=12(1+z),z=cos⁡ϑ\vec r_\phi=\vec e_z\ \Rightarrow\ F=\tfrac12(1+z),\qquad z=\cos\vartheta
    Haar measure is invariant under rotations of the Bloch sphere, so we may rotate r⃗ϕ\vec r_\phi to the north pole; zz is the height of r⃗ψ\vec r_\psi.
  3. Step 3:
    dA=sin⁡ϑ dϑ dα=∣dz∣ dα ⇒ Area{z1≤z≤z2}=2π(z2−z1)\mathrm dA=\sin\vartheta\,\mathrm d\vartheta\,\mathrm d\alpha=\lvert\mathrm dz\rvert\,\mathrm d\alpha\ \Rightarrow\ \mathrm{Area}\{z_1\le z\le z_2\}=2\pi(z_2-z_1)
    The area element of the unit sphere in spherical coordinates; with z=cos⁡ϑz=\cos\vartheta, dz=−sin⁡ϑ dϑ\mathrm dz=-\sin\vartheta\,\mathrm d\vartheta (Archimedes' hat-box theorem).
  4. Step 4:
    z∼U[−1,1] ⇒ F=12(1+z)∼U[0,1]z\sim\mathcal U[-1,1]\ \Rightarrow\ F=\tfrac12(1+z)\sim\mathcal U[0,1]
    The uniform distribution on the sphere gives every zone of equal height equal probability, and an affine map of a uniform variable is uniform: the same result as (d−1)(1−F)d−2(d-1)(1-F)^{d-2} at d=2d=2.

Result

d=2:PHaar(F)=1on [0,1]\boxed{d=2:\qquad P_{\mathrm{Haar}}(F)=1\quad\text{on }[0,1]}
Theorem

Second fidelity moment of the one-qubit HEA

For the one-qubit HEA ∣ψ(θ)⟩=∏lRZ(θl,1)RY(θl,0)∣0⟩\lvert\psi(\theta)\rangle=\prod_{l}R_Z(\theta_{l,1})R_Y(\theta_{l,0})\lvert0\rangle with LL layers and i.i.d. uniform angles on [0,2π)[0,2\pi), E[F]=12,E[F2]=13+16−L6.\mathbb E[F]=\frac12,\qquad\mathbb E[F^2]=\frac13+\frac{16^{-L}}{6}. The second moment exceeds the Haar value 13\frac13 for every LL, so the fidelity distribution is never exactly Haar, but the gap shrinks by a factor 1616 per layer. For L=1L=1, E[F2]=1132=0.34375\mathbb E[F^2]=\frac{11}{32}=0.34375.

Proof (11 steps)
  1. Step 1:
    F=12(1+r⃗⋅r⃗ ′),E[F2]=14(1+2 E[r⃗]⋅E[r⃗ ′]+E[(r⃗⋅r⃗ ′)2])F=\tfrac12\big(1+\vec r\cdot\vec r\,'\big),\qquad\mathbb E[F^2]=\tfrac14\Big(1+2\,\mathbb E[\vec r]\cdot\mathbb E[\vec r\,']+\mathbb E\big[(\vec r\cdot\vec r\,')^2\big]\Big)
    r⃗,r⃗ ′\vec r,\vec r\,' are the Bloch vectors of the two independent output states; expand the square and use independence for the cross term.
  2. Step 2:
    E[(r⃗⋅r⃗ ′)2]=∑a,bE[rarb] E[ra′rb′]=∥M∥F2,M=E[r⃗ r⃗ ⊤]\mathbb E\big[(\vec r\cdot\vec r\,')^2\big]=\sum_{a,b}\mathbb E[r_ar_b]\,\mathbb E[r'_ar'_b]=\lVert M\rVert_F^2,\qquad M=\mathbb E\big[\vec r\,\vec r^{\,\top}\big]
    (r⃗⋅r⃗ ′)2=∑a,brarbra′rb′(\vec r\cdot\vec r\,')^2=\sum_{a,b}r_ar_br'_ar'_b; independence factorises each term, and r⃗ ′\vec r\,' has the same law as r⃗\vec r.
  3. Step 3:
    RZ(β): (xyz)↦(xcos⁡β−ysin⁡βxsin⁡β+ycos⁡βz),RY(α) rotates (z,x) by αR_Z(\beta):\ \begin{pmatrix}x\\y\\z\end{pmatrix}\mapsto\begin{pmatrix}x\cos\beta-y\sin\beta\\x\sin\beta+y\cos\beta\\z\end{pmatrix},\qquad R_Y(\alpha)\ \text{rotates }(z,x)\text{ by }\alpha
    RP(α)R_P(\alpha) rotates the Bloch vector by the angle α\alpha about the axis PP (Bloch-sphere lessons of Stage 2).
  4. Step 4:
    Eβ[cos⁡β]=Eβ[sin⁡β]=Eβ[sin⁡βcos⁡β]=0,Eβ[cos⁡2β]=Eβ[sin⁡2β]=12\mathbb E_\beta\big[\cos\beta\big]=\mathbb E_\beta\big[\sin\beta\big]=\mathbb E_\beta\big[\sin\beta\cos\beta\big]=0,\qquad\mathbb E_\beta\big[\cos^2\beta\big]=\mathbb E_\beta\big[\sin^2\beta\big]=\tfrac12
    Averages over a full period for β∼U[0,2π)\beta\sim\mathcal U[0,2\pi); sin⁡βcos⁡β=12sin⁡2β\sin\beta\cos\beta=\frac12\sin2\beta and cos⁡2β=1+cos⁡2β2\cos^2\beta=\frac{1+\cos2\beta}2.
  5. Step 5:
    M=diag(mx,my,mz) → EβRZ  diag(mx+my2,mx+my2,mz),→ EαRY  diag(mx+mz2,my,mx+mz2)M=\mathrm{diag}(m_x,m_y,m_z)\ \xrightarrow{\ \mathbb E_\beta R_Z\ }\ \mathrm{diag}\Big(\tfrac{m_x+m_y}2,\tfrac{m_x+m_y}2,m_z\Big),\qquad\xrightarrow{\ \mathbb E_\alpha R_Y\ }\ \mathrm{diag}\Big(\tfrac{m_x+m_z}2,m_y,\tfrac{m_x+m_z}2\Big)
    For an independent uniform angle: E[(xcos⁡β−ysin⁡β)2]=mx+my2\mathbb E[(x\cos\beta-y\sin\beta)^2]=\frac{m_x+m_y}2, and the off-diagonal averages, e.g. E[(xcos⁡β−ysin⁡β)(xsin⁡β+ycos⁡β)]\mathbb E[(x\cos\beta-y\sin\beta)(x\sin\beta+y\cos\beta)] and E[(xcos⁡β−ysin⁡β)z]\mathbb E[(x\cos\beta-y\sin\beta)z], vanish by the previous line. So a diagonal MM stays diagonal.
  6. Step 6:
    E[r⃗]=Eβ[RZ(β)] Eα[RY(α)e⃗z]=Eβ[RZ(β)](Esin⁡α,0,Ecos⁡α)⊤=0 ⇒ E[F]=12\mathbb E[\vec r]=\mathbb E_\beta\big[R_Z(\beta)\big]\,\mathbb E_\alpha\big[R_Y(\alpha)\vec e_z\big]=\mathbb E_\beta\big[R_Z(\beta)\big](\mathbb E\sin\alpha,0,\mathbb E\cos\alpha)^\top=0\ \Rightarrow\ \mathbb E[F]=\tfrac12
    Already the first RYR_Y rotation with a uniform angle averages the Bloch vector to zero, and later rotations are linear. Then E[F]=12(1+E[r⃗]⋅E[r⃗ ′])=12\mathbb E[F]=\frac12(1+\mathbb E[\vec r]\cdot\mathbb E[\vec r\,'])=\frac12, the Haar value.
  7. Step 7:
    M0=diag(0,0,1),Ml=diag(pl,pl,ql),2pl+ql=E∣r⃗∣2=1,ql+1=pl+ql2=1+ql4M_0=\mathrm{diag}(0,0,1),\quad M_l=\mathrm{diag}(p_l,p_l,q_l),\quad2p_l+q_l=\mathbb E\lvert\vec r\rvert^2=1,\quad q_{l+1}=\frac{p_l+q_l}2=\frac{1+q_l}4
    Start from ∣0⟩\lvert0\rangle, r⃗=e⃗z\vec r=\vec e_z. One layer applies RYR_Y (mixing xx and zz: both become pl+ql2\frac{p_l+q_l}2) and then RZR_Z (mixing xx and yy), so the zz entry after the layer is pl+ql2\frac{p_l+q_l}2, and pl=1−ql2p_l=\frac{1-q_l}2 because pure states have ∣r⃗∣=1\lvert\vec r\rvert=1.
  8. Step 8:
    ql−13=14(ql−1−13) ⇒ qL=13+23 4−Lq_{l}-\tfrac13=\tfrac14\big(q_{l-1}-\tfrac13\big)\ \Rightarrow\ q_L=\tfrac13+\tfrac23\,4^{-L}
    Subtract the fixed point q∗=13q^\ast=\frac13 of q↦1+q4q\mapsto\frac{1+q}4 and iterate from q0=1q_0=1.
  9. Step 9:
    ∥ML∥F2=2pL2+qL2=(1−qL)22+qL2=13+32(qL−13)2=13+23 16−L\lVert M_L\rVert_F^2=2p_L^2+q_L^2=\frac{(1-q_L)^2}2+q_L^2=\frac13+\frac32\Big(q_L-\frac13\Big)^2=\frac13+\frac23\,16^{-L}
    Both sides equal 12−qL+32qL2\frac12-q_L+\frac32q_L^2; insert qL−13=23 4−Lq_L-\frac13=\frac23\,4^{-L}, so 32⋅49⋅16−L=23 16−L\frac32\cdot\frac49\cdot16^{-L}=\frac23\,16^{-L}.
  10. Step 10:
    E[F2]=14(1+13+23 16−L)=13+16−L6,L=1: 13+196=1132\mathbb E[F^2]=\tfrac14\Big(1+\tfrac13+\tfrac23\,16^{-L}\Big)=\frac13+\frac{16^{-L}}6,\qquad L=1:\ \frac13+\frac1{96}=\frac{11}{32}
    Insert into the first line with E[r⃗]=0\mathbb E[\vec r]=0. For Haar states M=13IM=\frac13I and ∥M∥F2=13\lVert M\rVert_F^2=\frac13, which reproduces E[F2]=13\mathbb E[F^2]=\frac13.
  11. Step 11:
    L=1:z=cos⁡θ0,0,Pr⁡(z≤ζ)=1−arccos⁡ζπ,fz(ζ)=1π1−ζ2≠12L=1:\quad z=\cos\theta_{0,0},\quad\Pr(z\le\zeta)=1-\frac{\arccos\zeta}{\pi},\quad f_z(\zeta)=\frac1{\pi\sqrt{1-\zeta^2}}\ne\frac12
    Where the excess comes from: uniform angles make the height of the Bloch vector arcsine-distributed (density largest near the poles ζ=±1\zeta=\pm1) instead of uniform, so very similar and very different pairs of states are over-represented.

∎

Derivation

The finite-sample floor of the KL estimate

  1. Step 1:
    P^j=Qj+εj,∑jεj=0\hat P_j=Q_j+\varepsilon_j,\qquad\sum_j\varepsilon_j=0
    Even a perfect Haar sampler produces a histogram that fluctuates around QQ; both vectors sum to 11.
  2. Step 2:
    (Q+ε)ln⁡(1+εQ)=(Q+ε)(εQ−ε22Q2+O(ε3))=ε+ε22Q+O(ε3)(Q+\varepsilon)\ln\Big(1+\frac\varepsilon Q\Big)=(Q+\varepsilon)\Big(\frac\varepsilon Q-\frac{\varepsilon^2}{2Q^2}+O(\varepsilon^3)\Big)=\varepsilon+\frac{\varepsilon^2}{2Q}+O(\varepsilon^3)
    Taylor expansion ln⁡(1+u)=u−u22+O(u3)\ln(1+u)=u-\frac{u^2}2+O(u^3) with u=ε/Qu=\varepsilon/Q.
  3. Step 3:
    DKL(P^ ∥ Q)=∑jεj22Qj+O(max⁡j∣εj∣3)D_{KL}(\hat P\,\Vert\,Q)=\sum_j\frac{\varepsilon_j^2}{2Q_j}+O\big(\max_j\lvert\varepsilon_j\rvert^3\big)
    Sum over the bins; the linear terms add up to ∑jεj=0\sum_j\varepsilon_j=0.
  4. Step 4:
    E[εj2]=Var[P^j]=Qj(1−Qj)N\mathbb E\big[\varepsilon_j^2\big]=\mathrm{Var}\big[\hat P_j\big]=\frac{Q_j(1-Q_j)}N
    The count of bin jj among NN independent samples is binomial with success probability QjQ_j.
  5. Step 5:
    E[DKL]≈∑j=0B−11−Qj2N=B−12N,B=20, N=2×104:  4.75×10−4\mathbb E\big[D_{KL}\big]\approx\sum_{j=0}^{B-1}\frac{1-Q_j}{2N}=\frac{B-1}{2N},\qquad B=20,\ N=2\times10^4:\ \ 4.75\times10^{-4}
    Leading order in 1/N1/N, using ∑jQj=1\sum_jQ_j=1. A KL estimate of order 5×10−45\times10^{-4} with these settings cannot be distinguished from Haar: only differences well above this floor are meaningful.

Result

ExprHaar sampler≈B−12N\boxed{\mathrm{Expr}_{\text{Haar sampler}}\approx\frac{B-1}{2N}}

The gallery samples random parameter vectors of the chosen ansatz (here the HEA with 33 qubits and 22 layers), draws the histogram of pairwise fidelities against the Haar bin probabilities QjQ_j for d=8d=8, and reports the KL divergence. For d=8d=8 the Haar density 7(1−F)67(1-F)^6 is concentrated near F=0F=0 with mean 18\frac18. Increase the number of layers and watch the histogram approach the Haar curve. With a single layer the circuit only produces the images under EE of product states; since EE is unitary, F=F0F1F2F=F_0F_1F_2 is then a product of three independent one-qubit HEA fidelities, so by the one-qubit moment theorem above E[F]=(12)3=18\mathbb E[F]=\big(\tfrac12\big)^3=\tfrac18 agrees with Haar while E[F2]=(1132)3≈0.0406\mathbb E[F^2]=\big(\tfrac{11}{32}\big)^3\approx0.0406 is far above the Haar value 28⋅9≈0.0278\tfrac2{8\cdot9}\approx0.0278: the one-layer histogram has a much heavier tail towards large fidelities.

Proposition

Real-amplitude circuits have a different fidelity law

If ∣ψ⟩\lvert\psi\rangle is uniformly distributed on the real unit sphere of Rd\mathbb R^d and ∣ϕ⟩\lvert\phi\rangle is a fixed real unit vector, then F=∣⟨ϕ∣ψ⟩∣2∼Beta(12,d−12)F=\lvert\langle\phi\vert\psi\rangle\rvert^2\sim\mathrm{Beta}\big(\frac12,\frac{d-1}2\big), with E[F]=1d\mathbb E[F]=\frac1d and E[F2]=3d(d+2)>2d(d+1)\mathbb E[F^2]=\frac3{d(d+2)}>\frac2{d(d+1)}. For d=4d=4 the density is 2πF−1/2(1−F)1/2\frac2\pi F^{-1/2}(1-F)^{1/2}, unbounded at F=0F=0. More generally, for every distribution of real unit vectors, E[F2]≥3d(d+2)\mathbb E[F^2]\ge\frac3{d(d+2)}. Hence an ansatz built from RYR_Y rotations and CNOTs only, whose states have real amplitudes, can never reproduce the Haar fidelity distribution, at any depth and for any distribution of its parameters.

Proof (10 steps)
  1. Step 1:
    ∣ψ⟩=x∥x∥,xk∼i.i.d.N(0,12),F=dx12x12+⋯+xd2\lvert\psi\rangle=\frac{x}{\lVert x\rVert},\qquad x_k\overset{\text{i.i.d.}}{\sim}\mathcal N(0,\tfrac12),\qquad F\overset{d}{=}\frac{x_1^2}{x_1^2+\cdots+x_d^2}
    The proof of the Gaussian lemma with orthogonal instead of unitary maps: the real Gaussian density π−d/2e−∥x∥2\pi^{-d/2}e^{-\lVert x\rVert^2} is invariant under orthogonal maps, which act transitively on the real sphere.
  2. Step 2:
    Pr⁡(xk2≤t)=2π∫0te−s2ds ⇒ fxk2(t)=t−1/2e−tπ=t1/2−1e−tΓ(12)\Pr(x_k^2\le t)=\frac2{\sqrt\pi}\int_0^{\sqrt t}e^{-s^2}\mathrm ds\ \Rightarrow\ f_{x_k^2}(t)=\frac{t^{-1/2}e^{-t}}{\sqrt\pi}=\frac{t^{1/2-1}e^{-t}}{\Gamma(\frac12)}
    Differentiate with the chain rule (ddtt=12t\frac{\mathrm d}{\mathrm dt}\sqrt t=\frac1{2\sqrt t}) and use Γ(12)=∫0∞t−1/2e−tdt=2∫0∞e−s2ds=π\Gamma(\frac12)=\int_0^\infty t^{-1/2}e^{-t}\mathrm dt=2\int_0^\infty e^{-s^2}\mathrm ds=\sqrt\pi: each xk2x_k^2 is Gamma(12,1)\mathrm{Gamma}(\frac12,1).
  3. Step 3:
    x22+⋯+xd2∼Gamma(d−12,1) ⇒ F∼Beta(12,d−12)x_2^2+\cdots+x_d^2\sim\mathrm{Gamma}\big(\tfrac{d-1}2,1\big)\ \Rightarrow\ F\sim\mathrm{Beta}\big(\tfrac12,\tfrac{d-1}2\big)
    The Gamma–Beta lemma, first for the sum (induction) and then for the ratio with a=12a=\frac12, b=d−12b=\frac{d-1}2.
  4. Step 4:
    E[Fm]=∫01ua+m−1(1−u)b−1du∫01ua−1(1−u)b−1du=Γ(a+m)Γ(a+b)Γ(a)Γ(a+b+m)\mathbb E[F^m]=\frac{\int_0^1u^{a+m-1}(1-u)^{b-1}\mathrm du}{\int_0^1u^{a-1}(1-u)^{b-1}\mathrm du}=\frac{\Gamma(a+m)\Gamma(a+b)}{\Gamma(a)\Gamma(a+b+m)}
    Moments of the Beta law from the Beta integral Γ(a)Γ(b)Γ(a+b)\frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)} evaluated in the lemma.
  5. Step 5:
    E[F]=aa+b=1d,E[F2]=a(a+1)(a+b)(a+b+1)=12⋅32d2(d2+1)=3d(d+2)\mathbb E[F]=\frac a{a+b}=\frac1d,\qquad\mathbb E[F^2]=\frac{a(a+1)}{(a+b)(a+b+1)}=\frac{\frac12\cdot\frac32}{\frac d2\big(\frac d2+1\big)}=\frac3{d(d+2)}
    Γ(z+1)=zΓ(z)\Gamma(z+1)=z\Gamma(z) with a=12a=\frac12, a+b=d2a+b=\frac d2.
  6. Step 6:
    3d(d+2)−2d(d+1)=3(d+1)−2(d+2)d(d+1)(d+2)=d−1d(d+1)(d+2)>0(d≥2)\frac3{d(d+2)}-\frac2{d(d+1)}=\frac{3(d+1)-2(d+2)}{d(d+1)(d+2)}=\frac{d-1}{d(d+1)(d+2)}>0\quad(d\ge2)
    The second moments differ, so the two laws differ. Combined with the bound of the next steps, no distribution of real states has the complex Haar second moment.
  7. Step 7:
    E[F2]=E[(ϕ⊤ψ)4]=∑a,b,c,eTabce2=∥T∥2,T=E[ψ⊗ψ⊗ψ⊗ψ]\mathbb E[F^2]=\mathbb E\big[(\phi^\top\psi)^4\big]=\sum_{a,b,c,e}T_{abce}^2=\lVert T\rVert^2,\qquad T=\mathbb E\big[\psi\otimes\psi\otimes\psi\otimes\psi\big]
    Any law of real unit vectors, with ϕ,ψ\phi,\psi independent and identically distributed: expand (ϕ⊤ψ)4=∑ϕaϕbϕcϕeψaψbψcψe(\phi^\top\psi)^4=\sum\phi_a\phi_b\phi_c\phi_e\psi_a\psi_b\psi_c\psi_e and factorise the expectation.
  8. Step 8:
    Gabce=δabδce+δacδbe+δaeδbc,⟨T,G⟩=3 E[∥ψ∥4]=3,∥G∥2=3d2+6dG_{abce}=\delta_{ab}\delta_{ce}+\delta_{ac}\delta_{be}+\delta_{ae}\delta_{bc},\qquad\langle T,G\rangle=3\,\mathbb E\big[\lVert\psi\rVert^4\big]=3,\qquad\lVert G\rVert^2=3d^2+6d
    Each of the three terms of GG contracts ψ⊗4\psi^{\otimes4} to (ψ⊤ψ)2=1(\psi^\top\psi)^2=1. In ∥G∥2\lVert G\rVert^2 each term has squared norm ∑δabδce=d2\sum\delta_{ab}\delta_{ce}=d^2, and each of the 66 ordered cross products, e.g. ∑δabδceδacδbe\sum\delta_{ab}\delta_{ce}\delta_{ac}\delta_{be}, forces a=b=c=ea=b=c=e and equals dd.
  9. Step 9:
    ∥T∥2≥⟨T,G⟩2∥G∥2=93d(d+2)=3d(d+2)\lVert T\rVert^2\ge\frac{\langle T,G\rangle^2}{\lVert G\rVert^2}=\frac9{3d(d+2)}=\frac3{d(d+2)}
    Cauchy–Schwarz inequality in the Euclidean space of 44-index tensors. The uniform law on the real sphere attains the bound, as computed above.
  10. Step 10:
    d=4:Γ(2)Γ(12)Γ(32)=1π⋅π2=2π,f(F)=2πF−1/2(1−F)1/2d=4:\quad\frac{\Gamma(2)}{\Gamma(\frac12)\Gamma(\frac32)}=\frac1{\sqrt\pi\cdot\frac{\sqrt\pi}2}=\frac2\pi,\qquad f(F)=\frac2\pi F^{-1/2}(1-F)^{1/2}
    Γ(32)=12Γ(12)\Gamma(\frac32)=\frac12\Gamma(\frac12). The density diverges at F=0F=0, while the complex Haar density 3(1−F)23(1-F)^2 equals 33 there.

∎

The second gallery uses the R_Y–CNOT ansatz on 22 qubits with one layer. All its gates have real matrices, so every output state has real amplitudes and, by the proposition, its fidelity histogram cannot converge to the Haar curve 3(1−F)23(1-F)^2 however many layers are added. Compare its KL divergence with the HEA at the same size: the RZR_Z rotations of the HEA are what give access to complex relative phases.

The TypeScript cell estimates the expressibility of the one-qubit HEA with the platform's simulator. It samples N=20 000N=20\,000 pairs for L=1,2,3L=1,2,3 layers and for exact Haar states generated from Gaussian amplitudes (the Box–Muller transform produces ∣gk∣2=−ln⁡u\lvert g_k\rvert^2=-\ln u, the Exp(1)\mathrm{Exp}(1) variable of the lemma), and prints the sample moments next to the exact value 13+16−L6\frac13+\frac{16^{-L}}6 and the KL divergence next to the floor B−12N\frac{B-1}{2N}.

TypeScriptrunnable
import { Statevector, mulberry32, runStatevector } from '@/lib/quantum-simulator';

// One-qubit hardware-efficient ansatz with L layers: R_Y(a_l) then R_Z(b_l) per layer (n = 1: no CNOT ladder).
function hea1(angles: number[]): Statevector {
  const ops = [];
  for (let k = 0; k < angles.length; k += 2) {
    ops.push({ kind: 'gate' as const, gate: 'RY' as const, targets: [0], params: [angles[k]] });
    ops.push({ kind: 'gate' as const, gate: 'RZ' as const, targets: [0], params: [angles[k + 1]] });
  }
  return runStatevector({ version: 1, numQubits: 1, numClbits: 0, ops });
}

const B = 20;      // histogram bins on [0, 1]
const N = 20000;   // sampled pairs of states per ansatz
const d = 2;       // Hilbert-space dimension of one qubit
const rng = mulberry32(44);

// Haar bin probabilities from the CDF 1 - (1 - F)^(d - 1); for d = 2 every bin has probability 1/B.
const haarBin = (j: number) => (1 - j / B) ** (d - 1) - (1 - (j + 1) / B) ** (d - 1);

// Haar-random qubit state g/|g| with Re g_k, Im g_k ~ N(0, 1/2): Box-Muller with |g_k|^2 = -ln u ~ Exp(1).
function haar1(): Statevector {
  const amp = () => {
    const r = Math.sqrt(-Math.log(1 - rng()));
    const phi = 2 * Math.PI * rng();
    return { re: r * Math.cos(phi), im: r * Math.sin(phi) };
  };
  return Statevector.fromAmplitudes([amp(), amp()], { normalize: true });
}

// Fidelity histogram of N independent pairs, its moments, and the expressibility KL(P_hat || P_Haar).
function expressibility(label: string, sample: () => Statevector, exactSecondMoment: number) {
  const counts = new Array<number>(B).fill(0);
  let m1 = 0;
  let m2 = 0;
  for (let s = 0; s < N; s++) {
    const F = sample().fidelity(sample());
    counts[Math.min(B - 1, Math.floor(F * B))]++;
    m1 += F / N;
    m2 += (F * F) / N;
  }
  let kl = 0;
  for (let j = 0; j < B; j++) {
    const p = counts[j] / N;
    if (p > 0) kl += p * Math.log(p / haarBin(j));
  }
  console.log(label + ': E[F] = ' + m1.toFixed(4) + '  E[F^2] = ' + m2.toFixed(4) + ' (exact ' + exactSecondMoment.toFixed(4) + ')  KL = ' + kl.toFixed(4));
}

for (const L of [1, 2, 3]) {
  const angles = () => Array.from({ length: 2 * L }, () => 2 * Math.PI * rng());
  expressibility('HEA n = 1, L = ' + L, () => hea1(angles()), 1 / 3 + 16 ** -L / 6);
}
expressibility('Haar (Gaussian)', haar1, 1 / 3);
console.log('finite-sample floor (B - 1)/(2N) = ' + ((B - 1) / (2 * N)).toFixed(6));

Expected output

HEA n = 1, L = 1: E[F] = 0.5014  E[F^2] = 0.3451 (exact 0.3438)  KL = 0.0103
HEA n = 1, L = 2: E[F] = 0.4987  E[F^2] = 0.3323 (exact 0.3340)  KL = 0.0003
HEA n = 1, L = 3: E[F] = 0.4981  E[F^2] = 0.3308 (exact 0.3334)  KL = 0.0004
Haar (Gaussian): E[F] = 0.5040  E[F^2] = 0.3381 (exact 0.3333)  KL = 0.0005
finite-sample floor (B - 1)/(2N) = 0.000475

The output agrees with the theory. For L=1L=1 the sample second moment is close to 1132\frac{11}{32} and the KL divergence is far above the floor: one layer covers the Bloch sphere, but with the wrong (arcsine) weights. Already for L=2L=2 the KL estimate is of the order of the sampling floor, and the Haar sampler itself lands at the floor. Expressibility is a property of the circuit and of the parameter distribution: drawing θ0,0\theta_{0,0} with density 12sin⁡ϑ\frac12\sin\vartheta on [0,π][0,\pi] instead of uniformly would make the one-layer one-qubit HEA exactly Haar.

Quick check

For two qubits (d=4d=4) and B=20B=20 bins, what is the Haar probability Q0=Pr⁡(F<0.05)Q_0=\Pr(F<0.05) of the first bin?

tolerance ±0.00001

Expressibility versus trainability, and hardware connectivity#

High expressibility sounds desirable: the target state should be reachable. But a circuit whose outputs are spread like Haar-random states also has a cost function that is almost the same everywhere. The following theorem makes this precise for the simplest cost, the expectation value of ZZ on one qubit.

Theorem

Cost values of Haar-random states concentrate exponentially

Let ∣ψ⟩\lvert\psi\rangle be Haar-random on nn qubits, d=2nd=2^n, and C(ψ)=⟨ψ∣Z0∣ψ⟩C(\psi)=\langle\psi\vert Z_0\vert\psi\rangle with Z0=Z⊗I⊗(n−1)Z_0=Z\otimes I^{\otimes(n-1)}. Then E[C]=0,Var[C]=1d+1=12n+1,Pr⁡(∣C∣≥ε)≤1ε2(2n+1).\mathbb E[C]=0,\qquad\mathrm{Var}[C]=\frac1{d+1}=\frac1{2^n+1},\qquad\Pr\big(\lvert C\rvert\ge\varepsilon\big)\le\frac1{\varepsilon^2(2^n+1)} .

Proof (7 steps)
  1. Step 1:
    C=∑k=0d−1sk pk,pk=∣ψk∣2,sk={+1,bit 0 of k is 0,−1,bit 0 of k is 1,∑ksk=0C=\sum_{k=0}^{d-1}s_k\,p_k,\qquad p_k=\lvert\psi_k\rvert^2,\qquad s_k=\begin{cases}+1,&\text{bit }0\text{ of }k\text{ is }0,\\-1,&\text{bit }0\text{ of }k\text{ is }1,\end{cases}\qquad\sum_ks_k=0
    Z0Z_0 is diagonal in the computational basis with eigenvalue (−1)k0(-1)^{k_0}; exactly half of the dd basis states have k0=0k_0=0.
  2. Step 2:
    pk=dEkE1+⋯+Ed∼Beta(1,d−1),E[pk]=1d,E[pk2]=2d(d+1)p_k\overset{d}{=}\frac{E_k}{E_1+\cdots+E_d}\sim\mathrm{Beta}(1,d-1),\qquad\mathbb E[p_k]=\frac1d,\qquad\mathbb E[p_k^2]=\frac2{d(d+1)}
    The Gaussian representation, applied to each coordinate (permuting coordinates is a unitary, so all pkp_k have the same law); moments from the Haar fidelity theorem.
  3. Step 3:
    ∑lpl=1 ⇒ E[pk2]+(d−1) E[pkpl]=E[pk∑lpl]=1d ⇒ E[pkpl]=1d(d+1) (k≠l)\sum_lp_l=1\ \Rightarrow\ \mathbb E[p_k^2]+(d-1)\,\mathbb E[p_kp_l]=\mathbb E\Big[p_k\sum_lp_l\Big]=\frac1d\ \Rightarrow\ \mathbb E[p_kp_l]=\frac1{d(d+1)}\ (k\ne l)
    All pairs k≠lk\ne l have the same joint law by permutation invariance; solve 1d−2d(d+1)=d−1d(d+1)\frac1d-\frac2{d(d+1)}=\frac{d-1}{d(d+1)} for the common value.
  4. Step 4:
    E[C]=∑ksk E[pk]=1d∑ksk=0\mathbb E[C]=\sum_ks_k\,\mathbb E[p_k]=\frac1d\sum_ks_k=0
    Linearity of expectation.
  5. Step 5:
    E[C2]=∑ksk2 E[pk2]+∑k≠lsksl E[pkpl]=d⋅2d(d+1)+(∑ksk)2−∑ksk2d(d+1)=2d+1−1d+1\mathbb E[C^2]=\sum_ks_k^2\,\mathbb E[p_k^2]+\sum_{k\ne l}s_ks_l\,\mathbb E[p_kp_l]=d\cdot\frac2{d(d+1)}+\frac{\big(\sum_ks_k\big)^2-\sum_ks_k^2}{d(d+1)}=\frac2{d+1}-\frac1{d+1}
    Expand the square; sk2=1s_k^2=1 and ∑k≠lsksl=(∑ksk)2−∑ksk2=0−d\sum_{k\ne l}s_ks_l=(\sum_ks_k)^2-\sum_ks_k^2=0-d.
  6. Step 6:
    Var[C]=E[C2]−E[C]2=1d+1=12n+1\mathrm{Var}[C]=\mathbb E[C^2]-\mathbb E[C]^2=\frac1{d+1}=\frac1{2^n+1}
    Combine the last two lines.
  7. Step 7:
    Pr⁡(∣C∣≥ε)≤Var[C]ε2=1ε2(2n+1),n=20, ε=0.01:  ≤104220+1≈9.54×10−3\Pr\big(\lvert C\rvert\ge\varepsilon\big)\le\frac{\mathrm{Var}[C]}{\varepsilon^2}=\frac1{\varepsilon^2(2^n+1)},\qquad n=20,\ \varepsilon=0.01:\ \ \le\frac{10^4}{2^{20}+1}\approx9.54\times10^{-3}
    Chebyshev's inequality. With 2020 qubits, fewer than one Haar-random state in a hundred has ∣⟨Z0⟩∣≥0.01\lvert\langle Z_0\rangle\rvert\ge0.01: over almost the entire state space the cost is flat at 00.

∎

From flat costs to vanishing gradients. If the states produced by an ansatz are distributed like Haar-random states, the theorem says that the cost takes essentially the same value almost everywhere in parameter space, and resolving its differences needs a number of measurement shots of order 1/Var[C]≈2n1/\mathrm{Var}[C]\approx2^n. McClean et al. (2018) showed that the same exponential suppression holds for the partial derivatives ∂θjC\partial_{\theta_j}C whenever the circuit parts before or after θj\theta_j are random enough (unitary 2-designs): the barren plateau phenomenon, derived in full in the lesson on barren plateaus. Holmes et al. (2022) turned this into a quantitative trade-off: the more expressive an ansatz, measured exactly by how close its distribution is to Haar, the smaller the variance of its gradients. A deep HEA with random initialisation is the textbook example of a circuit that is expressive and therefore hard to train. Shallow HEAs with local cost functions avoid the plateau (Cerezo et al. 2021): at depth O(log⁡n)O(\log n) the variance of the gradient of a cost built from few-qubit observables decays at most polynomially in nn. Initialising at θ=0\theta=0, where the circuit outputs ∣0n⟩\lvert0^n\rangle, is another way of starting far away from the Haar-like region.

Proposition

Rotations about X and Y from Z rotations and fixed gates

For every angle α\alpha: RX(α)=H RZ(α) HR_X(\alpha)=H\,R_Z(\alpha)\,H and RY(α)=S RX(α) S†=S H RZ(α) H S†R_Y(\alpha)=S\,R_X(\alpha)\,S^\dagger=S\,H\,R_Z(\alpha)\,H\,S^\dagger. Hence an HEA layer can be compiled into the parameter-free gates HH, SS, S†S^\dagger, CNOT and the rotations RZ(α)R_Z(\alpha), which superconducting processors implement as "virtual" gates by shifting the phase of later control pulses, at zero duration (McKay et al. 2017).

Proof (4 steps)
  1. Step 1:
    V e−iαP/2 V†=∑m≥0(−iα/2)mm!(VPV†)m=e−iα VPV†/2V\,e^{-i\alpha P/2}\,V^\dagger=\sum_{m\ge0}\frac{(-i\alpha/2)^m}{m!}\big(VPV^\dagger\big)^m=e^{-i\alpha\,VPV^\dagger/2}
    Conjugation by a unitary VV is multiplicative, VPmV†=(VPV†)mVP^mV^\dagger=(VPV^\dagger)^m, and continuous, so it passes through the exponential series.
  2. Step 2:
    HZH=12(111−1)(100−1)(111−1)=(0110)=X ⇒ HRZ(α)H=RX(α)HZH=\tfrac12\begin{pmatrix}1&1\\1&-1\end{pmatrix}\begin{pmatrix}1&0\\0&-1\end{pmatrix}\begin{pmatrix}1&1\\1&-1\end{pmatrix}=\begin{pmatrix}0&1\\1&0\end{pmatrix}=X\ \Rightarrow\ HR_Z(\alpha)H=R_X(\alpha)
    Matrix multiplication with H=H†H=H^\dagger.
  3. Step 3:
    SXS†=(100i)(0110)(100−i)=(0−ii0)=Y ⇒ SRX(α)S†=RY(α)SXS^\dagger=\begin{pmatrix}1&0\\0&i\end{pmatrix}\begin{pmatrix}0&1\\1&0\end{pmatrix}\begin{pmatrix}1&0\\0&-i\end{pmatrix}=\begin{pmatrix}0&-i\\i&0\end{pmatrix}=Y\ \Rightarrow\ SR_X(\alpha)S^\dagger=R_Y(\alpha)
    Matrix multiplication with S†=diag(1,−i)S^\dagger=\mathrm{diag}(1,-i).
  4. Step 4:
    RY(α)=S RX(α) S†=S H RZ(α) H S†R_Y(\alpha)=S\,R_X(\alpha)\,S^\dagger=S\,H\,R_Z(\alpha)\,H\,S^\dagger
    Substitute the first identity into the second. A layer RZ(β)RY(α)R_Z(\beta)R_Y(\alpha) thus costs two parameterised RZR_Z gates and four fixed gates per qubit.

∎

Proposition

SWAP overhead of non-adjacent CNOTs on a line

(a) SWAP=CNOT(a,b) CNOT(b,a) CNOT(a,b)\mathrm{SWAP}=\mathrm{CNOT}(a,b)\,\mathrm{CNOT}(b,a)\,\mathrm{CNOT}(a,b). (b) On a line of qubits coupled only to their neighbours, a CNOT between two qubits at distance k≥1k\ge1 implemented with nearest-neighbour SWAPs needs at least k−1k-1 SWAPs, and moving one qubit next to the other and back costs 6(k−1)+1=6k−56(k-1)+1=6k-5 CNOTs. (c) The ladder of the HEA consists of distance-11 CNOTs only, so its routing overhead is zero, whereas an all-to-all entangler on 44 qubits routed this way costs 3030 CNOTs instead of 66.

Proof (5 steps)
  1. Step 1:
    (xa,xb)→CNOT(a,b)(xa, xa⊕xb)→CNOT(b,a)(xb, xa⊕xb)→CNOT(a,b)(xb, xa)(x_a,x_b)\xrightarrow{\mathrm{CNOT}(a,b)}(x_a,\,x_a\oplus x_b)\xrightarrow{\mathrm{CNOT}(b,a)}(x_b,\,x_a\oplus x_b)\xrightarrow{\mathrm{CNOT}(a,b)}(x_b,\,x_a)
    Track the basis labels: xa⊕(xa⊕xb)=xbx_a\oplus(x_a\oplus x_b)=x_b and xb⊕(xa⊕xb)=xax_b\oplus(x_a\oplus x_b)=x_a.
  2. Step 2:
    both sides linear and equal on the basis ⇒ CNOT(a,b) CNOT(b,a) CNOT(a,b)=SWAP\text{both sides linear and equal on the basis}\ \Rightarrow\ \mathrm{CNOT}(a,b)\,\mathrm{CNOT}(b,a)\,\mathrm{CNOT}(a,b)=\mathrm{SWAP}
    Two linear maps that agree on a basis are equal (the lab of the SWAP lesson in Stage 2 builds exactly this circuit).
  3. Step 3:
    one neighbour SWAP changes ∣pos(a)−pos(b)∣ by at most 1 ⇒ #SWAP≥k−1\text{one neighbour SWAP changes }\lvert\mathrm{pos}(a)-\mathrm{pos}(b)\rvert\text{ by at most }1\ \Rightarrow\ \#\mathrm{SWAP}\ge k-1
    A SWAP of two neighbouring sites moves at most the two logical qubits on those sites, each by one position. If it exchanges aa and bb themselves their distance stays 11; otherwise at most one of them moves, by one. Reducing the distance from kk to 11 needs at least k−1k-1 SWAPs.
  4. Step 4:
    (k−1) SWAPs+CNOT+(k−1) SWAPs ⇒ 3(k−1)+1+3(k−1)=6k−5 CNOTs(k-1)\ \text{SWAPs}+\mathrm{CNOT}+(k-1)\ \text{SWAPs}\ \Rightarrow\ 3(k-1)+1+3(k-1)=6k-5\ \text{CNOTs}
    Move aa step by step towards bb, apply the now-local CNOT, and move aa back so that the qubit layout is restored; each SWAP costs three CNOTs by (a).
  5. Step 5:
    n=4:  3⋅(6⋅1−5)+2⋅(6⋅2−5)+1⋅(6⋅3−5)=3+14+13=30n=4:\ \ 3\cdot(6\cdot1-5)+2\cdot(6\cdot2-5)+1\cdot(6\cdot3-5)=3+14+13=30
    An all-to-all entangler has (42)=6\binom42=6 CNOTs: three pairs at distance 11, two at distance 22 and one at distance 33. The HEA ladder with its n−1=3n-1=3 CNOTs needs no routing at all.

∎

Quick check

On a line of qubits, how many CNOTs does a single CNOT between qubits at distance k=3k=3 cost when it is routed by nearest-neighbour SWAPs there and back?

tolerance ±0.5

The NumPy cell rebuilds the HEA from explicit Kronecker products: it reproduces the four lab fidelities (12\tfrac12 at zero, 11 and 00 for θ1,0,0=±π2\theta_{1,0,0}=\pm\tfrac\pi2, 14\tfrac14 for the rotation in the wrong layer), prints the prefix-parity permutation of the ladder and its order for n=2,…,6n=2,\dots,6, checks the normalisation, the mean, the second moment and the bin probabilities of the Haar density by Gauss–Legendre quadrature (exact for polynomials of degree below 4040), and iterates the averaged Bloch second-moment matrix of the one-qubit HEA layer by layer against 13+16−L6\frac13+\frac{16^{-L}}6.

Pythonrunnable
import numpy as np

I2 = np.eye(2)

def ry(t):
    c, s = np.cos(t / 2), np.sin(t / 2)
    return np.array([[c, -s], [s, c]], dtype=complex)

def rz(t):
    return np.diag([np.exp(-0.5j * t), np.exp(0.5j * t)])

def on(q, g, n):
    out = np.array([[1.0 + 0j]])
    for k in range(n):
        out = np.kron(out, g if k == q else I2)
    return out

def cnot(c, t, n):
    U = np.zeros((2 ** n, 2 ** n))
    for x in range(2 ** n):
        bits = [(x >> (n - 1 - k)) & 1 for k in range(n)]
        bits[t] ^= bits[c]
        U[int("".join(map(str, bits)), 2), x] = 1
    return U

def ladder(n):
    E = np.eye(2 ** n)
    for q in range(n - 1):
        E = cnot(q, q + 1, n) @ E
    return E

def hea_state(theta):  # theta[l, q, k]: layer l, qubit q, k = 0 (R_Y) or 1 (R_Z)
    L, n, _ = theta.shape
    psi = np.zeros(2 ** n, dtype=complex)
    psi[0] = 1
    for l in range(L):
        for q in range(n):
            psi = on(q, rz(theta[l, q, 1]) @ ry(theta[l, q, 0]), n) @ psi
        psi = ladder(n) @ psi
    return psi

# (1) The lab: fidelity with GHZ = (|000> + |111>)/sqrt(2) for the 3-qubit, 2-layer template.
ghz = np.zeros(8)
ghz[0] = ghz[7] = 1 / np.sqrt(2)
def lab(**angles):
    theta = np.zeros((2, 3, 2))
    for name, value in angles.items():
        _, l, q, k = name.split("_")
        theta[int(l), int(q), int(k)] = value
    return abs(np.vdot(ghz, hea_state(theta))) ** 2
for label, kw in [("all zero", {}), ("t_1_0_0 = +pi/2", {"t_1_0_0": np.pi / 2}),
                  ("t_1_0_0 = -pi/2", {"t_1_0_0": -np.pi / 2}), ("t_0_0_0 = +pi/2", {"t_0_0_0": np.pi / 2})]:
    print(f"{label:16s} F = {lab(**kw):.6f}")

# (2) The ladder is the prefix-parity permutation, and its order is the smallest power of two >= n.
E3 = ladder(3)
print("E:", " ".join(f"{x:03b}->{int(np.argmax(E3[:, x])):03b}" for x in range(8)))
orders = []
for n in range(2, 7):
    E, P, r = ladder(n), ladder(n), 1
    while not np.allclose(P, np.eye(2 ** n)):
        P, r = E @ P, r + 1
    orders.append(f"n={n}: {r}")
print("order of E:", ", ".join(orders))

# (3) Haar fidelity density (d - 1)(1 - F)^(d - 2): Gauss-Legendre quadrature is exact for these polynomials.
nodes, weights = np.polynomial.legendre.leggauss(20)
F, w = (nodes + 1) / 2, weights / 2
for d in (2, 4, 8):
    dens = (d - 1) * (1 - F) ** (d - 2)
    j = np.arange(20)
    Q = (1 - j / 20) ** (d - 1) - (1 - (j + 1) / 20) ** (d - 1)
    print(f"d = {d}: mass = {np.sum(w * dens):.10f}  E[F] = {np.sum(w * F * dens):.10f} (1/d = {1 / d:.10f})"
          f"  E[F^2] = {np.sum(w * F ** 2 * dens):.10f} (2/(d(d+1)) = {2 / (d * (d + 1)):.10f})  sum Q_j = {Q.sum():.10f}  Q_0 = {Q[0]:.6f}")

# (4) One-qubit HEA: average the Bloch second-moment matrix over uniform R_Y and R_Z angles, layer by layer.
def rot_avg(M, axis):  # 8 equispaced angles average cos^2, sin^2, sin cos exactly
    out = np.zeros((3, 3))
    for a in 2 * np.pi * np.arange(8) / 8:
        c, s = np.cos(a), np.sin(a)
        R = np.array([[c, 0, s], [0, 1, 0], [-s, 0, c]]) if axis == "y" else np.array([[c, -s, 0], [s, c, 0], [0, 0, 1]])
        out += R @ M @ R.T / 8
    return out
M = np.diag([0.0, 0.0, 1.0])
for L in range(1, 5):
    M = rot_avg(rot_avg(M, "y"), "z")
    p, q = M[0, 0], M[2, 2]
    print(f"L = {L}: M = diag({p:.6f}, {M[1, 1]:.6f}, {q:.6f})  E[F^2] = {(1 + np.sum(M ** 2)) / 4:.10f}"
          f"  1/3 + 16^-L/6 = {1 / 3 + 16.0 ** -L / 6:.10f}")

Expected output

all zero         F = 0.500000
t_1_0_0 = +pi/2  F = 1.000000
t_1_0_0 = -pi/2  F = 0.000000
t_0_0_0 = +pi/2  F = 0.250000
E: 000->000 001->001 010->011 011->010 100->111 101->110 110->100 111->101
order of E: n=2: 2, n=3: 4, n=4: 4, n=5: 8, n=6: 8
d = 2: mass = 1.0000000000  E[F] = 0.5000000000 (1/d = 0.5000000000)  E[F^2] = 0.3333333333 (2/(d(d+1)) = 0.3333333333)  sum Q_j = 1.0000000000  Q_0 = 0.050000
d = 4: mass = 1.0000000000  E[F] = 0.2500000000 (1/d = 0.2500000000)  E[F^2] = 0.1000000000 (2/(d(d+1)) = 0.1000000000)  sum Q_j = 1.0000000000  Q_0 = 0.142625
d = 8: mass = 1.0000000000  E[F] = 0.1250000000 (1/d = 0.1250000000)  E[F^2] = 0.0277777778 (2/(d(d+1)) = 0.0277777778)  sum Q_j = 1.0000000000  Q_0 = 0.301663
L = 1: M = diag(0.250000, 0.250000, 0.500000)  E[F^2] = 0.3437500000  1/3 + 16^-L/6 = 0.3437500000
L = 2: M = diag(0.312500, 0.312500, 0.375000)  E[F^2] = 0.3339843750  1/3 + 16^-L/6 = 0.3339843750
L = 3: M = diag(0.328125, 0.328125, 0.343750)  E[F^2] = 0.3333740234  1/3 + 16^-L/6 = 0.3333740234
L = 4: M = diag(0.332031, 0.332031, 0.335938)  E[F^2] = 0.3333358765  1/3 + 16^-L/6 = 0.3333358765

References

  1. A. Kandala, A. Mezzacapo, K. Temme, M. Takita, M. Brink, J. M. Chow and J. M. Gambetta, Hardware-efficient variational quantum eigensolver for small molecules and quantum magnets, Nature 549, 242–246 (2017).doi:10.1038/nature23879
  2. S. Sim, P. D. Johnson and A. Aspuru-Guzik, Expressibility and entangling capability of parameterized quantum circuits for hybrid quantum-classical algorithms, Advanced Quantum Technologies 2, 1900070 (2019).doi:10.1002/qute.201900070
  3. J. R. McClean, S. Boixo, V. N. Smelyanskiy, R. Babbush and H. Neven, Barren plateaus in quantum neural network training landscapes, Nature Communications 9, 4812 (2018).doi:10.1038/s41467-018-07090-4
  4. Z. Holmes, K. Sharma, M. Cerezo and P. J. Coles, Connecting ansatz expressibility to gradient magnitudes and barren plateaus, PRX Quantum 3, 010313 (2022).doi:10.1103/PRXQuantum.3.010313
  5. M. Cerezo, A. Sone, T. Volkoff, L. Cincio and P. J. Coles, Cost function dependent barren plateaus in shallow parametrized quantum circuits, Nature Communications 12, 1791 (2021).doi:10.1038/s41467-021-21728-w
  6. E. Grant, L. Wossnig, M. Ostaszewski and M. Benedetti, An initialization strategy for addressing barren plateaus in parametrized quantum circuits, Quantum 3, 214 (2019).doi:10.22331/q-2019-12-09-214
  7. D. C. McKay, C. J. Wood, S. Sheldon, J. M. Chow and J. M. Gambetta, Efficient Z gates for quantum computing, Physical Review A 96, 022330 (2017).doi:10.1103/PhysRevA.96.022330

Hands-on labs

  • Train a hardware-efficient ansatz to GHZ

    1 step · auto-graded

    The circuit is a fixed hardware-efficient ansatz on 33 qubits with 22 layers: in each layer every qubit gets RYR_Y then RZR_Z, followed by the CNOT ladder CNOT(0,1)\mathrm{CNOT}(0,1), CNOT(1,2)\mathrm{CNOT}(1,2). All 1212 angles start at 00, where the circuit outputs ∣000⟩\lvert000\rangle. Choose angle values so that the output is the GHZ state 12(∣000⟩+∣111⟩)\tfrac1{\sqrt2}(\lvert000\rangle+\lvert111\rangle), using the identity-at-zero and basis-state tracking arguments of the lesson.

    Open lab