Variational Quantum Simulation of Ultrastrong Light-Matter Coupling

We propose the simulation of quantum-optical systems in the ultrastrong-coupling regime using a variational quantum algorithm. More precisely, we introduce a short-depth variational form to prepare the groundstate of the multimode Dicke model on a quantum processor and present proof-of-principle results obtained via cloud access to an IBM device. We moreover provide an algorithm for characterizing the groundstate by Wigner state tomography. Our work is a first step towards digital quantum simulation of quantum-optical systems with potential applications to the spin-boson, Kondo and Jahn-Teller models.

We propose the simulation of quantum-optical systems in the ultrastrong-coupling regime using a variational quantum algorithm.More precisely, we introduce a short-depth variational form to prepare the groundstate of the multimode Dicke model on a quantum processor and present proofof-principle results obtained via cloud access to an IBM device.We moreover provide an algorithm for characterizing the groundstate by Wigner state tomography.Our work is a first step towards digital quantum simulation of quantum-optical systems with potential applications to the spinboson, Kondo and Jahn-Teller models.
Quantum simulation is one of the most prominent applications of quantum processors for solving problems in quantum physics and chemistry.Importantly, quantum simulation aims to circumvent the limited capabilities of classical computers to represent quantum states in exponentially large Hilbert spaces.Recently, a hybrid, quantum-classical simulation paradigm exploiting quantum variational principles has been introduced [1].Following this pioneering work, many other realizations of what is known as Variational Quantum Algorithm (VQA) have appeared in the literature [2][3][4][5].
VQAs have been shown to have some robustness against noise and thus appear appropriate for the current generation of Noise-Intermediate-Scale-Quantum (NISQ) processors [4,6,7].Although considerable effort has been devoted to solving proof-of-principle instances of problems in quantum chemistry [2][3][4][5] and optimization [8], the general applicability of this approach to other domains in physics is still a subject of debate and interest [9,10].Here, we use a VQA to simulate strongly interacting light-matter models.In particular, we focus on obtaining the groundstate of a set of two-level atoms coupled to electromagnetic modes, which is of fundamental interest and has practical applications for example for quantum-information processing and sensing [11][12][13][14][15].
The simplest case corresponds to that of a two-level atom coupled to a cavity mode and is described by the quantum Rabi Hamiltonian Here, ω q and ω c are the atomic and the electromagneticmode frequencies, σ µ (µ = x, y, z) the Pauli matrices and a (a † ) the annihilation (creation) operator for the oscillator, respectively.If the light-matter coupling constant, g, is small compared to the systems' frequencies, Eq. ( 1) reduces to the Jaynes-Cummings Hamiltonian [16].Under these conditions, the terms σ + a and σ − a † , where σ ± = (σ x ±iσ y )/2, lead to an exchange of a single excitation between the atom and the oscillator mode.Provided that g is greater than the decoherence rates of the atom and the cavity, this regime of light-matter interaction is referred to as strong coupling, and it is widely exploited for quantum-information processing purposes [17].
As g approaches a significant fraction of the bare atom and cavity frequencies, or becomes the largest energy scale in Eq. ( 1), the atom-cavity system enters the ultrastrong-(USC) and deep-strong-coupling (DSC) regimes, respectively [11,14,15,18,19].In these cases, the presence of the counter-rotating terms (σ + a † and σ − a) in Eq. ( 1) needs to be taken into account.Perturbation theory provides an accurate description for coupling strengths in the range of 10% − 30% of the system's frequencies, but has limited applicability beyond that regime [14].While an exact analytical solution in principle exists for Eq.(1) [20], larger systems involving multiple atoms and/or electromagnetic modes can only be handled numerically.
In the large-g limit, however, the mean cavity-mode occupation number and its quantum fluctuations are large and a sizable Fock space is required for numerical simulations.The total Hilbert-space dimension can thus quickly become unpractical for many-particle systems.This fact motivates the search for powerful analytical and numerical methods [11][12][13][21][22][23][24] and quantum-simulation algorithms [14,15,19,[25][26][27] for this problem.
We consider the generalization of Eq. ( 1) to N atoms and M electromagnetic modes, given by (2) where the constants {g ik } quantify the coupling strength between the i th atom (of frequency ω qi ) and the k th cavity mode (of frequency ω k ) referred below to as kmode.For M = 1, Eq. ( 2) reduces to the Dicke model, while the special case N = 1 corresponds to the multimode quantum Rabi model.Digital quantum simulation of such models requires the encoding of the bosonic modes into qubit registers.We choose to use a Single-

arXiv:1909.08640v1 [quant-ph] 18 Sep 2019
Excitation-Subspace (SES) encoding, in which the Fock space of a given k-mode is truncated to a maximum photon number n max k , and represented by a qubit register of size n max k + 1 [28][29][30][31].A mapping from the k-mode Fock space to the single-excitation subspace of the qubit register is then defined as where the tilde is used hereafter to indicate encoded states and operators.Importantly, under SES encoding, quadratic bosonic Hamiltonians lead to next-neighbor interactions at most [28].Indeed, the k-mode annihilation operator maps to a k → ãk = , where σ ± n k acts on the n k th qubit of the k-mode register.The 2local form of ãk relaxes connectivity requirements on the k-mode qubit register and thus leads to a reduced gate count.Other encodings can be found in Refs.[27,32,33].
Finding the groundstate |G of Eq. ( 2) by means of a VQA requires first to construct a proper variational form [1,6].That is, a unitary U (θ) parametrized by a real-valued vector θ, such that where noninteracting vacuum state, and θ * is obtained by classical minimization of the energy E(θ) = vac|U † (θ) HU (θ)| vac .Some intuition about a convenient choice of U (θ) can be gained from approximate disentangling transformations for Eq. ( 2) [22,23].We refer to such transformations indistinctly as polaron Ansätze.The simplest transformation is obtained for the case of N = 1, where it is useful to rotate H → H = P † HP by means of a qubitstate-dependent displacement of the k-modes where ω q is a renormalized frequency for the atom.As illustrated in Sect.IA of the Supplemental Material, the groundstate of H approaches the noninteracting groundstate of the atom-cavities system, |vac , in most coupling regimes.Therefore, the state P |vac approximates the groundstate |G in the laboratory frame.Exploiting this fact to prepare | G on a quantum computer requires compiling P from single-and two-qubit gates, for instance, using a Trotter decomposition.The need for reducing the Trotter error, however, can lead to quantum circuits of large depth.Moreover, this approach is sensitive to errors arising from imperfect qubit control and noise.As a way around this problem, we propose to leverage the structure of the polaron transformation to obtain a short-depth variational from.We do this by parameterizing the Trotter decomposition of P and letting the variational algorithm adjust the unitary such that the groundstate-Ansatz energy is minimized.The variational form has not only the purpose of discovering short-depth quantum circuits for synthesizing the USC groundstate, but also to potentially improving on the disentangling capabilities of Eq. ( 4).
We construct the variational form by choosing a convenient Trotter decomposition of P , first for the case of N = 1.We introduce two k-mode operators, X e k and X o k , which are defined such that P = , where {f k = g k /(ω k + ω q )} is a set of constants that will latter play the role of variational parameters.Although [ X e k , X o k ] = 0, X e k and X o k are respectively composed of commuting terms that act on even and odd sites of the k-mode qubit register (see the Supplemental Material, Sect.IB).The 2-local form of the encoded bosonic operators leads to an efficient implementation of the Trotter-expanded unitary where d k is the number of Trotter steps, that may vary with the k-mode index.As shown in Sect.IB of the Supplemental Material, the exponentials in this equation factorize exactly into a product of n max k controlledexchange gates acting on next-neighbor qubits of the k-mode register with the atom register being the control qubit.The implementation of Eq. ( 5) requires thus n max k × d k such gates per k-mode, adding to a total gate count of M k=1 n max k d k before quantum-circuit compilation.This number grows linearly with the number of kmodes, their Fock-space dimension and the order of the Trotter expansion (Trotter depth).Interestingly, since Eq. ( 5) parallelizes over the k-modes, its quantum-circuit depth does not scale with M .
For N > 1, the resulting variational form incorporates blocks of the form of Eq. (5) where the two-level-atom operator σ x → σ x i is now labeled by i ∈ [1, N ] and alternated among the respective qubit registers (see the Supplemental Material, Sect.IC).This observation leads to the more general expression where the coefficients f k → f s ik are variational parameters that depend on the Trotter step s ∈ [1, . . ., d ik ].Additionally, f s ik can also be made a function of the kmode photon number, such that f s ik → f s ik (n k ).As argued below, this trades shorter circuit depths for longer optimization runtime.
Important additional details apply, however, between the cases of N = 1 and N > 1.In particular, the case N > 1 requires Eq. ( 6) to be complemented by singlelayer short-depth variational form that acts on the atoms' registers.This extra step initializes the polaron variational circuit to the state | vac = converging to a 4% in the limit of large n max k (pink baseline) for g/ω = 0.8.This indicates that the minimum value of ∆en is reduced exponentially with n max k , reaching an absolute lower bound determined by the entanglement capabilities of the polaron Ansatz.Simulations do not include noise and are done using Qiskit [34].
is determined by an auxiliary optimization loop specified in Sect.IC of the Supplemental Material.
Fig. 1 shows the results for the (a) single-and (b) twomode Rabi Hamiltonian, and (c) single-and (d) twomode Dicke model for N = 2.The simulations assume the resonant case where atom and k-mode frequencies are set to ω qi = ω k ≡ ω and g ik ≡ g is swept in [0, ω].The resonance condition leads to strong entanglement between the atoms and cavity modes due to the energetically favorable exchange of excitations.To quantify the performance of the variational form, we define the error metric ∆ en = |(E vqe −E en )/E en |, accounting for the relative difference between the groundstate energy found by the VQA, E vqe , and the energy of the encoded groundstate, E en .An additional metric ∆ ex = |(E en −E ex )/E ex | quantifies the difference between E en and the numerically exact groundstate energy.We evaluate ∆ en and ∆ ex as a function of g/ω for circuits with Trotter depth d ik = d.The chosen Fock space truncation (see figure caption) leads to a small number of qubits while ensuring a relatively small ∆ ex .This choice seeks to reduce the quantum-hardware resources needed for simulation.
Results in panel (a) show that relative errors ∆ en below 1% are achieved by state-preparation circuits containing only 3 variational parameters (d = 3).A similar accuracy is obtained for circuits with d = 2 if additional pa-rameters dependent on the k-mode photon-number are incorporated (not shown).Remarkably, for d > 1, the energy of the variational Ansatz is significantly lower than H vac | n max k .The latter is the expectation value of Eq. ( 2) on the state P |vac | n max k within a truncated Fock space.This indicates that the variational algorithm can leverage the Trotter error to outperform the full polaron Ansatz under the same Fock-space restrictions and with very low circuit depth.Interestingly, we also find that the energy of the variational state falls below H vac | n max k →∞ in the full range of g/ω ∈ [0, 1] (not shown).The error metric ∆ ex remains below ∼ 2% in all the cases.
We observe a similar qualitative behavior for the twomode simulations in panel (b), although ∆ en increases to ∼ 2.5% for d = 4.The same accuracy is reached for circuits with d = 2 when variational parameters for each k-mode photon number are introduced (not shown).We find that the accuracy limit is both due to finite Fockspace truncation errors and the disentangling capabilities of the polaron Ansatz.Increasing the number of two-level atoms in the model, while keeping the number of qubits of the order of 10, leads to the results in panels (c)-(d) for which we find a maximum error of ∆ en 5% for d = 5 in the first case, and of ∆ en 8% for d = 4 in the second case.These results, however, are limited by Fock space truncation errors and can be improved by increasing the number of qubits in the simulations.It is worth noticing that, similarly to the case of N = 1, these variational circuits outperform the polaron Ansatz significantly for the same conditions.
The performance of the variational form may be improved further by means of simple modifications.For instance, a layer of a hardware-efficient (HE) gates [3] could be appended after each Trotter step, providing greater entangling capabilities for state preparation.Ideally, gates on such HE layers should conserve the number of excitations in the k-mode registers [30].Generalizations of Eq. ( 4) incorporating additional parameters are also a possibility [35].
As the number of qubits scales with ∼ (n max + 1) M , simulating the performance of the proposed VQA on a classical computer becomes quickly expensive.Moreover, circuits of larger depth and number of qubits could likely benefit from quantum devices tailored to compile the polaron Ansatz in fewer gates.An option is to engineer the required controlled-two-qubit gates directly on the quantum hardware.Sect.V of the Supplemental Material illustrates such special-purpose devices in the context of circuit QED.
The results of Fig. 1 suggest that the polaron variational form is a promising tool for investigating the USC groundstate in near-term quantum devices.For this reason, we implement the aforementioned strategy in currently available quantum hardware.Here, we use the IBM Q Poughkeepsie chip via the open-source framework Qiskit, taking advantage of the built-in SPSA opti-mizer [3,36] and the readout error mitigation techniques of Qiskit-Ignis [34].We use three qubits for the quantum simulation, two of them encoding the bosonic mode.The groundstate energies found this way, shown in Fig. 2 (star-shaped data points), are in good qualitative agreement with the theoretical estimations.
We find that the main limitations on the accuracy of the VQA are due to the level of noise in the quantum processor and to the capabilities of the SPSA optimizer given a finite number of optimization steps.To investigate the effect of the latter against the former, we perform the VQA with a desktop computer, assuming a larger number of optimization steps and the calibrated noise model of the quantum hardware.This produces a set of variational states with optimal parametrization according to the classical simulation.We then evaluate the energy expectation value of such states on the quantum processor, performing mitigation of readout errors.The result of this experiment (triangular-shaped data points) reach better accuracies than those obtained by means of the hybrid quantum-classical VQA.This suggests that noise processes on the quantum hardware prevent highaccuracy solutions to be reached in a reasonable number of optimization steps via cloud access, in the order of 150 SPSA trials.By controlling the level of noise in classical simulation, we also find that hybrid quantum-classical VQA solutions with ∆ en ∼ 1−2% for 150 SPSA trials are expected for noise levels one order of magnitude smaller than the present value.Note that in absence of noise, the number of optimizer steps required to reach numerical accuracy with respect to the reference value is very small in comparison, below 30 in the entire g/ω ∈ [0, 1] range.This allow us to conclude that the discrepancies encountered in the quantum-hardware runs are due the effect of noise and the limited optimizer calls rather than limitations of the proposed Ansatz.
Following this proof-of-principle demonstration, we present an alternative method for characterizing the prepared groundstate.This technique could be useful to probe entanglement metrics and to distinguish between nearly degenerate states.The latter situation occurs, for instance, within the groundstate manifold of the quantum Rabi model approaching the DSC regime.To this end, we introduce the joint Wigner function for a set of N qubits and M bosonic modes as generalizing the definition given in Ref. [37] for the case of N = M = 1.Here, {σ li i , l i ∈ [0, x, y, z]} are the Pauli matrices for the i th atom with Inversion of Eq. (7) gives the system's density matrix as Variational quantum-optics simulation of the Rabi model in resonance conditions on a quantum processor.Shown is the groundstate energy as a function of the coupling strength, both in units of ω.The cavity mode is encoded in a two-qubit register (n max k = 1).The light-blue bands enclose the range of results that are expected for 150 SPSA trials with levels of noise in the order of 0.1 and 1.0, relative to calibrated values (dotted lines).The black-dashed line is the encodedgroundstate energy.The star-shaped markers are the result of VQA runs for up to 150 SPSA trials on the quantum device.The pointing-up (pointing-down) triangular markers are the minimum (average) expectation value on quantum hardware of states that have been entirely optimized in the classical processor.The dispersion of such values is due to fluctuations in the level of noise of the quantum device between runs.Further details are provided in Sect.II of the Supplemental Material.
tegral is performed over d 2 α = M k=1 d 2 α k and the sum is extended to the 4 N possible values of l = (l 1 , . . ., l N ).This relation can be used for state reconstruction [37].
Expanding Eq. ( 7) in the Fock-state basis within the SES encoding we arrive to where Tr q is the trace operator over the atom registers and Ω ñ = ñ1 . . .ñM | D † (α)ρ D(α)|ñ 1 . . .ñM .W l (α) can be sampled by executing a quantum circuit that performs the necessary state-tomography gates.While these gates are simply single-qubit rotations for the atom registers, tomography gates correspond to the application of D(α) for the k-mode registers.Fortunately, the displacement operators can be easily implemented by a sequence of one-and two-qubit gates derived from a Trotter decomposition similar to that of the polaron transformation.Sect.III of the Supplemental Material includes further details.
We demonstrate this approach numerically for the case of a single atom and a cavity mode in resonance with g/ω = 1.Fig. 3 shows the reconstructed joint Wigner function W σ z (α) for an 8-qubit k-mode register.The result is compared to the numerically exact distribution, which is not affected by Trotter or Fock space truncation Trotter steps for D(α) errors.We observe a good qualitative agreement between the two distributions even for a few as 2 Trotter steps.This agreement improves as the number of Trotter steps used to implement D(α) is increased, although the discrepancy between the two distributions remains bounded from below due to finite-dimensional encoding errors.Finally, we discuss briefly the effect of common quantum error channels on the performance of the proposed VQA.It is worth highlighting that a SES encoding allows for damping errors of the form |ñ k → | 0 k to be detectable by joint-parity measurements of the k-mode register.This enables postselection of uncorrupted states, which can significantly reduce the impact of noise on the proposed variational algorithm.However, the downsides of using a SES encoding reside in two main points.First, this encoding trades shorter quantum-circuit depths for a relatively large qubit overhead compared to other possible encodings [32,33].On the flip side, this compromise might be leveraged by current quantum processors, which are mostly limited by decoherence rather than by the number of qubits [38].Second, noise channels that do not conserve the number of excitations in the k-mode qubit registers can become dominant for large qubit arrays (see Supplemental Material).This is a direct consequence of an exponential growth of the size of the comple-ment of SES with the number of qubits in the simulation.We also note that the observations above are generic to other proposals using SES encodings [28][29][30][31].
In conclusion, we introduced a short-depth and fewparameter variational form to study the interacting lightmatter groundstate of N atoms and M electromagnetic modes.We found that such a variational circuit can approximate the ultrastrong-coupling groundstate with very good accuracy.We implemented a proof-of-principle example on an IBM quantum processor, performing the mitigation of readout errors.Finally, we demonstrated the use of Wigner state-tomography to characterize the groundstate, and discussed the impact of noise on the variational algorithm.As the light-matter interaction Hamiltonian considered in this work is formally identical to the few-impurity spin-boson model, we envision applications to problems in condensed-matter physics for which the polaron transformation was originally introduced.The demonstration of quantum advantage by the variational approach introduced here is likely to require quantum-hardware of size and noise-resilience significantly beyond what is currently available.Our work is, however, a first step towards digital simulation of strongly interacting light-matter models with a quantum processor.In this section, we study the effect of the polaron transformation [Eq.( 4) of the main text] on the multimode Rabi Hamiltonian [N = 1 in Eq. ( 2) of the main text].We first consider a slightly more general unitary of the form where {f k } are real parameters to be determined below.By transforming H → H f = P † f HP f , where f = (f 1 , . . ., f M ), we arrive at Here, ω q = ω q e −2f •f is a renormalized frequency for the two-level atom, g k = g k − ω k f k are new light-matter coupling parameters and We have moreover defined the operator q f = exp(2σ x M k=1 f k a k ).The parameters f k in Eq. ( 2) can reduce the effective light-matter coupling strength g k and thus make the groundstate of such a Hamiltonian closer to that of the noninteracting case.This, however, holds if the higher-order corrections to σ z from the operator q † −f q f , which mixes the k-modes, can be made small at the same time.This leads to a compromise for the best value of f that can be resolved by optimization of such a parameter [1,2].More precisely, assuming that the groundstate of Eq. ( 2) results close to vacuum, minimizing the groundstate energy E f − ω q /2 leads to the optimal condition f k = g k /(ω k + ω q ) and to the implicit equation for the renormalized frequency of the atom.We note that the difficulty of solving the scalar equation Eq. ( 3) does not scale with the number of bosonic modes under consideration.This fact deserves to be highlighted as the optimal value of f is used to initialize the optimizer before executing the variational quantum algorithm presented in the main text and in more details in Sect.I B.
As discussed in Sect.I C, this is no longer true in the general case of N > 1, where the complexity of initialization scales exponentially with N .We can gain an understanding of how the polaron transformation works by analyzing the behavior of approximate solutions to Eq. ( 3) at the boundaries of the range 1 ≥ ω q /ω q ≥ 0. Indeed, for small coupling strengths g k , we expect ω q to differ only slightly from ω q and thus to be able to approximate ω q /ω q 1 − , where 1.In contrast, for large coupling strengths compared to the system frequencies, we expect ω q to vanish exponentially and thus ω q /ω q to hold.In the latter situation, the parameters of the polaron transformation approach the asymptotic scaling f k g k /ω k leading to g k → 0 (along with ω q → 0).This analysis indicates that the light-matter system effectively decouples in the polaron frame [see Eq. ( 2)], both in the strong and deep-strong coupling regimes.As demonstrated below, the disentangling capabilities of this transformation and the solution of Eq. (3) interpolate smoothly in the intermediate ultrastrong coupling regime, making this tool suitable for investigating the light-matter groundstate in a very broad range of parameters.
We confirm the intuition developed in the above paragraph by analyzing both the single-and two-mode quantum Rabi Hamiltonians in resonance conditions ω q = ω k ≡ ω, for g k ≡ g ∈ [0, ω].This is also the regime of parameters considered in the main text.Here, however, we are also interested in quantifying how well the state P |vac approximates the groundstate of the multimode Rabi Hamiltonian in all coupling regimes.Fig. 1 (a) compares the expectation value H vac to the exact groundstate energy of the single-mode (left) and twomode configurations (right).We observe that H vac follows closely the exact solution in the range g/ω ∈ [0, 2].This agreement is even improved for larger coupling strengths.Panels (b)-(d) show the effective atom frequency ω q and coupling strengths g k along with the f k parameters of the polaron transformation.Solid lines correspond to the numerical optimization of ω q , while dashed lines are obtained analytically from a series expansion of the right-hand side of Eq. (3) to second order (3) at the boundaries of 1 ≥ ω q /ωq ≥ 0. The asymptotic regime of f k g k /ω k is indicated in panel (c) (gray dashed lines).A similar line style in (d) indicates the coupling in the bare frame, g k .We observe that g k becomes exponentially small compared to g k as the bare coupling strength increases.This fact highlights the advantages of working in the polaron frame.
in the small parameter for both ω q /ω q 1 − and ω q /ω q regimes.As anticipated, we observe an exponential reduction of the ratio ω q /ω q along with the asymptotic tendency f k g k /ω k and g k → 0 as g increases.Although the energy of the state P |vac follows closely that of the groundstate in the laboratory frame in the DSC regime, we note that f k grows linearly with g k in such conditions.This unfavorable scaling would make the Trotter-decomposition-based approach employed in the main text rather inefficient, requiring quantum-circuits of larger depth to reach convergence.
The advantages of investigating the groundstate properties in the polaron frame are also emphasized by notic-ing that Eq. ( 2) can be recast into assuming g k /(ω k + ω q ) 1 and expanding q f to first order.Eq. ( 4) is reminiscent of the multimode Jaynes-Cummings model, with an additional term ∝ σ z a † k a k , that can either shift the frequency of the atom as a function of the number of photons in the k-modes (k = k ) or allow for the exchange of an excitation between two k-modes through the atom (k = k ).Provided that the effect of the latter interaction is only perturbative, the groundstate of Eq. ( 4) is the vacuum state.As shown in the main text, these conditions also lead to relatively small f k parameters, making it possible to construct a short-depth variational form by Trotter decomposition of the polaron transformation.

B. Quantum-circuit compilation of the polaron variational form
We now provide further details on the quantum circuit compilation of the polaron variational form designed for the multimode quantum Rabi model.Specifically, we seek to rewrite the unitary as a sequence of single-and two-qubit gates available on IBM's on-line quantum platform.The operator X e k in Eq. ( 5), which was introduced in the main text, is defined as X o k is defined analogously, with n k running over odds numbers in [0, n max k − 1].Since the Pauli products in Eq. ( 6) commute with each other, unitaries of the form exp(f s k σ x X e,o k /d k ) factorize into a sequence of controlled two-qubit gates that can be parallelized over the k-mode qubit registers.As shown in Fig. 2 (a), these two-qubit gates operate on pairs of next-neighbor qubits of the kmode registers and are controlled by the atom qubit.Panel (b) provides the compilation of the controlled twoqubit gates in elementary single-and two-qubit gates [3].
The atom and the k-mode registers are initialized to the states |0 and | 0 k = |10 . . .0 , respectively, which correspond to the noninteracting groundstates.This step if one wishes to introduce an extra degree of freedom, f s k (n k ), per controlled two-qubit gate.Importantly, the parameters of the variational form are initialized to the value specified by the polaron transformation in Eq. ( 4) of the main text.
is followed by d k sets of gates, where d k is the Trotter depth used to decompose the polaron unitary for mode k.In its simplest form, the polaron variational form has a single parameter per Trotter step.For a given s, the two-qubit-gate rotation angle is determined by such a parameter and the site index n k which introduces an additional scaling factor ∝ √ n k to this angle.The latter factor arises from the √ n k scaling of the matrix elements of the harmonic-oscillator ladder operators (a k and a † k ) in the photon-number basis.A version of the variational form with a larger number of parameters can easily be crafted by letting the two-qubit-gate rotation angle vary and potentially depart from the √ n k scaling.Moreover, an additional parameter can be introduced if the controlled two-qubit gates are implemented as in Fig. 2 (b), as this gate compilation uses two R z rotations that could be made independent from each other.In all cases, the initial value of the variational-form parameters must be set to that of the polaron transformation in Eq. ( 4) of the main text.An alternative option to the gate in Fig. 2 (b) is the implementation a three-qubit gate at a hardware level without the need of gate compilation.This possibility, investigated in depth in Sect.V, has a number of advantages with respect to the compiled version of such gate.In fact, a hardware-level implementation of the controlled two-qubit gates would lead to a reduction of the gate count, potentially enabling the simulation of systems of larger size.Moreover, the gate generator can be engineered to conserve the number of excitations in the k-mode registers and thus the SES encoding.Since this criteria is not guaranteed by the complied version of the gate in Fig. 2 (b), an error occurring half-way in such a sequence could severely impact the state fidelity.

C. Polaron transformation and variational form for the multimode Dicke Hamiltonian
In this section, we extend our approach to treat the multimode Dicke Hamiltonian [Eq.( 2) of the main text].To this end, the polaron transformation needs to be modified to include the N atoms as where the real parameters f ik , organized in the vectors f i = (f i1 , . . ., f iM ), for i = 1, . . ., N , now depend on both the atom and the k-mode indices [4].By transforming the Hamiltonian in Eq. ( 2) of the main text as H → H {fi} = P † {fi} HP {fi} , we find which generalizes Eq. ( 2).Here, we derive parameters analogous to those introduced for the quantum Rabi model, including renormalized atom frequencies ω qi = ω qi e −2fi•fi and light-matter coupling constants g ik = (g ik − ω k f ik ), along with the set of operators q fi = exp(2σ x i M k=1 f ik a k ) for i = 1, . . ., N .Unlike the N = 1 case, Eq. ( 8) also includes an effective two-body coupling ) between two atoms i and i , which is mediated by the k-modes.
Due to the latter interaction, the groundstate of Eq. ( 8) is not necessarily to that of the noninteracting case.Instead, we assume a groundstate Ansatz of the form M k=1 |ψ a | 0 k , where |ψ a is the groundstate of the effective spin Hamiltonian By minimizing the energy of H a , one finds a parameterization {f i } of the polaron transformation, that approximately decouples the atoms from the k-modes in Eq. ( 8) [4].In contrast to the case of the quantum Rabi model, the cost of finding a proper disentangling transformation now scales exponentially with system size.
Focusing on the limit of small N , Eq. ( 9) can still be handled on a classical processor which is used for diagonalizing H a and optimizing its groundstate energy as a function of {f i }.It is worth mentioning that quantum routines, such as variational eigensolvers and quantum annealing, might also be useful for this task.In particular, the latter method is attractive for large N due to the possibility of embedding Ising-type Hamiltonians with long-range interactions into physical models with bounded connectivity.
Regardless of how the above optimization is performed, a quantum circuit is needed to prepare the state |ψ a on the atom registers.This initialization step is followed by the application of the polaron variational form ) where σ x i is the Pauli-X operator for the i th atom, and d ik is the Trotter order of the polaron unitary involving this qubit and the k-mode labeled by k.Furthermore, f s ik are the parameters of the polaron variational form, adding to a total of N i=1 M k=1 d ik parameters, which scales linearly with the number of atoms.Additional parameters can be introduced by allowing f s ik → f s ik (n k ) to depend on the k-mode photon-number index.Fig. 3 shows a schematic of the variational form including the initialization step.We assume that |ψ a can be synthesized by a set of d a hardware-efficient layers acting on the atom registers [5].Although this is not scalable to large N , the choice of a hardware-efficient approach is motivated by the lack of structure in Eq. ( 9).As previously shown in Fig. 2, the polaron variational form contains sets of controlled two-qubit gates acting on the k-mode registers.In the present case, the control qubit is swept across the atom registers, while the number of Trotter steps d ik may vary form one set to the other.Finally, we note that the compilation of the controlled two-qubit gates in one-and two-qubit gates is the same as in Fig. 2 (b).

