Phase diagram reconstruction of the Bose-Hubbard model with a Restricted Boltzmann Machine wavefunction
Vladimir Vargas-Calderón, Herbert Vinck-Posada, Fabio A. González
PPhase diagram reconstruction of the Bose-Hubbard model with a RestrictedBoltzmann Machine wavefunction
Vladimir Vargas-Calder´on ∗ and Herbert Vinck-Posada Grupo de Superconductividad y Nanotecnolog´ıa, Departamento de F´ısica,Universidad Nacional de Colombia, Bogot´a, Colombia
Fabio A. Gonz´alez
Computing Systems and Industrial Eng. Department,Universidad Nacional de Colombia, Bogot´a, Colombia (Dated: April 28, 2020)Recently, the use of neural quantum states for describing the ground state of many- and few-body problems has been gaining popularity because of their high expressivity and ability to handleintractably large Hilbert spaces. In particular, methods based on variational Monte Carlo haveproven to be successful in describing the physics of bosonic systems such as the Bose-Hubbard model.However, this technique has not been systematically tested on the parameter space of the Bose-Hubbard model, particularly at the boundary between the Mott insulator and superfluid phases. Inthis work, we evaluate the capabilities of variational Monte Carlo with a trial wavefunction given bya Restricted Boltzmann Machine to reproduce the quantum ground state of the Bose-Hubbard modelon several points of its parameter space. To benchmark the technique, we compare its results tothe ground state found through exact diagonalization for small one-dimensional chains. In general,we find that the learned ground state correctly estimates many observables, reproducing to a highdegree the phase diagram for the first Mott lobe and part of the second one. However, we findthat the technique is challenged whenever the system transitions between excitation manifolds, asthe ground state is not learned correctly at these boundaries. Nonetheless, we propose a methodto discard noisy probabilities learned in the ground state, which improves the quality of the resultsproduced by the method.
I. INTRODUCTION
The dimension of a Hilbert space that describes thepossible states of a many-body system scales exponen-tially with the number of one-body states and withthe number of particles. In most cases, this featurecomes as an important practical difficulty for physi-cists to study many-body problems, because approximatetechniques have to be implemented in order to performsimulations (e.g., dynamical mean-field theory , densitymatrix renormalization group (DMRG) ). One of themost studied quantum objects is the ground state of amany-body system for a number of reasons: particlestend to occupy the lowest energy states first, as dic-tated by the aufbau principle; also excited states inheritthe ground state structure; and most prominently, it isan object that encodes changes in observables that ex-hibit quantum phase transitions. Recently, the proposalof variational wavefunctions with neural networks andmachine learning-inspired wavefunctions for the groundstate and its optimization through variational MonteCarlo (VMC) has proven to be successful in approxi-mating with high fidelity the ground state of several con-densed matter many- and few-body systems. Some exam-ples are atomic and molecular systems , the transverse-field Ising model with quenching , the Heisenbergmodel and its anti-ferromagnetic version , the quan-tum harmonic oscillator in electric field , the Hubbardmodel , the J - J Heisenberg model and its antiferro-magnetic version , the XXZ model , and the dipolar Bose-Hubbard system .The success in approximating the ground state comesfrom two sources. The first one is technical, where thehigh expressivity and capabilities of neural networks toapproximate arbitrary functions are exploited, as wellas the thoroughly studied optimization methods in ma-chine learning . The second one is physical: despite theHilbert space being exponentially large with respect toone-body states and number of particles, only a small setof those states are needed to describe the ground state .Here the definition of small can vary, as it will be seen inthis work.Even though reported results indicate that the quan-tum ground state can be represented through neural net-work wavefunctions, most of these studies focus on theconvergence of the energy, and few of them use the tech-nique to characterize quantum correlations near a quan-tum phase transition to study the phase space of a Hamil-tonian comprehensively. Indeed, the convergence of theenergy comes quicker than the convergence of the state,as any first-order error in the variational state leads toonly a second-order error in its corresponding energy .Therefore, focusing on energy convergence might be mis-leading when asserting the power of VMC with neu-ral network trial wavefunctions. In contrast to previ-ous works, we intensively and systematically test the ca-pabilities of VMC to reconstruct, through a RestrictedBoltzmann Machine (RBM) wavefunction, the quantumground state of a one-dimensional Bose-Hubbard systemthroughout the superfluid (SF) and Mott insulator (MI) a r X i v : . [ c ond - m a t . d i s - nn ] A p r phases. We chose the Bose-Hubbard model because ofthe known numerical difficulties near the MI-SF bound-ary, and also because the problem involves bosons, whoseHilbert space is larger than fermions, thus increasing thecomputational difficulty in the study. Consequently, wecontrast exact diagonalization solutions with the onesprovided by the VMC-RBM method in systems of 5 and8 sites for a fine mesh in the Bose-Hubbard parame-ter space that exhibit the first and part of the secondMott lobes. It is only due to the sweeping of the Bose-Hubbard parameter space that we are able to see thatpoints near the MI-SF boundary are challenging for theVMC-RBM method. However, we also find that someof the differences between the learned ground state andthe ground state found through exact diagonalization aredue to noise that can be cleaned. Therefore, in this work,we also contribute with a state cleaning technique thatremoves noise that we identify comes from over-fitting.Nonetheless, we point out that there have been previ-ous efforts to examine the learning capability of differ-ent neural network wavefunctions for the ground state ofthe Bose-Hubbard Hamiltonian. However, as mentionedbefore, these studies do not extensively test this recentVMC framework in the Bose-Hubbard phase space. Forinstance, a permutation symmetric RBM implemented inNetKet has been used to study the energy and parti-cle density convergence at two points in the SF and MIphases, finding a difficulty in the convergence to the nu-merically exact particle density in the SF phase . How-ever, a qualitatively good location of the boundary be-tween the SF and MI phases for the region enclosing thefirst two Mott lobes was achieved using the particle den-sity as a discriminator . RBMs have also been used toshow how for a fixed number of bosons in the system, thelearned ground state is able to replicate the numericallyexact order parameter of the quantum phase transition .However, fixing the number of bosons greatly reduces thesize of the Hilbert space that is to be sampled with VMC,implying that a more complete sampling of the basisstates is possible, making the technique more exact. Also,the Bose-Hubbard model has been considered as a toymodel for trying full-forward neural networks (FFNNs)and convolutional neural networks (CNNs) wavefunc-tions. For instance, the ground state of one and two-dimensional finite lattices with parabolic confinement po-tential for a fixed number of bosons was approximatedthrough an FFNN ansatz . In the case of no confine-ment and periodic boundary conditions (which introducedisplacement symmetry), both FFNNs and RBMs wereused in large 1D lattices of up to 40 sites . CNNs havealso been introduced and compared with FFNN wave-functions to approximate this ground state . A CNN isalso proposed as a means to introduce the ground stateat finite temperature in a 1D lattice for larger on-siteboson interaction than hopping interaction . Despitethe important results derived from these works, they donot focus on the ability of the neural network ansatz toreproduce the ground state near the quantum phase tran- sition boundary. Nevertheless, unsupervised learning hasbeen previously used to classify states as belonging to theMott insulator or SF phases in the Bose-Hubbard modelwith previously obtained states or physical quantities, forseveral values of the order parameters .We believe that the ability of VMC technique to repro-duce the ground state near the quantum phase transitionhas to be extensively tested, which is why we pay close at-tention to the description of the variational ground statenear the SF-MI phase for a small number of sites to com-pare it with the numerically exact ground state. This pa-per is organized as follows. Section II exposes the mainfeatures of the Bose-Hubbard Hamiltonian and describesthe SF-MI transition in the system. In this section, anoverview of the VMC technique is given, discussing thedifficulties that must be faced and overcome in order tolearn the ground state of the Bose-Hubbard model. Insection III, the main results are given, including energyconvergence, the overlap of the approximated and exactground states, the reproduction of the phase diagram viatwo different order parameters, as well as tomographiesindicating relevant Fock states for the ground state. Fi-nally, in section IV conclusions of this work are given. II. MODEL AND VARIATIONAL MONTECARLOA. Bose-Hubbard Hamiltonian
The Bose-Hubbard Hamiltonian describes the inter-actions between bosons that can occupy sites in a d -dimensional lattice. These interactions are characterizedby a hopping energy t , an on-site interaction U and achemical potential µ , so that the grand-canonical Hamil-tonian reads ( (cid:126) = 1) ˆ H = − t (cid:88) (cid:104) ij (cid:105) (ˆ a † i ˆ a j + H.c.) + U (cid:88) i ˆ n i (ˆ n i − − µ (cid:88) i ˆ n i , (1)where ˆ a i is the annihilation operator at site i , andˆ n i = ˆ a † i ˆ a i is the number operator. The notation (cid:104) ij (cid:105) indicates that the sum runs over pairs of neighbor sitesin the lattice. The Bose-Hubbard model is able to re-produce experimental results in Josephson-junction net-works and in lattices of ultra-cold atoms . Thelatter offers precise control of the lattice parameters .Theoretically, a lot of attention has been devoted bothto understand the quantum phases of the system (theones that arise from eq. (1) and from the disordered orextended Bose-Hubbard model with longer range interac-tions) , as well to calculate the quantum phase transi-tion boundaries , whose precision has improved overthe years with better calculation techniques and comput-ing power, revealing features such as the re-entrance phe-nomenon , where for particular values of the chemicalpotential, the system switches between the MI phase tothe SF phase, and back to the MI phase before definitelyentering the SF phase after an increase of t/U .For simplicity, we will restrict our analyses to the d = 1 case, where only two quantum phases are possi-ble. When the on-site interaction energy is much largerthan the hopping energy, the latter becomes negligible,and the Hamiltonian is written as the sum of indepen-dent Hamiltonians for each site U ˆ n i (ˆ n i − − µ ˆ n i , whichcan be immediately diagonalized by the number basis.The corresponding eigenenergies are U n i ( n i − − µn i ,which reach minimum values for fixed µ and U at n i = max { , (cid:100) µ/U (cid:101)} (note that all sites are equiva-lent). Moreover, in this regime, the expected varianceof the local number operator is 0, i.e. (cid:104) ˆ n i (cid:105) − (cid:104) ˆ n i (cid:105) = 0.This regime characterizes the MI phase. On the otherhand, when the hopping energy is much larger thanthe on-site interaction energy, the latter becomes neg-ligible. Thus, the Hamiltonian can also be written asthe sum of independent Hamiltonians, but in momen-tum space, where ˜ a k = N − / (cid:80) Nj =1 ˆ a j e − ix j p k / (cid:126) is theannihilation boson operator in momentum representa-tion. Here, x j = c × j , where c is the lattice constant,and p k = 2 πk (cid:126) / ( N × c ). Each independent Hamiltonianin momentum space ( (cid:80) k ( − t cos(2 πk/N ) − µ )˜ a † k ˜ a k ) haseigenenergies − t cos(2 πk/N ) − µ , which reach their min-imum when all bosons condense with 0 momenta. Notethat the energies are independent of the states’ occupa-tion, meaning that the ground state is degenerate forany number of particles. This regime is known as theSF phase, characterized by a delocalized wave function(formally described by algebraic decaying spatial corre-lations , which is why the SF phase in 1D is not atrue Bose-Einstein condensate) which has a non-zero ex-pected variance of the local number operator. In fact,the probability distribution for the local occupation isPoissonian, meaning that (cid:104) n i (cid:105) − (cid:104) n i (cid:105) = (cid:104) n i (cid:105) . For afixed chemical potential and an infinite number of sites,there exists a continuous phase transition from the MIphase to the SF phase as U decreases with respect to t , apart from a range in µ in the first Mott lobe, wherere-entrance exists. B. Variational Monte Carlo with RestrictedBoltzmann Machine Wave Function
In this section, a review of the method introducedby Carleo and Troyer is given, where a VMC setup isused to find the ground state, formally written as a trialwavefunction given by an RBM. In the Fock space basis,the ground wave function of the Bose-Hubbard systemcan be written as | Ψ (cid:105) = ∞ (cid:88) n =0 ,n =0 ,... Ψ( n , n , . . . ) | n , n , . . . (cid:105) , (2)where n i is the occupation number of the i -th site,and | Ψ( n , n , . . . ) | are the probability amplitudes cor-responding to the Fock states | n , n , . . . (cid:105) . In order to map the wave function to a computer, both the num-ber of sites and the number of possible particles in eachsite have to be truncated. We will refer to the numberof sites as N and to the maximum number of particlesin each site as M −
1. The coefficients Ψ( n , n , . . . , n N )are approximated by an RBM. RBMs are generative neu-ral networks, formally described by a bipartite undi-rected graph such as the one shown in fig. 1, wherethere is a layer of visible neurons denoted by v thatare used to input real data, and a layer of hidden neu-rons denoted by h that are used as latent variables ofthe model . In particular, the wavefunction coefficientstake the form Ψ( n , n , . . . , n N ) ≈ ψ θ ( n , . . . , n N ) = (cid:80) h e − E RBM ( v ( n ) , h ) , where E RBM ( v ( n ) , h ) is the energyof the RBM , and θ are the variational parameters ofthe RBM. As a short-hand notation, an occupation con-figuration is denoted as n , and it is inputted to the visiblelayer of the RBM. However, the configuration first needsto be one-hot encoded as follows: each occupation n i isencoded into an M -component vector whose j -th com-ponent is δ j,n i , j = 0 , , . . . , M −
1; then, the vectors forevery occupation are concatenated into v ( n ). Moreover,if the N H hidden neurons are restricted to binary values1 or −
1, the approximated wave function coefficients canbe written as ψ θ ( n ) = e (cid:80) j a j v j ( n ) N H (cid:89) (cid:96) =1 b (cid:96) + (cid:88) j W (cid:96),j v j ( n ) , (3)where a j , b (cid:96) and W are the complex-valued visible bias,hidden bias and connection matrix of the RBM, respec-tively, which are to be learned to approximate the groundstate. FIG. 1. Illustration of the used RBM, where each site occupa-tion is one-hot encoded into M visible neurons, depicted withdifferent colors for different sites in the top layer. There are N H hidden neurons in the bottom layer connected with thevisible neurons through weights W (cid:96),j that connect the (cid:96) -thhidden neuron with the j -th visibile neuron. Thus, the approximation of the wave function coeffi-cients is done through the adjustment of the parameters θ : { a j , b (cid:96) , W (cid:96),j } that minimize the energy (cid:104) ψ θ | ˆ H | ψ θ (cid:105) .At each step of the minimization, a set M of configura-tions n is sampled from | ψ θ ( n ) | using the Metropolis-Hastings algorithm, so that the energy can be efficientlyestimated as (cid:104) ψ θ | ˆ H | ψ θ (cid:105) ≈ |M| (cid:88) n ∈M (cid:88) n (cid:48) (cid:104) n | ˆ H | n (cid:48) (cid:105) ψ θ ( n (cid:48) ) ψ θ ( n ) . (4)More explicitly, the following steps are carried out togenerate the sample M . At the first iteration, a state n is randomly proposed. Then, at the i -th iteration:1. Under some updating rule, propose a new state n (cid:48) i from n i .2. With probability min { , | ψ θ ( n (cid:48) i ) /ψ θ ( n i ) | } acceptthe state n (cid:48) i , i.e. n i +1 ← n (cid:48) i . If it is not accepted,then n i +1 ← n i .To build the sample M a total of 1000 iterations areperformed. Then, through either stochastic gradient de-scent or stochastic reconfiguration , the energy in eq. (4)can be minimized, producing a new set of parameters θ , as explained by Carleo and Troyer . Both the sam-pling and minimization of the energy with respect to theRBM parameters are repeated iteratively until the RBMparameters converge, resembling an Expectation Maxi-mization algorithm . A schematic representation of thevariational Monte Carlo technique is shown in fig. 2. FIG. 2. Representation of the variational Monte Carlo tech-nique. With randomly initialized parameters θ , a set of statesfrom the Hilbert space H is sampled. By minimizing the en-ergy defined in eq. (4), the parameters θ are updated. Thesetwo steps are repeated n times with the objective of samplingthe states (in the occupation basis) that are relevant for theground state, depicted by a red blot, with the correct proba-bility distribution. Once the RBM has been trained, the state can be re-constructed by sampling states from | ψ θ ( n ) | . An im-portant issue that immediately arises is that for an un-known target probability distribution | Ψ( n ) | , the sam-pling can be too small or too large. If it is small, impor-tant information about the ground state might not be taken into account, whereas if it is large, noisy probabil-ity from other states can be taken into account. How-ever, the small or large nature of the sampling is relativeto the number of states that contribute significantly tothe ground state. Critically, if the number of states thatcontribute significantly to the ground state is minimal,then the sampling might never be able to visit thosestates in the Hilbert space (especially for intractablylarge Hilbert spaces). This is the case of the Mott in-sulator phase, where the ground state corresponds to | n = m, n = m, . . . (cid:105) for the m -th Mott lobe. III. RESULTS
We swept several values of chemical potential and hop-ping energy corresponding to the first and part of thesecond Mott lobe in the t/U – µ/U space, with U = 1,performing 12000 sampling and optimization steps withNetKet , where 1000 Metropolis-Hastings steps weredone for each sampling step. This was done for threedifferent scenarios: for 5 sites, we used 8 and 20 hiddenneurons, and for 8 sites, we used 11 hidden neurons. Inall cases, the maximum number of bosons allowed persite was 4. A. Energy convergence
Since the VMC minimizes the energy, we check thatthe energy computed with eq. (4) converges by measur-ing its variance for the last 500 sampling-optimizationsteps, as well as by measuring the absolute error of theenergy when compared to the exact ground state energyfound through Lanczsos diagonalization. For the afore-mentioned three scenarios, the variance for the last 500sampling-optimization steps is shown in fig. 3(a)-(c), andthe corresponding absolute errors with respect to the ex-act ground state energy are shown in fig. 3(d)-(f). It isseen that in the Hamiltonian parameter space, the ma-jority of the energies have low-variance, showing conver-gence towards a value that is in excellent agreement withthe exact ground-state energies.Since the Hilbert spaces are small enough to computeall the probability amplitude coefficients for every statein the Fock basis, we can also compute any observable ˆ O as (cid:104) ψ θ | ˆ O | ψ θ (cid:105) = (cid:80) n ∈H | ψ θ ( n ) | (cid:104) n | ˆ O | n (cid:105) (cid:80) n ∈H | ψ θ ( n ) | . (5)In particular, the absolute error of the energy computedthrough eq. (5) is shown in fig. 4. It is clear that in thecase of 5 sites and 20 neurons, a very large region bothin the SF and MI phases show a strong disagreementbetween the energy calculated through eqs. (4) and (5),showing that even though the energy converged, the statedid not. Another recurrent pattern in the energies calcu-lated through eqs. (4) and (5) is an arc of high absolute FIG. 3. Energy variance and absolute error of the last 500sampling-optimization steps. Dashed lines show the phaseboundaries computed for 128 sites with DMRG by Ejima et al. . In (a), the two first Mott lobes and the SF region arelabeled explicitly. (a)-(c) show the variance for the last 500sampling-optimization steps for the cases of 5 sites and 8 hid-den neurons, 8 sites and 11 hidden neurons, and 5 sites and 20hidden neurons, respectively. (d)-(f) show the correspondingabsolute errors between the average energy value for the last500 sampling-optimization steps and the exact ground stateenergy. errors formed in the left-most side of the Mott lobes,which we will address later. FIG. 4. Absolute error between the RBM state expected en-ergy computed through eq. (5) and the exact ground stateenergy for (a) 5 sites and 8 hidden neurons, (b) 8 sites and 11hidden neurons and (c) 5 sites and 20 hidden neurons. Dashedlines are the MI-SF boundaries as in fig. 3.
B. Overlap
The RBM state non-convergence for the case of 5 sitesand 20 hidden neurons is further confirmed when we mea-sure the overlap between the exact ground state | ψ exact (cid:105) and the RBM ground state | ψ θ (cid:105) , shown in fig. 5(c). It isnow clear why the expected energy with respect to thecomplete RBM state shown in fig. 4(c) presents large er-rors when compared to the expected energy with respectto the exact ground state: it appears that the RBM hasnot learned the ground state in the bottom-right regionof the plot, which is a region that covers part of theMott lobe, as well as part of the SF phase region. Infact, since the minimization is carried out with eq. (4),the RBM only learns relative probability amplitudes be-tween the sampled states, and due to a large number of hidden neurons, the RBM overfits. Moreover, in the caseof 5 sites and 8 hidden neurons, and the case of 8 sitesand 11 hidden neurons, the RBM finds difficulty in learn-ing the ground state in the limit between the MI and theSF phase as shown in fig. 5(a) and (b). The difficultyin learning those states, and in general, in treating theground state near the MI-SF boundary comes from theKosterlitz-Thouless-like quantum phase transition in 1Dsystems , where an exponentially small Mott gap ex-ists . We see once again that there are arcs of low over-lap points formed in the left-most side of the Mott lobes.Within the SF phase, there are also lines of low overlap,which appear because of the finite size effects of the siteschain. Note that there are as many of these fictitiousboundaries as there are sites in the periodic chain understudy. Nevertheless, when comparing the 5 sites cases, itis seen that for values of µ/U > , cf. fig. 5(c)) thanwhen only 8 hidden neurons are used (see fig. 5(a)). FIG. 5. Overlap |(cid:104) ψ exact | ψ θ (cid:105)| between the exact and RBMground states for (a) 5 sites and 8 hidden neurons, (b) 8 sitesand 11 neurons, and (c) 5 sites and 20 neurons. Dashed linesare the MI-SF boundaries as in fig. 3. Since the advantage of VMC over exact diagonalizationcomes for intractably large Hilbert spaces, it is not alwayspossible to compute the probability amplitudes for all ofthe Fock states basis. In such a case, the RBM state canbe built as | ˜ ψ θ (cid:105) = (cid:80) n ∈M (cid:48) ψ θ ( n ) | n (cid:105) (cid:112)(cid:80) n ∈M (cid:48) | ψ θ ( n ) | , (6)where M (cid:48) is a sample of Fock states sampled with theMetropolis-Hastings algorithm, with an acceptance prob-ability of min { , p GC ( n i +1 ) /p GC ( n i ) } , where p GC ( n ) = Z − exp (cid:16) − (cid:104) n | ˆ H | n (cid:105) (cid:17) is the probability associated withthe grand canonical ensemble (note that the term µ (cid:104) ˆ N (cid:105) has already been introduced in the Hamiltonian). Thisstrategy was used to generate a sample M (cid:48) of up to 2048Fock states yielding a state | ˜ ψ θ (cid:105) for every point in thephase diagram. The overlap between the exact groundstate | ψ exact (cid:105) and the sampled RBM ground state | ˜ ψ θ (cid:105) isshown in fig. 6 for the three studied scenarios. Compar-ing fig. 5 with fig. 6, it is seen that the retrieved sampledstate “cleans” the RBM state by removing Fock statesfor which the probability amplitudes were badly learned,and leaves the Fock states that are relevant to the ac-tual ground state. Despite this cleaning, the larger theHilbert space, the more states have to be sampled to con-sider all relevant Fock states, as it is seen that 2048 Fockstates are insufficient to capture the ground state in thecase of 8 sites shown in fig. 6(b); however, note that thelow overlap in the arc from fig. 5(b) almost completelydisappears after the cleaning, cf. fig. 6(b). An importantphenomenon is that the overlap of the sampled RBMstate diminishes as t gets larger. This happens becauselarger values of t/U imply larger delocalization of theground wave function, which involves more Fock states.Moreover, the Hilbert space size for 8 sites consists of5 = 390625 Fock states, which is why the sampling for2048 results in poor representations of the ground state,especially in the SF phase. Sampling more Fock stateseventually reconstructs the exact ground state with veryhigh overlap, except in the boundary between the MI andSF phases (data not shown). FIG. 6. Overlap | (cid:104) ψ exact | ˜ ψ θ (cid:105)| between the exact and sam-pled RBM ground states for (a) 5 sites and 8 hidden neurons,(b) 8 sites and 11 hidden neurons and (c) 5 sites and 20 hid-den neurons, for a maximum of 2048 states sampled from theHilbert space. Dashed lines are the MI-SF boundaries as infig. 3. C. Order parameter
The phase diagram of the Bose-Hubbard model canbe reconstructed by measuring quantities in the sam-pled RBM ground state that exhibit the phase transition.As mentioned before, we choose the variance of the lo-cal number operator, in particular, of the first site. Infig. 7(a) and (b), the order parameter measured with theexact ground state is shown for 5 and 8 sites, respectively.It is seen that in the Mott insulator phase, the varianceof the number of bosons in the first lattice site is nearto 0, but not exactly 0 because of finite size effects. Onthe other hand, fig. 7(d) and (e) show the order param-eter for the sampled RBM ground state with 2048 Fockstates for the cases of 8 and 20 hidden neurons for 5 sites,which show excellent agreement with their exact counter-part fig. 7(a). However, at the boundary between the MIand SF phases, there are absolute errors that could be ashigh as 0.15, which come from the difficulty of learningthe ground state near the phase transition boundary. Inspite of these differences at the boundaries, it is clear thatthe learned RBM ground state mimics the re-entrancefound in finite 1D chains . Finally, fig. 7(c) shows theorder parameter for the exact RBM ground state for 8 sites and 11 hidden neurons. Figure 7(f) also shows theorder parameter for 8 sites and 11 hidden neurons butwith a different sampling limit. Instead of fixing a num-ber of Fock states to be sampled (2048 and 4096 wereinsufficient, data not shown), we build the sample M (cid:48) byaccepting states with the Metropolis-Hastings algorithmuntil there is a run of 400 consecutive state proposals thatdo not raise a new accepted state into M (cid:48) . In this case,the state is cleaned (note that the arc of badly learnedorder parameter almost completely disappears). A clearindicator of the number of Fock states to represent theground state emerges: for very low values of t , within thefirst Mott lobe, only one sampled state ( | , , , , (cid:105) ) isneeded to reproduce the exact ground state; in the SFphase, up to 30000 states are sampled before hitting a400 streak of no newly accepted states into the samplethat constitutes the RBM ground state. FIG. 7. Order parameter Var(ˆ n ) for the ground state ob-tained through exact diagonalization for 5 sites (a), and 8sites (b); for the sampled RBM ground state with 2048 Fockstates for 5 sites with 8 hidden neurons (d) and 20 hiddenneurons (e); and for the exact RBM ground state for 8 siteswith 11 neurons (c) and the sampled RBM ground state (f).The white line is the 0.001 contour line, and the yellow onecorresponds to the 0.01 contour line in all plots. Dashed linesare the MI-SF boundaries as in fig. 3. Other quantities can be used as an order parameter,which more explicitly relate to quantum correlations suchas entanglement. For instance, a partial trace carried outover all degrees of freedom except the first site yields areduced density matrix ρ ( t, µ ), from which the linearentropy S ( t, µ ) = 1 − Tr (cid:8) ρ ( t, µ ) (cid:9) can be measured (thevon Neumann entanglement entropy can also be used,e.g., ). Figure 8 shows the absolute errors between theexact and the sampled RBM ground states for the threestudied scenarios, where the errors in the MI-SF bound-ary become large. D. Tomography
In all the studied scenarios, there are problems forlearning the ground state at the MI-SF boundaries,as well as the mini-plateaus boundaries within the SF
FIG. 8. Absolute error for the linear entropy at the first siteof the chain between the exact ground state and the sampledRBM ground state for (a) 5 sites, 8 hidden neurons and 2048Fock states, (b) 8 sites, 11 hidden neurons and 4096 Fockstates, and (c) 5 sites, 20 hidden neurons, and 2048 Fockstates. Dashed lines are the MI-SF boundaries as in fig. 3. phase. In order to understand the differences betweenthe learned ground state and the one obtained throughexact diagonalization, we performed a study of the com-position of the ground states for 5 sites and 8 hiddenneurons. For that reason, we examined the probabil-ity amplitudes of the ground state (both RBM learnedand exact) at 11 different points in the t - µ space, fixing t = 0 .
1, and with µ at the middle and border of eachplateau in the MI and SF phases. We also examined theRBM and exact ground states for very small t at themiddle of the first Mott lobe, as indicated by the blackdot in fig. 9, where the ground state for the RBM and forexact diagonalization was | , , , , (cid:105) with an associatedprobability amplitude of 99.99%, as expected.Note that the Bose-Hubbard chain is invari-ant (up to a phase) under displacements and in-versions, i.e. there are displacement and in-version operators that act as follows on Fockstates: ˆ T | n , . . . , n N − , n N (cid:105) = e iφ | n N , n , . . . , n N − (cid:105) for displacement, and ˆ I | n , n , . . . , n N − , n N (cid:105) = e iϕ | n N , n N − , . . . , n , n (cid:105) for inversion. Therefore, in or-der to perform a tomography, we must take into accountthat all Fock states that belong to the same rung definedby displacement and inversion operations are equivalent.Now, each Fock state can be brought to a canonical Fockstate through a consecutive application of displacementand inversion operators onto the original Fock state. Thiscanonical Fock state is selected as the lexicographicallysmallest one, after every possible application of displace-ment and inversion operators, as in . In fig. 9, we showthe probability amplitude distribution of the Fock statemanifold rungs for the exact and the RBM ground states.More explicitly, the exact ground state is | Ψ (cid:105) = (cid:88) n Ψ( n ) | n (cid:105) = (cid:88) i (cid:88) n ∈ R i Ψ( n ) | n (cid:105) , (7)where i indexes the rungs, and R i is the i -th rung. Theprobability amplitude corresponding to a rung is, there-fore, the sum of the probability amplitudes of all of itsFock states, and the bars from fig. 9 show those rungsprobability amplitudes. Reading the plot from rightto left, i.e. starting with the smallest value of µ , wesee that within the first Mott lobe, both the RBM and the exact ground states show very similar distributionsover three rungs of the one-filling manifold: | , , , , (cid:105) , | , , , , (cid:105) and | , , , , (cid:105) . Increasing µ up to the firstMI-SF boundary (slightly within the SF phase), we seethat the probability amplitude distribution starts to dif-fer between the exact and the RBM ground states, andthe RBM struggles to identify if the ground state is nowin an excitation manifold above the one-filling manifold,represented by the rung | , , , , (cid:105) , or in the MI phasewhich is represented by the rung | , , , , (cid:105) . Increasingagain µ in order to be at the middle of the first plateauwithin the SF phase shows that both the RBM and theexact ground states have similar probability amplitudedistributions, even though the RBM assigns a probabilityto other rungs (not shown in the plot because their contri-bution is less than 0.01). It is clear that this first plateauis mostly represented by the rung | , , , , (cid:105) which is inthe 6-th excitation manifold. If we continue to increase µ we see that near the boundaries between the MI-SFphases and between the plateaus within the SF phase,the RBM and the exact ground states exhibit differencesin the probability amplitude distributions. On the con-trary, in the middle of those plateaus, the probability am-plitude distributions from both RBM and exact groundstates are very similar (which is also seen in fig. 5(a)).Moreover, each time the tomography moves onto a newplateau (with a higher value of µ ), the excitation mani-fold increases by one, up to the 10-th excitation manifold,which corresponds to the two-filling manifold at the sec-ond Mott lobe, characterized by the rung | , , , , (cid:105) aswell as the rungs | , , , , (cid:105) and | , , , , (cid:105) . IV. CONCLUSIONS
In this work, we systematically tested the capabilitiesof VMC with a trial ground wavefunction given by anRBM on the one-dimensional Bose-Hubbard model. Themotivation for the technique comes from the possibilityof incorporating it into the toolbox of quantum physicsto tackle theoretical problems that are difficult to studynumerically due to the intractably large Hilbert spaces.Thus, it is first needed to intensively test the techniqueto reproduce known results, and it is also needed to theo-retically explain why the technique works (this is a chal-lenging open question which involves the question of whyneural networks work well). Only if the community iden-tifies the strengths and weaknesses of the technique isit possible to use it to explore problems of interest thatinvolve vast Hilbert spaces. Regarding the model un-der study, we repeatedly found differences between theexact ground state found through the exact diagonaliza-tion of the Bose-Hubbard Hamiltonian and the learnedground state. We did so in one-dimensional chains of 5and 8 sites. In order to better learn the ground state, thevariational trial wavefunction was enriched in the case of5 sites by increasing the number of hidden units in theRBM’s hidden layer. Although results improved near the
FIG. 9. Quantum tomography of the probability amplitude of Fock state manifold rungs for the ground state found throughexact optimization and through VMC with an RBM wave function ansatz. In the lower part of the figure, a 3D plot of theorder parameter for the exact ground state is shown for the case of 5 sites. Each bar color represents a manifold rung depictedthrough its canonical Fock state as the lexicographically smallest one (e.g. the rung {| , , , , (cid:105) , | , , , , (cid:105) , . . . , | , , , , (cid:105)} is represented by | , , , , (cid:105) ). A black dot in the middle of the first Mott lobe indicates that a tomography was made thereas well. MI-SF boundary, other regions of the Bose-Hubbard pa-rameter space suffered from badly learned ground states.For that reason, we proposed a sampling technique thatcleans the ground state, getting rid of contributions thatarise because of over-fitting. However, the VMC-RBMmethod to find the ground state was proven to yield goodresults not only in computing the energy (with low errorsand clear signs of convergence in most of the parame-ter space) but also in computing other observables whichexplicitly involve quantities related to quantum correla-tions, such as the linear entropy of one of the chain’s sites.In particular, we reconstructed with excellent accuracythe phase diagram of the Bose-Hubbard model using thelocal occupation variance as an order parameter for mostof the parameter space. Nevertheless, we also carried outquantum state tomographies to understand the compo-sition of the ground states at the interfaces between theMI-SF phases and also within the excitation manifoldplateaus formed for finite chains in the SF phase, whichrevealed significant differences between the Fock statesdistributions of the exact ground state and the learnedground state.Accordingly, we must raise the attention that the VMCtechnique is consistently challenged near the quantumphase transition of the Bose-Hubbard model. We do not doubt that the high-expressivity of neural quantumstates are capable of improving those results under moreextensive sampling, or more complex variational trialwavefunctions apart from RBMs; however, we report thatsome questions have to be answered with precision beforeconfidently using VMC as a method to accurately explorethe physics of a many-body problem without the supportof any other numerical method. In the first place, we ob-served that not every point of the phase diagram requiresthe same amount of computational work to learn theground state. Most notably, the ground state in the MIphase is mostly explained by only one Fock state, whereasthe number rapidly grows for ground states in the SFphase. Therefore, automated ways of stopping samplingwithin the sampling steps should be taken into consider-ation. In particular, when using the Metropolis-Hastingsalgorithm, we found that if a Markov chain of a certainlength was formed with no new sampled states, we couldstop sampling, yielding a high-quality ground state. Sec-ondly, the VMC technique showed badly learned groundstates along arcs formed within the Mott lobes. Whythese form remains unanswered, as they are not relatedto the structure of the MI-SF boundaries. However, afterwe clean the learned ground state, these arcs almost com-pletely disappear. Thus, it might be the case that theyare related to over-fitting, which traduces in noisy proba-bility from Fock states that are not relevant to the groundstate. Finally, concerning the first question, the numberof optimization steps required to learn the ground statehas to be better understood, as it is only clear when the energy converges, but not the state . ACKNOWLEDGMENTS
We want to thank Andr´es Urquijo and J. P. Restrepo-Cuartas for useful discussions at the early stage of thisstudy. ∗ [email protected] A. Georges, G. Kotliar, W. Krauth, and M. J. Rozenberg,Dynamical mean-field theory of strongly correlated fermionsystems and the limit of infinite dimensions, Rev. Mod.Phys. , 13 (1996). U. Schollw¨ock, The density-matrix renormalization group,Rev. Mod. Phys. , 259 (2005). G. Carleo and M. Troyer, Solving the quan-tum many-body problem with artificial neu-ral networks, Science , 602 (2017),https://science.sciencemag.org/content/355/6325/602.full.pdf. J. Han, L. Zhang, and W. E, Solving many-electronSchr¨odinger equation using deep neural networks, Jour-nal of Computational Physics , 108929 (2019),arXiv:1807.07014. K. Choo, A. Mezzacapo, and G. Carleo, Fermionic neural-network states for ab-initio electronic structure (2019),arXiv:1909.12852 [physics.comp-ph]. D. Pfau, J. S. Spencer, A. G. de G. Matthews, andW. M. C. Foulkes, Ab-initio solution of the many-electronschr¨odinger equation with deep neural networks (2019),arXiv:1909.02487 [physics.chem-ph]. J. Hermann, Z. Sch¨atzle, and F. No´e, Deep neural net-work solution of the electronic schr¨odinger equation (2019),arXiv:1909.08423 [physics.comp-ph]. J. Kessler, F. Calcavecchia, and T. D. K¨uhne, Artificialneural networks as trial wave functions for quantum montecarlo (2019), arXiv:1904.10251 [physics.comp-ph]. M. Ruggeri, S. Moroni, and M. Holzmann, Nonlinear Net-work Description for Many-Body Quantum Systems inContinuous Space, Physical Review Letters , 205302(2018). G. Carleo, Y. Nomura, and M. Imada, Constructing exactrepresentations of quantum many-body systems with deepneural networks, Nature communications , 1 (2018). O. Sharir, Y. Levine, N. Wies, G. Carleo, and A. Shashua,Deep autoregressive models for the efficient variationalsimulation of many-body quantum systems, Phys. Rev.Lett. , 020503 (2020). M. Hibat-Allah, M. Ganahl, L. E. Hayward, R. G. Melko,and J. Carrasquilla, Recurrent neural network wavefunc-tions (2020), arXiv:2002.02973 [cond-mat.dis-nn]. S. Czischek, M. G¨arttner, and T. Gasenzer, Quenchesnear Ising quantum criticality as a challenge for artifi-cial neural networks, Physical Review B , 1 (2018),arXiv:1803.08321. Y. Nomura, A. S. Darmawan, Y. Yamaji, and M. Imada,Restricted Boltzmann machine learning for solvingstrongly correlated quantum systems, Physical Review B , 1 (2017), arXiv:1709.06475. P. Teng, Machine-learning quantum mechanics: Solving quantum mechanics problems using radial basis functionnetworks, Phys. Rev. E , 033305 (2018). A. Szab´o and C. Castelnovo, Neural network wave func-tions and the sign problem (2020), arXiv:2002.04613 [cond-mat.str-el]. K. Choo, T. Neupert, and G. Carleo, Two-dimensionalfrustrated J1-J2 model studied with neural networkquantum states, Physical Review B , 10.1103/Phys-RevB.100.125124 (2019), arXiv:1903.06713. X. Liang, W. Y. Liu, P. Z. Lin, G. C. Guo, Y. S. Zhang, andL. He, Solving frustrated quantum many-particle modelswith convolutional neural networks, Physical Review B ,1 (2018), arXiv:1807.09422. A. Haim, R. Kueng, and G. Refael, Variational-correlationsapproach to quantum many-body problems (2020),arXiv:2001.06510 [cond-mat.str-el]. P. Rosson, M. Kiffner, J. Mur-Petit, and D. Jaksch, Char-acterizing the phase diagram of finite-size dipolar Bose-Hubbard systems, Physical Review A , 13616 (2020),arXiv:1909.09099. X. Gao and L. M. Duan, Efficient representation of quan-tum many-body states with deep neural networks, NatureCommunications , 1 (2017), arXiv:1701.05039. S. Sun, Z. Cao, H. Zhu, and J. Zhao, A survey of optimiza-tion methods from a machine learning perspective (2019),arXiv:1906.06821 [cs.LG]. L. E. Ballentine,
Quantum mechanics: a modern devel-opment (World Scientific Publishing Company, 1998) pp.293–294. G. Carleo, K. Choo, D. Hofmann, J. E. Smith, T. Wester-hout, F. Alet, E. J. Davis, S. Efthymiou, I. Glasser, S.-H.Lin, et al. , Netket: A machine learning toolkit for many-body quantum systems, SoftwareX , 100311 (2019). K. McBrian, G. Carleo, and E. Khatami, Ground statephase diagram of the one-dimensional bose-hubbard modelfrom restricted boltzmann machines, in
Journal of Physics:Conference Series , Vol. 1290 (IOP Publishing, 2019) p.012005. Z. Li,
Quantum simulation of quantum phases of matterwith interacting photons , Ph.D. thesis, National Universityof Singapore (2019). H. Saito, Solving the Bose–Hubbard model with machinelearning, Journal of the Physical Society of Japan ,10.7566/JPSJ.86.093001 (2017), arXiv:1707.09723. K. Choo, G. Carleo, N. Regnault, and T. Neupert, Sym-metries and Many-Body Excitations with Neural-NetworkQuantum States, Physical Review Letters , 167204(2018), arXiv:1807.03325. H. Saito and M. Kato, Machine learning technique to findquantum many-body ground states of bosons on a lattice,Journal of the Physical Society of Japan , 1 (2018), arXiv:1709.05468. N. Irikura and H. Saito, Neural-network quantum statesat finite temperature, , 1 (2019), arXiv:1911.02774. Y. H. Liu and E. P. Van Nieuwenburg, Discriminative Co-operative Networks for Detecting Phase Transitions, Phys-ical Review Letters , 176401 (2018), arXiv:1706.08111. P. Huembeli, A. Dauphin, and P. Wittek, Identifying quan-tum phase transitions with adversarial neural networks,Physical Review B , 15 (2018), arXiv:1710.08382. P. Broecker, F. F. Assaad, and S. Trebst, Quantum phaserecognition via unsupervised machine learning, arXiv e-prints , arXiv:1707.00663 (2017), arXiv:1707.00663 [cond-mat.str-el]. J. P´erez D´ıaz,
Detection of quantum phase transitions viamachine learning algorithms , Master’s thesis, UniversitatPolit`ecnica de Catalunya (2019). M. P. Fisher, P. B. Weichman, G. Grinstein, and D. S.Fisher, Boson localization and the superfluid-insulatortransition, Physical Review B , 546 (1989). A. van Oudenaarden and J. E. Mooij, One-DimensionalMott Insulator Formed by Quantum Vortices in JosephsonJunction Arrays, Physical Review Letters , 4947 (1996). R. Baltin and K. H. Wagenblast, Quantum phase transi-tions for bosons in one dimension, Europhysics Letters ,7 (1997), arXiv:9705261 [cond-mat]. C. Bruder, R. Fazio, and G. Schon, The bose-hubbardmodel: from josephson junction arrays to optical lattices,Annalen der Physik , 566 (2005). M. Greiner, O. Mandel, T. Rom, A. Altmeyer, A. Widera,T. W. H¨ansch, and I. Bloch, Quantum phase transitionfrom a superfluid to a Mott insulator in an ultracold gas ofatoms, Physica B: Condensed Matter , 11 (2003). D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, andP. Zoller, Cold bosonic atoms in optical lattices, PhysicalReview Letters , 3108 (1998). I. Bloch, J. Dalibard, and W. Zwerger, Many-body physicswith ultracold gases, Reviews of Modern Physics , 885(2008), arXiv:0704.3011. T. St¨oferle, H. Moritz, C. Schori, M. K¨ohl, andT. Esslinger, Transition from a strongly interacting 1d su-perfluid to a mott insulator, Phys. Rev. Lett. , 130403(2004). I. B. Spielman, W. D. Phillips, and J. V. Porto, Mott-insulator transition in a two-dimensional atomic bose gas,Phys. Rev. Lett. , 080404 (2007). N. Gemelke, X. Zhang, C.-L. Hung, and C. Chin, In situobservation of incompressible mott-insulating domains inultracold atomic gases, Nature , 995 (2009). G. Grynberg, P. Horak, and C. Mennerat-Robilliard, Spa-tial diffusion of atoms cooled in a speckle field, EurophysicsLetters (EPL) , 424 (2000). C. Kollath, U. Schollw¨ock, J. von Delft, and W. Zwerger,Spatial correlations of trapped one-dimensional bosons inan optical lattice, Phys. Rev. A , 031601 (2004). P. Lugan, D. Cl´ement, P. Bouyer, A. Aspect, M. Lewen-stein, and L. Sanchez-Palencia, Ultracold bose gases in 1ddisorder: From lifshits glass to bose-einstein condensate,Phys. Rev. Lett. , 170403 (2007). J. K. Freericks and H. Monien, Strong-coupling expansionsfor the pure and disordered Bose-Hubbard model, PhysicalReview B , 2691 (1996). S. Ejima, H. Fehske, F. Gebhard, K. Zu M¨unster, M. Knap,E. Arrigoni, and W. Von Der Linden, Characterizationof Mott-insulating and superfluid phases in the one- dimensional Bose-Hubbard model, Physical Review A -Atomic, Molecular, and Optical Physics , 1 (2012),arXiv:1203.1120. M. A. Cazalilla, R. Citro, T. Giamarchi, E. Orignac, andM. Rigol, One dimensional bosons: From condensed mat-ter systems to ultracold gases, Reviews of Modern Physics , 1405 (2011), arXiv:1101.5337. N. Teichmann, D. Hinrichs, M. Holthaus, and A. Eckardt,Process-chain approach to the Bose-Hubbard model:Ground-state properties and phase diagram, Physical Re-view B - Condensed Matter and Materials Physics , 1(2009). G. G. Batrouni, R. T. Scalettar, G. T. Zimanyi, andA. P. Kampf, Supersolids in the bose-hubbard hamilto-nian, Phys. Rev. Lett. , 2527 (1995). W. Krauth, N. Trivedi, and D. Ceperley, Superfluid-insulator transition in disordered boson systems, Phys.Rev. Lett. , 2307 (1991). W. Krauth and N. Trivedi, Mott and superfluid transitionsin a strongly interacting lattice boson system, EurophysicsLetters (EPL) , 627 (1991). P. Niyaz, R. T. Scalettar, C. Y. Fong, and G. G. Ba-trouni, Ground-state phase diagram of an interacting bosemodel with near-neighbor repulsion, Phys. Rev. B , 7143(1991). T. K¨uhner and H. Monien, Phases of the one-dimensionalBose-Hubbard model, Physical Review B - CondensedMatter and Materials Physics , R14741 (1998),arXiv:9712307 [cond-mat]. S. Ejima, H. Fehske, and F. Gebhard, Dynamic proper-ties of the one-dimensional Bose-Hubbard model, Epl ,10.1209/0295-5075/93/30002 (2011), arXiv:1102.2028. T. D. K¨uhner and S. R. White, One-dimensional Bose-Hubbard model with nearest-neighbor interaction, Physi-cal Review B - Condensed Matter and Materials Physics , 12474 (2000), arXiv:9906019 [cond-mat]. G. G. Batrouni and R. T. Scalettar, World-line quantummonte carlo algorithm for a one-dimensional bose model,Phys. Rev. B , 9051 (1992). V. A. Kashurnikov, A. V. Krasavin, and B. V. Svis-tunov, Mott-insulator-superfluid-liquid transition in a one-dimensional bosonic hubbard model: Quantum montecarlo method, Journal of Experimental and TheoreticalPhysics Letters , 99 (1996). G. G. Batrouni, R. T. Scalettar, and G. T. Zimanyi, Quan-tum critical phenomena in one-dimensional bose systems,Phys. Rev. Lett. , 1765 (1990). N. Elstner and H. Monien, Dynamics and thermodynamicsof the bose-hubbard model, Phys. Rev. B , 12184 (1999). W. Koller and N. Dupuis, Variational cluster perturbationtheory for bose–hubbard models, Journal of Physics: Con-densed Matter , 9525 (2006). M. Pino, J. Prior, and S. R. Clark, Capturing the re-entrant behavior of one-dimensional Bose-Hubbard model,Physica Status Solidi (B) Basic Research , 51 (2013). P. Smolensky, Information processing in dynamical sys-tems: Foundations of harmony theory, in
Parallel Dis-tributed Processing: Explorations in the Microstructure ofCognition, Vol. 1: Foundations (MIT Press, Cambridge,MA, USA, 1986) p. 194–281. S. Sorella, M. Casula, and D. Rocca, Weak bind-ing between two aromatic rings: Feeling the van derwaals attraction by quantum monte carlo methods,The Journal of Chemical Physics , 014105 (2007), https://doi.org/10.1063/1.2746035. A. P. Dempster, N. M. Laird, and D. B. Rubin, Maxi-mum likelihood from incomplete data via the em algo-rithm, Journal of the Royal Statistical Society. Series B(Methodological) , 1 (1977). T. Giamarchi,
Quantum physics in one dimension , Vol. 121(Clarendon press, 2003). If instead the sample is built from the distribu-tion | ψ n | , then the acceptance probability becomesmin { , | ψ θ ( n i +1 ) /ψ θ ( n i ) | } and the state build througheq. (5) inherits the RBM state problems (data not shown). P. Buonsante and A. Vezzani, Ground-state fidelity andbipartite entanglement in the bose-hubbard model, Phys.Rev. Lett. , 110601 (2007). S.-h. Park, C. Park, and M.-C. Cha, Critical point of theone-dimensional boson hubbard model, Journal of the Ko-rean Physical Society , 1553 (2004). In this regard, NetKet team has included the Gelman-Rubin statistics .73