Boson Sampling

Boson sampling is the sharpest answer physics has yet given to a deceptively simple question: what is the least a quantum device can do and still be impossible to simulate on any classical computer? Scott Aaronson and Alex Arkhipov's 2011 answer strips a quantum computer down to its barest optical bones. Send $n$ indistinguishable single photons into an $m$ -mode interferometer built only from beam splitters and phase shifters — no nonlinearity, no entangling gates, no measurement-dependent feedforward — and simply record which detectors click. The machine computes no function and answers no yes/no question. It only draws samples from the distribution of output photon patterns. Yet that pale, passive act is, under widely believed conjectures, beyond the reach of every classical algorithm. The reason lives in one algebraic accident. Because photons are bosons, their many-particle amplitudes are built by symmetrizing over all the ways indistinguishable particles can be routed to detectors, and that symmetric sum is exactly the permanent of a submatrix of the interferometer's unitary. The permanent looks like the determinant with every minus sign erased — and that erased sign is the difference between a quantity Gaussian elimination handles in $O(n^3)$ steps and one Valiant proved is $\#\mathrm{P}$ -hard, as hard as counting the satisfying assignments of any Boolean formula. Fermions antisymmetrize and give determinants, so free-fermion optics is classically trivial; bosons symmetrize and give permanents, so free-boson optics is classically forbidding. This article follows that single sign from the Hong–Ou–Mandel dip through the Aaronson–Arkhipov collapse theorem — an efficient classical simulator would collapse the polynomial hierarchy — out to Jiuzhang and Borealis, the photonic machines that turned this complexity gadget into a claim of quantum advantage.

Key takeaways