II. PERFORMANCE OF THE VQA ON A QUANTUM PROCESSOR
A. VQA simulations with different hardware-noise levels Mitigating the effect of noise is one of the greatest challenges for near term quantum computers.In order to first quantify this effect, it is useful to investigate the performance of the VQA for a variable noise strength.One way to modify the effective level of noise that a quantum algorithm is subject to, is to perform the quantum gates necessary for the computation having made these artificially slower.This enhances the effect of any decoherence channel and thus leads to an increased noise strength.Although we do not have low-level access to the pulses applied on the quantum hardware, like in Ref. [6], we can simulate the effect of a variable noise strength on the VQA by modifying the error model accordingly.This strategy also allow us to simulate a noise level below the calibrated values for the quantum device in use, which are provided by Qiskit [7].
We modify the noise level in simulation defining a noise factor, η noise , such that where T 1 and T error probabilities and r readout is the readout error probability, respectively.The quantities in Eq. ( 11) which are labeled as "device" correspond to calibrated values for the day that the runs were executed.Fig. 4 shows the result of the simulations with variable noise strength.We use the Simultaneous-Perturbation-Stochastic-Approximation (SPSA) method as classical optimizer, with parameters α = 0.602 and γ = 0.101 defined in Ref. [5].We perform 25 calibration steps to compute the parameters a and c and another 100 SPSA trials for the actual optimization procedure.From the simulations, we conclude that VQA results with acceptable accuracy with respect to the exact groundstate energy would be attainable with noise levels that are one order of magnitude smaller than the actual ones.We note that extrapolation schemes to the zero-noise regime, like the ones discussed in Refs.[6,8], can potentially help to further mitigate the effects of noise.

