Four algorithms extend classical simulation of photonic quantum circuits

Jean Senellart, Stephen Wein and Pierre-Emmanuel Emeriau, at Quandela, present four sampling algorithms in a preprint submitted to arXiv on November 8, 2026. The supplied record does not establish peer review. This is mathematical and computational research, not an experiment with a new photonic device.

The question matters for assessing quantum advantage: a device’s output must be compared with what classical computers can reproduce. Here, sampling means generating random detection patterns with the probabilities predicted by the circuit, rather than calculating every possible outcome.

Building an output one photon at a time

The starting point is Clifford and Clifford’s algorithm for single photons passing through a linear interferometer, a network that mixes optical modes. It constructs an output sequentially while adding input photons in a uniformly random order.

The authors identify a mathematical consistency condition that lets this strategy work more broadly. Their “set-chain lemma” relates non-negative weights for a partially constructed output to weights for smaller input sets. When that identity holds, the algorithm can choose the next output mode using one weight evaluation per candidate mode.

That does not make each evaluation cheap. The method still involves permanent calculations, whose cost grows exponentially in the baseline algorithm. The advance is avoiding additional expensive sums that a more direct calculation would require.

Four extensions, with different costs

Sources in a superposition of vacuum and one photon. These sources are not equivalent to a classical lottery deciding which emitters fire: coherence changes the output probabilities. The algorithm first draws the total photon number, then uses a rectangular version of the permanent to retain that coherence. Its exponential cost scaling matches that of single photons under input loss at the same emission probability, although extra work remains within each step.

For 34 sources with an emission probability of 0.8, the measured runtime was about 0.23 times the lossless 34-photon baseline. Against the more relevant baseline of single photons with matching loss, however, it was about 4.1 times slower. Fewer detected photons explain the reduction against the lossless case; coherence is not computationally free.

Entangled photon trains. A sequential emitter with an internal memory can produce a photonic matrix-product state. The method compresses the contribution of emissions outside the active set into a sum of squared amplitudes, rather than enumerating every omitted label combination.

On the tested trains, extending to 20 emissions, the number of terms fell to one or two near the final, most expensive steps. Whether that favorable behavior persists for longer trains remains open. If the term count grows exponentially, another existing sampling approach could eventually become faster.

Adaptive linear optics. Here, intermediate photon counts determine later optical settings. The algorithm draws those counts, initializes a compatible partial chain and rejects continuations that would put additional photons into already measured modes.

For 38 photons and twice as many optical modes, a cost model estimated runtime factors of 1.03, 1.07 and 1.10 for one, two and three adaptive modes. These are model estimates, not measured end-to-end timings. Separate measurements found approximately 2.1 and 4.1 chains started per accepted sample for one and two adaptive modes; early rejection limits wasted work.

A non-linear phase element. The straightforward chain identity fails here. Instead, the authors express the circuit as a coherent sum of linear circuits, sample from those branches and apply a rejection rule that restores the required interference.

For one Kerr phase element at a phase parameter of π/2, the theoretical acceptance probability is one-half. Combining this with measured branch-evaluation costs gives a modeled overhead of about 4.2 times the single-photon sampler. This is not universal: the response’s structure matters, and adding multiple sites can multiply the expense sharply.

What the checks establish

The authors tested 19 small instances using 40,000 samples each, comparing their output distributions with independently coded probability calculations. All discrepancies were within 1.7 standard deviations of the variation expected from finite sampling. Selected mathematical results were also machine-checked in Lean, but that formalization does not cover every extension or certify the GPU implementation.

“Exact” describes the ideal algorithms in exact arithmetic. The implementation uses floating-point arithmetic, so agreement within sampling noise is supporting evidence, not proof of numerically exact output.

The benchmarks used one NVIDIA H100 and specified random-interferometer geometries. Comparisons with competing methods relied on cost models, not timed implementations. Noise coverage is restricted too: general photon distinguishability, inefficient adaptive detectors and loss between non-linear sites are not covered.

The result broadens classical simulation tools. It does not make arbitrary photonic circuits easy or rule out quantum advantage.