Tensor-Network Methods for Classical Circuit Simulation

A pure state of nn qubits, ∣ψ⟩=∑s1…snψs1…sn∣s1…sn⟩|\psi\rangle=\sum_{s_1\dots s_n}\psi_{s_1\dots s_n}|s_1\dots s_n\rangle, is literally a tensor with nn indices, each of dimension 22. Written out in full it takes 2n2^n complex amplitudes — the wall that every classical simulation of quantum computing slams into. Tensor-network methods are the discovery that this wall is an illusion for the states that actually matter. Ground states of local gapped Hamiltonians, the outputs of shallow circuits, weakly entangled dynamics: all of them occupy a vanishingly thin, highly structured sliver of Hilbert space, and inside that sliver the giant tensor factorizes into a network of small tensors joined by contracted bond indices. The cost of storing and manipulating the state is then set not by 2n2^n but by how much entanglement has to flow across those bonds. This deep-tier synthesis takes the parameterized-circuit and expressibility machinery of aq-qci-parameterized-quantum-circuits-ansatze and aq-qci-expressibility-trainability-tradeoffs and turns it into a concrete question of classical resources. Three ideas organize everything. First, the workhorse ansatz — the Matrix Product State (MPS) — represents a one-dimensional state as a product of matrices of size at most χ×χ\chi\times\chi, and a single theorem fixes χ\chi: the bond dimension across any cut equals the Schmidt rank there, so entanglement is the exact currency of classical cost. Second, to simulate a circuit you either evolve an MPS gate by gate, re-compressing by SVD (the TEBD algorithm), or treat the whole circuit as one network and contract it in a near-optimal order whose cost is exponential in the graph's treewidth. Third, truncation makes the approximation certifiable: the squared weight of the discarded Schmidt values is exactly the simulation error. The whole subject reduces "is this circuit classically hard?" to "how fast does the entanglement grow?"

Key takeaways
  • A state on nn qubits is an order-nn tensor with 2n2^n amplitudes; a tensor network factorizes it into a contracted graph of small tensors whose cost scales with the network's connectivity, not with 2n2^n.
  • Matrix Product States represent 1D states with bond dimension χ\chi; they exactly capture any state whose bipartite entanglement entropy is at most log⁡2χ\log_2\chi across every cut, and evolve one gate in time O(χ3)O(\chi^3).
  • Circuit simulation contracts the network built from the input ∣0⟩⊗n|0\rangle^{\otimes n}, the gate tensors, and the final observable; the dominant cost is governed by the contraction-order-dependent treewidth of the circuit graph.
  • Truncation is the decisive lever: after each gate you keep only the largest χ\chi singular values, and the discarded weight ϵ=∑k>χλk2\epsilon=\sum_{k>\chi}\lambda_k^2 is exactly the squared simulation error — a tunable, certified approximation.
  • Area-law 1D states, shallow and low-degree circuits, and Clifford circuits are efficiently simulable; generic deep or 2D circuits drive χ∼2n/2\chi\sim 2^{n/2}, which is precisely the hardness foundation of random-circuit-sampling supremacy claims.
  • Key idea: entanglement, measured by the Schmidt rank across each cut, is the true cost currency — a tensor network simulates a state in time polynomial in nn and in the bond dimension χ\chi, and χ\chi need only equal the Schmidt rank, so low-entanglement states are classically cheap and circuit simulation reduces to bounding entanglement growth.
The core insight

Entanglement, measured by the Schmidt rank across each cut, is the true cost currency: a tensor network simulates a quantum state in time polynomial in nn and polynomial in the bond dimension χ\chi, and χ\chi need only equal the Schmidt rank. Low-entanglement (area-law) states are therefore classically efficient, and circuit simulation reduces to bounding the growth of entanglement.

Where this lives today

Tensor networks are simultaneously the leading classical competitor to quantum hardware and the sharpest theoretical instrument for judging quantum advantage. On the competitor side, DMRG remains the gold-standard baseline for one- and quasi-one-dimensional quantum chemistry and materials, routinely reaching machine-precision ground-state energies with bond dimensions of a few thousand — which is exactly why any Variational Quantum Eigensolver claim must be benchmarked against it, and why a genuine quantum advantage requires problems where DMRG's χ\chi provably explodes, typically strongly correlated 2D systems. The most public collision came with random-circuit sampling: Google's 2019 Sycamore experiment on 53 qubits at depth 20 estimated 10,000 years for the best classical supercomputer, and within a few years tensor-network contraction with optimized orderings and index slicing (Pan and Zhang; Gray and Kourtis; Huang and collaborators) had cut that estimate to a matter of days on a GPU cluster, by attacking the very treewidth of equation (7). The same tools now serve as a pre-screen for variational ansätze — before committing a Parameterized Quantum Circuit to hardware, one simulates it with a bounded-χ\chi MPS to see whether its entanglement stays low (classically tractable, hence advantage-free) or grows fast enough to signal both quantum hardness and the onset of barren plateaus, tying this article directly to the expressibility–trainability trade-off. And the networks are models in their own right: MPS and Tree Tensor Networks are trained by sweeping optimization as classical machine-learning models, with the bond dimension playing the role of model capacity. The open frontier is the boundary itself — every improvement in contraction heuristics narrows the region where quantum devices are provably hard to simulate, and locating that shrinking boundary precisely is one of the live questions of the field.