B. Simulations on the quantum processor
For the quantum-hardware runs we use three qubits out of the twenty available on the IBM Q Poughkeepsie chip.The connectivity map of the device is shown in Fig. 5.The average error rates recorded throughout our experiments were (5.25 ± 0.212) × 10 −3 and (3.75 ± 0.364) × 10 −2 for single-quit gates and CNOT gates, respectively.These error rates are highly depended on the actual date on which the experiment took place.Mitigation of readout errors is done with the standardized methods provided by Qiskit-Ignis.There, a measurement calibration matrix is used to identify readout errors by preparing 2 N basis states and estimating the probability distribution of such state, with N the num-FIG.5. Connectivity map of IBM Q Poughkeepsie quantum processor.Qubits 14, 18 and 19 (highlighted) were used for our experiment.The choice of qubits was depended on the respective single and two quit gate errors.ber of qubits in simulation.The probability distribution of an unknown state can then be corrected based on these estimates.The calibration matrix was updated after every run or every 120 minutes of wall-clock time for the VQA runs.
Beyond coherent and incoherent errors on quantum hardware, the main limitation to greater accuracy has been found to be the classical optimizer.Indeed, the SPSA method fails to acquire the expected solution in several cases, potentially getting stuck into local minima.We confirm this hypothesis indirectly by performing the optimization of the variational Ansatz in simulation assuming the calibrated noise model for the device.We then compute the expectation value of such variational states on quantum hardware.This additional experiment can reach significantly better accuracy than those obtained by running the optimization over quantumhardware energy estimations, as shown in Fig. 2 of the main text.We therefore expect new and more powerful optimization algorithms to enable higher accuracy VQA results.

