Finite-Size Scaling on a Digital Quantum Simulator

Finite-Size Scaling on a Digital Quantum Simulator Using Quantum Restricted Boltzmann Machine The critical point and the critical exponents for a phase transition can be determined using the Finite-Size Scaling (FSS) analysis. This method assumes that the phase transition occurs only in the infinite size limit. However, there has been a lot of interest recently in quantum phase transitions occurring in finite size systems such as a single two-level system interacting with a single bosonic mode e.g., in the Quantum Rabi Model (QRM). Since these phase transitions occur at a finite system size, the traditional FSS method is rendered inapplicable for these cases. For cases like this, we propose an alternative FSS method in which the truncation of the system is done in the Hilbert space instead of the physical space. This approach has previously been used to calculate the critical parameters for stability and symmetry breaking of electronic structure configurations of atomic and molecular systems. We calculate the critical point for the quantum phase transition of the QRM using this approach. We also provide a protocol to implement this method on a digital quantum simulator using the Quantum Restricted Boltzmann Machine algorithm. Our work opens up a new direction in the study of quantum phase transitions on quantum devices. Keywords: finite-size scaling, quantum phase transition, quantum simulator, quantum restricted Boltzmann machine, quantum rabi model 1 INTRODUCTION A phase transition occurs whenever the thermodynamic functions of a system become nonanalytic e.g. as a liquid changes into a gas, the density of the system changes discontinuously. If the phase transition occurs at a finite temperature T ≠ 0, the transition is called a classical phase transition (CPT) as it is dominated by thermal fluctuations. On the other hand, if the transition occurs by tuning some parameter in the system’s Hamiltonian as T → 0, it is called a quantum phase transition (QPT) since it is dominated by quantum fluctuations. A CPT appears only when the system is infinite i.e., in the thermodynamic limit [1]. On the other hand, a QPT doesn’t necessarily require the thermodynamic limit. Recently there has been a lot of interest in QPTs occurring in finite size light-matter interaction systems [2–7]. Quantum Rabi Model (QRM) describes the interaction of a two-level system with a bosonic field mode (see Eq. 1 for the Hamiltonian.) This model has gained a lot of significance in the study of ultrastrong light-matter coupling regimes where the so-called counterrotating terms can not be ignored [8]. Quantum Rabi Model has been shown to exhibit a QPT [2]. Namely, when the energy separation of the two levels in the system Ω becomes infinitely large compared to the frequency of the bosonic mode ω0, the ground state of the Hamiltonian undergoes a phase transition from a normal phase to a superradiant phase as the light-matter coupling exceeds the critical value. Moreover, the ground state of the JaynesCummings model (JCM) which can be obtained from the QRM by performing the rotating-wave approximation has also been shown to exhibit the normal-superradiant phase transition [3]. Later on, a more general anisotropic QRM in which the rotating and counter-rotating terms can have different coupling strengths was also considered [4]. The QRM and JCM are limiting cases of this model. It was shown that the ground state for this more general case also undergoes the normal-superradiant phase transition. The phase transition in QRM has also been demonstrated experimentally using a 171Yb+ ion in a Paul trap [7]. This experimental demonstration of a phase transition in a single two-level system has incited a lot of interest since this opens up an avenue for studying critical phenomena in controlled, small quantum systems. In CPTs and some QPTs (which require N → ∞), a finitesize scaling (FSS) analysis can be done to extract the critical point and the critical exponents of the transition [1, 9]. While this procedure is inapplicable to the QPTs discussed above since these phase transitions occur at a finite system size, the phase transitions in these paradigmatic light-matter interaction models occur only in the limit Ω/ω0 → ∞ and FSS analysis can be done in Ω/ω0 [2–4] instead. In this paper, however, we propose a different approach to study such phase transitions. We apply the FSS in Hilbert space method [10–15] to the QPT in Quantum Rabi Model. In this approach, the truncation of the system is done not in the physical space but in the Hilbert space. The set of basis states spanning the infinite dimensional Hilbert space is truncated to a finite set and the scaling ansatz is employed in terms of the size of this set. This approach has previously been developed and applied to a single particle in Yukawa potential [11, 13] and the problem of finding electronic structure critical parameters for atomic and molecular systems [10, 12, 14–16]. In recent years, digital and analog quantum simulators have emerged as a promising platform for the simulation of quantum phenomena. Quantum simulators have already been used to study phase transitions using the method of partition function zeros [17] and the Kibble-Zurek mechanism [18, 19]. In this paper, we present a protocol to implement the finite-size scaling method on a digital quantum simulator. We use the Quantum Restricted Boltzmann Machine (QRBM) algorithm to find the critical point of the Quantum Rabi model. This paper is organized as follows. In Section 2, we explain the theory of Quantum Rabi Model, Finite-Size Scaling and the Quantum Restricted Boltzmann Machine. In Section 3, we present our results obtained using the exact diagonalization method and QRBM. Finally in Section 4, we discuss our results and future prospects of studying quantum phase transitions on quantum devices. 2 THEORY 2.1 Quantum Rabi Model The QRM describes a two-level system interacting with a bosonic field mode. The Hamiltonian is [2], HRabi
Ω σ z + ω0 a† a − λσ x
a + a†  2 (1) where we’ve chosen Z = 1. Here, σz and σx are the Pauli Z and X matrices respectively, Ω is the energy separation between the two levels in the system, ω0 is the frequency of the bosonic mode and λ is the system-environment coupling strength. The parity operator Π
† eiπ(a a+|↑〉〈↑|) commutes with HRabi. So, HRabi√has a Z2 symmetry. This model has a critical point at g
2λ/ ω0 Ω
gc
1 in the limit Ω/ω0 → ∞ [2]. Ω/ω0 → ∞ is analogous to the thermodynamic limit for this case, and in experiments where Ω/ω0 has to be finite, we’ll observe finite-size effects like in any other phase transition [2]. For g < 1, the system is in the normal phase and the ground state is |ϕ0np (g)〉
S[rnp (g)]|0〉|↓〉 where S[x]
exp[x2 (a†2 − a2 )] and rnp (g)
− 14 ln(1 − g2 ). The rescaled ground state energy and photon number are eG (g)
ω0 ω0 † Ω 〈HRabi 〉
−ω0 /2 and nG (g)
Ω 〈a a〉
0 respectively. For g > 1, the system is in a superradiant phase and the ground state is two-fold degenerate, |ϕ0sp (g)〉
D[± αg ]S[rsp (g)]|0〉|↓± 〉 here rsp (g)
− 14 ln(1 − g−4 ) and D[α]
exp[α(a† − a)]. |↓± 〉 is the 2λα negative eigenvalue eigenstate of 2g1 2 σ z ± g2 Ωg σ x where   Ω αg
4g2 ω0 (g4 − 1). The rescaled ground state energy and photon number are eG (g)
ωΩ0 〈HRabi 〉
−ω0 (g2 + g−2 )/4 and nG (g)
ωΩ0 〈a† a〉
(g2 − g−2 )/4 respectively. As shown in Figures 1A,B, d2eG/dg2 is discontinuous at g = gc = 1, indicating a continuous phase transition and nG
ωΩ0 〈a† a〉 is an order parameter for this phase transition. In the normal phase, nG is zero whereas in the superradiant phase, Z2 symmetry is spontaneously broken and nG becomes non-zero. We can also write effective low-energy Hamiltonians in both the normal and the superradiant phases. For g < 1, HRabi can be reduced to the following effective Hamiltonian [2], Hnp
ω0 a† a − ω0 g2 Ω 2
a + a†  − . 4 2 (2) The system’s degrees of freedom have been removed by projecting to |↓〉〈↓|, since this is a low energy description. Similarly, for g > 1 the effective Hamiltonian can be written as [2], Hsp
ω0 a† a − ω0 Ω 2
a + a†  −
g2 + g−2 , 4g4 2 (3) where this time around the Hamiltonian has been projected along |↓± 〉〈↓± |. In Section 3, we’ll use Hnp and Hsp to find the critical point of the model. 2.2 Finite-Size Scaling The FSS method is widely used to determine the critical points and the critical exponents in phase transitions [1]. To demonstrate the method, consider that we have an infinite 2d system that undergoes a classical phase transition at a critical temperature T = Tc [9]. Suppose Q is a quantity that becomes singular at T = Tc with some power law behavior FIGURE 1 | Phase Transition in Quantum Rabi Model. (A) The rescaled ground state energy eG /ω0
〈HRabi 〉/Ω and (d2eG/dg2)/ω0 as functions of g. The discontinuity in (d2eG/dg2)/ω0 at g = gc = 1 indicates a countinuous phase transition. (B) The order parameter nG
ωΩ0 〈a† a〉 as a function of g. nG becomes non-zero when the Z2 symmetry is spontaneously broken at g > gc = 1. Q∞ (T) ~|T − Tc |−ω . (4) We can also think of this system as an infinite collection of infinite stripes, where the stripes are infinitely extended along one direction and stacked along the perpendicular direction. Now suppose there are only an N number of stripes. If N is finite, Q should be regular at T = Tc since finite systems cannot have non-analyticities at T ≠ 0. The singularity at T = T c should appear only when N → ∞. The finite size scaling hypothesis assumes the existence of a scaling function FQ such that QN (T) ≃ Q∞ (T)FQ
Nξ ∞ (T), (5) where QN is the observable Q for a system with N stripes and Q∞ corresponds to the system in the thermodynamic limit. ξ∞ is the correlation length for the infinite system. Eq. 5 is valid when N is large. The correlation length also diverges as a power law near the critical point, ξ ∞ (T) ~|T − Tc |−] . (6) Substituting Eqs 4, 6 in Eq. 5, QN (T) ≃ |T − Tc |−ω FQ (N|T − Tc |] ). (7) Since QN(T) should be regular at T = Tc, the scaling function should cancel the divergence due to |T − Tc|−ω. Therefore, the scaling function should be of the form FQ(x) ~ xω/] as x → 0. We should then have, QN (Tc ) ~ Nω/] . (8) If we define a function ΔQ(T; N, N′) such that ΔQ
T; N, N′
log
QN (T)QN′ (T) , log
NN′ (9) then the value of this function at T = Tc, ΔQ(Tc; N, N′) ≃ ω/] is independent of N and N′. Therefore, for three different values N, N′ and N′′, the curves ΔQ(T; N, N′) and ΔQ(T; N′, N′′) will intersect at the critical point T = Tc. This is how we can locate the critical point using the finite size scaling hypothesis. We can also find the critical exponents ω and ]. Noting from Eq. 4 that zQ∞ (T) ~|T − Tc |−(ω+1) . zT (10) Therefore, we should have ΔzQ/zT(Tc; N, N′) ≃ (ω + 1)/]. Define a new function Γω(T; N, N′) such that Γω
T; N, N′
ΔQ
T; N, N′ . ΔzQ/zT
T; N, N′ − ΔQ
T; N, N′ (11) The value of this function at the critical point Γω(Tc; N, N′) ≃ ω is independent of N and N′ and gives us the critical exponent ω. Then ] can be determined using ]≃ ω . ΔQ
Tc ; N, N′ (12) As we’ve already stated in the Introduction, this method cannot be used for the kinds of phase transitions we are interested in which occur at a finite system size. However, for such cases we can consider an extension of the approach discussed above [10–16]. In this extended approach, instead of truncating the system in the physical space, the system is truncated in the Hilbert space [16]. The FSS ansatz looks exactly the same except that N now represents the size of the set of basis states which spans the truncated Hilbert space [16]. Moreover, the temperature T will be replaced by the parameter g which is being tuned across the critical point. This approach has been shown by Kais and co-workers to work in the case of a particle in Yukawa potential [11, 13] and the calculation of electronic structure critical parameters for atomic and molecular systems [10, 12, 14–16]. 2.3 Quantum Restricted Boltzmann Machine Solving quantum many-body problems accurately has been a taxing numerical problem since the size of the wavefunction scales exponentially. The idea of taking advantage of the aspects of Machine Learning (ML) related to dimensionality reduction and feature extraction to capture the most relevant information came from the work by Carleo and Troyer [20], which introduced the idea of representing the many-body wavefunction in terms of an Artificial Neural Network (ANN) to solve for the ground states and time evolution for spin models, with a Restricted Boltzmann Machine (RBM) as the chosen architecture for this ANN. More recently, the critical behavior of the quantum Hall plateau transition based on wavefunctions has been studied in a 2D disordered electron system with the usage of a Convolutional Neural Network (CNN) [21]. However, we focus on using an RBM architecture in this work. An RBM consists of a visible layer and a hidden layer with each neuron in the visible layer connected to all neurons in the hidden layer but the neurons within a layer are not connected to each other. The quantum state is ψ expanded in the basis |x〉:  ψ 〉
 ψ (x)|x〉 (13) The Neural Network Quantum State (NQS) [20] describes the wavefunction ψ(x) to be written as ψ(x; θ), where θ represents the parameters of the RBM. ψ(x; θ) is now written in terms of the probability distribution that is obtained from the RBM as follows: 1 ψ (x; θ) ∝  e2 i ai σ zi +j bj hj +ij wij σ zi hj (14) {h} where, σ zi is the Pauli z operator at ith site, σ zi and hj take values { + 1, − 1}, θ = {ai, bj, wij} are the trainable bias and weight parameters of the RBM. Using stochastic optimization, the energy E(θ) is minimized. This work was extended to obtain the ground states of the Bose-Hubbard model [22] and for the application of quantum state tomography [23]. With the rapid developments in the domains of ML and Quantum Computing (QC), the appetite for integrating ideas in both of these areas has been growing considerably. The last decade has seen a surge in the application of classical ML for quantum matter, wherein these methods have been adopted to benchmark, estimate and study the properties of quantum matter [24–27], with recently showing provable classification efficiency in classifying quantum states of matter [28]. RBM based ansätzes have been shown to capture entanglement transitions [29] and using an RBM with local sparse connectivity achieves higher accuracy compared to its dense counterpart when applied to disordered quantum Ising chains [30]. The protocols and algorithms related to ML implementable on a quantum system so-called Quantum machine Learning [31] is expected to have the potential of changing the course of fundamental scientific research [32] along with industrial pursuit. In lieu of today’s Noisy Intermediate Scale Quantum (NISQ) devices, the ideas which utilize both classical and quantum resources, such that the part of the problem which has an exponential scaling is implemented on the quantum platform while the rest are dealt with classically, are being carefully investigated for various applications. Such algorithms are known as classical-quantum hybrid algorithms. In the work by Xia and Kais [33], a modified RBM with three layers was introduced, the third layer to account for the sign of the wavefunction, to solve for the ground state energies of molecules (see Figure 2). Now, the parametrized wavefunction ψ(x; θ) is written as a function of P(x) along with a sign function s(x): {h} P(x)