The physical picture

Picture the 2n2^n-amplitude tensor as a brick wall, and imagine slicing it with a knife at any point along its length. The number of independent "channels" — Schmidt coefficients — that connect the two halves is the only thing that measures how entangled the state is across that cut. A tensor network is nothing but the wiring diagram that carries exactly those channels and no more: the bond dimension χ\chi is the number of wires crossing the cut. If a state needs only a few wires anywhere along the wall, it is cheap to store; if every cut is saturated with 2n/22^{n/2} wires, you are back to the full 2n2^n. The first figure makes this trade concrete with real numbers: an MPS chain replaces one monstrous tensor by nn tiny ones, and the storage collapses from petabytes to megabytes.

Figure
Factorizing the order-n amplitude tensor into a Matrix Product State chain. The single tensor on the left carries 2ⁿ complex amplitudes, one dangling physical leg per qubit; the chain on the right stores n small tensors joined by bond indices of size χ, costing n·d·χ² numbers. For n = 50 and χ = 64 that is 1.13×10¹⁵ amplitudes, about 18 petabytes, versus 4.10×10⁵ numbers, about 6.6 megabytes — a compression of 2.7×10⁹
One giant tensor (2ⁿ numbers) factorizes into a chain of small tensors (n·d·χ²) ψ s₁ s₂ … sₙ (n legs) 2ⁿ amplitudes n=50: 1.13×10¹⁵ ≈ 18 PB SVD χ χ s₁ sₙ n·d·χ² numbers n=50, χ=64: 4.10×10⁵ ≈ 6.6 MB compression 2ⁿ / (n·d·χ²) ≈ 2.7×10⁹ when the entanglement is bounded

The reason this works — and the reason it must sometimes fail — is the Schmidt decomposition. At every bipartition A∣BA|B the state can be written ∣ψ⟩=∑kλk∣αk⟩A∣βk⟩B|\psi\rangle=\sum_k\lambda_k|\alpha_k\rangle_A|\beta_k\rangle_B with at most χ=min⁡(2∣A∣,2∣B∣)\chi=\min(2^{|A|},2^{|B|}) nonzero terms, and an MPS is simply these local Schmidt decompositions glued together site by site, the bond index running over the Schmidt vectors. A two-qubit gate is a small order-4 tensor dropped onto two legs; applying it can at most multiply the Schmidt ranks on either side of the cut, so entanglement grows multiplicatively with circuit depth in the worst case. That multiplicative growth is the whole story of when simulation stays cheap and when it explodes. The second figure is the quantity a practitioner actually watches: the bond dimension a 1D brickwork circuit demands as its depth increases, and the depth at which it overruns any classical memory budget.

Figure
Bond dimension needed to represent a generic 1D brickwork circuit versus its depth, on a logarithmic axis. Each layer can roughly square the local Schmidt rank, so the light-cone bound is χ = 2^(d/2): χ = 32 at depth 10, 1024 at depth 20, and about 10⁹ at depth 60. A large-workstation MPS ceiling of χ = 2¹⁶ ≈ 6.5×10⁴ is overrun near depth 32 — the tractability cliff beyond which exact MPS simulation fails
The punchline: χ = 2^(d/2) overruns a χ = 2¹⁶ budget near depth 32 — the tractability cliff 10⁰ 10² 10⁴ 10⁶ 10⁸ 0 10 20 30 40 50 60 circuit depth d bond dimension χ workstation ceiling χ = 2¹⁶ ≈ 6.5×10⁴ d=10: χ=32 d=20: χ=1024 d=32: χ=65536 (cliff) d=60: χ≈1.1×10⁹

The core equations

The starting point is the definition that makes the whole subject possible: a quantum state is a tensor. Reading the amplitudes ψs1…sn\psi_{s_1\dots s_n} as the entries of an nn-index array is not an analogy but an identity, and every tensor-network method is a statement about how that array factorizes.

∣ψ⟩=∑s1,…,sn=01ψs1s2⋯sn ∣s1⟩⊗⋯⊗∣sn⟩(1)|\psi\rangle=\sum_{s_1,\dots,s_n=0}^{1}\psi_{s_1 s_2\cdots s_n}\,|s_1\rangle\otimes\cdots\otimes|s_n\rangle \qquad (1)

where: ψs1⋯sn\psi_{s_1\cdots s_n} — order-nn amplitude tensor, 2n2^n complex entries; si∈{0,1}s_i\in\{0,1\} — physical index of qubit ii; nn — number of qubits. Read it as: the exponential cost of quantum states lives entirely in the number of indices, not in any single number being large. The one-dimensional factorization replaces that single array with a chain of matrices. Each site carries a small tensor A[i]siA^{[i]s_i} that is a matrix in the bond indices for every value of the physical index, and the whole amplitude is recovered by matrix multiplication along the chain.