III. SAMPLING THE JOINT WIGNER FUNCTION
We now provide details about the sampling of the joint Wigner distribution in Eq. ( 8) of the main text.We begin by noticing that measurements of the Pauli products σ l1 1 . . .σ l N N can be done in the computational basis, provided that a set of single-qubit gates {R li } are executed on the atom register prior to qubit readout.Taking this into consideration, we now focus on the case where the Pauli string σ l1 1 . . .σ l N N contains only σ z and/or identity operators and no prior rotation of ρ is needed.
Eq. ( 8) of the main text makes use of the probability distribution of the displaced density matrix D † (α)ρ D(α).In order to sample such a distribution, D † (α) needs first to be compiled into single-and twoqubit gates.Since this joint-displacement operator is a product of single-mode displacements of the form D † (α k ), we only provide the quantum-circuit compilation for the latter unitary.To this end, it is convenient to introduce the real (α R k ) and imaginary (α I k ) parts of the displacement parameter α k = α R k + iα I k , and to expand the displacement operator as where " " indicates an equivalence up to a global phase.Making use of the site operators of the k-mode registers, we find Splitting the k-mode registers into even-and odd-index qubit subsets, the exponentiation of the operators in Eq. ( 12) can be implemented by a Trotter expansion, as it was done for the polaron variational form.Fig. 6 summarizes the procedure for sampling the joint Wigner function with a set of tomography gates applied on the ultrastrong-coupling groundstate synthesized by the polaron variational form in panel (a).A first set of two-qubit gates, compiled in panel (b), implements a displacement operator along the imaginary-α k axis with a Trotter order d I .This is followed by a similar set of gates, compiled in panel (c), implementing a displacement along the real-α k axis with Trotter depth d R .Additionally, single-qubit tomography gates are applied on the atom registers.The circuit is terminated by readout of both the atom and k-mode registers.An histogram of counts (d) is constructed by repeating this procedure for a fixed (l, α) pair.The joint Wigner function can then be computed from this histogram approximating the trace operator in Eq. ( 8) of the main text by N i=1 βq i c(q 1 , . . ., q N ; ñ1 , . . ., ñM ), (13) where c(q 1 , . . ., q N ; ñ1 , . . ., ñM ) are the normalized counts for the basis vector |q 1 , . . ., q N ; ñ1 , . . ., ñM in which q i ∈ [0, 1] is the state of the i th atom, and β qi ∈ [0, 1] accounts for the presence of a σ z i operator before qubit readout.More precisely, such a parameter is set according to the rules β qi = 1 − q i if l i = z, and β qi = 0 if l i = 0. We note that the approximate relation in Eq. ( 13) can be replaced by an exact equivalence in the limit of large counts.
The reconstruction error scales with the amplitude of the displacement parameters |α k | 2 , although it can be reduced by increasing the Trotter order in the implementation of the unitaries D † (α k ).The finite Fock-state truncation of the encoded bosonic modes sets an upper bound to the accuracy of the reconstructed joint Wigner distribution, as demonstrated in Fig. 3 of the main text.Given that the quantum circuit corresponding to such operators is appended to that of the polaron variational form, this procedure would ultimately be limited by the strength of the noise in the quantum processor.However, we find that Trotter depths as small as 2 are enough to demonstrate the qualitative features of the joint Wigner distribution.
This section discuses the so-called memory error of a small qubit register encoding a bosonic mode.We consider both phase-and bit-flip error channels acting on a copy of the maximally entangled state |ψ ∼ n max k n k =0 |n k , in absence of logical gates.Specifically, we compute the state fidelity where ρ = |ψ ψ| ⊗ |ψ ψ| and ρ ∼ Here, E C are multiqubit error channels obtained by composition of single-qubit ones, E q,C .The latter have the general form E q,C (•) = i , where {E i } are the Kraus operator for the channel C. Denoting the error probability with r, we define the phase-flip channel by x correspond to the bit-flip channel.Fig. 7 shows the result for the state fidelity assuming a k-mode register containing up to a maximum of 7 qubits.As anticipated in the main text, we find that bitflip errors are the most relevant as the size of the k-mode register is increased.This can be understood intuitively by looking at the complement SES of the SES subspace used for the encoding.Since the number of basis vectors in SES, and thus the dimension of this subspace, grows exponentially with n max k , a noise operator breaking the SES symmetry could significantly affect the state fidelity in the limit of large n max k .On the other hand, it is worth noticing that this might not necessarily limit the performance of near-term algorithms requiring only a small number of qubits.Alternatives for scaling-up to larger devices include the use of qubits with naturally long T 1 times, or a different encoding for the bosonic modes [10,11].Future work will investigate the performance crossover of the various possible encodings as the variational circuit is scaled up.Finally, we note that state fidelity, although standard, is a strong metric to evaluate the performance of our variational algorithm, and provides only a qualitative estimation of the impact on the energy of the variational ansatz.

