Computing Backends

To run a simulation, Perceval is integrating different computing backends implemented from state of the art algorithms. Each of these backends have different specificities and features that we describe in that section.

Features

Sampling or Weak Simulation

Sampling is the simulation task closest to the actual running of a physical circuit. Given known input states, the sampling will produce output states one at a time as it would be observed by ideal detectors. Sampling is considered as a weak simulation of a circuit, since it does not give explicit output distribution, nor the nature of the mixture states generated by the circuit but a mere observation of individual outputs.

Sampling has been studied in length (see [Clifford and Clifford, 2018], [Clifford and Clifford, 2020]), in the context of Boson Sampling, as a classical computing challenge pushing further the limit of the quantum supremacy.

Strong Simulation

Strong(er) Simulation (see [Heurtel et al., 2022]) of a circuit should provide an access to complete distribution and nature of the output state. Compared to Sampling, circuit designer are interested in the actual probabilistic distribution of the outputs and their exact characteristic.

In particular, for a \(m\) port circuit, one would like to know the exact expected probability of detecting photons on a given port, and not a mere estimation based on sampling observations. We also differentiate here the ability of getting the probability (or probability amplitude) of a single output, and the possibilities of getting probabilities of all the different outputs for a give input state.

Also, even for a simple circuit showing an equiprobable probability of detecting \(\ket{0,1}\) and \(\ket{1,0}\), we would want to know if the output is \(\frac{1}{\sqrt 2}(\ket{0,1}+\ket{1,0})\) or \(\frac{1}{\sqrt 2}(\ket{0,1}-\ket{1,0})\) which are very distinct states.

Finally, a fine-grained simulation would need not only to give output state probability but also probability amplitude. Indeed, probability amplitude is required for further evolution of the output states, but also analysis of polarization for circuit with polarization support, etc.

Strongest Simulation

Beyond simulation of perfect circuit describes by unitary matrix, goal of Perceval is also to model non unitary phenomenon like loss of photons, noise, time delays, and more. Ideal simulators should take these phenomenon into accounts.

The Backends

Perceval has 6 different built-in back-ends with the support of optimized C++ library.

Comparison Table

Features Name

CliffordClifford2017

SLOS

Naive

NaiveApprox

Stepper

MPS

Sampling Efficiency

\(\mathrm{O}(n2^n+poly(m,n))\)

\(\mathrm{O}(mC_n^{n+m-1})\)

N/A [1]

N/A [1]

N/A [1]

N/A [1]

Single output Efficiency

N/A

N/A

\(\mathrm{O}(n2^n)\)

\(\mathrm{O}(n)\)

\(\mathrm{o}(N_cC_n^{n+m-1})\)

\(\mathrm{o}(N_cC_n^{n+m-1})\)

Full Distribution Efficiency

N/A

\(\mathrm{O}(nC_n^{n+m-1})\)

\(\mathrm{O}(n2^nC_n^{n+m-1})\)

\(\mathrm{O}(nC_n^{n+m-1})\)

\(\mathrm{o}(N_cC_n^{n+m-1})\)

\(\mathrm{o}(N_cC_n^{n+m-1})\)

Probability Amplitude

No

Yes

Yes

Yes

Yes

Yes

Support Symbolic Computation

No

Yes

No

No

Yes

No

Support of Time-Circuit

No

No

No

No

Yes

No

Practical Limits

\(n\approx30\)

\(n,m<20\)

\(n\approx30\)

where:

  • \(n\) is the number of photons

  • \(m\) is the number of modes

  • \(N_c\) is the number of elementary circuits

CliffordClifford2017

This backend is the implementation of the algorithm introduced in [Clifford and Clifford, 2018]. The algorithm, applied to Boson Sampling, aims to produce provably correct random samples from a particular quantum mechanical distribution. Its time and space complexity are respectively \(\mathrm{n2^n+mn^2}\) and \(\mathrm{m}\) (in addition to matrix storing). The algorithm has been implemented in C++, and uses an adapted Glynn algorithm [Glynn, 2010] to efficiently compute \(n\) simultaneous sub-permanents.

Recently, the same authors have proposed a faster algorithm in [Clifford and Clifford, 2020] with an average time complexity of \(\mathrm{n\rho_\theta^n}\) for a number of modes \(m=\theta n\) which is linear in the number of photons \(n\), where:

\[\rho_\theta = \frac{(2\theta+1)^{2\theta+1}}{(4\theta)^{\theta}(\theta+1)^{\theta+1}}\]

For example, if we were to work with dual rail path encoding (ignoring for now the number of auxiliary modes required), we would typically work with \(\theta=2\), and the average performance is then \(\mathrm{n(\frac{5^5}{8^23^3})^n} \approx \mathrm{n1.8^n}\).

SLOS

The Strong Linear Optical Simulation SLOS algorithm developed by a subset of the present authors is introduced in [Heurtel et al., 2022]. It unfolds the full computation path in memory, leading to a remarkable time complexity of \(\mathrm{nC_n^{n+m-1}}\) for computing the full distribution. The current implementation also allows restrictive sets of outputs, with average computing time in \(\mathrm{n\rho_\theta^n}\) for single output computation. As discussed in [Heurtel et al., 2022], Boson Sampling with SLOS is possible with the time complexity of [Clifford and Clifford, 2020], though it has not yet been implemented in the current version of Perceval.

The tradeoff in this approach is a huge memory usage in \(\mathrm{nC^{n+m-1}_n}\) that limits usage on personal computers to circuits with \(\approx 20\) photons and to \(\approx 24\) photons on super-computers.

Naive

This backend implements direct permanent calculation and is therefore suited for single output probability computation with small memory cost. Both Ryser’s [Ryser, 1963] and Glynn’s [Glynn, 2010] algorithms have been implemented. Extra-care has been taken on the implementation of these algorithms, with usage of different optimisation techniques including native multithreading and SIMD vectorisation primitives. Benchmarking of these algorithms and comparison with the implementation present in the The Walrus library is provided in following figure:

_images/performance-permanent.png

Comparison of the average time [2] to calculate a permanent of an \(n\times n\) Haar random matrix. The processor is a 32 core, 3.1GHz Intel Haswell. For The Walrus, version 0.19 is used and installed from pypi. The Ryser implementation is run on 4 or 32 threads. The Glynn implementation is run on a single thread. What is interesting to note is that all implementations have convergence to the theoretical performance but the factor between optimised and less optimised implementation still makes a perceptible time difference for the end-user.

NaiveApprox

This backend does the same computations as Naive, but uses Gurvits estimate to compute the permanent. Aside of usual probability() and prob_amplitude() methods, it offers a 99% confidence interval on the probability, or a 99% sure error bound on the amplitude. A better accuracy can be obtained with a higher iteration count.

With this approximated backend, you can achieve a few probability estimates for high photon counts. Example given:

>>> from perceval.utils import BasicState, Matrix
>>> from perceval.backends import NaiveApproxBackend
>>> from perceval.components import Unitary
>>>
>>> circuit_size = 60
>>> n_photons = 30
>>> backend = NaiveApproxBackend(100_000_000)
>>> backend.set_circuit(Unitary(Matrix.random_unitary(size)))
>>> input_state = BasicState([1]*n_photons + [0]*(size-n_photons))
>>> backend.set_input_state(input_state)
>>> interval = backend.probability_confidence_interval(BasicState([1]*n_photons + [0]*(size-n_photons)))
>>> print(f"Probability in {interval}")
Probability in [6.051670221391749e-20, 1.5297683283662674e-19]

MPS

Matrix Product State (MPS) is based on a type of tensor network simulation, which gives an approximation of the output states [Schollwöck, 2011], [Oh et al., 2021]. As the Stepper, MPS backend does the computation on each component of the circuits one-by-one, and not on the whole unitary, but has the unique feature of performing approximate state evolution. The states are represented by tensors, which are then updated at each component. These tensors can be seen as a big set of matrices, and the approximation is done by choosing the dimension of these matrices, called the bond dimension.

Stepper

This simulator takes a totally different approach. Without computing the circuit’s overall unitary matrix first, it applies the unitary matrix associated with the components in each layer of the circuit one-by-one, simulating the evolution of the statevector. The complexity of this backend is therefore proportional to the number of components. It enables simple debugging of circuits by exposing intermediate states.

The Stepper formerly a backend is a simulator since Perceval 0.9.

Footnotes