ψs1⋯sn=∑a1,…,an−1Aa1[1]s1Aa1a2[2]s2⋯Aan−1[n]sn(2)\psi_{s_1\cdots s_n}=\sum_{a_1,\dots,a_{n-1}} A^{[1]s_1}_{a_1}A^{[2]s_2}_{a_1 a_2}\cdots A^{[n]s_n}_{a_{n-1}} \qquad (2)

where: A[i]siA^{[i]s_i} — site tensor at ii; aia_i — bond (virtual) index of dimension ≤χ\le\chi; χ\chi — bond dimension, the maximal Schmidt rank over any cut. Read it as: storage falls from 2n2^n to n d χ2n\,d\,\chi^2 with d=2d=2, polynomial in nn whenever χ\chi stays bounded. Everything about cost traces back to the Schmidt decomposition at a chosen bipartition. It is the microscope that counts, cut by cut, how many terms are genuinely needed — and that count is precisely the bond dimension the MPS must carry across that cut.

∣ψ⟩=∑k=1χAλk ∣αk⟩A⊗∣βk⟩B,∑kλk2=1(3)|\psi\rangle=\sum_{k=1}^{\chi_A}\lambda_k\,|\alpha_k\rangle_A\otimes|\beta_k\rangle_B,\qquad \sum_k\lambda_k^2=1 \qquad (3)

where: λk≥0\lambda_k\ge 0 — Schmidt coefficients (singular values), ordered λ1≥λ2≥⋯\lambda_1\ge\lambda_2\ge\cdots; ∣αk⟩A,∣βk⟩B|\alpha_k\rangle_A,|\beta_k\rangle_B — orthonormal Schmidt vectors on each side; χA\chi_A — Schmidt rank, the number of nonzero λk\lambda_k and the required bond dimension. Read it as: the required χ\chi is a property of the state at each cut, not a free parameter of the method. The bridge from Schmidt rank to physics is the entanglement entropy. The eigenvalues of the reduced density operator are the squared Schmidt coefficients, and an entropy over χA\chi_A outcomes can never exceed log⁡2χA\log_2\chi_A — so bounded entanglement directly caps the bond dimension.

S(ρA)=−∑kλk2log⁡2λk2 ≤ log⁡2χA(4)S(\rho_A)=-\sum_k \lambda_k^2\log_2\lambda_k^2\ \le\ \log_2\chi_A \qquad (4)

where: ρA=TrB∣ψ⟩⟨ψ∣\rho_A=\mathrm{Tr}_B|\psi\rangle\langle\psi| — reduced density operator on AA; S(ρA)S(\rho_A) — von Neumann bipartite entanglement entropy; χA\chi_A — Schmidt rank across the cut. Read it as: an area-law state with S=O(1)S=O(1) needs only χ=O(1)\chi=O(1), while a volume-law state with S∼n/2S\sim n/2 forces χ∼2n/2\chi\sim 2^{n/2}. Approximation enters through a single, honest knob. Keeping only the largest χ′\chi' Schmidt values produces the best low-rank approximation, and the squared weight of everything you threw away is the squared error of the approximation — a number you can compute and report.

∣ψ~⟩=∑k=1χ′λk∣αk⟩∣βk⟩,ϵ=∑k>χ′λk2=∥∣ψ⟩−∣ψ~⟩∥2(5)|\tilde\psi\rangle=\sum_{k=1}^{\chi'}\lambda_k|\alpha_k\rangle|\beta_k\rangle,\qquad \epsilon=\sum_{k>\chi'}\lambda_k^2=\big\||\psi\rangle-|\tilde\psi\rangle\big\|^2 \qquad (5)

where: χ′<χ\chi'<\chi — truncated bond dimension actually retained; ϵ\epsilon — discarded weight, the certified squared approximation error; λk\lambda_k — singular values, of which the smallest are dropped. Read it as: no simulation result is meaningful without quoting ϵ\epsilon (or the fidelity 1−ϵ1-\epsilon) alongside it. To run a circuit as an MPS you apply one gate at a time. A two-qubit gate merges two neighbouring site tensors into a single object Θ\Theta, which is then re-split by SVD and truncated back to bond dimension χ\chi — the elementary step of the TEBD algorithm.

Θai−1ai+1sisi+1=∑ai,ti,ti+1Utiti+1sisi+1 Aai−1ai[i]ti Aaiai+1[i+1]ti+1(6)\Theta^{s_i s_{i+1}}_{a_{i-1}a_{i+1}}=\sum_{a_i,t_i,t_{i+1}}U^{s_i s_{i+1}}_{t_i t_{i+1}}\,A^{[i]t_i}_{a_{i-1}a_i}\,A^{[i+1]t_{i+1}}_{a_i a_{i+1}} \qquad (6)