V. CIRCUIT-QED IMPLEMENTATION OF THE CONTROLLED-EXCHANGE GATES
With the purpose of reducing the gate count of the polaron variational form, we now present a superconducting-qubit implementation of a controlledexchange gate.We stress, however, that the proposed approach could be leveraged by any other quantumhardware platform with native interactions similar to those found in a standard circuit-QED setup.Below, we provide an ideal implementation of the gate interaction and then suggest a superconducting circuit that approaches the ideal scheme.
We first consider the case of a single atom and k-mode The legend applies for both left and right plots.We observe a significant decrease of the state fidelity as the size of the qubit register is increased.As discussed in the main text, bit-flip errors are expected to be dominant as the register size is scaled up to simulate a larger number of modes with also greater Fock-state truncation.The inset shows the state fidelity as a function of n max k for an error probability of 10 −1 .
registers.The frequency of the physical qubit corresponding to the atom is denoted by ω ν , while the frequencies of the qubits belonging to the k-mode register are denoted by ω n k , with n k ∈ [0, 1, . . ., n max k ].Note that these frequencies are not related to the parameters of the problem that one wishes to simulate.With the purpose of engineering a controlled two-qubit gate, we assume the two qubits of the k-mode register-labeled by µ ∈ [n k , n k + 1]-to be independently coupled to the atom qubit with a time-dependent interaction strength.This situation is described by a 3-qubit Hamiltonian of the form where we take, in particular, Ω µ = Ω 0 µ +2ε µ sin[(ω µ −ω ν + δ µ )t+φ µ ].Here, Ω 0 µ is an always-on interaction strength, 2ε µ is the modulation amplitude, δ µ is a frequency detuning with respect to the µ-ν transition and φ µ is a relative phase.Counter-rotating terms of the form σ + µ σ + ν and its Hermitian conjugate have been omitted after a RWA.
To make the three-qubit-gate interaction explicit, we now perform a standard time-dependent Schrieffer-Wolff transformation with generator conceived to remove first-order interaction terms by the condition H 0 + [H 0 , S] + H int − i Ṡ = 0, where 1, we expand the transformed Hamiltonian up to second order in the interaction strength, and move to a frame rotating at frequencies ω ν for the atom qubit and ω µ + δ µ for the k-mode qubits, where the modulated interaction is resonant.Setting the phase of the drives as φ n k = 0 and φ n k +1 = −π/2, and performing a second RWA, we find the effective Hamiltonian where drive parameters have been chosen to satisfy −δ µ + (Ω 0 µ ) 2 /(ω µ − ω ν ) + ε 2 µ /δ µ = 0, with δ µ = 2(ω µ − ω ν ) + δ µ .The drive condition removes terms ∝ σ z µ from Eq. ( 16) and makes the three-qubit interaction resonant in the current frame.Moreover, ξ n k has been defined as an effective exchange-interaction rate between the two neighboring qubits of the k-mode register, that is mediated by the atom qubit and given by Additionally, we derive a shift to the frequency of the atom qubit given by δω ν = −(Ω 0 µ ) 2 /(ω µ − ω ν ) − ε 2 µ /δ µ , due to the interaction with the two other qubits and the presence of the drive.
Evolution under the Hamiltonian in Eq. ( 16) generates the desired controlled-exchange gate operation exp )] which is key to the polaron variational form.Due to the presence of the term ∝ δω ν in Eq. ( 16), an unintentional R z rotation on the atom qubit needs to be corrected for by applying an additional single-qubit gate.Modulation amplitudes ε µ /2π of the order of 10 MHz and typical values of δ µ /2π of the order of the GHz lead to controlled-exchange rates ξ n k /2π in the range 0.1 − 0.5 MHz.Despite this number being small compared to standard rates of one-and two-qubit gates in superconducting-qubit architectures, counting with a direct implementation of the three-qubit gate still provides a significant advantage with respect to its compiled counterpart in Fig. 2. In fact, the proposed gate is designed to conserve the excitation number of the k-mode registers and thus the SES encoding.Furthermore, while the gate time of the direct implementation is proportional to the desired rotation angle, the compiled version of the gate has an approximately fixed gate time determined by the number of CNOT gates in the circuit.This important difference would be leveraged further as the Trotter order of the polaron variational form is increased, making the controlled rotations closer to the identity.
Having presented an ideal model for the controlledexchange gate, we now elaborate on a possible superconducting-circuit implementation of the Hamiltonian in Eq. (14).In particular, we consider an architecture made of transmon qubits and tunable couplers (see Fig. 8), similar to that studied in Ref. [12].Using couplers to mediate parametric interactions allows us FIG. 8. Variational quantum-optics (VQO) superconducting processor.(a) Schematic of a controlled-exchange gate between three transmon qubits.The qubit in red plays the role of the atom register, controlling the switch on and off of an exchange interaction between two neighboring transmons belonging to the k-mode register.(b) Superconducting-qubit implementation of the concept in (a).A tunable-coupler (light green) is introduced to mediate the interaction between the atom and k-mode registers.(c) Device for the VQO simulation of the ultrastrong interaction between N atoms (in red) and two bosonic modes (light blue, left and right).A superconducting resonator acting as a quantum bus is require to enable long-range interactions between the atoms and the coupler modes.Moreover, the bus mode enables dispersive two-qubit gates between the atom qubits if made frequencytunable [13], which are required for state preparation in the case of N > 1.
to remove the need for frequency tunability of the qubit modes resulting in greater coherence times.Standard circuit quantization of the unit-cell device in Fig. 8 (b), followed by a two-level and rotating-wave approximations leads to the Hamiltonian where ω b c [Φ ext ] denotes the bare frequency of the tunable coupler, ω b ν is the bare frequency of the atom qubit, {ω b µ } are the bare frequencies of two neighboring qubits in the k-mode register and {g ν , g µ } are the respective coupling strengths between such qubits and the coupler.We assume that the coupler frequency can be tuned and modulated by a external magnetic flux Φ ext through the coupler's SQUID loop.Following [12], we perform the adiabatic elimination of the coupler mode by means of a Schrieffer-Wolff transformation in the dispersive regime Assuming that the coupler mode remains in its groundstate at all times, we derive an effective Hamiltonian of the form H eff = H ideal + H err , where H ideal is the ideal interaction model given in Eq. ( 14) with frequency parameters H err is a spurious off-resonant term coupling directly the two qubits of the k-mode register.Implementation of the interaction model in Eq. ( 16) from H eff requires a twotone modulation of Φ ext (t) at frequencies ω µ −ω ν +δ µ for µ ∈ [n k , n k +1].We observe that the effect of H err can be exactly canceled by tuning the qubits to the destructiveinterference condition ∆ However, this also leads to a small interaction strength for the controlled-exchange operation.A better alternative is to consider the two k-mode qubits being coupled to an additional ancillary mode whose frequency is chosen to counteract the effect of H err .Moreover, if the k-mode qubits are properly detuned the residual interaction only leads to a frequency renormalization of the drive condition above and to an off-resonant controlled-exchange interaction between ν and n k (n k + 1) via n k + 1 (n k ) which can be dropped by means of a RWA.As anticipated, we find that for typical circuit parameters and without optimization, the gate-interaction rate ξ n k can reach values in the range of 0.1 − 0.5 MHz, assuming a modulation amplitude between 25 − 50% of Ω 0 µ [12].Further improvements on the speed of the gate might be enabled by optimization of the proposed circuit, the use of other possible coupling schemes implementing Eq. ( 16), or optimal control techniques [14].
The proposed implementation may be scaled-up to a larger number qubits, as shown schematically in Fig. 8 (c) for the case of N atoms and two k-modes.A cavity bus mode is used to enable long-range interactions between the tunable couplers and the qubits playing the role of atoms registers.Moreover, the bus mode allows for the implementation of two-qubit gates between the atom qubits, which are necessary to initialize the polaron variational form for N > 1.We note that the controlled-exchange gates can be parallelized over even and odd qubits of the k-mode registers.Finally, we stress that scaling-up to a larger number of qubits entails issues that are beyond the scope of the present work and require to be examined in greater detail.The analysis of this section, however, provides a path forward towards the implementation of variational-quantum optics algorithms on special-purpose hardware.

