Variational Simulation of Schwinger's Hamiltonian with Polarisation Qubits
O. V. Borzenkova, G. I. Struchalin, A. S. Kardashin, V. V. Krasnikov, N. N. Skryabin, S. S. Straupe, S. P. Kulik, J. D. Biamonte
VVariational Simulation of Schwinger’s Hamiltonian with Polarisation Qubits
O. V. Borzenkova, ∗ G. I. Struchalin, A. S. Kardashin, V. V. Krasnikov, N. N. Skryabin, S. S. Straupe, S. P. Kulik, † and J. D. Biamonte ‡ Skolkovo Institute of Science and Technology, 3 Nobel Street, Moscow 121205, Russian Federation Quantum Technology Centre and Faculty of Physics, M.V. Lomonosov Moscow State University,1 Leninskie Gory Street, Moscow 119991, Russian Federation (Dated: September 22, 2020)The numerical emulation of quantum physics and quantum chemistry often involves an intractablenumber of degrees of freedom and admit no known approximations in a general form. In practice,representing quantum-mechanical states using available numerical methods become exponentiallymore challenging with increasing system size. Recently quantum algorithms implemented as varia-tional models, have been proposed to accelerate such simulations. Here we study the effect of noiseon the quantum phase transition in the Schwinger model, within a variational framework. Theexperiments are built using a free space optical scheme to realize a pair of polarization qubits andenable any two-qubit state to be experimentally prepared up to machine tolerance. We specificallyexploit the possibility to engineer noise and decoherence for polarization qubits to explore the limitsof variational algorithms for NISQ architectures in identifying and quantifying quantum phase tran-sitions with noisy qubits. We find that noise does not impede the detection of the phase transitionpoint in a large range of noise levels.
I. INTRODUCTION
The numerical emulation of quantum systems under-pins a wide assortment of science and engineering andtouches on fields ranging from statistical and quantumphysics to biology and even to life- [1, 2] and behavioral-sciences [3–5]. A physical simulator bootstraps onephysical system to emulate the properties of another.While the time and memory required in the simulation ofphysical systems, particularly strongly correlated many-body quantum systems, using traditional computers of-ten scales exponentially in the system size, the same isnot always true for the physics-based quantum simulator.Indeed, Richard Feynman first speculated that instead ofviewing the simulation of quantum systems using classi-cal computers as a no-go zone due to its apparent com-putational difficulty, Feynman argued [6] that physicalsystems themselves naturally posses the computationalcapacity to be harnessed and used.Variational approaches to optimization and simulationof eigenstates [7–12] have been used recently to port ideasfrom machine learning [13] to enhance algorithms withquantum processors [13–15]. These approaches rely onan iterative quantum-to-classical variational procedure.Proven to be a universal model of quantum computa-tion in [16]—where the ansatz circuits are proven to beuniversal in [17]—the variational approach to quantumcomputation arose naturally as the pathway between astatic simulator and a fully programmable gate-basedquantum information processor. The variational modelof quantum computation is the algorithmic workhorse of ∗ [email protected] † https://quantum.msu.ru ‡ https://quantum.skoltech.ru the current (NISQ: Noisy Intermediate-Scale Quantum)technology era.Up-to-date experimental works realize variational algo-rithms on different quantum hardware including super-conducting qubits [10, 11], trapped atoms [8, 9, 12] andphotonic quantum processors [7, 18]. The most commonapplication of quantum variational algorithms includesquantum chemistry problems. Despite the fact that themain purpose of the algorithm is finding eigenvalues andeigenvectors, [19] show that variational techniques canalso find excited states, and various other proposals fur-ther expand the limits of applicability [20].The variational quantum eigensolver (VQE) [7] per-forms classical optimization to minimize a mean Hamilto-nian expected value—found by quantum hardware. Thepurpose of this algorithm is to determine the eigenvaluesof a particular Hamiltonian, which describes a physicalsystem, for example, the interaction of spins or electronicsystems [8, 9]. A classical computer initially sets a vec-tor of parameters θ = { θ i } for i ∈ N and an experimen-tal setup prepares a quantum state | ψ ( θ ) (cid:105) parametrizedby these control parameters. After that, the state ismeasured, and the evaluation of the mean Hamiltonianvalue occurs. The parameters θ are adjusted to find theground-state energy: E min ( θ ) = min θ (cid:104) ψ ( θ ) | H | ψ ( θ ) (cid:105) (1)Therefore, the problem consists of using classical opti-mization algorithms to select optimal parameters θ cor-responding to the (ideally) minimal value of energy.Here we report an experimental implementation ofVQE in a photonic system. We target the exploration ofa quantum phase transition in the Schwinger model. Wespecifically exploit the possibility to engineer noise anddecoherence [21] for polarization qubits to explore thelimits of variational algorithms for NISQ architectures in a r X i v : . [ qu a n t - ph ] S e p HWP1QWP1
Dichroic mirror
Filter L D H W P ( ) Q W P ( ) C o i n c i d e n c e c i r c u i t Wollaston prism
Polarization-entangled photons
HWP2 L C V R FIG. 1. Experimental setup implementing the variational quantum eigensolver algorithm using a pair of polarization qubits.Half-wave plates (HWP3, HWP4) and quarter-wave plates (QWP2, QWP3) in each channel prepare the desired variational state | ψ ( θ ) (cid:105) . Wollaston prisms implement projective measurements in the computational basis (horizontal and vertical polarization).HWP1 and QWP1 control the polarization of the pump beam, while HWP2 rotates the polarization by 90 ◦ to ensure the samepolarization of the pump for the two pumping directions for the nonlinear crystal (PPKTP) inside the Sagnac interferometer.Liquid crystal variable retarder (LCVR) is used to arificially introduce dephasing noise on the first qubit, when necessary. identifying and quantifying quantum phase transitionswith noisy qubits. II. ENCODING THE ANSATZ STATE INPOLARIZATION QUBITS
Ansatz state preparation is strongly connected withthe particular experimental realization, because elementswith tunable parameters used in the setup determine theparametrization of the ansatz. Therefore, we start fromthe experimental scheme to clarify the origin of the cho-sen ansatz.We implement a VQE algorithm using polarization-encoded qubits. The experimental setup scheme is shownin Fig. 1. The initial state preparation is carried out bya two-photon source based on spontaneous parametricdown-conversion process (SPDC) in the Sagnac interfer-ometer [22]. A 405-nm laser diode beam is divided by apolarization beam-splitter (PBS), which makes it possi-ble to pump a 30-mm long periodically poled KTP non-linear crystal (PPKTP) in two opposite directions. As aresult of a type-II SPDC, pairs of signal and idler pho-tons with orthogonal polarizations are generated in bothdirections. Then each photon pair is divided on the PBSand sent to different arms of the scheme. Thus, at theoutput of the two-photon source, we have the following entangled state: | ψ in (cid:105) = α ( θ , θ ) | HV (cid:105) + β ( θ , θ ) | V H (cid:105) , (2)where the coefficients α and β depend on the angularpositions θ and θ of waveplates QWP1 and HWP1,which are placed in the pump beam. By rotating QWP1and HWP1, we can alter the degree of entanglement ofthe initial state. The photon pairs are coupled to single-mode fibers and transferred to the measurement part ofthe setup. Motorized quarter-wave (QWP2, QWP3) andhalf-wave (HWP3, HWP4) plates are placed in each armafter the single-mode fiber channel, allowing to obtainany polarization state at the output. Finally, the Wol-laston prism spatially separates the vertical and the hor-izontal polarizations to detect the prepared states usingsingle-photon detectors in each of the arms. According tothe measurement results, the classical algorithm transfersthe new parameter values to the motorized plates untilthe optimal set of parameters is obtained.We should note that estimation of a single mean valueof a Hamiltonian requires projective measurements inseveral bases, while the Wollaston prism projects onlyonto | H (cid:105) and | V (cid:105) states. To change the basis one may usean additional pair of QWP and HWP retarders, mountedjust before the Wollaston prism. However, we chose amore economic setup, where the local unitary transfor-mation of the initial state and the transformation of the | (cid:105) U QWP ( θ ) U HWP ( θ ) • U QWP ( θ ) U HWP ( θ ) | (cid:105) X X U
QWP ( θ ) U HWP ( θ )FIG. 2. Schematic of the VQE algorithm. The initial state isprepared by three single-qubit gates and a Controlled-X gate. U QWP ( θ ) and U HWP ( θ ) are used to control the initial statein experiment. Other four single-qubit transformations serveboth for the ansatz state preparation and the measurementbasis change. measurement basis are compiled together, e. g., (cid:104) H | BU HWP ( θ ) U QWP ( θ ) = (cid:104) H | U HWP ( θ (cid:48) ) U QWP ( θ (cid:48) ) , (3)here U ( θ ) is a transformation of a corresponding wave-plate with an axis angle θ and B is a unitary matrix thatchanges the basis. New angles θ (cid:48) are calculated auto-matically in our algorithm to perform measurements indesired bases.By mapping experimental optical elements to the gatemodel we arrive at the ansatz preparation circuit withsix tunable parameters θ i , i = 1 , . . . , θ i physically correspond to thewaveplates’ rotation angles. A general waveplate with aphase shift δ and an axis position θ performs the trans-formation U ( δ, θ ): U ( δ, θ ) = V ( θ ) D ( δ ) V † ( δ ) , (4) V ( θ ) = θ ) − iσ y sin( θ ) , D ( δ ) = e iδ | (cid:105)(cid:104) | . A controlled-X gate corresponds to the SPDC of thepump photon in the nonlinear crystal.Taking into account ansatz preparation scheme, ourVQE algorithm implementation consists of the four mainsteps:1. SPDC source emits the initial entangled state | ψ in (cid:105) (2).2. Once the initial state has been prepared, a localunitary transformation U ⊗ U is applied to getthe probe state | ψ ( θ ) (cid:105) : | ψ ( θ ) (cid:105) = U ( θ , θ ) ⊗ U ( θ , θ ) | ψ in ( θ , θ ) (cid:105) . (5)Unitaries U and U are composed of the waveplatetransformations: U = U HWP ( θ ) U QWP ( θ ), U = U HWP ( θ ) U QWP ( θ ).3. The cost function E ( θ ) = (cid:104) ψ ( θ ) | H | ψ ( θ ) (cid:105) is cal-culated by summing up measurement results withcoefficients depending on the problem Hamiltonian.Since usually Hamiltonian is expressed as a linearcombination of Pauli observables and our setup al-lows only projective measurements, we first shoulddecompose the Hamiltonian as a linear combina-tion of projectors onto eigenbases of Pauli matri-ces. Change of basis is carried out according to therule (3). 4. The value of E ( θ ) is minimized as a function ofthe parameters θ using a classical optimizer rou-tine. In particular we use simultaneous perturba-tion stochastic approximation (SPSA) algorithm.The reader is referred to Appendix A for details. III. THE SCHWINGER MODEL
The Schwinger model describes interactions betweenDirac fermions via photons in a two-dimensional space.In Ref. [12], the authors map the model to the lat-tice model of an electron-positron array. The SchwingerHamiltonian exhibits a quantum phase transition: thesignature of which (in finite dimensions) allows us to de-termine new features in VQE behavior and clarify itsrobustness to noise.The Schwinger Hamiltonian H N describes electron-positron pair creation and annihilation, their interactionand takes into account the particle mass: H N = w N − (cid:88) j =1 [ σ + j σ − j +1 + H.c. ] + m N (cid:88) j =1 ( − j σ zj − g N (cid:88) j =1 L j . (6)It consists of the three terms: the first one is responsi-ble for the interaction of an electron and a positron, thesecond depends on bare mass m of the particles, and thethird stands for the energy of the electric field. We as-sume the coefficients w = g = 1 and only consider the de-pendence of the Hamiltonian ground energy on the baremass. The operators in the third term are given by L j = (cid:15) − j (cid:88) l =1 [ σ zl + ( − l ] , (7)where we set the background electric field parameter (cid:15) to zero.The problem Hamiltonian can be encoded in the mul-tiqubit system by using its decomposition into Paulistrings: P α = σ α ⊗ σ α ⊗ . . . ⊗ σ α N N with single-qubitPauli operators σ α i i ∈ { I, σ xi , σ yi , σ zi } as H = (cid:88) a h α P α , (8)where N denotes the number of qubits and h α ∈ R arereal coefficients. In further consideration, we will use thisrepresentation. We carried out numerical simulationsand experiments for the case of two qubits, for whichthe Schwinger Hamiltonian takes the form H = σ x σ x +2 σ y σ y − σ z + 12 σ z σ z + m σ z − σ z ) . (9)The quantum phase transition manifests itself in thebehavior of the order parameter (cid:104) O (cid:105) = 12 N ( N − (cid:88) j>i (cid:104) (1+( − i σ zi )(1+( − j σ zj ) (cid:105) . (10) - - - - - - -
10 m E i g e n v a l u e (a) - - O r d e r p a r a m e t e r (b) FIG. 3. Eigenvalue (a) and order parameter (b) versusbare mass m . Solid lines—analytical solution, cyan points—simulation, red points—experimental result. For polarization-encoded pair of qubits the order param-eter is simply a projector onto | V H (cid:105) state: (cid:104) O (cid:105) = 14 (cid:104) (1 − σ z )(1 + σ z ) (cid:105) = (cid:104)| V H (cid:105)(cid:104)
V H |(cid:105) . (11)Two-qubit Schwinger Hamiltonian has four non-degenerate eigenvalues E , . . . , E . Two intermediateeigenvalues, E = 2 and E = 1, are constant and donot depend on the mass m . The largest and the smallesteigenvalues, E , = 1 / ± (cid:112) m + m + 17 /
4, vary with m in a symmetric manner.We are interested in the ground energy of the Hamilto-nian that corresponds to the minimal eigenvalue E min ≡ E . The graph of its dependence on m is depicted inFig. 3a. Also Fig. 3b shows the order parameter versusmass. The solid lines correspond to the exact analyt-ical solutions, dots represent the results of simulationsand experiment. A phase transition signifies itself in therapid change of the order parameter from one to zero andit is expected near the point m = − /
2, where (cid:104) O (cid:105) = 1 / m = − /
2, we found a discrepancy between analyti-cal solutions and VQE simulations. The Hamiltonian H ( m = − /
2) has the ground energy E min = − / | (cid:105) − | (cid:105) ) / √ | Ψ − (cid:105) . Adistinguishing feature of the singlet state is its invari-ance under local unitary rotations U ∈ SU(2): | Ψ − (cid:105) = ( U ⊗ U ) | Ψ − (cid:105) . Therefore, the target function E ( θ ) re-mains constant on some parameter manifold. Note thatthis plateau does not change with the mass m , because ∀ m : (cid:104) Ψ − | H ( m ) | Ψ − (cid:105) = − / m = − / m (cid:54) = − / E min < − /
2, while residing near the plateau with E = − /
2. So the landscape of E ( θ ) in the puncturedneighbourhood of m = − / m is far away from the phase transition point, the plateaudoes not strongly influence the results, because E min ismuch lower than E = − / | Ψ − (cid:105) being theHamiltonian eigenvector for m = − /
2. A more generalview on the cause of the convergence problem is that itappears any time, when the ansatz is general enough toperform arbitrary local unitary transformations and theHamiltonian ground state is close to some Bell state (notnecessarily | Ψ − (cid:105) ). Indeed, all Bell states are equivalentunder local transformations, so we can find a local mapthat brings a Bell state | ψ (cid:105) to a singlet one | Ψ − (cid:105) : | Ψ − (cid:105) = ( W ⊗ W ) | ψ (cid:105) , (12)where W , are some single-qubit unitary matrices. Con-sequently, an arbitrary Bell state | ψ (cid:105) is invariant underthe following transformation: ∀ U ∈ SU(2) : | ψ (cid:105) = ( W † U W ⊗ W † U W ) | ψ (cid:105) . (13)If ansatz circuit is general enough to prepare differenttransformations of the form (13), then the plateau in thelandscape of E ( θ ) appears. Therefore, when the Hamil-tonian ground state is close to the Bell state, the nearbyplateau will create flat valley landscape.The simplest opportunity to get around poor optimizerconvergence is by a correct choice of the initial point.We gathered statistics for 10 random initial points θ for m = − /
2, 0, 1 /
2, and 10 and found that near thephase transition the algorithm sticks to the plateau muchfrequently than to the proper minimum (see Appendix Bfor details).
IV. EFFECTS OF NOISE
Compared to other types of quantum computers, pho-ton circuits have low intrinsic noise levels. This meansthat we can add noise to the system in a controlled man-ner and get the dependencies of the parameters of in-terest on the noise level. We took advantage of this toevaluate the effect of noise on the phase transition that ϵ - - - - - - -
10 m S i m u l a t e d e i g e n v a l u e (a) - - - - - - -
10 m E i g e n v a l u e (b) ϵ - - O r d e r p a r a m e t e r (c) - - O r d e r p a r a m e t e r (d) FIG. 4. Noise simulations for dephasing of both qubits (a), (c) and one qubit (b), (d). Figs. (a) and (b) show the minimaleigenvalue dependence on m and Figs. (c) and (d)—the dependence of the order parameter. Red lines correspond to noiselesssimulations, the color blur corresponds to the increase of noise strength (cid:15) from 0 . . we observed without the noise. We expect that as thedegree of dephasing increases, the phase transition willblur until it disappears completely. This will allow us toestimate the acceptable noise level in the system imple-menting VQE to identify quantum phase transitions.The origin of the noise model used is connected withour experimental implementation. We artificially intro-duce noise to the system with liquid crystal variable re-tarders (LCVR) that allow us to change the phase of thespecific polarization component of the light field. If thephase shift δ varies during the data acquisition time, thenthis leads to effective decoherence of the system state.The noise channel E ( ρ ) is thus the transformation (4)averaged over δ taken from some interval (depending onthe noise strength). The explicit action of the noise chan-nel is ρ (cid:48) = E ( ρ ) = (cid:88) j =1 E j ρE † j , (14) E i = V ( θ ) D i ( δ ) V † ( θ ) ,D ( δ ) = (cid:114) − (cid:15) (cid:18) e iδ
00 1 (cid:19) , D ( δ ) = (cid:114) (cid:15) (cid:18) e iδ − (cid:19) , where E j are the Krauss operators, θ is a LCVR axisangle, and δ is a mean retardance. Noise strength iscontrolled by the parameter (cid:15) , 0 ≤ (cid:15) ≤
1. We set θ = π/ (cid:15) ranging from 0 to 1 with a 0 . m . As expected, the presence of noise in thesystem prevents the algorithm from converging to the ex-act eigenvalue, and noise escalation leads to convergencedeterioration (Fig. 4).Finding appropriate eigenvalue becomes challengingfor the case of simultaneous dephasing in both channels,and full dephasing ( (cid:15) = 1) leads to degeneracy—the al-gorithm converges to 1 for any m . The phase transitionin the order parameter blurs with increasing noise anddisappears for (cid:15) = 1. Full dephasing makes the order pa-rameter constant and equal 1 / m . In the case ofa single noise channel the phase transition remains visi-ble even with (cid:15) = 1, while the maximum value of orderparameter is halved. V. CONCLUSION
Quantum phase transitions as metal-insulator transi-tion and transition between quantum Hall liquid states,can be predicted and inquired by quantum algorithms.As we experimentally demonstrated, noise does not im-pede the detection of the phase transition point in a largerange of noise levels. Only completely dephasing chan-nels acting on both qubits prevent finding it in our model.This result demonstrates the noise-tolerance of VQE notonly from speed and quality of convergence perspectivebut also from a practical point of view of determiningthe parameters of the Hamiltonian corresponding to aquantum phase transition.We observe slow VQE convergence near the phase tran-sition point and connect this behavior with the Hamil-tonian ground state’s closeness to the two-qubit singletstate. It seems to be a common effect for a combina-tion of sufficiently general ansatz circuits and Hamilto-nians, where the ground state exhibits additional sym-metry. This hypothesis should be verified in future re-search. Possible approaches to circumvent poor conver-gence may include quantum approximate optimizationalgorithm (QAOA) [15], because it uses specific ansatzadjusted for the target Hamiltonian.
ACKNOWLEDGMENTS
The Skoltech team acknowledges support from the re-search project,
Leading Research Center on QuantumComputing (agreement No. 014/20). The MSU team ac-knowledges financial support from the Russian Founda-tion for Basic Research (RFBR Project No. 19-32-80043and RFBR Project No. 19-52-80034) and support underthe Russian National Technological Initiative via MSUQuantum Technology Centre.
Competing interests:
The authors declare no competing interests.
Data andcode availability:
The data that supports this studyare available within the article. The code for generat-ing the data will be made available on GitHub after thispaper is published.
Appendix A: Classical optimizer
The target function under minimization E ( θ ) is themean Hamiltonian value, but in the experiment, onlyrandom samples of E ( θ ) obtained by repetitive measure-ments are available. So experimental VQE is a stochasticoptimization problem. We use a simultaneous perturba-tion stochastic approximation (SPSA) algorithm [24] as aclassical optimizer in our VQE implementation. It is use-ful for high-dimensional problems, where the gradient ofthe objective function is not directly available, becauseSPSA requires only two function evaluations per itera-tion for any number of parameters in the optimizationproblem.Single SPSA iteration proceeds as follows:1. Generate a random vector ∆ with elements being ± g : g = E ( θ + b ∆ ) − E ( θ − b ∆ )2 b ∆ . (A1)3. Move to the new point θ (cid:48) : θ (cid:48) = θ − a g . (A2)Scalar variables a and b are called meta parameters . Theparameter a describes the iteration step and b definesfinite difference to calculate the gradient. They changewith the number of iterations k according to schedule: a ( k ) = a − a f k . + a f , b ( k ) = b − b f k . + b f . (A3)Usually final values a f and b f are set to zero to en-sure convergence in the limit k → ∞ . However, we usenonzero a f and b f to track a slow drift of the experimen-tally prepared probe state | ψ ( θ ) (cid:105) over time [25]. The driftoccurs mainly due to the instability of polarization trans-formation in optical fibers connecting the SPDC sourceand measurement part of the setup.Moreover, we find out influence of mass parame-ter m on VQE convergence—closeness to phase transi-tion makes it slowly. So we adjust meta parameters foreach m as a ,f ( m ) = ¯ a ,f . m + 1 , b ,f ( m ) = ¯ b ,f . m + 1 (A4)In our simulations and the experiment we used ¯ b = 0 . b f = 0 .
002 and tried different ¯ a and ¯ a f to find trade-offbetween the number of iterations and accuracy. For ¯ a =0 .
01 and ¯ a f = 0 .
003 convergence is slowly, especially inthe experiment. After different simulations we chose ¯ a =0 .
05 and ¯ a f = 0 . Appendix B: Initial point statistics
We carried out numerical simulations of the VQE algo-rithm for m = 0 , − / , ,
10 to investigate how the choiceof an initial point θ affects convergence and explorethe set of obtained solutions. Recall that the SchwingerHamiltonian H (9) undergo phase transition of the or-der parameter at m = − /
2, so points m = 0 and m = 1are nearby and symmetric w. r. t. phase transition and m = 10 is an example of a distant point. To collectstatistics, we execute the VQE algorithm 10 times foreach m starting from freshly generated random initialpoints θ . The points are distributed uniformly in a six-dimensional hypercube with the side length equal to π ,which coincides with the period of the target function E ( θ ).Each VQE run results in the final point θ , the energylevel E ( θ ) (eigenvalue), and the order parameter (cid:104) O (cid:105) .Fig. 5 shows histograms of (cid:104) O (cid:105) for m = 0 and m = 1. N u m b e r o f r un s FIG. 5. Order parameter histograms for m = 0 (darker, left)and m = − VQEruns. Red vertical lines present analytical solutions. Thehistograms are nearly a reflection of each other around (cid:104) O (cid:105) =1 / As one can see, there is a sharp peak near a wrong value (cid:104) O (cid:105) = 1 / (cid:104) O (cid:105) ≈ .
38 and (cid:104) O (cid:105) ≈ .
62 for m = 0 and m = −
1, respectively. As it was said in themain text, (cid:104) O (cid:105) = 1 / m . As expected, the histogram for m = 10 is unimoal and centered around the exact eigen-value for the corresponding Hamiltonian. When m ap-proaches phase transition point m = − /
2, the secondpeak emerges around E = − /
2, which is precisely theHamiltonian eigenvalue for m = − /
2. This erroneouspeak is small for m = 1, but it becomes even higher thanthe true one for m = −
1. Closeness to phase transitionpoint changes convergence statistics dramatically—lessthan 1% of simulations reveal proper values for m = − θ . First, we bringall found points θ to a hypercube that corresponds tothe target function period. After that, for graphic pur-poses, we decreased the dimensionality of obtained solu-tions θ from six to three using principal component anal-ysis (PCA). PCA finds a lower-dimensional hyperplane inthe original space, which has a minimum average squareddistance from the points to the hyperplane. Then pointsare projected to the approximating hyperplane. PCAhelps to keep the real structure of the initial space andfind any clusters of points with the same values.Fig. 7 presents obtained PCA projections for m = − .
5, 0, and 1 in two views, which will be called “top”and “side” for convenience. Color shows target functionvalues E ( θ ). Blue points correspond to erroneous eigen-values for the given m . The overall structure is similarfor different m , especially for − . m = 0.For m = 1, the fraction of good solutions increases, blueareas slowly disappear. This suggests that the algorithm - - - - - - - - N u m b e r o f r un s (a) m = 10 - - - - - - - N u m b e r o f r un s (b) m = 1 - - - - - - - - N u m b e r o f r un s (c) m = 0 FIG. 6. Histograms of eigenvalues that are found during 10 VQE runs for m = 10 (a), m = 1 (b), and m = 0 (c). Redvertical lines show analytical solutions. converges to the desired point with higher probability,which is in perfect agreement with the histogram in Fig.6b. There are two classes of proper solutions: pointsfrom the first class are located “inside” the erroneousplateau and for the second lie “outside”. However, itappears difficult to isolate areas of initial points θ thatcan guarantee finding the true minimum or lead to oneor another class of solutions, and additional research isrequired. - - - - - - (a) Top, m = − . - - - - - - (b) Side, m = − . - - - - - - - (c) Top, m = 0. - - - - - - - (d) Side, m = 0. - - - - - - (e) Top, m = 1. - - - - - - (f) Side, m = 1. FIG. 7. Principal component analysis of VQE solution points θ for m = − . m = 0 (c, d), and m = 1 (e, f). Leftfigures correspond to “top” view and right ones to “side” projection. Color shows target function values E ( θ ).[1] N. Lambert, Y.-N. Chen, Y.-C. Cheng, C.-M. Li, G.-Y.Chen, and F. Nori, Nature Physics , 10 (2013). [2] F. Neukart, G. Compostella, C. Seidel, D. Von Dollen, S. Yarkoni, and B. Parney, Frontiers in ICT , 29 (2017).[3] T. Werlang, C. Trippe, G. Ribeiro, and G. Rigolin, Phys-ical review letters , 095702 (2010).[4] J. Abadie, B. P. Abbott, R. Abbott, T. D. Abbott,M. Abernathy, C. Adams, R. Adhikari, C. Affeldt,B. Allen, G. Allen, et al. , Nature Physics , 962 (2011).[5] M. Lubasch, J. Joo, P. Moinier, M. Kiffner, andD. Jaksch, Physical Review A , 010301 (2020).[6] R. P. Feynman, Foundations of Physics , 507 (1986).[7] A. Peruzzo, J. McClean, P. Shadbolt, M.-H. Yung, X.-Q.Zhou, P. J. Love, A. Aspuru-Guzik, and J. L. Obrien,Nature communications , 4213 (2014).[8] M.-H. Yung, J. Casanova, A. Mezzacapo, J. Mcclean,L. Lamata, A. Aspuru-Guzik, and E. Solano, ScientificReports , 3589 (2014).[9] Y. Shen, X. Zhang, S. Zhang, J.-N. Zhang, M.-H. Yung,and K. Kim, Physical Review A , 020501 (2017).[10] P. J. OMalley, R. Babbush, I. D. Kivlichan, J. Romero,J. R. McClean, R. Barends, J. Kelly, P. Roushan,A. Tranter, N. Ding, et al. , Physical Review X , 031007(2016).[11] A. Kandala, A. Mezzacapo, K. Temme, M. Takita,M. Brink, J. M. Chow, and J. M. Gambetta, Nature , 242 (2017).[12] C. Kokail, C. Maier, R. van Bijnen, T. Brydges, M. K.Joshi, P. Jurcevic, C. A. Muschik, P. Silvi, R. Blatt, C. F.Roos, et al. , Nature , 355 (2019).[13] J. Biamonte, P. Wittek, N. Pancotti, P. Rebentrost,N. Wiebe, and S. Lloyd, Nature , 195 (2017), arXiv:1611.09347 [quant-ph].[14] V. Akshay, H. Philathong, M. Morales, and J. Biamonte,Physical Review Letters (2020), 10.1103/phys-revlett.124.090504.[15] E. Farhi, J. Goldstone, and S. Gutmann, “A quantumapproximate optimization algorithm,” (2014), unpub-lished, arXiv:1411.4028.[16] J. Biamonte, “Universal variational quantum computa-tion,” (2019), arXiv:1903.04500 [quant-ph].[17] M. E. Morales, J. Biamonte, and Z. Zimbor´as, arXivpreprint arXiv:1909.03123 (2019).[18] J. Carolan, M. Mohseni, J. P. Olson, M. Prabhu,C. Chen, D. Bunandar, M. Y. Niu, N. C. Harris, F. N.Wong, M. Hochberg, et al. , Nature Physics , 322(2020).[19] J. R. McClean, S. Boixo, V. N. Smelyanskiy, R. Bab-bush, and H. Neven, Nature Communications (2018),10.1038/s41467-018-07090-4.[20] O. Higgott, D. Wang, and S. Brierley, Quantum , 156(2019).[21] A. Pechen, Physical Review A (2011), 10.1103/phys-reva.84.042106.[22] A. Fedrizzi, T. Herbst, A. Poppe, T. Jennewein, andA. Zeilinger, Optics Express , 15377 (2007).[23] H. H. Rosenbrock, The Computer Journal , 175 (1960).[24] J. C. Spall, IEEE Transactions on Automatic Control ,332 (1992).[25] O. Granichin and N. Amelina, IEEE Transactions on Au-tomatic Control60