where: Utiti+1sisi+1U^{s_i s_{i+1}}_{t_i t_{i+1}} — order-4 gate tensor on sites i,i+1i,i+1; Θ\Theta — merged two-site tensor, re-split by SVD and truncated; ai−1,ai+1a_{i-1},a_{i+1} — outer bonds preserved during the local update. Read it as: the entire cost of one gate is the cost of one small SVD, which is O(χ3)O(\chi^3). The complementary paradigm treats the whole circuit as one network and asks only for a scalar — an amplitude or an expectation value — by summing over every shared index in a clever order. The cost of the best order is controlled by a single graph quantity, the treewidth.

contraction cost=2 O(tw(G))⋅poly(∣G∣)(7)\text{contraction cost}=2^{\,O(\mathrm{tw}(G))}\cdot\mathrm{poly}(|G|) \qquad (7)

where: GG — the circuit graph (vertices are tensors, edges are shared indices); tw(G)\mathrm{tw}(G) — its treewidth; for a depth-dd 1D circuit tw=O(d)\mathrm{tw}=O(d), for a 2D L×LL\times L grid tw=O(L)\mathrm{tw}=O(L). Read it as: simulation is fixed-parameter tractable in "how tangled the circuit is," and treewidth is the precise measure of that tangle.

Derivation AWhy the bond dimension equals the Schmidt rank

The goal is to prove the theorem that makes entanglement the exact currency of classical cost: the minimal bond dimension of an MPS at any cut is neither more nor less than the Schmidt rank of the state across that cut. This is what licenses the whole enterprise — it says an MPS is not a lossy heuristic but an exact representation whose size is pinned by a physical property of the state, so a state with little entanglement is provably cheap and a state with a lot is provably expensive. We place the MPS in mixed canonical form about the bond between sites ii and i+1i+1 and read off the Schmidt decomposition directly.

Step 1. Bring every tensor left of the cut into left-canonical gauge, ∑s(A[j]s)†A[j]s=1\sum_s (A^{[j]s})^\dagger A^{[j]s}=\mathbb{1} for j≤ij\le i, by a sweep of QR or SVD decompositions. Why: this makes the left tensors isometries, so the vectors they define on the bond index aia_i are automatically orthonormal — they are the Schmidt vectors ∣αk⟩A|\alpha_k\rangle_A.

Step 2. Symmetrically bring every tensor right of the cut into right-canonical gauge, ∑sA[j]s(A[j]s)†=1\sum_s A^{[j]s}(A^{[j]s})^\dagger=\mathbb{1} for j>ij>i. Why: now the right tensors are isometries too, giving orthonormal Schmidt vectors ∣βk⟩B|\beta_k\rangle_B indexed by the same bond aia_i, so both halves are expressed in orthonormal bases labelled by one index.

Step 3. With both gauges fixed, the diagonal bond matrix Λ\Lambda sitting on the cut carries the remaining weights, and the state reads ∣ψ⟩=∑ai=1χiΛai∣αai⟩A∣βai⟩B|\psi\rangle=\sum_{a_i=1}^{\chi_i}\Lambda_{a_i}|\alpha_{a_i}\rangle_A|\beta_{a_i}\rangle_B. Why: this is literally the Schmidt decomposition of equation (3) with λai=Λai\lambda_{a_i}=\Lambda_{a_i}, so the bond index and the Schmidt index are the same object.

Step 4. Compare ranks. Any MPS with bond dimension χi\chi_i produces a Schmidt form of rank at most χi\chi_i, so χi≥rank ρAi\chi_i\ge\mathrm{rank}\,\rho_{A_i}; the canonical form just built attains equality, giving χimin⁡=rank ρAi=#{λk>0}\chi_i^{\min}=\mathrm{rank}\,\rho_{A_i}=\#\{\lambda_k>0\}. Why: no gauge choice can represent more Schmidt terms than the bond can hold, and the canonical construction shows none are wasted — the minimum is exactly the Schmidt rank.

Go deeper · Derivation B · 1 minOptimality of singular-value truncation and the certified error
Derivation BOptimality of singular-value truncation and the certified error

The aim is to show that throwing away the smallest Schmidt values is not merely convenient but provably optimal, and that the price you pay is exactly the discarded weight — a quantity you can print next to your answer as a warranty. Without this result, truncation would be an uncontrolled approximation and every MPS circuit simulation would be a guess; with it, TEBD becomes a method whose error is known in advance and accumulates predictably. We take the merged two-site tensor Θ\Theta from equation (6), read off its SVD, and ask for the best rank-χ′\chi' replacement in Frobenius norm.

Step 1. Write Θ=∑kλk ukvk†\Theta=\sum_k\lambda_k\,u_k v_k^\dagger in its singular-value decomposition and note the Frobenius norm is unitarily invariant, so ∥Θ−Θ~∥F2=Tr[(Θ−Θ~)†(Θ−Θ~)]\|\Theta-\tilde\Theta\|_F^2=\mathrm{Tr}[(\Theta-\tilde\Theta)^\dagger(\Theta-\tilde\Theta)] may be evaluated in the singular basis. Why: working in the SVD basis turns a matrix-approximation problem into a plain sum over singular values, where the optimum is transparent.