MFIG. 1 .
FIG. 1. Groundstate-energy estimation for the single-and two-mode version of the Rabi and Dicke models in resonance conditions.Fock-space truncation in panels (a)-(d) is set to n max k = 3, 3, 5 and 4 corresponding to 5, 9, 8 and 12 qubits, respectively.The legend is shared between all panels, although (a) and (d) display results only for di < 4 and d ik < 5, respectively.We show the error metrics ∆en (left scale) and ∆ex (light-blue dashed line and right scale) defined in the main text, along with H vac|n max k (triangular markers).The inset in panel (b) shows H vac|n max k

FIG. 1 .
FIG.1.Performance of polaron transformation in the strong-, ultrastrong-and deep-strong-coupling regimes for the singleand two-mode quantum Rabi models in resonance conditions.(a) Comparison between H vac and the exact groundstate energy.The present scale emphasizes the qualitative agreement between these two quantities, although it does not capture relatively small deviations which can be appreciated in the main text.Panels (b)-(d) show the polaron transformation parameters (black solid lines) obtained by numerical optimization.Dashed lines are analytical approximations of Eq. (3) at the boundaries of 1 ≥ ω q /ωq ≥ 0. The asymptotic regime of f k g k /ω k is indicated in panel (c) (gray dashed lines).A similar line style in (d) indicates the coupling in the bare frame, g k .We observe that g k becomes exponentially small compared to g k as the bare coupling strength increases.This fact highlights the advantages of working in the polaron frame.

