Introduction to GBS

Gaussian Boson Sampling (GBS) is a special-purpose model of photonic quantum computation, first introduced in Ref. [13]. In its most general form, GBS consists of preparing a multi-mode Gaussian state and measuring it in the Fock basis. It differs from universal photonic circuits only because it does not employ non-Gaussian gates and it restricts measurements to the Fock basis. In general, the output distribution of a GBS device cannot be simulated in polynomial time with classical computers [2][13]. Applications of GBS aim to harness the computational power unique to GBS to perform useful tasks.

GBS offers great versatility in the scope of problems that it can encode. This led to the appearance of several GBS algorithms with applications to quantum chemistry [45], graph optimization [44][49], molecular docking [40], graph similarity [50] [51][52], and point processes [48]. The GBS applications layer is built with the goal of providing users with the capability to implement these GBS algorithms using only a few lines of code. Programming the GBS device, generating samples, and classical post-processing of the outputs are taken care of automatically by built-in functions.

The GBS distribution

A general pure Gaussian state can be prepared from a vacuum state by a sequence of single-mode squeezing, multimode linear interferometry, and single-mode displacements. It was shown in Ref. [13] that for a Gaussian state with zero mean – which can be prepared using only squeezing followed by linear interferometry – the probability \(\Pr(S)\) of observing an output \(S=(s_1, s_2, \ldots, s_m)\), where \(s_i\) denotes the number of photons detected in the \(i\)-th mode, is given by

\[\Pr(S) = \frac{1}{\sqrt{\text{det}(\mathbf{Q})}}\frac{\text{Haf}(\mathbf{\mathcal{A}) }}{s_1!s_2!\cdots s_m!},\]


\[\begin{split}\begin{align} \mathbf{Q}&:=\mathbf{\Sigma} +I/2\\ \mathbf{\mathcal{A}} &:= \mathbf{X} \left(I- \mathbf{Q}^{-1}\right),\\ \mathbf{X} &:= \left[\begin{smallmatrix} 0 & I \\ I & 0 \end{smallmatrix} \right], \end{align}\end{split}\]

and \(\mathbf{\Sigma}\) is the covariance matrix of the Gaussian state. The matrix function \(\text{Haf}(\cdot)\) is the hafnian. When the state is pure, the matrix \(\mathbf{\mathcal{A}}\) can be written as \(\mathbf{\mathcal{A}}=\mathbf{A}\oplus \mathbf{A}^*\) and the distribution becomes

\[\Pr(S) = \frac{1}{\sqrt{\text{det}(\mathbf{Q})} } \frac{|\text{Haf}(\mathbf{A}_S)|^2}{ s_1!\ldots s_m!},\]

where \(\mathbf{A}\) is an arbitrary symmetric matrix with eigenvalues bounded between \(-1\) and \(1\). Therefore, besides encoding Gaussian states into a GBS device, it is also possible to encode symmetric matrices \(\mathbf{A}\). Notably, these may include adjacency matrices of graphs, and kernel matrices.

Programming a GBS device

In GBS without displacements, indicating gate parameters is equivalent to specifying the symmetric matrix \(\mathbf{A}\). Employing the Takagi-Autonne decomposition, we can write

\[\mathbf{A} = \mathbf{U}\, \text{diag}(\lambda_1, \lambda_2, \ldots, \lambda_m)\, \mathbf{U}^T,\]

where \(\mathbf{U}\) is a unitary matrix. The matrix \(\mathbf{U}\) is precisely the unitary operation that specifies the linear interferometer of a GBS device. The values \(0\leq \lambda_i < 1\) uniquely determine the squeezing parameters \(r_i\) via the relation \(\tanh (r_i) = \lambda_i\), as well as the mean photon number \(\bar{n}\) of the distribution from the expression

\[\bar{n} = \sum_{i=1}^M \frac{\lambda_i^2}{1-\lambda_i^2}.\]

It is possible to encode an arbitrary symmetric matrix \(\mathbf{A}\) into a GBS device by rescaling the matrix with a parameter \(c>0\) so that \(cA\) satisfies \(0\leq \lambda_i < 1\) as in the above decomposition. The parameter \(c\) controls the squeezing parameters \(r_i\) and the mean photon number. Overall, a GBS device can be programmed as follows:

  1. Compute the Takagi-Autonne decomposition of \(\mathbf{A}\) to determine the unitary \(\mathbf{U}\) and the values \(\lambda_1, \lambda_2, \ldots, \lambda_m\).

  2. Program the linear interferometer according to the unitary \(\mathbf{U}\).

  3. Solve for the constant \(c>0\) such that \(\bar{n} = \sum_{i=1}^M \frac{(c\lambda_i)^2}{1-(c\lambda_i)^2}\).

  4. Program the squeezing parameter \(r_i\) of the squeezing gate \(S(r_i)\) acting on the \(i\)-th mode as \(r_i=\tanh^{-1}(c\lambda_i)\).

The GBS device then samples from the distribution

\[\Pr(S) \propto c^{k} \frac{|\text{Haf}(\mathbf{A}_S)|^2}{s_1!\ldots s_m!},\]

with \(k = \sum_{i}s_{i}\).

The GBS applications layer includes functions for sampling from GBS devices that are programmed in this manner. It also includes a function for sampling more general Gaussian states, which are useful for applications to quantum chemistry. These functions can be found in the sample module.

GBS algorithms work by choosing a clever way of encoding problems into a GBS device and generating many samples, which are then be post-processed by classical techniques. The applications layer contains a dedicated module for each of the known GBS algorithms.