Step 2. Invoke the Eckart–Young–Mirsky theorem: among all matrices of rank at most χ′\chi', the minimizer keeps the χ′\chi' largest singular values, Θ~⋆=∑k=1χ′λkukvk†\tilde\Theta^\star=\sum_{k=1}^{\chi'}\lambda_k u_k v_k^\dagger, and every other rank-χ′\chi' matrix does at least as badly. Why: this is the classical linear-algebra fact that makes "keep the biggest, drop the smallest" the correct and not merely a plausible rule.

Step 3. Evaluate the residual: min⁡rank Θ~≤χ′∥Θ−Θ~∥F2=∑k>χ′λk2=ϵ\min_{\mathrm{rank}\,\tilde\Theta\le\chi'}\|\Theta-\tilde\Theta\|_F^2=\sum_{k>\chi'}\lambda_k^2=\epsilon, and because ∑kλk2=∥∣ψ⟩∥2=1\sum_k\lambda_k^2=\||\psi\rangle\|^2=1 the discarded weight is directly the squared state-distance error. Why: the error you make is nothing but the tail of the Schmidt spectrum you deleted, so it is computed for free during the SVD you were already doing.

Step 4. Propagate over the circuit: renormalizing after each truncation, the total fidelity is bounded below by the product of per-gate factors, F≳∏gates(1−ϵg)F\gtrsim\prod_{\text{gates}}(1-\epsilon_g). Why: multiplying the small per-gate infidelities gives a running, certifiable accuracy budget for the whole simulation, which is exactly what you report as its warranty.

Go deeper · Derivation C · 2 minThe O(n χ³) cost of one TEBD sweep
Derivation CThe O(n χ³) cost of one TEBD sweep

The purpose is to pin down, operation by operation, why a one-dimensional circuit simulates in time polynomial in the qubit number and cubic in the bond dimension, so that "efficient for area-law states" is a precise flop count rather than a slogan. This is the number a practitioner budgets against: it tells you that doubling χ\chi costs eight times the work, and that the qubit number enters only linearly, which is why 1D DMRG and TEBD scale to hundreds of sites. We count a single two-site update with physical dimension d=2d=2 and bond dimension χ\chi, then multiply out over a full brick-wall layer.

Step 1. Form the two-site tensor by contracting two χ×d×χ\chi\times d\times\chi site tensors over their shared bond, Θ(asi),(si+1b)=∑aiA[i]A[i+1]\Theta_{(a s_i),(s_{i+1} b)}=\sum_{a_i}A^{[i]}A^{[i+1]}, at a cost of O(d2χ3)O(d^2\chi^3). Why: this contraction is a matrix product over the internal bond, and matrix products over a dimension-χ\chi index carry the cubic factor that will dominate.

Step 2. Apply the gate by multiplying Θ\Theta against the d2×d2d^2\times d^2 gate matrix on its physical indices, costing O(d4χ2)O(d^4\chi^2). Why: the gate touches only the two physical legs, so it is cheaper than the SVD by a factor of χ\chi and is subleading for d=2d=2.

Step 3. Split the result back into two site tensors by the SVD Θ′=UΣV†\Theta'=U\Sigma V^\dagger, keep the largest χ\chi singular values, and reabsorb Σ\Sigma to restore canonical form; the SVD of the (dχ)×(dχ)(d\chi)\times(d\chi) matrix dominates at O(d3χ3)O(d^3\chi^3). Why: the SVD is the expensive step, and its cubic scaling in χ\chi sets the cost of the whole update — everything else is lower order.

Step 4. Multiply out one brick-wall layer of O(n)O(n) gates and a depth-dd circuit of dd layers: Cgate=O(χ3)C_{\text{gate}}=O(\chi^3), Csweep=O(nχ3)C_{\text{sweep}}=O(n\chi^3), and Ccircuit=O(ndχ3)C_{\text{circuit}}=O(nd\chi^3), valid only while χ\chi stays bounded. Why: the linear count of gates per layer and layers per circuit turns the per-gate cube into the canonical O(ndχ3)O(nd\chi^3) — polynomial precisely when the entanglement, and hence χ\chi, does not blow up. *As a concrete check: for χ=64\chi=64, χ3=2.62×105\chi^3=2.62\times10^5, so one n=50n=50 sweep costs ≈1.3×107\approx 1.3\times10^7 multiply-adds — confirming the linear-in-nn, cubic-in-χ\chi scaling.*

Limiting cases

Whether a circuit is classically simulable is decided by a single question — how the entanglement entropy, and hence the bond dimension, scales with system size nn and depth dd:

CaseResultMeaning
Area-law / gapped 1D ground statesS(ρA)=O(1) ⇒ χ=O(1)S(\rho_A)=O(1)\ \Rightarrow\ \chi=O(1)Constant bond dimension suffices; DMRG and TEBD simulate these in time O(nχ3)O(n\chi^3) — the regime where tensor networks dominate outright.
Shallow / low-treewidth circuitscost=2O(tw(G))⋅poly(n)\text{cost}=2^{O(\mathrm{tw}(G))}\cdot\mathrm{poly}(n)Markov–Shi: contraction is fixed-parameter tractable in the circuit-graph treewidth; a depth-dd circuit on a path has tw=O(d)\mathrm{tw}=O(d).
Clifford circuitscost=O(poly(n))\text{cost}=O(\mathrm{poly}(n)) regardless of SSGottesman–Knill: the stabilizer tableau simulates even maximally entangled Clifford circuits efficiently — entanglement alone is not sufficient for hardness; non-stabilizer "magic" is also required.
Generic deep / 2D supremacy circuitsS(ρA)∼volume ⇒ χ∼2n/2S(\rho_A)\sim\text{volume}\ \Rightarrow\ \chi\sim 2^{n/2}Bond dimension grows exponentially and exact PEPS contraction is #P-hard; this is the hardness foundation of random-circuit-sampling experiments.

Worked examples

1 · Truncation error and the fidelity budget (fully worked)

After a gate you compute the normalized Schmidt spectrum of squared singular values {0.80, 0.15, 0.04, 0.01}\{0.80,\,0.15,\,0.04,\,0.01\} across a cut. You truncate to bond dimension χ′=2\chi'=2. Find the discarded weight, the fidelity of this single step, and then the largest per-gate error you can tolerate if a 10410^4-gate circuit must finish with fidelity at least 0.90.9.

Step 1. The squared singular values already sum to 0.80+0.15+0.04+0.01=10.80+0.15+0.04+0.01=1, so the state is normalized and the λk2\lambda_k^2 are genuine Schmidt probabilities. Truncating to χ′=2\chi'=2 keeps the two largest, {0.80,0.15}\{0.80,0.15\}, and drops {0.04,0.01}\{0.04,0.01\}.

Step 2. By equation (5) the discarded weight is the sum of the dropped squared values, ϵ=0.04+0.01=0.05\epsilon=0.04+0.01=\mathbf{0.05}. The retained weight is 1−ϵ=0.951-\epsilon=0.95.

Step 3. Renormalizing the truncated state and overlapping it with the original gives fidelity Fstep=∣⟨ψ∣ψ~⟩∣2=1−ϵ=0.95F_{\text{step}}=|\langle\psi|\tilde\psi\rangle|^2=1-\epsilon=\mathbf{0.95} — the discarded weight is the step infidelity, exactly as the Eckart–Young derivation promised. A 5% error in one step is enormous; a real simulation must keep ϵ\epsilon orders of magnitude smaller.

Step 4. Over G=104G=10^4 gates the fidelities multiply, F≳(1−ϵg)GF\gtrsim(1-\epsilon_g)^{G}. Requiring F≥0.9F\ge 0.9 gives ϵg≤−ln⁡(0.9)/G=0.1054/104≈1.05×10−5\epsilon_g\le -\ln(0.9)/G=0.1054/10^4\approx\mathbf{1.05\times10^{-5}}. Each gate may discard at most about 10−510^{-5} of the weight; at the illustrative ϵg=10−4\epsilon_g=10^{-4} the same circuit would finish at F≈(1−10−4)104≈e−1=0.37F\approx(1-10^{-4})^{10^4}\approx e^{-1}=0.37 — a stark demonstration that truncation error compounds and must be budgeted from the outset.

2 · Exponential compression versus the volume-law wall

Compare the classical cost of storing a 5050-qubit state exactly, as an MPS with bond dimension χ=64\chi=64, and as an MPS forced to represent a volume-law state. Use complex128 (16 bytes per amplitude).

Step 1. Exact storage of the full amplitude tensor needs 250=1.13×10152^{50}=\mathbf{1.13\times10^{15}} complex numbers, which at 16 bytes each is 1.8×10161.8\times10^{16} bytes =18 petabytes=\mathbf{18\ \text{petabytes}} — beyond any single machine.

Step 2. The MPS at χ=64\chi=64 stores n d χ2=50×2×642=4.10×105n\,d\,\chi^2=50\times2\times64^2=\mathbf{4.10\times10^{5}} complex numbers, about 6.66.6 megabytes — a compression of 250/(ndχ2)≈2.7×1092^{50}/(n d\chi^2)\approx\mathbf{2.7\times10^{9}}. This is the payoff of bounded entanglement: a half-billion-fold shrink that turns an impossible problem into a laptop one.

Step 3. By equation (4), that χ=64\chi=64 can only hold entanglement up to S≤log⁡264=6S\le\log_2 64=\mathbf{6} bits across each cut. A generic deep circuit drives the half-chain entropy toward the volume-law maximum S∼n/2=25S\sim n/2=25 bits, which by equation (4) demands χ∼225≈3.4×107\chi\sim 2^{25}\approx\mathbf{3.4\times10^{7}}.

Step 4. At that bond dimension the MPS would store ndχ2≈50×2×(3.4×107)2≈1.2×1017n d\chi^2\approx 50\times2\times(3.4\times10^7)^2\approx\mathbf{1.2\times10^{17}} numbers — larger than the 2502^{50} it was meant to compress, plus a TEBD cost scaling as χ3≈4×1022\chi^3\approx4\times10^{22} per gate. The method that beat 2n2^n by nine orders of magnitude for a low-entanglement state offers no advantage at all once the state is volume-law — the compression lives and dies with the entanglement.

Physics of the portrait

The Matrix Product State is where a computational trick from condensed-matter physics turned into the sharpest classical lens on quantum computing, and its history is a relay of ideas passed between the two fields.

Origins · 1 minmatrix product state and how this was discovered
MP
matrix product state
the χ-bond factorization of a many-body wavefunction, from White's DMRG in 1992 to Vidal's TEBD in 2003 · classical simulation of quantum states and circuits

The Matrix Product State is not a person but the central object of this subject, and its lineage runs through both condensed-matter physics and quantum information. Its computational birth was Steven White's Density-Matrix Renormalization Group in 1992, a variational method that reached astonishing precision for one-dimensional quantum lattices without anyone yet naming the state it produced. In 1995 Stellan Östlund and Stefan Rommer recognized that the DMRG fixed point was exactly a product of matrices — the MPS — retrofitting a transparent structure onto an empirical success. The quantum-computing chapter opened in 2003, when Guifré Vidal proved that a circuit generating only limited entanglement can be simulated classically in polynomial time by evolving an MPS and re-truncating by SVD, the TEBD algorithm, recasting entanglement as the resource that separates easy circuits from hard ones. A year later Frank Verstraete and Ignacio Cirac generalized the picture to two dimensions with Projected Entangled Pair States and rebuilt the whole formalism in the language of quantum information, and in 2007 Vidal added the Multiscale Entanglement Renormalization Ansatz for critical systems. The circuit-simulation viewpoint was completed in 2008 when Igor Markov and Yaoyun Shi showed that simulating a quantum circuit is contracting a tensor network, with a cost controlled by the graph's treewidth — the result behind every modern attempt to spoof a quantum-supremacy experiment on classical hardware.

History
1992
— White introduces the Density-Matrix Renormalization Group for 1D quantum lattices
1995
— Östlund and Rommer identify the DMRG fixed point as a matrix product state
2003
— Vidal's TEBD efficiently simulates slightly entangled quantum circuits
2004
— Verstraete and Cirac introduce PEPS and recast MPS in quantum-information language
2008
— Markov and Shi cast circuit simulation as tensor contraction with a treewidth bound

Going further

Texts and reviews. Nielsen and Chuang, Quantum Computation and Quantum Information (10th anniv. ed.), Ch. 2 (tensor products, Schmidt decomposition) · Preskill, Lecture Notes on Quantum Computation, Ch. 4 (entanglement, Schmidt rank) · Watrous, The Theory of Quantum Information, Ch. 2 (Schmidt decomposition, operator-Schmidt) · Kitaev, Shen, Vyalyi, Classical and Quantum Computation (treewidth and simulation complexity) · Schollwöck, "The density-matrix renormalization group in the age of matrix product states," Ann. Phys. 326, 96 (2011) · Vidal, "Efficient classical simulation of slightly entangled quantum computations," PRL 91, 147902 (2003) · Verstraete, Cirac, Murg, "Matrix product states, projected entangled pair states, and variational renormalization group methods," Adv. Phys. 57, 143 (2008) · Markov and Shi, "Simulating quantum computation by contracting tensor networks," SIAM J. Comput. 38, 963 (2008) · Orús, "A practical introduction to tensor networks," Ann. Phys. 349, 117 (2014). The key results used throughout are: the entropy bound S≤log⁡2χS\le\log_2\chi follows from concavity of the Shannon entropy; the discarded-weight error formula ϵ=∑k>χ′λk2\epsilon=\sum_{k>\chi'}\lambda_k^2 follows from Eckart–Young–Mirsky optimality and the normalisation ∑kλk2=1\sum_k\lambda_k^2=1; the TEBD cost is O(nχ3)O(n\chi^3); and the contraction cost is 2O(tw)poly(∣G∣)2^{O(\mathrm{tw})}\mathrm{poly}(|G|). The worked numbers are exact: 250=1.126×10152^{50}=1.126\times10^{15}; MPS storage 50⋅2⋅642=4.096×10550\cdot2\cdot64^2=4.096\times10^{5}, compression 2.7×1092.7\times10^{9}; truncation ϵ=0.05\epsilon=0.05, Fstep=0.95F_{\text{step}}=0.95; fidelity budget ϵg≤−ln⁡(0.9)/104≈1.05×10−5\epsilon_g\le-\ln(0.9)/10^4\approx1.05\times10^{-5} and (1−10−4)104≈0.37(1-10^{-4})^{10^4}\approx0.37; light-cone χ=2d/2\chi=2^{d/2} crossing χ=216\chi=2^{16} at d=32d=32. The Sycamore figures are from Google (2019) and the subsequent classical contraction results (Pan and Zhang; Gray and Kourtis).

Prerequisites. Parameterized Quantum Circuits and Ansatz Design · Variational Quantum Eigensolver · Expressibility, Entangling Capability, and Trainability Trade-offs. Leads to. Barren Plateaus in Variational Landscapes · Quantum Approximate Optimization Algorithm · Quantum Machine-Learning Models.

Ask the tutor: Ask me to walk the Schmidt-rank–bond-dimension theorem line by line, to work out the bond dimension and TEBD cost for a circuit of your size, or to trace how a tensor-network contraction spoofed a quantum-supremacy experiment.

Is a Matrix Product State an approximation or an exact representation?

Both, depending on the bond dimension you allow. Given a large enough χ\chi — up to 2n/22^{n/2} at the middle cut — an MPS represents any state exactly, because the bond dimension can always be made equal to the Schmidt rank. The approximation only enters when you cap χ\chi below what the state demands and truncate the Schmidt spectrum, at which point the error is the certified discarded weight ϵ\epsilon. So an MPS is an exact reparameterization that becomes an approximation exactly when, and by exactly the amount that, you refuse to pay for the entanglement present. For area-law states the exact χ\chi is already small, so no approximation is needed at all.

Why does entanglement, rather than sheer state complexity, set the cost?

Because the bond dimension across a cut is the Schmidt rank across that cut, and the Schmidt rank is a direct measure of how entangled the two halves are. A state can have all 2n2^n amplitudes nonzero and still be a product state with Schmidt rank one — utterly cheap — while a state with far fewer nonzero amplitudes can be maximally entangled and expensive. The MPS is blind to everything about the state except how many Schmidt channels cross each cut, so it is entanglement, not amplitude count or any other notion of complexity, that it charges you for. This is why the whole subject reframes "is this circuit hard?" as "how much entanglement does it build?"

If tensor networks can spoof supremacy experiments, was quantum advantage disproved?

No — it was sharpened. The classical contraction results did not simulate the quantum device faithfully; they exploited that the specific circuits were shallow enough, and the fidelity low enough, that a cleverly ordered tensor contraction could reproduce the sampling statistics within the same low fidelity the hardware achieved. What they demonstrated is that the particular depth-20, 53-qubit instance was not as far beyond classical reach as first estimated, which pushes experimentalists toward deeper circuits and larger, higher-treewidth hardware graphs. The contest is a moving boundary: each side forces the other to be more precise about exactly where classical simulation becomes intractable, and that boundary is the real scientific output.

Why is two dimensions so much harder than one, even at fixed bond dimension?

Because contracting the network is easy only when the graph is nearly a path. A 1D MPS is a chain with treewidth O(1)O(1), so its norm and expectation values contract in O(nχ3)O(n\chi^3) by sweeping transfer matrices. A 2D PEPS lives on a grid whose treewidth scales as its linear size LL, and by equation (7) the contraction cost is 2O(L)2^{O(L)}; worse, exact contraction of a PEPS is #P-hard in general, so no ordering rescues it. Raising χ\chi does not help because the obstruction is the connectivity of the network, not the size of the tensors — the extra spatial dimension gives entanglement too many paths to spread, and there is no one-dimensional order in which to sum the indices cheaply.

What exactly does the treewidth measure, and why does contraction cost depend on it?

Treewidth is a graph invariant that quantifies how far a graph is from being a tree: a tree has treewidth one, a path has treewidth one, and a densely connected grid has treewidth growing with its size. When you contract a tensor network you sum over shared indices in some order, and each step creates an intermediate tensor whose number of open legs is bounded by the width of the tree decomposition you are following. A width-ww decomposition guarantees no intermediate ever has more than w+1w+1 legs, so each summation costs at most 2w+12^{w+1} and the whole contraction costs 2O(w)poly(∣G∣)2^{O(w)}\mathrm{poly}(|G|). Treewidth is therefore the precise, order-independent measure of "how tangled" a circuit is, and it is what makes shallow, geometrically local circuits tractable and deep all-to-all circuits hard.

How do I know my truncated simulation is trustworthy?

Watch the discarded weight. After every SVD truncation you get ϵ\epsilon for free — it is the sum of the squared singular values you dropped — and the running product ∏g(1−ϵg)\prod_g(1-\epsilon_g) lower-bounds the fidelity of your whole simulation against the exact untruncated result. If that product stays close to one, your answer is trustworthy to that fidelity; if it collapses, your bond dimension is too small for the entanglement the circuit is generating and the result is unreliable no matter how clean it looks. A simulation reported without its discarded weight or final fidelity is uninterpretable, because the same tidy-looking output could be exact or could be dominated by truncation. Convergence in χ\chi — checking that the answer stops moving as you raise the bond dimension — is the companion diagnostic.

Before this

Where this leads