k-mode register atom FIG. 2 .
FIG. 2. Quantum-circuit compilation of the polaron variational form [Eq. (5)] for one of the k-mode registers.(a) Schematic representation of the variational form as a sequence of controlled two-qubit gates.Such gates correspond to the unitaries exp(f s k σ x X e k /d k ) and exp(f s k σ x X o k /d k ), performed at each Trotter step labeled by s.(b) Compilation of the controlled two-qubit gate in one-and two-qubit gates available on the IBM quantum hardware [3].Here, H denotes the Hadamard gate, Y = Rx(π/2), and Rz, R z are rotations around the Z axis using the conventional notation in which Rµ(θ) = exp(−iθσ µ /2) for µ = x, y, z.If the parameters f s k are taken to be independent of the k-mode photon number, then Rz = Rz(f s k √ n k + 1/d k ) and R z = R † z , where n k labels the sites within the k-mode register.If, in contrast,f s k → f s k (n k ) is allowed to vary from site to site, then Rz = Rz[f s k (n k ) √ n k + 1/d k ] and R z = R † z -or R z = Rz[f s k (n k ) √ n k + 1/d k ]if one wishes to introduce an extra degree of freedom, f s k (n k ), per controlled two-qubit gate.Importantly, the parameters of the variational form are initialized to the value specified by the polaron transformation in Eq. (4) of the main text.