Boson sampling is a deliberately non-universal model: $n$ identical photons pass through an $m$ -mode passive linear-optical interferometer described by an $m\times m$ unitary $U$ , and one samples the output detection pattern — nothing more.
The transition amplitude between input and output Fock states is the permanent of an $n\times n$ submatrix of $U$ ; output probabilities are squared moduli of permanents divided by occupation factorials.
Computing the permanent is $\#\mathrm{P}$ -hard (Valiant, 1979) while the determinant is $O(n^3)$ -easy; the two differ only by the sign $\mathrm{sgn}(\sigma)$ , so free bosons are hard for exactly the reason free fermions are easy.
Aaronson and Arkhipov proved that an efficient exact classical simulator would place a $\#\mathrm{P}$ -hard problem in $\mathrm{BPP}^{\mathrm{NP}}$ and collapse the polynomial hierarchy to $\mathrm{PH}=\Sigma_3^p$ — a consequence essentially nobody believes.
The approximate-sampling hardness rests on two conjectures — the Permanent-of-Gaussians Conjecture and the Permanent Anti-Concentration Conjecture — made relevant by taking $U$ Haar-random and $m\gtrsim n^2$ , so submatrices look i.i.d. complex Gaussian.
Gaussian boson sampling swaps single photons for squeezed light, easing sources and turning permanents into hafnians (still $\#\mathrm{P}$ -hard); loss and distinguishability are the experimental adversaries that can restore classical simulability.
Key idea: indistinguishable bosons through a linear interferometer interfere so that amplitudes are matrix permanents, and because the permanent is $\#\mathrm{P}$ -hard while the determinant of free fermions is easy, sampling photon patterns is classically intractable unless the polynomial hierarchy collapses.

The core insight

Indistinguishable bosons passing through a linear interferometer interfere so that transition amplitudes are matrix permanents; because the permanent is $\#\mathrm{P}$ -hard while the determinant (free fermions) is easy, sampling photon patterns is classically intractable unless the polynomial hierarchy collapses.

Where this lives today

Boson sampling has become the photonic front line of the quantum-advantage debate. The 2020 USTC Jiuzhang experiment (Pan Jianwei, Lu Chaoyang and collaborators) used squeezed-light sources, a 100-mode interferometer, and threshold detectors to detect up to about 76 photons, and estimated that reproducing its output would take the then-fastest supercomputer billions of years; Jiuzhang 2.0 and 3.0 pushed the photon numbers and detection rates further. In 2022 Xanadu's Borealis demonstrated a programmable, time-multiplexed Gaussian boson sampler with 216 squeezed modes, adding the reconfigurability earlier demonstrations lacked. But the story is a moving target: classical algorithms have repeatedly clawed back ground by exploiting the very imperfections — photon loss, partial distinguishability, finite squeezing — that real devices suffer, with tensor-network and Markov-chain "spoofers" (for instance the loss-exploiting algorithms of Oh, Jiang, and collaborators) matching or beating some advantage claims on the relevant benchmarks. The open theoretical questions are exactly the load-bearing conjectures: whether the Permanent-of-Gaussians Conjecture is true, whether anti-concentration holds, and precisely where the loss threshold sits that separates a genuinely hard sampler from a classically simulable one. On the applications side, because Gaussian-boson-sampling probabilities are hafnians that count weighted perfect matchings, groups have proposed GBS devices as special-purpose heuristics for molecular vibronic spectra, dense-subgraph finding, and graph-similarity kernels — useful to watch, but without the advantage proofs that the bare sampling task enjoys.

The physical picture

Start with two photons and one beam splitter, the atom of the whole subject. Two identical photons arrive at the two input ports of a balanced 50:50 splitter at the same instant. Classically you would expect them to sort independently, giving a coincidence (one photon out each side) half the time. Quantum mechanically the two amplitudes for a coincidence — both transmit, or both reflect — are equal in size and opposite in sign, so they cancel exactly. The photons refuse to leave separately; they always bunch into the same output port. This is the Hong–Ou–Mandel effect, and the cancellation is nothing other than a $2\times2$ permanent evaluating to zero. Every richer boson-sampling instance is this dip scaled up: for $n$ photons in $m$ modes, the amplitude for a given output pattern is an $n\times n$ permanent, a coherent sum over all $n!$ ways indistinguishable particles can be assigned to detectors.

Figure

Hong–Ou–Mandel bunching as a 2×2 permanent: two identical photons enter a balanced 50:50 beam splitter in state 1,1, and the coincidence amplitude equals Per of the beam-splitter matrix = (1/√2)(−1/√2) + (1/√2)(1/√2) = 0, so the photons never exit separately — they always bunch, giving P(2,0) = P(0,2) = ½ and P(1,1) = 0

The single sign that makes the coincidence cancel is not incidental — it is the entire complexity story in miniature. Antisymmetrize the multi-particle wavefunction, as fermions demand, and the routing sum acquires a $\mathrm{sgn}(\sigma)$ factor, turning it into a determinant that Gaussian elimination evaluates in $O(n^3)$ time; free-fermion (matchgate) optics is therefore classically simulable. Symmetrize, as bosons demand, and every sign becomes $+1$ , turning the sum into a permanent whose best known exact algorithm (Ryser's) still costs $O(n\,2^n)$ and whose exact evaluation Valiant proved $\#\mathrm{P}$ -hard. Bosons interfere in a way no efficient bookkeeping can untangle. The second figure sets the two branches side by side.

Figure

The determinant–permanent fork drawn from one interferometer: n photons through an m-mode unitary U give output probabilities equal to the squared modulus of Per(U from S,T). Fermionic antisymmetry would give det(A) = Σ sgn(σ) ∏ A, computable by Gaussian elimination in O(n³) — free fermions are easy. Bosonic symmetry gives Per(A) = Σ ∏ A, whose best exact algorithm is Ryser's O(n·2ⁿ) and which is #P-hard — the one missing sign is the whole gap between tractable and classically intractable

Why sample rather than compute? A boson sampler cannot hand you the value of any single permanent: each output probability is exponentially small, so reading one off to useful precision would demand exponentially many repetitions. The machine's power is purely distributional — it draws from a distribution defined by exponentially many permanents at once. The hardness proof turns exactly this feature into a weapon, feeding the hypothetical classical sampler's existence into Stockmeyer's approximate-counting theorem. And the hardness is robust precisely because the interferometer is taken Haar-random and large, $m\sim n^2$ : in that dilute regime an $n\times n$ submatrix of $U$ behaves like an i.i.d. complex Gaussian matrix, so the operative statement concerns permanents of generic Gaussian instances, not hand-crafted worst cases.

The core equations

The central formula of the model reads off a probability from a permanent. A passive interferometer sends each input photon into a coherent superposition of output modes, and the amplitude to find a specified output occupation is the coherent sum over all assignments of photons to detectors — the permanent of the submatrix of $U$ picked out by the input and output patterns. The occupation factorials normalize for the bosonic bunching of repeated modes.

P(S) = \frac{\big|\mathrm{Per}(U_{S,T})\big|^2}{s_1!\cdots s_m!\, t_1!\cdots t_m!} \qquad (1)

where: $P(S)$ — probability of the output occupation pattern $S=(s_1,\dots,s_m)$ ; $U$ — the $m\times m$ unitary of the linear-optical network; $U_{S,T}$ — the $n\times n$ matrix formed by taking $t_j$ copies of column $j$ and $s_i$ copies of row $i$ of $U$ ; $s_i,t_i$ — output/input occupation of mode $i$ , with $\sum_i s_i=\sum_i t_i=n$ ; $\mathrm{Per}(A)=\sum_{\sigma\in\mathcal S_n}\prod_i A_{i\sigma(i)}$ . Reading: the whole physics is compressed into which submatrix $U_{S,T}$ each detection pattern selects — the interferometer realizes a chosen $U$ , and the photons compute the permanent. The permanent and determinant are twins separated by a single factor. Written as sums over the symmetric group, they run over the same $n!$ products of matrix entries; the determinant weights each by the sign of the permutation, the permanent weights each by $+1$ .

\mathrm{Per}(A) = \sum_{\sigma\in\mathcal S_n} \prod_{i=1}^n A_{i,\sigma(i)}, \qquad \det(A) = \sum_{\sigma\in\mathcal S_n} \mathrm{sgn}(\sigma)\prod_{i=1}^n A_{i,\sigma(i)} \qquad (2)

where: $\mathcal S_n$ — the group of permutations of $n$ symbols; $\mathrm{sgn}(\sigma)=\pm1$ — the parity of $\sigma$ ; $A$ — an $n\times n$ matrix. Reading: the alternating sign lets row operations simplify a determinant without changing its value; strip it and no such reduction survives. Valiant (1979) proved the permanent of a $\{0,1\}$ matrix is $\#\mathrm{P}$ -complete — it counts perfect matchings of a bipartite graph, encoding the solution count of any NP problem. The best general exact algorithm for the permanent is Ryser's inclusion–exclusion formula. It does not beat the exponential wall, but it replaces the naive $n\cdot n!$ of summing permutations with a sum over the $2^n$ subsets of columns.

\mathrm{Per}(A) = (-1)^n \sum_{S\subseteq\{1,\dots,n\}} (-1)^{|S|} \prod_{i=1}^n \sum_{j\in S} A_{ij} \qquad (3)

where: $S$ — a subset of column indices; $|S|$ — its cardinality; the outer sign $(-1)^n$ and inner $(-1)^{|S|}$ implement inclusion–exclusion. Reading: the cost is $O(n\,2^n)$ arithmetic operations. This exponential is the yardstick every classical boson-sampling simulator is measured against — and the reason a few dozen photons already strain the largest supercomputers. For the approximate-hardness argument the relevant object is not a worst-case matrix but a random Gaussian one, and a clean exact identity anchors the analysis: the mean squared modulus of the permanent of an i.i.d. standard complex Gaussian matrix is exactly $n!$ .

\mathbb E\big[\,|\mathrm{Per}(X)|^2\,\big] = n!\,, \qquad X_{ij}\ \text{i.i.d. } \mathcal N(0,1)_{\mathbb C} \qquad (4)

where: $\mathbb E$ — expectation over the Gaussian entries; $\mathcal N(0,1)_{\mathbb C}$ — a standard complex Gaussian, $\mathbb E[X_{ij}\overline{X_{kl}}]=\delta_{ik}\delta_{jl}$ . Reading: $n!$ is the natural scale of a Gaussian permanent, so the Permanent-of-Gaussians Conjecture asks for additive approximation to $\pm\,\epsilon\,n!$ ; that this typical size is large — not concentrated near zero — is the content of the anti-concentration conjecture that the proof also needs. The payoff is the collapse chain. An efficient classical sampler, fed into Stockmeyer's counting theorem, would let an $\mathrm{NP}$ -oracle machine estimate an individual output probability — hence a $\#\mathrm{P}$ -hard permanent — and Toda's theorem then drags the entire polynomial hierarchy down to its third level.

\mathcal A\in\mathrm{FBPP} \;\Longrightarrow\; \#\mathrm{P}\subseteq\mathrm{BPP}^{\mathrm{NP}} \;\xrightarrow{\ \text{Toda}\ }\; \mathrm{PH}=\Sigma_3^p \qquad (5)

where: $\mathcal A$ — a hypothetical exact classical sampler in randomized polynomial time $\mathrm{FBPP}$ ; $\mathrm{BPP}^{\mathrm{NP}}$ — randomized polynomial time with an NP oracle; $\mathrm{PH}$ — the polynomial hierarchy; $\Sigma_3^p$ — its third level. Reading: the conclusion $\mathrm{PH}=\Sigma_3^p$ is disbelieved as strongly as $\mathrm P\neq\mathrm{NP}$ , so the contrapositive is read as evidence that no efficient exact classical sampler exists — a hardness statement resting on standard complexity conjectures rather than on any physical assumption. Finally, the experimentally dominant variant replaces single-photon Fock inputs with squeezed vacuum, keeping the state Gaussian through the network and converting permanents into hafnians of a matrix built from the squeezing and $U$ .

P_{\mathrm{GBS}}(S) = \frac{1}{\sqrt{\det(\sigma_Q)}}\,\frac{\mathrm{Haf}(A_S)}{s_1!\cdots s_m!} \qquad (6)

where: $\mathrm{Haf}(A)=\sum_{\mathrm{PM}}\prod_{(i,j)}A_{ij}$ — the hafnian, summed over perfect matchings of $2k$ indices (the symmetric-matrix analogue of the permanent); $A_S$ — the submatrix selected by the click pattern $S$ ; $\sigma_Q$ — the $Q$ -function covariance matrix normalizing the Gaussian state. Reading: the hafnian is also $\#\mathrm{P}$ -hard, so the hardness survives; and because the hafnian counts weighted perfect matchings of a graph, GBS devices double as heuristic samplers for graph and chemistry problems.

Derivation AWhy the output amplitude is a permanent

The goal is to derive equation (1) from first principles — the second-quantized action of a passive interferometer on Fock states — and to see the permanent emerge not as an analogy but as a forced consequence of bosonic commutation. We take $n$ photons entering in distinct modes (a collision-free input pattern $T$ ) and ask for the amplitude to detect them in a collision-free output pattern $S$ , then contract creation and annihilation operators until the sum over routings appears literally as $\mathrm{Per}(U_{S,T})$ .

Step 1. A passive network implements a unitary $\varphi(U)$ acting on creation operators by the interferometer matrix, $\varphi(U)\,\hat a_j^\dagger\,\varphi(U)^\dagger=\sum_k U_{kj}\,\hat b_k^\dagger$ . Why: beam splitters and phase shifters are quadratic in the mode operators, so each input mode fans out into a fixed superposition of output modes read from the columns of $U$ .

Step 2. Build the input as a product of creation operators, $|T\rangle=\prod_{j\in T}\hat a_j^\dagger\,|0\rangle$ , and evolve it to $\varphi(U)|T\rangle=\prod_{j\in T}\big(\sum_k U_{kj}\hat b_k^\dagger\big)|0\rangle$ . Why: since $\varphi(U)|0\rangle=|0\rangle$ , evolving the state just evolves each creation operator.

Step 3. Project onto the collision-free output $\langle S|=\langle 0|\prod_{i\in S}\hat b_i$ , one annihilation per occupied mode. Why: reading a detector pattern applies the matching annihilation operators, one clean $\hat b_i$ per clicked detector.

Step 4. Contract. The surviving terms pair each output annihilation with a distinct input factor — a bijection $\sigma$ — giving $\langle 0|\prod_i\hat b_i\prod_j(\sum_k U_{kj}\hat b_k^\dagger)|0\rangle=\sum_{\sigma\in\mathcal S_n}\prod_i U_{i,\sigma(i)}$ . Why: the vacuum expectation forces every creation to be annihilated once, and those pairings are the symmetric group — with no sign, because bosonic operators commute.

Step 5. The unit-weight sum over $\mathcal S_n$ is the permanent, $\langle S|\varphi(U)|T\rangle=\mathrm{Per}(U_{S,T})$ ; squaring with occupation factorials gives $P(S)=|\mathrm{Per}(U_{S,T})|^2/(\prod_i s_i!\prod_i t_i!)$ . Why: anticommuting (fermionic) operators would have inserted $\mathrm{sgn}(\sigma)$ and produced a determinant — the statistics pick the twin. *Verified: reproduces the standard Aaronson–Arkhipov (2013) result; the factorial normalization and sign-free contraction are the load-bearing steps.*

Go deeper · Derivation B · 1 minHong–Ou–Mandel bunching from a two-by-two permanent

Derivation BHong–Ou–Mandel bunching from a two-by-two permanent

The goal is to compute the simplest non-trivial boson-sampling instance completely — two photons on one balanced beam splitter — and to show that the famous coincidence cancellation is arithmetic on a $2\times2$ permanent, with the vanishing traceable to the very sign-freedom that makes the general problem hard. We fix the beam-splitter unitary, form the coincidence amplitude as a permanent, evaluate it, and then use probability conservation to pin the bunching probability at one.

Step 1. Take a balanced 50:50 beam splitter in the Hadamard convention, $U=\tfrac1{\sqrt2}\begin{pmatrix}1&1\\1&-1\end{pmatrix}$ , with one photon in each input mode, so $T=S=(1,1)$ . Why: with both input and output modes singly occupied, the selected submatrix $U_{S,T}$ is the full $2\times2$ matrix, and no occupation factorials complicate the amplitude.

Step 2. Write the coincidence amplitude as the permanent of that matrix, $\langle 1,1|\varphi(U)|1,1\rangle=\mathrm{Per}\begin{pmatrix}\tfrac1{\sqrt2}&\tfrac1{\sqrt2}\\\tfrac1{\sqrt2}&-\tfrac1{\sqrt2}\end{pmatrix}$ . Why: by the previous derivation the amplitude for one-in-one-out is exactly the permanent of the full interferometer matrix — this is boson sampling with $n=2$ , $m=2$ .

Step 3. Evaluate: the permanent adds the two diagonal-and-off-diagonal products with a plus sign, $\big(\tfrac1{\sqrt2}\big)\big(-\tfrac1{\sqrt2}\big)+\big(\tfrac1{\sqrt2}\big)\big(\tfrac1{\sqrt2}\big)=-\tfrac12+\tfrac12=0$ . Why: the two indistinguishable routing amplitudes are equal and opposite, so they interfere destructively — this exact cancellation is the Hong–Ou–Mandel dip, and it would read $-1$ (no cancellation of this kind) had we taken the determinant.

Step 4. Conclude by conservation. Since the total output probability is one and the coincidence channel $|1,1\rangle$ has amplitude zero, the two bunched outcomes $|2,0\rangle$ and $|0,2\rangle$ must share unit probability; the $\sqrt{2!}$ normalization of a doubly occupied mode splits it evenly, $P(2,0)=P(0,2)=\tfrac12$ . Why: the photons leave together with certainty — a purely quantum, statistics-driven correlation with no classical analogue, and the seed from which all boson-sampling interference grows. *Verified: matches the textbook HOM result; recomputed as $\mathrm{Per}=0$ and $P(\text{coincidence})=0$ , $P(\text{bunch})=1$ .*

Go deeper · Derivation C · 2 minHow an exact classical sampler collapses the polynomial hierarchy

Derivation CHow an exact classical sampler collapses the polynomial hierarchy

The goal is a self-contained sketch of the Aaronson–Arkhipov exact-case theorem: if some classical randomized algorithm samples exactly from the boson-sampling distribution in polynomial time, then the polynomial hierarchy collapses to its third level. The chain runs from embedding a hard permanent as an output probability, through Stockmeyer's approximate counting, to Toda's theorem — using no physical assumptions, only standard complexity theory.

Step 1. Encode a hard instance. Choose an interferometer and input so that a designated output $S^\star$ has probability $\Pr_D[S^\star]=|\mathrm{Per}(X)|^2/N$ for an essentially arbitrary integer matrix $X$ and a known normalization $N$ . Why: because output probabilities are squared permanents, the sampler's distribution secretly contains the value of any permanent we care to embed — the machine is a permanent oracle in disguise, if only we could read one probability.

Step 2. Apply Stockmeyer's theorem. Given the code of a classical sampler $\mathcal A$ for $D$ , an $\mathrm{FBPP}^{\mathrm{NP}}$ machine can multiplicatively approximate any individual probability $\Pr_D[S^\star]$ to within $(1\pm1/\mathrm{poly})$ . Why: Stockmeyer's approximate counting estimates the acceptance probability of any efficiently describable randomized process using an NP oracle — the existence of the classical circuit, not the physical device, is what this step consumes.

Step 3. Read off the permanent. Combining Steps 1 and 2, $|\mathrm{Per}(X)|^2$ becomes multiplicatively approximable in $\mathrm{BPP}^{\mathrm{NP}}$ ; for nonnegative integer permanents a good multiplicative estimate pins the exact value, which is $\#\mathrm{P}$ -hard. Why: a $\#\mathrm{P}$ -hard quantity has now been computed inside $\mathrm{BPP}^{\mathrm{NP}}$ , a class sitting low in the hierarchy — an apparent contradiction waiting for Toda to make it formal.

Step 4. Invoke Toda's theorem, $\mathrm{PH}\subseteq\mathrm{P}^{\#\mathrm{P}}$ . If $\#\mathrm{P}$ problems can be solved in $\mathrm{BPP}^{\mathrm{NP}}$ , then all of $\mathrm{PH}$ collapses into $\mathrm{BPP}^{\mathrm{NP}}\subseteq\Sigma_3^p$ . Why: Toda proved the entire polynomial hierarchy reduces to a single call to a counting oracle, so trivializing counting trivializes the hierarchy above the third level.

Step 5. Conclude. The polynomial hierarchy would collapse to $\mathrm{PH}=\Sigma_3^p$ ; since this is overwhelmingly believed false, no efficient exact classical sampler for boson sampling should exist. Why: this is an unconditional-in-spirit hardness result — it needs only the belief that $\mathrm{PH}$ is infinite; the harder approximate-sampling theorem additionally invokes the Permanent-of-Gaussians and anti-concentration conjectures. *Verified: this is the exact-case argument of Aaronson–Arkhipov (2013, Theorem 1); the logical chain Stockmeyer + Toda $\Rightarrow \mathrm{PH}=\Sigma_3^p$ is the standard formulation.*

Limiting cases

The model's character swings sharply with the particle statistics, the photon-to-mode ratio, and the coherence of the photons — each limit either preserves or destroys the classical hardness.

Case	Result	Meaning
Free fermions	$\langle f\|U\|f\rangle=\det(U_{S,T})$	Antisymmetric statistics give a determinant; matchgate circuits simulate classically in $O(n^3)$ — the easy complementary case that makes the boson result surprising.
Free bosons	$\langle s\|\varphi(U)\|s\rangle=\mathrm{Per}(U_{S,T})/\sqrt{\mu_S\mu_T}$	Symmetric statistics give a permanent; sampling is classically hard under standard conjectures.
Dilute regime	$m\gtrsim n^2$	Outputs are collision-free with high probability (each mode at most one photon), and $n\times n$ submatrices of Haar-random $U$ look Gaussian — the regime where the hardness proof applies.
Distinguishable photons	$P_{\mathrm{dist}}(S)=\mathrm{Per}\big(\|U_{S,T}\|^2\big)$	Losing coherence replaces the complex permanent by a permanent of a nonnegative matrix, approximable in randomized polynomial time (Jerrum–Sinclair–Vigoda) — the quantum advantage evaporates.
Heavy loss	effective $n_{\mathrm{eff}}\ll n$	Enough photon loss lets classical algorithms match the device by simulating only the few surviving photons; a loss threshold separates hard from easy.

Worked examples

1 · The permanent and determinant of the same matrix (fully worked)

Take the integer matrix $A=\begin{pmatrix}1&2&3\\4&5&6\\7&8&9\end{pmatrix}$ . Compute both its permanent and its determinant to see, concretely, how one sign flip changes everything. (The matrix is not a unitary submatrix — the point is the arithmetic that a real interferometer would run on complex entries.)

Step 1. List the six permutations of three columns and their signed products. The permanent and determinant sum the same six products $\prod_i A_{i,\sigma(i)}$ ; only the determinant attaches $\mathrm{sgn}(\sigma)$ . The products are $1\cdot5\cdot9=45$ , $1\cdot6\cdot8=48$ , $2\cdot4\cdot9=72$ , $2\cdot6\cdot7=84$ , $3\cdot4\cdot8=96$ , $3\cdot5\cdot7=105$ .

Step 2. Add with all-plus signs for the permanent: $\mathrm{Per}(A)=45+48+72+84+96+105=\mathbf{450}$ .

Step 3. Add with permutation signs for the determinant — the three even permutations $(45,84,96)$ enter with $+$ , the three odd ones $(48,72,105)$ with $-$ : $\det(A)=45-48-72+84+96-105=\mathbf{0}$ . The matrix is singular, yet its permanent is a robust $450$ .

Step 4. Read the lesson. The determinant's cancellation is what lets Gaussian elimination shortcut it; the permanent has no such structure, and for a genuine $30\times30$ or $50\times50$ interferometer submatrix there is no known shortcut past Ryser's $O(n\,2^n)$ . Per(A) = 450 while det(A) = 0 — the identical six products, differing only in sign, give a tractable quantity and an intractable one.

2 · How many photons make a classical supercomputer give up

A classical simulator evaluates a boson-sampling probability with Ryser's formula at cost $\approx n\,2^n$ operations. Estimate the cost at $n=30$ , $n=50$ , and $n=76$ photons, and connect the last figure to a real experiment.

Step 1. At $n=30$ : $2^{30}\approx1.07\times10^{9}$ , so $n\,2^n\approx 30\times1.07\times10^{9}\approx\mathbf{3.2\times10^{10}}$ operations — a fraction of a second on a laptop. Thirty photons are firmly classical.

Step 2. At $n=50$ : $2^{50}\approx1.13\times10^{15}$ , giving $n\,2^n\approx50\times1.13\times10^{15}\approx\mathbf{5.6\times10^{16}}$ operations — hours to days on a large cluster, and this is per probability; sampling requires estimating many. Fifty photons sit at the ragged edge.

Step 3. At $n=76$ : $2^{76}\approx7.56\times10^{22}$ , so $n\,2^n\approx76\times7.56\times10^{22}\approx\mathbf{5.7\times10^{24}}$ operations. At an optimistic $10^{18}$ operations per second (exascale), a single permanent would take $\sim\!5.7\times10^{6}$ seconds — about two months — and no exact sampler can afford one such evaluation per sample.

Step 4. Anchor it. The 2020 USTC Jiuzhang Gaussian-boson-sampling experiment detected up to $\sim76$ photons across 100 modes, precisely the regime where the exact classical cost explodes; this scaling is the quantitative backbone of the quantum-advantage claim. At n = 76 the exact classical cost (~5.7×10²⁴ operations) is beyond any existing computer — while the photonic device samples in microseconds.

Physics of the portrait

Boson sampling is a theorist's invention that became an experimentalist's flagship. It was born not in a laboratory but in a complexity-theory argument aimed squarely at a physical question: is nature, at the level of a few dozen non-interacting photons, already beyond classical reach?

Origins · 1 minScott Aaronson and how this was discovered

Scott Aaronson

b. 1981, quantum complexity theorist · co-author of the 2011 boson-sampling proposal and its hardness theorems

Scott Aaronson, with his then-student Alex Arkhipov at MIT, introduced boson sampling in "The Computational Complexity of Linear Optics" (STOC 2011; Theory of Computing 2013), and in doing so reshaped what "quantum advantage" could mean. The intellectual backdrop was a tension. Knill, Laflamme and Milburn had shown in 2001 that linear optics with adaptive measurement and feedforward is universal for quantum computation — but that adaptivity is experimentally brutal. Aaronson and Arkhipov asked the opposite question: how little can you keep and still be classically intractable? Their answer discards adaptivity entirely, keeping only Fock-state inputs, a passive Haar-random interferometer, and photon counting, and rests the hardness on the oldest hard object in the counting world — Leslie Valiant's 1979 proof that the permanent is $\#\mathrm{P}$ -complete. The masterstroke was to make hardness robust to noise by arguing about permanents of Gaussian matrices via the Permanent-of-Gaussians and anti-concentration conjectures, so that even an approximate, error-prone sampler would collapse the polynomial hierarchy. Aaronson, a prolific expositor whose blog Shtetl-Optimized shaped a generation's intuition about quantum computing, framed the result deliberately as evidence, not engineering: boson sampling is useless as a computer and precisely for that reason a clean probe of the Extended Church–Turing Thesis. Within two years small photonic demonstrations appeared; within a decade Jiuzhang and Borealis pushed it to the advantage frontier — and to the running battle with classical spoofing algorithms that continues today.

History

1979

— Valiant proves that computing the permanent is #P-complete

1987

— Hong, Ou and Mandel observe two-photon bunching, the HOM dip

2011

— Aaronson and Arkhipov propose boson sampling and prove the collapse theorems

2013

— first small photonic boson-sampling experiments demonstrate three- and four-photon interference

2020

— USTC's Jiuzhang claims quantum advantage via Gaussian boson sampling with ~76 detected photons

Going further

References. Aaronson and Arkhipov, "The Computational Complexity of Linear Optics," Theory of Computing 9 (2013) 143 (STOC 2011) — the founding paper and both hardness theorems · Valiant, "The complexity of computing the permanent," Theoretical Computer Science 8 (1979) 189 — $\#\mathrm{P}$ -completeness of the permanent · Nielsen and Chuang, Quantum Computation and Quantum Information, Ch. 4 and 7 — linear-optics formalism and complexity classes · Kitaev, Shen and Vyalyi, Classical and Quantum Computation — counting classes and the polynomial hierarchy · Preskill, Lecture Notes on Quantum Computation, Ch. 6 — context and the Extended Church–Turing Thesis · Hamilton et al., "Gaussian Boson Sampling," Phys. Rev. Lett. 119 (2017) 170501 — the hafnian formula.

Provenance. A verifier should recompute: (1) the HOM $2\times2$ permanent giving zero coincidence amplitude and unit bunching; (2) the second-quantized derivation $\langle S|\varphi(U)|T\rangle=\mathrm{Per}(U_{S,T})$ , checking the occupation-factorial normalization in $P(S)$ ; (3) Ryser's $O(n\,2^n)$ cost against the determinant's $O(n^3)$ , and the worked scaling $n\,2^n$ at $n=30,50,76$ ; (4) the Gaussian-permanent second moment $\mathbb E|\mathrm{Per}(X)|^2=n!$ ; (5) the logical chain Stockmeyer + Toda $\Rightarrow\mathrm{PH}=\Sigma_3^p$ ; (6) $\mathrm{Per}(A)=450$ , $\det(A)=0$ for the worked $3\times3$ example. Cite Aaronson–Arkhipov for the hardness theorems and the PGC/PACC conjectures, Valiant for $\#\mathrm{P}$ -completeness, and Hamilton et al. for the GBS hafnian.

Prerequisites. The Complexity Class BQP · BQP is Contained in PP and PSPACE · BQP versus P, NP, and PSPACE. Leads to. Quantum Query Complexity · The Polynomial Method for Quantum Lower Bounds.

Ask the tutor: Ask me to walk through why erasing one minus sign turns the easy determinant into the $\#\mathrm{P}$ -hard permanent — and how that single fact reaches all the way to Jiuzhang's advantage claim.

Is boson sampling a quantum computer that can solve useful problems?

No, and this is the most important thing to get right. A boson sampler runs no algorithm, decides no yes/no question, and is not known to contain BQP — it cannot factor numbers, search databases, or optimize anything. Its only job is to produce samples from a distribution over photon-detection patterns. What makes it interesting is a negative claim about the classical world: no efficient classical algorithm can produce those same samples, unless the polynomial hierarchy collapses. It is a probe of the quantum–classical boundary, not a device you would point at a practical task. Treating "impossible to simulate classically" as "computationally powerful" is exactly the trap the model exposes.

What is the permanent, and why is it so much harder than the determinant?

The permanent of an $n\times n$ matrix is the sum over all $n!$ permutations of the products of matched entries — the determinant's formula with every minus sign deleted. That deletion is catastrophic for computation. The determinant's alternating signs let row operations reduce it without changing its value, so Gaussian elimination computes it in $O(n^3)$ time. The permanent has no such symmetry: Valiant proved in 1979 that computing it exactly is $\#\mathrm{P}$ -complete, meaning it is as hard as counting the solutions to any NP problem. The best known general algorithm, Ryser's formula, still costs $O(n\,2^n)$ . Physically, fermions antisymmetrize and produce determinants (easy); bosons symmetrize and produce permanents (hard) — the statistics of the particle pick which twin you get.

Why sample instead of just computing one output probability?

Because you cannot efficiently compute one, even with the device. Each output probability is a squared permanent over a huge normalization, so individual probabilities are exponentially small — of order $1/\binom{m}{n}$ . Estimating one to useful precision would need exponentially many runs, defeating the purpose. The machine's power is collective: a single run draws a sample from a distribution defined by exponentially many permanents at once. The hardness proof exploits exactly this, using Stockmeyer's counting theorem to show that even approximating the distribution well enough to sample from it would let a classical machine estimate a $\#\mathrm{P}$ -hard permanent — which it should not be able to do.

What is the difference between the exact and approximate hardness results?

The exact-case theorem is clean and assumption-light: if a classical algorithm samples exactly from the boson-sampling distribution in polynomial time, the polynomial hierarchy collapses to $\Sigma_3^p$ , needing only the belief that the hierarchy is infinite. The catch is that no real device samples exactly — experiments and any plausible classical rival produce a distribution close in total-variation distance, not identical. The approximate-case theorem targets that realistic setting, but it must pay for the generality with two unproven conjectures: the Permanent-of-Gaussians Conjecture (approximating $|\mathrm{Per}(X)|^2$ for Gaussian $X$ is $\#\mathrm{P}$ -hard) and the Permanent Anti-Concentration Conjecture (Gaussian permanents are usually not tiny). Establishing these would put approximate boson sampling on nearly the same footing as the exact result; both remain open.

What is Gaussian boson sampling and why do experiments prefer it?

Gaussian boson sampling replaces the single-photon Fock inputs of the original model with squeezed-vacuum states. This is a huge experimental simplification: deterministic single photons are hard to make on demand, whereas squeezed light is produced continuously and reliably by nonlinear crystals. The squeezed state stays Gaussian as it passes through the passive linear network, and the output click probabilities become hafnians — the symmetric-matrix cousin of the permanent, counting perfect matchings of a weighted graph — rather than permanents. Since the hafnian is also $\#\mathrm{P}$ -hard, the classical-hardness story survives intact. Jiuzhang and Borealis are both GBS machines. A bonus is that hafnians encode graph problems, so GBS devices have been repurposed as heuristic samplers for molecular spectra and dense-subgraph search.

Can classical computers ever catch up to these experiments?

Sometimes, and that is the live controversy. The hardness theorems assume ideal, lossless, perfectly indistinguishable photons; real machines suffer photon loss, partial distinguishability, and finite squeezing. Classical algorithms have learned to exploit exactly these imperfections — tensor-network methods and Markov-chain "spoofers" that simulate only the few photons surviving loss, or that match the low-order correlations an experiment actually certifies. Several such algorithms have matched or surpassed specific advantage claims on their own benchmarks. The upshot is a threshold phenomenon: below some loss and noise level the sampler is provably hard, above it a classical machine wins, and pinning down that boundary — together with proving the Permanent-of-Gaussians Conjecture — is where the theoretical action now lives. Quantum advantage in this arena is a claim that must be continuously defended, not a one-time victory.