x′ a σz + i i i e e {h} b h + j j j a σ z′ + i i i ij b h + j j j wij σ zi hj ij wij σ z′ h i j ⎢ ⎝c +  di σ i ⎞ ⎠⎤⎥⎦ ⎣⎛ s(x)
tanh⎡ (15) (16) i The wavefunction ansatz in terms of the RBM can be expressed as [33]:   ψ 〉
 P(x)s(x)|x〉 (17) x A quantum circuit comprising of a single-qubit (Ry) and multi-qubit y-rotation gates (C1 − C2 − Ry) is employed, to sample the Gibbs distribution. The utilization of Ry gates caters to the bias parameter of visible and hidden layers part of the distribution, while C1 − C2 − Ry gates tend to the weights part of the distribution. In the work by Sureshbabu et al. [34], the implementation of such a circuit on IBM-Q devices was shown, wherein a new ancillary qubit is introduced to store the value corresponding to every C1 − C2 − Ry gate (Figure 3). The term n denotes the number of visible qubits and m denotes the number of hidden units. In this formalism, the number of ancillary qubits required is n × m. Starting all the qubits from a |0〉, the Ry and C1 − C2 − Ry rotations are performed, and a measurement is performed on all the qubits. If all the ancillary qubits are in |1〉, then the sampling is deemed successful and the states corresponding to the first m + n qubits provide the distribution P(x). The joint probability distribution defined over the parameters of the circuit θ = {a, b, w} and a set of y = {σz, h} is given by: P