FIG. 3 .
FIG.3.Quantum circuit compilation of the variational form for the multimode Dicke problem for 3 atoms and a single kmode register.The variational form is initialized on the state M k=1 |ψa | 0 k by means of hardware-efficient layers acting on the atom registers.The hardware-efficient variational form is followed by layer of controlled two-qubit gates identical to those in 2 (b).Here, the control qubit is swept over the atom registers, and the number of Trotter steps may vary with the atom and k-mode indices.

FIG. 6 .
FIG. 6.Sampling the joint Wigner distribution corresponding to the atom k-mode registers.(a) Quantumcircuit compilation of the displacement operators applied on the k-mode registers (here shown for a single mode), along with the single-qubit gates needed for tomography of the atom registers.These gates are executed on the ultrastrongcoupling groundstate synthesized by the polaron variational form.The k-mode displacement unitary is split in real and imaginary components, implemented by gate sequences with Trotter depth d R and d I , respectively, which may in general be different.The imaginary component is depicted first, and makes use of the two-qubit gates in (b) where Rz = Rz(αI k √ ñk + 1/d I ) and R z = R † z .What follows is the implementation of the real component of the displacement operator, which executes the two-qubit gates in (c) whereRz = Rz(α R k √ ñk + 1/d R ) and R z = R † z .Qubit readout is performed at the end of the quantum circuit in (a), leading to the histogram of counts in (d) after several repetitions of the experiment.This allows for reconstruction of the joint Wigner function as described in the main text.

10 − 3 10 −2 10 −FIG. 7 .
FIG.7.Effect of phase-flip (a) and bit-flip (b) error channels on the fidelity of a maximally entangled state as a function of the total error probability.We consider a k-mode qubit register of size n max k + 1 ≤ 7. The legend applies for both left and right plots.We observe a significant decrease of the state fidelity as the size of the qubit register is increased.As discussed in the main text, bit-flip errors are expected to be dominant as the register size is scaled up to simulate a larger number of modes with also greater Fock-state truncation.The inset shows the state fidelity as a function of n max
FIG.4.Simulation of the VQA performance under variable noise strength.ηnoise = 1 corresponds to the realistic noise model provided by Qiskit for the IBM Q Poughkeepsie device.The shared area indicates the range of results which are possible with levels of noise ηnoise = 0.1 − 2. Extrapolation to zero-noise indicates that results with acceptable accuracy are only attainable with levels of noise one order of magnitude smaller than those in the current generation of devices.
2 are single-qubit relaxation and dephasing times, r 1q−g and r 2q−g are single-and two-qubit gate