y, θ
a σz+ i i i e {y} b h + j j j a σ z′ + i i i e ij wij σ zi hj b h + j j j ij wij σ z′ h i j (18) The probability of successful sampling can be improved by rewriting the distribution P(y, θ) as Q(y, θ) and setting |w | k
max(1, 2ij )[33, 35]: 1 Q
y, θ
ek  a σz + i i i 1 k {y} e  b h + j j j a σ z′ + i i i ij b h + j j j wij σ zi hj  ij wij σ z′ h i j (19) Firstly, the QRBM is implemented classically, i.e., the quantum circuit is simulated on a classical computer. This execution caters to the ideal results that can be obtained through the QRBM algorithm. Then, the quantum circuit is implemented on the Digital Quantum Simulator, the qasm FIGURE 2 | Restricted Boltzmann Machine architecture. The first layer is the visible layer with bias parameters denoted by ai. The second layer is the hidden layer with bias parameters denoted by bj. The third layer is the sign layer with bias parameters denoted by c. The weights associated with the connections between the visible neurons and the hidden neurons are designated by wij. The weights associated with the connections between the visible neurons and the neuron of the sign layer are designated by di. FIGURE 3 | The quantum circuit to sample the Gibbs distribution. n is the number of qubits belonging to the visible layer and m is the number of qubits belonging to the hidden layer. There are m × n ancillary qubits. simulation backend. This simulator is part of the highperformance simulators from IBM-Q. The circuit is realized using IBM’s Quantum Information Software Toolkit titled Qiskit [36]. Though no noise model was utilized, as a result of finite sampling, statistical fluctuations in the values of probabilities in observing the circuit in the measurement basis, are present in the obtained results. Having obtained the distribution Q(y, θ), the probabilities are raised to the power of k, to get P(y, θ). Following this, the sign function is computed classically, thereby calculating |ψ〉. Then, the expectation value for the Hamiltonian H [〈Ψ|H|Ψ〉] is computed to get the energy, which is minimized using gradient descent to obtain the ground state eigenenergy of H. The resource requirements demanded by this algorithm are quadratic. The number of qubits required are (m + n) to encode the visible and hidden nodes, and (m × n) to account for the ancillary qubits. Hence, the number of qubits scales as O(mn). The number of Ry gates required are (m + n) and the number of C1 − C2 − Ry gates required are (m × n). In addition, each C1 − C2 − Ry gate requires 6n X-gates to account for all the states spanned by the control qubits. Therefore, the number of gates required also scales as O(mn). Obtaining the ground states or minimum eigenvalues of a given matrix using exact diagonalization has a complexity of ≈ j3, with j being the dimension of the column space for the given matrix [37]. 3 RESULTS 3.1 Exact Diagonalization In this section, we demonstrate the calculation of the critical point of the Quantum Rabi model using the Finite-Size Scaling method. As discussed before, the phase transition in QRM occurs only in the limit Ω/ω0 → ∞. This limit is not straightforward to implement in HRabi given in Eq. 1. Instead, we have considered the effective low-energy Hamiltonians Hnp and Hsp given in Eqs 2, 3 respectively. In Hnp and Hsp, Ω is involved only in a constant term which can be removed from the Hamiltonians and the limit Ω/ω0 → ∞ can then be easily imposed. In Hnp and Hsp, the degrees of freedom of the two-level system have been traced out and the only degrees of freedom we have are those of the bosonic mode. Let’s first consider the normal phase Hamiltonian Hnp. The Hilbert space for this Hamiltonian is spanned by the familiar harmonic oscillator number states {|0〉, |1〉, |2〉, . . .}. We can truncate the full Hilbert space to an N-dimensional Hilbert space spanned by {|0〉, |1〉, . . . , |N − 1〉} to apply the finite-size scaling analysis. In this restricted Hilbert (N) space, of H np can √  be found by using a|m〉
√ the matrix form † m|m − 1〉 and a |m〉
m + 1|m + 1〉. Once we have the matrix form, we can then use the exact diagonalization (N) method to find the ground state of H(N) np with energy Enp . Consider the scaling law for the ground state energy in the vicinity of the critical point g = gc, E
g ~|g − gc |α . (20) Here E is the ground state energy. We slightly modify the formula in Eq. 9 to take into account the difference in the signs of the exponents in Eqs 4, 20. The new formula with Q = E is, (N′) logE(npN)
gEnp
g ΔHnp
g; N, N′
, log
N′N (21) We plot the curves ΔHnp (g; N, N + 2) for N = 8, 10, . . ., 30 in Figure 4A. We then plot the intersection points g(N) np of the curves ΔHnp (g; N − 4, N − 2) and ΔHnp (g; N − 2, N) as a function of N as shown in Figure 4B. To find the limit of g(N) np as N → ∞, we used the Bulirsch-Stoer algorithm [44, 45]. The limit was (np)
0.999996. calculated to be g(N) np → 0.999996. So gc In a similar way, we then consider Hsp. The curves ΔHsp (g; N, N + 2) are plotted in Figure 4C for N = 8, 10, . . ., 30 and the intersection points g(N) sp are plotted in Figure 4D as a function of N. In this case, the extrapolation to N → ∞ gives the
0.999987. Both the calculated values of gc(np) critical value g(sp) c and g(sp) are very close to the exact value gc = 1. c 3.2 Quantum Restricted Boltzmann Machine Now we illustrate the implementation of the FSS method using the QRBM algorithm. The results are shown in Figure 5. Figure 5A,C show the results for Hnp and Hsp using the classical implementation of the algorithm respectively. Whereas, Figure 5B,D correspond to the results for Hnp and Hsp when the algorithm is implemented using the qasm simulator from IBM-Q respectively. The QRBM algorithm is run for N = 8, 10, 12, 14, 16. For the case of N = 8, the number of qubits associated with the visible nodes equals 3, the number of qubits associated with the hidden nodes equals 3, and 9 ancillary qubits were used. The quantum circuit consists of 6 Ry gates associated with the bias parameters, 9 C1 − C2 − Ry gates associated with the weights. Since, each C1 − C2 − Ry gate requires 6 X-gates, a total of 54 Xgates were used. For the case of N = 10,..,16, the number of qubits associated with the visible nodes equals 4, the number of qubits associated with the hidden nodes equals 4, and 16 ancillary qubits were used. The quantum circuit consists of 8 Ry gates associated with the bias parameters, 16 C1 − C2 − Ry gates associated with the weights. Since, each C1 − C2 − Ry gate requires 6 X-gates, a total of 96 X-gates were used. Starting from random initialization, all parameters are updated via gradient descent. A learning rate of 0.01 was chosen and the algorithm is run for around 30,000 iterations. In order to assist with the convergence to the minimum eigenenergies, warm starting is employed. The method of warm starting is essentially initializing the parameters of the current point with the parameters of a previously converged point of calculation, which helps in avoiding the convergence to a local minima. The black curves plotted in the insets in Figure 5 represent the deviation of the QRBM results (black dashed curves) from the exact diagonalization results (blue solid curves). They were calculated using the average of the quantity |Δ(ED)(g) − Δ(QRBM)(g)/Δ(ED)(g)|× 100 over all the four curves. An enlarged version of the error plots is shown in Figure 6 can be found in the Supporting Information section. For each case the overall error close to g = 1.000 is not more than ~ 5% which implies convergence to the right result. Moreover, for the case of Hsp, we notice that the error is very small for the classical implementation i.e., ~ < 1% throughout the range of the graph. The critical point using Hnp was found to be gc(np)
1.008 for both the classical and qasm implementations. Similarly, the FIGURE 4 | Finite-Size Scaling for Quantum Rabi model. We used N = 8, 10, . . ., 32. (A) Graphs of ΔHnp (g; 8, 10), ΔHnp (g; 10, 12), . . . , ΔHnp (g; 30, 32) as a function (np) (N) (N) (N)
0.999996. (C) of g. (B) Intersection points g(N) np where ΔHnp (gnp ; N − 4, N − 2)
ΔHnp (gnp ; N − 2, N), as a function of 1/N. As N → ∞, gnp → 0.999996. So, gc (N) (N) where Δ (g ; N − 4, N − 2)
Δ (g Graphs of ΔHsp (g; 8, 10), ΔHsp (g; 10, 12), . . . , ΔHsp (g; 30, 32) as a function of g. (D) Intersection points g(N) Hsp Hsp sp sp sp ; N − 2, N), as a (sp)
0.999987. function of 1/N. As N → ∞, g(N) sp → 0.999987. So, gc FIGURE 5 | QRBM Implementation of FSS for QRM. The light blue line represents results obtained from exact diagonalization and dashed black line represents QRBM results. (A) Classical implementation of QRBM corresponding to normal phase, graphs of ΔHnp (g; 8, 10), ΔHnp (g; 10, 12), . . . , ΔHnp (g; 14, 16) as a function of g. (B) QRBM in both the cases implemented on qasm simulator corresponding to normal phase, graphs of ΔHnp (g; 8, 10), ΔHnp (g; 10, 12), . . . , ΔHnp (g; 14, 16) as a function of g. The g(np) c is calculated to be 1.008. (C) Classical implementation of QRBM corresponding to superradiant phase, graphs of ΔHsp (g; 8, 10), ΔHsp (g; 10, 12), . . . , ΔHsp (g; 14, 16) as a function of g. (D) QRBM implemented on qasm simulator corresponding to superradiant phase, graphs of ΔHsp (g; 8, 10), ΔHsp (g; 10, 12), . . . , ΔHsp (g; 14, 16) as a function in both the cases is calculated to be 0.996. The inset plots display the mean percentage error between the exact diagonalization results and QRBM results. of g. The g(sp) c 4 DISCUSSION AND OUTLOOK FIGURE 6 | Error plots from the insets of Figure 5. (A–D) are the insets shown in Figures 5A–D respectively. critical point for the case of Hsp was found to be g(sp)
0.996 for c both the classical and qasm implementations. Here we notice that although, the convergence for the data obtained from both the classical and qasm implementations turns out to be the same for both Hnp and Hsp, such a perfect match appears to be somewhat coincidental. Here, again the Bulirsch-Stoer algorithm [44, 45] which sets the criteria used to deduce these convergence results. The convergence plots are shown in Figure 7 have been added to the Supporting Information section. In this paper, we have used the Finite-Size Scaling in Hilbert Space approach to calculate the critical point of the Quantum Rabi Model. We used the low-energy effective Hamiltonians for both the normal and superradiant phases respectively to show that the critical point is g c ≈ 1. The original FSS approach in which the truncation is done in the physical space has been widely used to calculate critical points and critical exponents since its inception. However, that approach was not applicable to Quantum Phase Transitions which occur at a finite system size. With the rise in interest in QPTs occurring in these finite size systems, our approach provides a natural extension of the original FSS method to study such phase transitions. To our knowledge, this is the first time this approach has been used to study a QPT in a light-matter interaction system. We have also provided a recipe for the implementation of this method on a universal quantum computer using the Quantum Restricted Boltzmann Machine algorithm. It was shown that results obtained from the classical gate simulation match those obtained from the IBM-Q’s qasm simulator. Such an implementation scales quadratically while the exact diagonalization scales cubically in the best case and exponentially in the worst case. Looking forward, we are interested in applying this approach to other QPTs such as the QPT in anisotropic QRM. We would also like to use our method to calculate the critical exponents in addition to the critical points in these phase transitions. It would also be interesting to see if this approach can be used to predict any new phase transition for some other non-integrable model. FIGURE 7 | Convergence diagrams for results in Figure 5. (A–D) correspond to convergence results for data in Figures 5A–D respectively. The same procedure was used as the one shown in Figures 4B,D. Another very promising research direction is to implement the FSS method for phase transitions in classically intractable manybody models such as exotic electronic and magnetic systems. These include general quantum materials, for example where Coulomb potential leads to a gapped spectrum in energy, including in direct band-gap semiconductors in the thermodynamic limit. Conventionally speaking, it might be necessary to resort to the original finite-size scaling in the physical space approach for these systems since they exhibit criticality only in the limit N → ∞. However, the ground state of an appropriately truncated Hamiltonian could be deduced using the QRBM algorithm as shown in the paper towards efficient implementation on a digital quantum simulator. A simile can also be drawn between a many-body bulk gap separating a continuum of excited states from the ground state manifold to the gapped Rabi model discussed in this paper. Such an approach can be useful in emergent topological systems, such as in Weyl semimetals, 1-D Kitaev spin chains, quantum spin liquids, and others, on which there is a tremendous explosion of interest [38–43]. Topological phase transitions are devoid of any conventional order parameter and a quantum solution deriving from the approach outlined in this paper can help us bypass resource and scaling limitations of DMRG and exact REFERENCES 1. Goldenfeld N. Lectures on Phase Transitions and the Renormalization Group; Frontiers in Physics. Addison-Wesley Publishing Company (1992). 2. Hwang M-J, Puebla R, Plenio MB. Quantum Phase Transition and Universal Dynamics in the Rabi Model. Phys Rev Lett (2015) 115:180404. doi:10.1103/ physrevlett.115.180404 3. Hwang M-J, Plenio MB. Quantum Phase Transition in the Finite JaynesCummings Lattice Systems. Phys Rev Lett (2016) 117:123602. doi:10.1103/ physrevlett.117.123602 4. Liu M, Chesi S, Ying Z-J, Chen X, Luo H-G, Lin H-Q. Universal Scaling and Critical Exponents of the Anisotropic Quantum Rabi Model. Phys Rev Lett (2017) 119:220601. doi:10.1103/physrevlett.119.220601 5. Hwang M-J, Rabl P, Plenio MB. Dissipative Phase Transition in the Open Quantum Rabi Model. Phys Rev A (2018) 97:013825. doi:10.1103/physreva.97. 013825 6. Xie Y-F, Chen X-Y, Dong X-F, Chen Q-H. First-order and Continuous Quantum Phase Transitions in the Anisotropic Quantum Rabi-Stark Model. Phys Rev A (2020) 101:053803. doi:10.1103/physreva.101.053803 7. Cai M-L, Liu Z-D, Zhao W-D, Wu Y-K, Mei Q-X, Jiang Y, et al. Observation of a Quantum Phase Transition in the Quantum Rabi Model with a Single Trapped Ion. Nat Commun (2021) 12:1126. doi:10.1038/s41467-021-21425-8 8. Forn-Díaz P, Lamata L, Rico E, Kono J, Solano E. Ultrastrong Coupling Regimes of Light-Matter Interaction. Rev Mod Phys (2019) 91:025005. doi:10. 1103/revmodphys.91.025005 9. Derrida B, Seze LD. In Application of the Phenomenological Renormalization to Percolation and Lattice Animals in Dimension 2. In: JL CARDY, editor. Current Physics–Sources and Comments, 2. Elsevier (1988). p. 275–83. doi:10. 1016/b978-0-444-87109-1.50024-8 10. Neirotti JP, Serra P, Kais S. Electronic Structure Critical Parameters from Finite-Size Scaling. Phys Rev Lett (1997) 79:3142–5. doi:10.1103/physrevlett. 79.3142 11. Serra P, Neirotti JP, Kais S. Finite-size Scaling Approach for the Schrödinger Equation. Phys Rev A (1998) 57:R1481–R1484. doi:10.1103/physreva.57.r1481 12. Serra P, Neirotti JP, Kais S. Electronic Structure Critical Parameters for the Lithium Isoelectronic Series. Phys Rev Lett (1998) 80:5293–6. doi:10.1103/ physrevlett.80.5293 diagonalization approaches to calculate the critical point and the critical exponents. DATA AVAILABILITY STATEMENT The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author. AUTHOR CONTRIBUTIONS SK and AB designed the research problem, BK and SS performed the calculations, all authors discussed the results and wrote the paper. ACKNOWLEDGMENTS We acknowledge funding by the US Department of Energy, Office of Science, National Quantum Information Science Research Centers, Quantum Science Center. 13. Serra P, Neirotti JP, Kais S. Finite Size Scaling in Quantum Mechanics. J Phys Chem A (1998) 102:9518–22. doi:10.1021/jp9820572 14. Kais S, Serra P. Quantum Critical Phenomena and Stability of Atomic and Molecular Ions. Int Rev Phys Chem (2000) 19:97–121. doi:10.1080/ 014423500229873 15. Shi Q, Kais S. Finite Size Scaling for Critical Parameters of Simple Diatomic Molecules. Mol Phys (2000) 98:1485–93. doi:10.1080/00268970009483354 16. Kais S, Serra P. Finite-Size Scaling for Atomic and Molecular Systems. Adv Chem Phys (2003) 125:1–99. doi:10.1002/0471428027.ch1 17. Francis A, Zhu D, Huerta Alderete C, Johri S, Xiao X, Freericks JK, et al. Manybody Thermodynamics on Quantum Computers via Partition Function Zeros. Sci Adv (2021) 7:eabf2447. doi:10.1126/sciadv.abf2447 18. Keesling A, Omran A, Levine H, Bernien H, Pichler H, Choi S, et al. Quantum Kibble–Zurek Mechanism and Critical Dynamics on a Programmable Rydberg Simulator. Nature (2019) 568:207–11. doi:10. 1038/s41586-019-1070-1 19. Dupont M, Moore JE. Quantum Criticality Using a Superconducting Quantum Processor (2021). arXiv:2109.10909. 20. Carleo G, Troyer M. Solving the Quantum many-body Problem with Artificial Neural Networks. Science (2017) 355:602–6. doi:10.1126/science.aag2302 21. Li Z, Luo M, Wan X. Extracting Critical Exponents by Finite-Size Scaling with Convolutional Neural Networks. Phys Rev B (2019) 99:075418. doi:10.1103/ physrevb.99.075418 22. Saito H. Solving the Bose-Hubbard Model with Machine Learning. J Phys Soc Jpn (2017) 86:093001. doi:10.7566/jpsj.86.093001 23. Torlai G, Mazzola G, Carrasquilla J, Troyer M, Melko R, Carleo G. Neuralnetwork Quantum State Tomography. Nat Phys (2018) 14:447–50. doi:10. 1038/s41567-018-0048-5 24. Carrasquilla J, Melko RG. Machine Learning Phases of Matter. Nat Phys (2017) 13:431–4. doi:10.1038/nphys4035 25. Deng D-L, Li X, Das Sarma S. Machine Learning Topological States. Phys Rev B (2017) 96:195145. doi:10.1103/physrevb.96.195145 26. Butler KT, Davies DW, Cartwright H, Isayev O, Walsh A. Machine Learning for Molecular and Materials Science. Nature (2018) 559:547–55. doi:10.1103/ physrevb.96.195145 27. Carleo G, Cirac I, Cranmer K, Daudet L, Schuld M, Tishby N, et al. Machine Learning and the Physical Sciences. Rev Mod Phys (2019) 91:045002. doi:10. 1103/revmodphys.91.045002 28. Huang H-Y, Kueng R, Torlai G, Albert VV, Preskill J. Provably Efficient Machine Learning for Quantum many-body Problems (2021). arXiv preprint arXiv:2106.12627. 29. Medina R, Vasseur R, Serbyn M. Entanglement Transitions from Restricted Boltzmann Machines. Phys Rev B (2021) 104:104205. doi:10.1103/physrevb. 104.104205 30. Pilati S, Pieri P. Simulating Disordered Quantum Ising Chains via Dense and Sparse Restricted Boltzmann Machines. Phys Rev E (2020) 101:063308. doi:10. 1103/PhysRevE.101.063308 31. Biamonte J, Wittek P, Pancotti N, Rebentrost P, Wiebe N, Lloyd S. Quantum Machine Learning. Nature (2017) 549:195–202. doi:10.1038/ nature23474 32. Sajjan M, Li J, Selvarajan R, Sureshbabu SH, Kale SS, Gupta R, et al. Quantum Computing Enhanced Machine Learning for Physico-Chemical Applications (2021). arXiv preprint arXiv:2111.00851. 33. Xia R, Kais S. Quantum Machine Learning for Electronic Structure Calculations. Nat Commun (2018) 9:4195–6. doi:10.1038/s41467-01806598-z 34. Sureshbabu SH, Sajjan M, Oh S, Kais S. Implementation of Quantum Machine Learning for Electronic Structure Calculations of Periodic Systems on Quantum Computing Devices. J Chem Inf Model (2021) 61(6):2667–2674. doi:10.1021/acs.jcim.1c00294 35. Sajjan M, Sureshbabu SH, Kais S. Quantum Machine-Learning for Eigenstate Filtration in Two-Dimensional Materials. J Am Chem Soc (2021) 143: 18426–45. doi:10.1021/jacs.1c06246 36. Aleksandrowicz G, Alexander T, Barkoutsos P, Bello L, Ben-Haim Y, Bucher D, et al. Qiskit: An Open-Source Framework for Quantum Computing. Zenodo (2019) 16:11. doi:10.5281/zenodo.2562111 37. Harris CR, Millman KJ, van der Walt S, Gommers R, Virtanen P, Cournapeau D, et al. Array Programming with NumPy. Nature (2020) 585:357–62. doi:10.1038/s41586-020-2649-2 38. Satzinger KJ, Liu YJ, Smith A, Knapp C, Newman M, Jones C, et al. Realizing Topologically Ordered States on a Quantum Processor. Science (2021) 374: 1237–41. doi:10.1126/science.abi8378 39. Kitaev A. Anyons in an Exactly Solved Model and beyond. Ann Phys (2006) 321:2–111. January Special Issue. doi:10.1016/j.aop.2005.10.005 40. Liu Y-J, Shtengel K, Smith A, Pollmann F. Methods for Simulating String-Net States and Anyons on a Digital Quantum Computer. arXiv:2110.02020 2021. 41. Xiao X, Freericks JK, Kemper AF. Determining Quantum Phase Diagrams of Topological Kitaev-Inspired Models on NISQ Quantum Hardware. Quantum (2021) 5:553. doi:10.22331/q-2021-09-28-553 42. Wen X-G. Colloquium: Zoo of Quantum-Topological Phases of Matter. Rev Mod Phys (2017) 89:041004. doi:10.1103/revmodphys.89.041004 43. Hasan MZ, Kane CL. Colloquium: Topological Insulators. Rev Mod Phys (2010) 82:3045–67. doi:10.1103/revmodphys.82.3045 44. Bulirsch R, Stoer J. Numerical Treatment of Ordinary Differential Equations by Extrapolation Methods. Numer Math (1966) 8:1–13. doi:10.1007/ bf02165234 45. Henkel M, Schutz G. Finite-lattice Extrapolation Algorithms. J Phys A: Math Gen (1988) 21:2617–33. doi:10.1088/0305-4470/21/11/019 APPENDIX A: using the following rules Bulirsch-Stoer Algorithm For hN = 1/N where N = 0, 1, 2, . . ., the Bulirsch-Stoer algorithm can be used to find the limit of a function T(hN) as N → ∞44,45. For demonstration, consider that we only have T(hN) for N = 0, 1, 2, 3, then the following rows are computed successively, T(−1N)
0 T(0N)
T(hN ) 0 1 2 3 T (1) 0 T 1(0) T (2) 0 T 1(1) T (0) 2 T 1(2) T (1) 2 T 3(0) T 0(3) (23) N) N+1) T(m≥1
T(m−1 ω T (0) 0 (22) −1 N) hN T(N+1) − T(m−1 N+1) N) +
T(m−1 − T(m−1   1 − (m−1  − 1 N+1) (N+1) hN+m Tm−1 − Tm−2 (24) where ω is a free parameter determined by minimizing (0) (i+1) − T(i) ε(i) m |. The final answer is T3 . m
|Tm