Instantons in the Hofstadter butterfly: difference equation, resurgence and quantum mirror curves

We study the Harper-Hofstadter Hamiltonian and its corresponding non-perturbative butterfly spectrum. The problem is algebraically solvable whenever the magnetic flux is a rational multiple of 2π. For such values of the magnetic flux, the theory allows a formulation with two Bloch or θ-angles. We treat the problem by the path integral formulation, and show that the spectrum receives instanton corrections. Instantons as well as their one loop fluctuation determinants are found explicitly and the finding is matched with the numerical band width of the butterfly spectrum. We extend the analysis to all 2-instanton sectors with different θ-angle dependence to leading order and show consistency with numerics. We further argue that the instanton-anti-instanton contributions are ambiguous and cancel the ambiguity of the perturbation series, as they should. We hint at the possibility of exact 2-instanton solutions responsible for such contributions via Picard-Lefschetz theory. We also present a powerful way to compute the perturbative fluctuations around the 1-instanton saddle as well as the instanton-anti-instanton ambiguity by using the topological string formulation.


Introduction
The spectral problem for electrons on a two-dimensional square lattice in a uniform magnetic field was originally considered by Harper in 1955 [1], where an elegant difference equation was derived. More than 20 years later in 1976 Hofstadter derived a recursive equation which allowed him to plot the spectrum as a function of the magnetic field, now known as the Hofstadter butterfly [2]. Due to the magnetic effect, the electron spectrum shows a rich structure. Recently, a novel link between a two-dimensional electron lattice system and a Calabi-Yau geometry was found in [3]. It was pointed out in [3] that this Hofstadter's spectral problem is related to another spectral problem appearing in the mirror JHEP01(2019)079 geometry of the toric Calabi-Yau manifold known as local F 0 [4]. The interesting point of this relation is that the magnetic effect is interpreted as a kind of quantum deformations of the Calabi-Yau geometry. One can probe quantum Calabi-Yau geometry by the 2d electron lattice system in the magnetic field. The correspondence was generalized to the triangular lattice and another Calabi-Yau manifold [5].
In the present paper, our goal is a more quantitative understanding of this relation as well as the non-perturbative and resurgent structure of the spectrum. We here focus on the band structure of the Harper-Hofstadter problem in the weak magnetic limit. In this regime, we can treat the magnetic flux perturbatively. The perturbative expansion of the energy spectrum can explain the position (the center) of the band for each Landau level. However, it does not explain the width of bands because the band width is non-perturbative in the weak magnetic flux limit. Such non-perturbative corrections are caused by quantum mechanical tunneling effects. We will demonstrate that the non-perturbative band width is explained by instanton effects in the path integral formalism. This was observed long ago in [6] (see also [7] for the WKB approach to the problem). However here we will focus on the resurgent properties intimately related to these instantons, or more correctly to the multi-instanton contributions which we discuss in some details.
Technically, we have a very efficient way to compute the perturbative expansion of the energy spectrum around the trivial saddle [8,9], but this efficient way is not applicable for the computation of semiclassical expansion around the other nontrivial saddles. To our knowledge, there are no systematic ways to compute the semiclassical expansions around the instanton saddles in the Harper-Hofstadter model. We employ several approaches to extract this information. One is a brute force numerical approach, which we use as a check. The second is a path-integral approach, where we find the exact saddle of the path-integral action and the one-loop fluctuation. We push this computation to the 2-instanton sector and find matching results to the numerics. The instanton analysis is performed only to the leading-order in perturbation theory, and is not easily extended to perturbative corrections around the instanton saddles.
To extract corrections around instanton saddles, we employ a rather unconventional approach. We use the connection with a toric Calabi-Yau threefold, local F 0 , and find that the non-perturbative band width is captured by the free energy of the refined topological string on this geometry. Using this remarkable connection, we can efficiently compute the semiclassical fluctuation around the 1-instanton saddle by using the string theory technique, called the refined holomorphic anomaly equations [10][11][12]. Our approach here is conceptually very similar to the previous works [13][14][15] on certain quantum mechanical systems. 1 We would like to emphasize that here we have a realistic electron system where string theory techniques can be applied.
The structure of the rest of the paper is as follows. In section 2 we quickly review the eigenvalue problem of the Harper-Hofstadter model and its exact solutions when the magnetic flux φ is 2π times a rational number. We argue that in the latter case there are two Bloch's angles which can be turned on, while only one of them can be turned on

JHEP01(2019)079
if the magnetic flux has a generic value, as in the case of a trans-series solution of the Harper-Hofstadter model. In section 3, we make a trans-series ansatz for the energy in the small φ limit. We then compute the leading order contribution in the 1-instanton sector for the ground state energy by a path integral calculation, and find that it agrees with the numerical results. In section 4, we perform further path integral calculations in the 2-instanton sector, and compare with the numerical results. The imaginary part of the instanton-anti-instanton sector is extracted numerically using the well-known relation to the large order growth of the perturbative energy. Inspired by [15][16][17], we also find in section 5 the fluctuations in the 1-instanton and instanton-anti-instanton sector can be computed from topological string on local F 0 . Finally we conclude and list some open problems in section 6.

The Harper-Hofstadter problem
To prepare for the other sections, we quickly review in this section the classic results on the Harper-Hofstadter model [1,2], including the formulation of its eigenvalue problem, and the exact solutions when the magnetic flux is 2π times a rational number. We make the careful distinction that there are two Bloch's angles in this case while only one of them can be turned on if the value of the magnetic flux is generic.

The eigenvalue problem of the Harper-Hofstadter equation
The Harper-Hofstadter model describes an electron in a two dimensional lattice potential with a uniform magnetic flux in the perpendicular direction. Let the lattice spacing be a, and suppose the electron momentum has components k x and k y in the two directions. The energy of the electron before turning on the magnetic flux is, up to a normalization E = − 1 2 e ikxa + e −ikxa + e ikya + e −ikya + 2 . (2.1) We have chosen for later convenience a particular normalization so that the energy vanishes for zero electron momentum. In this convention, the energy forms a single band 0 ≤ E ≤ 4. After we turn on the magnetic flux, quantum mechanically we get the Hamiltonian operator by replacing the momentum k by the operator 2 π := p − A. Notice that p is the canonical momentum. Upon the gauge transformation A → A + ∇Λ, the Hamiltonian is only invariant up to a canonical transformation p → p + ∇Λ. Under such a canonical transformation, the state of the Hilbert space transforms as |Ψ → e iΛ(x,y) |Ψ . Notice that the momentum π generally depends on the coordinates. Indeed this is reflected in the fact that the commutator [π x , π y ] = iF xy (x, y) (2.2) where F xy (x, y) = ∂ x A y − ∂ y A x is the xy component of the field-strength tensor of A, i.e. the magnetic field through the xy-plane at the point (x, y). Henceforth, we consider the case where the magnetic field is uniform: F xy (x, y) = B.

JHEP01(2019)079
Replacing 3 x = π x a, y = π y a, we have that the lattice Hamiltonian becomes H = − 1 2 e ix + e −ix + e iy + e −iy + 2 , (2. 3) with the commutation relation [x, y] = iφ, (2.4) where φ = Ba 2 is the flux of the magnetic field through the plaquette. We will also use the exponentiated notation T x = e ix , T y = e −iy (2.5) with the commutation relation T x T y = e iφ T y T x (2.6) so that the Hamiltonian can be written as We regard x and y as the canonical operators, and can now look at eigenstates |ψ in the x-representation, i.e. define ψ(x) = x|ψ where |x is an eigenstate of x with eigenvalue x, so that which is just the difference equation.

Symmetries and θ-angles
The Hamiltonian (2.3) clearly commutes with the symmetry operators 4 each of which generates a group Z. The labelling above is becausẽ But we generally haveT Since for a generic value of φ ∈ R the operators commute up to a phase, we can say that the physical symmetry group Z × Z acts projectively. Let us first choose that φ = 2π/Q , Q ∈ Z .
(2.12) 3 Despite the notation, x and y are not the original coordinates of the system, but are proportional to the magnetic translation operators. 4 Similar operators also play an important role in the context of quantum mechanics associated with toric Calabi-Yau threefolds [3].

JHEP01(2019)079
In that case the two operators commute, and the symmetry Z × Z is no longer acting projectively. Now we can project to simultaneous eigenstates of the operatorsT x andT y , i.e. we can demand thatT The angles θ x and θ y are Bloch's angles for the x and y translations. Notice however that they can only be defined in this way if 2π/φ ∈ Z.
Next, we consider more general case that φ/(2π) = P/Q ∈ Q , (2.14) where P, Q are coprime integers. Then we have that Clearly, if P = 1, the generatorsT x ,T y must be supplemented by the generator 5 I P = e i 2π P , and the Z × Z must be centrally extended by Z P .
What about θ-angles? In this case we have that [(T x ) P , (T y ) P ] = 0, and we can define θ x , θ y angles by the simultaneous eigenstate of (T x ) P and (T y ) P . Alternatively, in this case we also have [(T x ), (T y ) P ] = 0, so we could equally define the two θ-angles as eigenstates of these two operators. Finally if P = n 2 is a perfect square, we have that [(T x ) n , (T y ) n ] = 0 and we can define θ-angels accordingly as well. In most cases we will only consider P = 1. Then, all these definitions of θ-angles coincide, and we are back to the scenario (2.12).
Finally if φ/2π is irrational, thenT where α/(2π) = −2π/φ is irrational as well. The additional generator I α = e iα generates the group Z, so the Z × Z is centrally extended by Z. In this case we are allowed only one θ-angle, which we can get as an eigenstate of eitherT x orT y but not both simultaneously.

Exact solutions for rational magnetic flux
It is well-known that the eigenvalue problem (2.8) can be solved exactly if the rationality condition (2.14) is satisfied [2]. Let us set where P, Q are two coprime integers and Q > 0. The underlying reason of the exact solvability is that in the case of (2.17) we can project onto simultaneous eigenstates of the powersT P x andT P y , as these two operators commute. This will allow, as we shall see, for a finite-dimensional representation of the operators T x and T y , in which the Hamiltonian (2.7) is written, and give us an algebraic equation for the eigenvalue problem. Note that in this case T x , T y are also shift operators, as because there always exists an integer k such that e − 2πiQ P k = IP .

JHEP01(2019)079
Recall that in this case we can define θ-angles as eigenvalues ofT P x = T Q y andT P y = T Q x . Now let us for the moment choose θ x = θ y = 0 mod 2π, i.e.
(T (0) In other words we impose periodic boundary conditions on physical states under the shift x → x − 2πP and y → y − 2πP . The algebra (2.6), which now reads has a finite dimensional representation in terms of the clock and shift matrices where q = e iφ = e 2πiP Q . Note also that (T (0) as it should. Now let us introduce the twisted boundary condition through the replacement while the algebra (2.6) is intact. Alternatively, the twisted boundary condition is equivalent to a deformation of the Hamiltonian. Using the notation k x = θ x /Q, k y = θ y /Q, we can write the Hamiltonian operator depending on k x and k y as (2.24) so that the characteristic equation det(H − E I Q×Q ) = 0 is given by  Figure 1. The Hofstadter butterfly plots energy levels E kx,ky (N, φ) with 0 ≤ k x , k y ≤ 2π/Q against magnetic flux φ ∈ 2πQ for the Harper-Hofstadter model. We take φ/2π to be P/Q for any coprime pairs of positive integers such that P ≤ Q and Q ≤ 30. with As in [18], it is straightforward to check that Using the symmetry under the mapping (k x , k y ) → (k y , −k x , ), one finds that the equation (2.25) can be simplified to F P/Q (E, 0, 0) + 4 = 2(cos(θ x ) + cos(θ y )) , (2.28) with the Bloch's angels θ x = Qk x , θ y = Qk y . It is then a simple job to get eigen-energy E by solving (2.28). We notice that the equation (2.28) depends on the value of P only through the polynomial F P/Q (E, 0, 0). Note (2.28) indicates that the minimal ranges for the Bloch's angles θ x , θ y are as they should. By varying the values of θ x , θ y , the eigen-energies E(θ x , θ y ) form bands. The two edges of a energy band correspond to (θ x , θ y ) = (0, 0), (π, π). If we turn off one Bloch's angle, the energy band width is reduced to its one half. We reproduce in figure 1 the famous plot of the Hofstadter butterfly, which is a plot of the energy bands as a function of the magnetic flux φ when φ/2π is rational.

Why trans-series expansion?
We are interested in the energy spectrum of the Harper-Hofstadter model in the weak flux limit φ → 0. As discussed in section 2, with generic values of φ, we should use the Hamiltonian operator (2.3) for the twisted boundary condition with only one Bloch's angle. Throughout this paper, however, we consider the weak flux limit with the specific form then we can introduce two distinct Bloch's angles θ x and θ y simultaneously. We want to understand the spectral behavior in the limit (3.1). To do so, it is useful to treat φ as a continuous parameter even in the specific case (3.1). Since the Hamiltonian is a Laurent polynomial of e ix and e iy , we can use the Mathematica package BenderWu [8,9] to compute its perturbative energy. 6 The first few orders are as follows where N is the Landau level of the eigen-energy. We note the agreement with earlier studies [19,20]. The perturbative energy (3.2) or even its Borel resummation cannot be the full answer. First of all, the higher order terms of the perturbative series have the same sign, and thus its Borel transform of the perturbative series has poles on the positive axis, leading to ambiguity in the Borel resummation. This ambiguity is an indication that the energy receives non-perturbative corrections. We discuss the ambiguity in detail in section 4.4. Second, the perturbative series clearly does not depend on Bloch's angles, thus itself alone cannot explain the energy bands. As a result, the band spectrum should have the transseries expansion, with the explicit dependence of θ x and θ y in instanton sectors. The trans-series expansion of the spectrum should take the following form: The leading perturbative contribution is given by (3.2). The k-instanton sector is exponentially suppressed by a factor e −kA/φ with a constant A. 7 Our goal in this paper is to reveal this trans-series structure. In particular, we will show explicit forms for a few instanton sectors. In the subsequent subsections, we will first compute the leading (1-loop) order contribution to the 1-instanton for the ground state energy by an honest path integral computation, and then compare them with the 6 The Hamiltonians considered in [9] consist of operators e x and e y with [x, y] = i . To translate it into our case here, one has to identify = −φ. See subsection 5 for detail. We have however updated the BenderWu package version 2.2 with a function BWDifferenceArray, which allows of mixed inclusion of terms e x , e p , e ix , e ip . The package is available on Wolfram Package site. 7 More precisely, beyond the one-instanton order, in general logarithmic corrections of the form log φ also appear. Therefore, the full trans-series expansion consists of three kinds of trans-monomials: φ, e −A/φ and log φ.

JHEP01(2019)079
prediction from the numerical analysis. We further argue that the quantum fluctuations in the one-instanton sector for any energy level can be read off from the topological string theory on local F 0 . In the next section, we will investigate the 2-instanton sector.

Path integral in one-instanton sector
The problem of instantons in the Harper/Hofstadter problem was first discussed in [6] where the authors computed the one-instanton and its one-loop determinant numerically.
Here we will re-derive these instanton solutions and compute analytically the one-loop fluctuations in the instanton sectors of the ground state energy.
To begin with, let's reproduce the Hamiltonian operator (2.3) for convenience The cosine potential has infinitely many degenerate vacua located at x = 2πn x , y = 2πn y , n x , n y ∈ Z . (3.5) Classically we have complete freedom of whether to identify different vacua as physically equivalent. This is not possible quantum mechanically for generic values of φ as we shall see.
Treating φ as the Planck constant, the above Hamiltonian can be associated with the Euclidean path integral with boundary conditions for x and y to be specified momentarily. The partition function above is related to the eigen-energies E(N ) with levels N = 0, 1, 2, . . . of the Hamiltonian H by so that the ground state energy E(0) can be obtained through the Euclidean path integral in the large β limit. Before we continue we should emphasize that the action of the above path integral is similar to that of the phase-space quantum mechanical system where x is identified with a coordinate, and y is identified with a momentum. The difference is that here we do not have a purely Gaussian dependence on the "momentum" y. For this reason we cannot integrate it out. Still one may hope to analyze the problem semi-classically. But there are several issues here. Firstly the semi-classics of path-integrals is to this day not a completely understood subject, but it has become clear recently that the correct interpretation of it is via the Picard-Lefschetz (PL) theory [21][22][23][24][25][26][27][28][29][30]. The PL theory analysis is by far not a straightforward matter, and requires the identification of saddles which contribute in the semi-classical expansion. As we shall see all such saddles of the action above will be on complex x, y trajectories. We do not a priori know whether such saddles should contribute.
To determine it we should compute the so-called intersection number of the co-thimble (we JHEP01(2019)079 refer the reader to the cited literature for details). This is a difficult task way beyond our current understanding. We will find some instanton solutions and argue that they must contribute on physical grounds. We will check quantitatively their contribution against numerics and find exact agreement.
Secondly it is not clear whether a continuum limit of the above path-integral exists. The path-integral is typically obtained by slicing the Boltzmann weight into N pieces, and inserting a complete set of states in between. This amounts to a lattice discretization of the path-integral, with a lattice spacing = β/N . Upon integration over the momentum, the resulting path-integral has a Gaussian suppression factors e −(...
. As we take the continuum limit → 0 the path of x is forced to be smoother and smoother. No such smoothness seems to be justified in the continuum-limit of the phase-space path integral above. Still as we shall see the semiclassical analysis passes many non-trivial checks against the numerical brute-force calculation.
The boundary conditions of the path integral can be made strictly periodic. This amounts to saying that values of coordinates (x, y) and (x + 2πn x , y + 2πn y ) are physically distinct for any n x , n y ∈ Z. In this case the above Lagrangian has a shift symmetry which takes x → x + 2πn x and y → y + 2πn y , with n x,y ∈ Z. Now let us consider the values of x and x + 2π to be physically equivalent. In other words we are gauging the shift symmetry of the scenario above, projecting the full Hilbert space down to eigenstates of a shift symmetry operator. Without a θ x -term, the projection will be to singlets of the shift operator. Gauging the symmetry amounts to saying that the boundary conditions must be relaxed to include periodicity of x(t) up to a 2π shift, i.e.
where m x is to be summed over. The integers m x can be viewed as holonomies of the Z-valued gauge field which we have to sum over in order to project to a subspace of singlets under the shift symmetry x → x + 2π.
Notice however that after gauging the x-shift symmetry, shifting y to y + 2πn y we get an additional phase in the partition function The above is only x-shift symmetry is gauged) and that y → y + 2π is a global symmetry we must have that 8 φ = 2π/Q. This is of course evident from the Hilbert space picture, but it is satisfying to see it in the path-integral. Incidentally we can say that there is a 't Hooft anomaly between the two (Z) x and (Z) y shift symmetries, so that the system must break at least one of the two to saturate the anomaly.
Since we are assuming that φ = 2π/Q, we can insert the two θ-angles by introducing the terms θ yẏ 2π and θ xẋ 2π . The path integral can be treated by the saddle-point approximation if φ is small. The main contribution comes from the perturbative saddle for which x = y = 0 at any time t. This solution does not break the translational symmetry on the time-circle,

JHEP01(2019)079
and all its modes are Gaussian. The perturbative partition function can be expanded in powers of φ using the Feynman diagrams. The result will be the perturbative partition function which we denote as Z 0 . In turn this is related to the perturbative energies as follows where E pert (N ) is the perturbative energy at level N . On the other hand, the contributions of the partition function can be classified by their topological winding number, i.e.
where Z 0 is the expansion around the trivial saddle point (i.e. the perturbative vacuum), and it is responsible for perturbative contributions E pert (0) while Z n =0 come from different instanton sectors (n counts the instanton number). The constant C above may be UV divergent, and may be removed by the appropriate definition of the path integral measure. Further all Z n -s are UV divergent. However all the UV divergences are the same, and soẐ n is UV finite. The constant C therefore factorizes, and is of no physical consequence as it cancels in the observables. The dilute instanton gas approximation makes now the following assumption: the multi-instanton contributions factorize to 1-instanton contributions. Sô Summing over n we simply have Now theẐ ±1 is given bŷ where K is the measure of the 1-instanton configurations, including the perturbative corrections, and θ is the relevant θ-angle coupling to the instantons. 9 Therefore the 1-instanton correction to the ground state energy is given by To get this correction we need to compute K. 9 In the Harper-Hofstadter problem we will have two types of instantons which tunnel in x-and ydirections respectively. So we may have two θ-angles: θ = θx or θ = θy coupling to the tunneling events x → x + 2π and y → y + 2π. Recall that these θ angles can only be defined when 2π/φ ∈ Z.

JHEP01(2019)079
Let us first consider the partition function Z(β, θ) in the trivial vacuum given by x = y = 0 by expanding in x and y up to quadratic terms and performing the Gaussian integral to get Now we consider the 1-instanton sector. For this purpose, we need to solve for the 1-instanton configuration. The equations of motion for the partition function (3.6) is We solve these equations in the appendix A to give the 1-instanton solution where t 0 is a free parameter interpreted as the center of the instanton. Note that x 1 (t) starts from 0 in t = −∞ and reaches 2π in t = +∞, and thus it indeed has topological charge 1, while y 1 (t) is always imaginary and its imaginary value reaches the maximum cos −1 (3) at t = t 0 . This means that we are considering the instanton tunneling in the x-direction. We call it an x-instanton. We plot x 1 (t) and −iy 1 (t) in figure 2. The profile of an anti-instanton is obtained by simply the time-reversal transformation. 10 We also notice that the Hamiltonian is constant as it should be, with the help of the e.o.m. (3.18a), which will be of use later. There exists in fact another type of 1-instanton due to the fact that the Hamiltonian function is also periodic in y. In the example of the Harper-Hofstadter model, one can easily find the new instanton due to the symmetry of the theory under the map Applying this map to the instanton solution (3.19), we get a new instanton solution with the x-and y-profiles exchanged (up to a minus sign). We call it a y-instanton, since it has a non-trivial topological charge in the y-direction, but a trivial topological charge in the x-direction. This instanton does not couple to θ x . Instead it couples to the θ y -angle. Let us compute the action of the 1-instanton configuration (3.19), in the limit β → ∞.

The action of the instanton is computed analytically in appendix A and it reads
where C is the Catalan's constant. Now we compute the one-loop partition function in the 1-instanton sector, by performing the expansion dt cos x 1 · δx 2 + cos y 1 · δy 2 − 2i δẋδy .
(3.24) We can first integrate out δy. However notice that in doing so we will get a nontrivial factor in front of the path-integral, because the coefficient of δy 2 is not a constant. To avoid this, let us first replace δỹ = √ cos y 1 δy and δx = δx/ √ cos y 1 . 12 Notice that this replacement keeps the measure invariant i.e. D(δx)D(δỹ) = D(δx)D(δy). Upon integrating out the δỹ, we get 12 Note that cos y1 > 0, because y1 is purely imaginary on the instanton trajectory.

JHEP01(2019)079
The operatorÕ has a zero mode given by ψ 0 (t) = N −1ẋ1 (t) √ cos y 1 , as can be checked. Here is the normalization factor. So the above expression of the one loop weight of the instanton cannot be correct. The zero mode originates from the time-translation symmetry of the theory. In other words, field fluctuations which only change the location of the instanton do not change the action, and the modes in this direction must be treated exactly (i.e. beyond the Gaussian approximation).
To find the measure of the instanton we must first separate out the zero mode, which we denote by t 0 . We will get that where the prime indicates that the zero mode has been excluded from the determinant. The µ above is the measure of the instanton moduli t 0 (or the moduli space metric). It is given by (see appendix B) so that the one-loop instanton contribution to the partition function is given by The contribution is of course divergent, as the functional determinant is infinite in the continuum. We therefore normalize it with respect to the perturbative partition function. The normalized 1-instanton partition function is given bŷ Comparing with (3.14), we find that the prefactor K entering formula (3.15) is given by 13 As we show in the appendix C, the ratio of determinants is given by .
(3.33) Note that in obtaining the above result, we have used Dirichlet boundary conditions for the space of function acted on by the operators O and O 0 . Since we will only be interested in 13 A possible minus sign can be absorbed into the θ angle.

JHEP01(2019)079
the limit β → ∞, the boundary conditions will not matter. However if one is interested in computing the instanton contributions to higher energy levels, a computation with periodic boundary conditions is necessary.
Further since we only care about the limit β → ∞, we can make convenient approximations. We notice that the 1-instanton configuration (A.6), (A.7) has the following asymptotic formẋ where Besides, the integrand of the integral over t is small when t is close to ±β/2, so the integral over t is saturated away from them. So regarding the two integrals over t and t , only the − β 2 t β 2 region is important. The two integrals can be approximated by Pulling these two integrals out of the integral of t, the latter becomes (ẋ 1 / √ cos y 1 , Apply (3.34) in the remaining part of the determinant evaluation, we find in we get that the factor K in (3.32) is given by (3.38) The anti-instanton partition function is the same but with the opposite topological charge. Therefore using (3.15), the leading order 1-instanton correction to the ground state energy given by x-instanton coupled to θ x is Since we have two kinds of instantons coupled to θ x and θ y respectively, the full 1-instanton correction is finally given by We will see that it indeed agrees with the numerical results given in (3.43).

Comparison with numerical analysis
In the previous subsection, we gave the path integral analysis in the Harper-Hofstadter model. The spectrum of the Harper-Hofstadter Hamiltonian turned out to receive nonperturbative corrections in the weak flux limit. The eigen-energy takes the form of the trans-series expansion (3.3). As was seen in the previous subsection, the one-instanton sector consists of x-and y-instantons and their anti-instantons. From the isotropy, it should take the form as where A is the instanton action given by (3.22), is an unknown coefficient, and P 1-inst fluc is the perturbative fluctuation around the one-instanton saddle. We have computed N (1) (0, φ) in the previous subsection, but it is not easy to compute it for excited states in the same way. In this subsection, we predict N (1) (N, φ) and P 1-inst fluc from the numerical analysis. To extract one-instanton contribution, the best way is probably to investigate the width of the energy bands. As is clear from the trans-series form of the eigen-energies (3.3) with (3.41), the band width is controlled, at the leading order, by the one-instanton sector. From the angle dependence of (3.41), one can easily see that two band edges correspond to (θ x , θ y ) = (0, 0) and (θ x , θ y ) = (π, π). Therefore the band width is given by where · · · represents the higher instanton contributions, which are irrelevant here. We compute this band width for various values of φ = 2π/Q and N using the exact formula (2.28), and fix unknown parameters in the formula above by numerical fitting. The strategy is the same as the one used in [31]. We refer the reader to this work for details. We first confirm that the exponential decay of the band width in φ → 0 is actually explained by the instanton action A = 8C. By the numerical fitting, we also find the explicit form of the coefficient N (1) (N, φ): For the lowest Landau level N = 0, this is indeed in agreement with the path integral result (3.40). Finally, we make a comment that with numerical calculation we can also go beyond the leading order contribution. With the help of the Richardson transformation, as explained in [31], we find the fluctuation in the 1-instanton sector to be While checking this result is very hard with instanton calculus of path-integrals, our topological string theory analysis of section 5 essentially obtains the same result.

Two-instanton calculation
Now we wish to go beyond the dilute instanton gas approximation, and compute the contributions of the two-instanton sector to the leading order in semi-classics. Recall that we have two types of instantons, which we will call I x and I y , where I x is a tunneling event in x, i.e. it takes x → x + 2π, while I y is a tunneling event in y → y + 2π.
We will consider all kinds of two-instanton events, ranging from "pure" correlations Before computing their interactions, we should stress that the contribution of such events has long been subject to debates. Particularly tricky is the instanton-anti-instanton contribution [IĪ], which is a priori ill-defined. This is because when instanton and anti-instanton are close to each other the configuration is indistinguishable from the perturbative vacuum, and it is not clear how such configurations should be taken into account (see [29] for an incomplete list of references on the topic).
If we naively superpose the well-separated instanton and anti-instanton, where we label their separation by τ , the action will be an increasing function of τ . Such a configuration spends most of the time in one of the vacua (say x = 0) and then tunnels to the other vacuum (x = 2π), lingering there for the time τ , and then returns back to the original (x = 0) vacuum. The action of such a configuration is approximately (as we will demonstrate in the next subsection) where the exponential contribution is the "classical" interaction 14 of the instanton-antiinstanton pair. The contribution of such a class of configurations to the partition function would then be 15 where φ is the coupling constant, t 0 is the "center of mass" location of the pair, and K is the one-loop measure of the individual (anti-)instantons. The integral over t 0 will simply produce one power of β, while the rest of the expression will be related to the IĪ contribution to the energy. The integral over their separation is, however, an awkward operation. As we shall see in the next subsection (4.18), the interaction constant B is negative, so the integral is saturated by its lower limit τ ∼ 0, where the approximations of the above expression are invalid, and where the notion of the instanton-anti-instanton is ill defined.
14 The term "classical" is used to reflect the 1/φ dependence of the interaction, but it is a bit of a misnomer, because an instanton-anti-instanton event is in fact a large-quantum fluctuation, and is in no way classical. 15 The subtracted unity is to control the IR divergence due to the uncorrelated instantons. Since uncorrelated instantons have already been taken into account by the instanton gas approximation it should be subtracted here to avoid double counting.

JHEP01(2019)079
Bogomolny [32] and Zinn-Justin [33,34] argued long time ago that the ill-defined IĪ amplitude is connected with the ambiguity of the Borel sum of the perturbation theory. They correctly argued that the definition of the IĪ amplitude must be ambiguous in the same way that the perturbation theory is. A prescription which is now dubbed the Bogomolny-Zinn-Justin (BZJ) prescription, is to take the coupling φ to be negative, so that the above integral is saturated away from τ ∼ log(1/φ) 1, where the approximations are valid. The above integral over τ is then performed to produce a correction to the energy where γ E is the Euler's constant, and Γ(•, •) the incomplete gamma function. The last term is exponentially small when φ < 0 so it is normally dropped. Further the expression is ambiguous if we now send φ from negative to positive values in the upper or lower complex half-plane, because of the appearance of the log. Moreover the ambiguity is exactly canceled by the ambiguity in the Borel sum of the perturbation theory. This was one of the great successes of resurgence in quantum mechanics and our understanding of its relationship with path-integrals. The BZJ prescription, however revolutionary, causes some unease. Perhaps the most uncomfortable one is that it requires dropping a factor which is exponentially small when φ < 0, but becomes exponentially large when the correct limit φ > 0 is taken. In recent years it became increasingly evident that at the heart of the correct interpretation of the BZJ result is the Picard-Lefschetz theory -a generalization of the steepest decent method to multi-integral (or indeed path-integral) cases. In fact it was only recently that a resolution of this puzzle was proposed by the interpretation of the instanton-anti-instanton pair as a saddle point at infinity [29], which establishes a concrete method for a systematic calculation of the semi-classical expansion in path integrals. The procedure is roughly as follows (We refer the reader to [29] for details.): 1) Consider an instanton-anti-instanton configuration for the case of finite time β.
2) Note that if the instanton and the anti-instanton are at opposite ends of a temporal circle, the configuration becomes a saddle point. Since the action can be decreased by bringing the pair closer together, the saddle point in question is "unstable".
3) Treat the saddle point with Picard-Lefschetz theory, i.e. instead of integrating over a cycle of real instanton-anti-instanton separation, replace the cycle with the Lefschetz thimble integral (i.e. the "steepest decent cycle"), along which the action is monotonically increasing.
4) Note that the imaginary part of the thimble integral is ambiguous depending on whether Im φ is greater or smaller than zero, and that the ambiguity cancels the Borel sum ambiguity of the path-integral, while the real part is identical to the BZJ result above, provided that we drop the incomplete-gamma term.
In particular the ambiguity, which comes from the imaginary part, is given by

JHEP01(2019)079
where we used our result (3.38). We would like to point out that the ambiguity does not contain the interaction term for the ground-state energy, i.e. it is independent of the constant B which parametrizes the instanton-anti-instanton interactions. This is in fact clarified by the thimble integration procedure in [29], summarized above. The ambiguity comes from the vicinity of the critical point at infinity, which, for a finite temporal extent, is the instanton-anti-instanton pair at opposing ends of the temporal circle. Since the saddle is "unstable" with regards to the perturbations in the real field space, the proper thimble integration will force us to integrate along the direction of imaginary separation, 16 inducing an imaginary factor in the result. This is the ambiguity, and in this case it is saturated in the vicinity of the IĪ saddle. When we take β → ∞, this vicinity of the IĪ saddle moves to infinity, where the instanton and anti-instanton are decorrelated, and all dependence on the interactions vanishes.
On the other hand, the real part is given by which depends explicitly on the instanton interaction parameter B, and which compared to corresponding terms in (4.5) has already changed sign inside the logarithm. We leave the computation of the interaction parameter to the next subsection.
Let us now consider other two-instanton events with nonzero topological charges. As opposed to quantum mechanics these in addition to the pure types (second line in (4.1)) include also the mixed types (4.2). It is straightforward to repeat the same analysis as for the instanton-instanton events, and it yields always more or less the same results. However, a crucial difference is that, as we shall see in the next subsection, the interaction constant B is positive if the two instantons are identical. In this case we have no need to change the sign of the coupling in the logarithm and the ambiguity is absent. The energy correction contains the real part only, which is (4.7) with B replaced by −B. Furthermore the mixed type as in (4.2) also must be present on physical grounds, because we expect terms of the kind cos(θ x ) cos(θ y ) to be present in the ground-state energy θ-dependence. Quite unexpectedly their interactions are all found to be purely imaginary, and thus the individual corresponding energy correction is complex. Nevertheless, when all eight mixed 2-instanton events in (4.2) are considered, the total energy correction is real.
In fact, according to the logic of [29] all these contributions should correspond to exact saddles of the QM problem on a compact S 1 time. What is especially interesting is that while the instanton-instanton and instanton-anti-instanton events of the same type (4.1) have their counterpart in quantum mechanics and are thus not very surprising by analogy, 17 the mixed type (4.2) are a different matter. Yet as we shall see their naive BZJ amplitudes 16 The contour IĪ separation parameter τ along the thimble eventually bends and becomes parallel to the real axis in the complex τ -plane, which gives the real contribution. 17 Such saddles can be thought of as the motions in a periodic inverted-potential which either oscillate with a period β between two peeks of the inverted potential (i.e. between two classical vacua of the potential)a saddle that corresponds to an instanton-anti-instanton pair -or roll with the "energy" slightly higher than the peak of the inverted potential so that in precisely one period of the imaginary time β the particle winds twice -a saddle corresponding to an instanton-instanton event.

JHEP01(2019)079
are in an extremely good agreement with the numerics, so we are inclined to believe that such saddles must also exist.

Instanton interactions
This section is devoted to the determination of instanton-interaction constants B in various instanton events, and the corresponding correction to the ground state energy. To start with, we would like to write down a general formula with any given choice of correlations. The ansatz goes as follows: let's consider a superposition of two instanton events, where x α and x β (y α and y β ) are either instanton or anti-instanton in the x (y) direction. We shift the solutions to separate two events (x α , y α ) and (x β , y β ) by τ 0. We further make the assumption that the "tunneling" of (x α , y α ) takes place when t 0, while for (x β , y β ) it happens when t 0. As a consequence, it shows that x α (0) or y α (0) differs from x α (+∞) or y α (+∞) by exponentially small terms, while x β (0) or y β (0) differs from x β (−∞) or x β (−∞) by exponentially small terms.
The two-instanton action can be split into two parts, For ease of notation, we denote the two terms on the r.h.s. as S − and S + respectively. Let's first concentrate on S − . Our assumption implies that δx = x β − x β | −∞ and δy = y β − y β | −∞ can be treated as small perturbations in this region. Let us Taylor expand S − up to the second order (4.10) Here we have used the fact that both x β (−∞) and y β (−∞) must be integer multiples of 2π. It can be shown that sum of the two terms 0 −∞ dt( 1 2 ∂ 2 x H| xα δx 2 + 1 2 ∂ 2 y H| yα δy 2 ) always has the order o(e −2τ ), while we will look for the interactions of order o(e −τ ). Indeed the reader may wonder why we even bothered to expand up to a quadratic order in δx and δy in the first place. The reason is that the terms δyδẋ is special because it is a total derivative, and will contribute a finite amount to the o(e −τ ) order. Then recall our normalization and the equations of motion (3.18a)(3.18b) as well as the conservation law (3.20), we can further simplify it to be If we choose the normalization carefully such that y β | −∞ = 0 and −i ∞ −∞ dtẋ α y α = A, we can obtain, at the leading order, a very neat formula (4.12)

JHEP01(2019)079
The same story, with the only difference that δx = x α − x α | +∞ and δy = y α − y α | +∞ are treated as small perturbations, goes for S + and yields given the choice of renormalization y α | ∞ = 0 and −i ∞ −∞ dtẋ β y β = A (again we observe that ∞ 0 dt( 1 2 ∂ 2 x H| x β (δx ) 2 + 1 2 ∂ 2 y H| y β (δy ) 2 ) has no contribution at this order). Summing up S + and S − we get the following approximation of two-instanton action (4.14) Notice that 2A is already accounted for by dilute instanton gas approximation, while the remaining part will yield the exponential contribution predicted in (4.3). Now it is time to consider concrete examples and plug in instanton solutions. First of all, it is convenient to recall the asymptotic behaviors of instanton solutions Let us go over all the 2-instanton events listed in (4.1), (4.2).
• "Pure" instanton-anti-instanton. A superposition can be chosen as to be negative.
• "Pure" instanton-instanton. The four events in the second line of (4.1) also have the same interaction, so we only need to compute [I x I x ]. A superposition can be chosen as which is positive.

JHEP01(2019)079
• "Mixed" events. We work out explicitly the [I x I y ] pair. A superposition satisfying the constraints above can be chosen as

(4.22)
Plugging into our general formula (4.14), we get On the other hand, if we consider the [I y I x ] correlation, we need to shift our solutions   With all the formulas in hand, we are able to determine various contributions to the ground state energy due to various 2-instanton events. We assume that up to 2-instantons, the ground state energy of the Harper-Hofstadter model has the most general transseries form cos θ x cos θ y , (4.27) which respects the symmetry θ x → θ y as well as θ x → −θ x , θ y → −θ y . The 1-instanton correction has already been discussed in section 3. We will check the various 2-instanton corrections in this section. Let's first look at the E IĪ 0 term. From the discussion above, we know that this term has both real and imaginary parts, given by (4.7) and (4.6) respectively for an individual 2-instanton event. Reading off the instanton interaction constant B IxĪx from (4.17), and summing up all four events in the first line of (4.1), we find the real correction is while the imaginary correction, i.e. the ambiguity is

JHEP01(2019)079
Next, E II 0 receives contribution from both I x I x andĪ xĪx events, which is the same as the sum of I y I y andĪ yĪy . Since the instanton-interaction B IxIx (4.21) is negative, the energy correction is real, and we find Finally, all the eight mixed events listed in (4.2) contribute to E IImix 0 . Although each individual event has imaginary instanton-interaction, as one sees in (4.23), (4.25), (4.26), and thus gives complex correction to the ground state energy, one can check that the imaginary contributions cancel and the total contribution of all the eight events is real. It reads

Numerical studies of two-instanton sector
In this subsection, we perform a numerical study of the various 2-instanton corrections to the ground state energy, and compare them with the predictions computed in the last subsection. And we will find perfect agreement. We confine ourselves to the real parts of the corrections, and leave the study of the imaginary part (ambiguity) to the next subsection. The trans-series (4.27) of the ground state energy already gives us a hint as how to extract 2-instanton corrections numerically. We have thus all the 2-instanton contributions can be obtained from the following linear contributions (4.33) Since the r.h.s. can be computed exactly when φ = 2π/Q for Q ∈ N, these simple linear formulas allow us to easily compute 2-instanton corrections indicated on the l.h.s. up to at least 3-instanton corrections for a sequence of Q up to very large Q, with very good accuracy for large Q. Note that here E IĪ 0 actually refers to only its real part; the imaginary value cannot be computed in this way as it cancels in physical observables. To check the imaginary part of the 2-instanton sector we can analyze the perturbation series, and match its lateral Borel sum with the ambiguity from the instantons. This will be discussed in the next section. We also notice that since E 0 (π/2, π/2) = E 0 (0, π) (cf. (2.28)), we have We comment that in using (4.33) we will make at most an error exponentially suppressed by a one-instanton factor. We demonstrate this explicitly by comparing with the results of Fourier transformation in appendix D Let us thus focus on E IĪ 0 and E II 0 . For a sequence of Q = 2π/φ, we expect improving agreement with the path integral predictions (4.28), (4.30) as Q increases. Note that from numerics, we only obtain the combination E pert 0 + E IĪ 0 , and we have to remove E pert 0 by hand by subtracting the Borel-Padé sum of the perturbative ground state energy. Poles in the Borel plane, which are responsible for the ambiguity (4.29), are dealt with by Cauchy principal value integration. This additional complication limits the range of Q we can push for. The agreement between the numerical results and the path integral predictions is excellent, as demonstrated in the matching digits plots figure 3.

Large order growth and ambiguity of energy
According to the resurgence theory, the large order growth of the perturbative energy expansion is controlled by the ambiguity (imaginary part) of energy, which receives contributions from instanton sectors with topological charge zero (see for instance [35]). The first such sector is the instanton-anti-instanton sector [IĪ] including all four events listed in the first line of (4.1). The imaginary energy correction from this sector is is the Stoke's constant, related to the ambiguity of the lateral Borel resummation of the perturbative expansion, and b N is the leading exponent of φ in the instantonanti-instanton sector. Let us denote the perturbative expansion by The resurgent analysis then suggests the following relation (4.36) We will use this relation to compute numerically the imaginary part of E IĪ . We start with the ground state with N = 0. We compute a In this process, it is convenient to use the Richardson transformation to accelerate the convergence (see for instance [36] for details). It is easy to check that these numerical estimations reproduce the exact values so that in the leading order, we have which agrees with the path integral calculation (4.29). As in the 1-instanton sector, once the analytic values of S (0) and A are fixed, numerically we can go beyond the leading order and further extract the values of a , · · · (4.41)

JHEP01(2019)079
These coefficients should give the perturbative fluctuation around the instanton-antiinstanton saddle. We repeat the same computation for higher energy levels. Observing the general structure (4.36), we find that In addition, we extract the coefficients a (4.43) From these data, we could construct the [IĪ] contribution to the imaginary part of the eigen-energy where the function A(N, φ) is nothing else but the "non-perturbative" A-function appearing in the Zinn-Justin-Jentschura exact quantization conditions [37,38] in conventional quantum mechanics. In our example, the first few terms of A(N, φ) read where ν = N + 1/2.

JHEP01(2019)079 5 Instanton fluctuation from topological string
Here we reveal an interesting connection between the fluctuation parts P 1-inst fluc , P IĪ fluc and topological string theory.
Before our analysis, we would like to remind the reader that the Harper-Hofstadter model is closely related to a Calabi-Yau threefold called the canonical bundle of F 0 , also known as local F 0 in string theory community, as first pointed out in [3]. According to mirror symmetry, all the Gromov-Witten invariants of local F 0 are encoded in an algebraic curve, called mirror curve, whose equation reads 18 e x + e −x + e y + e −y = u . One can obtain another QM model by promoting x and y without the rotation, i.e., one considers the Hamiltonian We choose a normalization of H F 0 slightly different from that in the literature to match the normalization of (2.3) we use in this paper. Motivated by topological string considerations [39,40], this QM model has been thoroughly studied, both its spectrum [4,[41][42][43][44][45][46] and its wave functions [47][48][49] (see also [50,51]). This has led to exciting development of nonperturbative completion of topological string theory and topological string / spectral theory duality [4,[52][53][54][55][56][57][58], which in turn inspired a new procedure to solve non-perturbatively QM models [15][16][17], as well as the discovery of a new class of exactly solvable deformed QM models [59]. We would like to point out that on the one hand, the Hamiltonian (5.2) and that of the Harper-Hofstadter model are rather different in nature. The former is confining and has a discrete spectrum, while the Harper-Hofstadter Hamiltonian is periodic and thus has a rich band structure. On the other hand, the spectra of the two Hamiltonians are closely related in the semi-classical regime. In fact, the perturbative eigen-energies of H F 0 was computed in [15], also using the BenderWu package [8,9], and it is easy to check that they are related to the perturbative eigen-energies of H(0, 0) by the map We will see in later sections that many results [15] also apply for the Harper-Hofstadter model as well with appropriate modification.

JHEP01(2019)079
The large order growth of the perturbative energy of H F 0 has been analyzed in detail in [15], and it is incorporated in the leading non-perturbative correction 19 to the perturbative series. It is revealed in [15] that this non-perturbative correction can be obtained from the refined free energies in the Nekrasov-Shatashvili limit of topological string theory on the Calabi-Yau threefold local F 0 . We will demonstrate that we can obtain the 1-instanton correction (and the instanton-anti-instanton correction) of the Harper-Hofstadter model from their data by applying the map → −φ. This is not obvious at first glance because the 1-instanton correction here is the half order of the non-perturbative correction in [15]. This is a consequence of the fact that the 1-instanton sector and the instanton-anti-instanton sector are closely interrelated, as suggested in [31].
Let us quickly review the results of [15] concerning the spectrum of H F 0 . The perturbative eigen-energy can be computed also by using the BenderWu package [8,9], and the first few terms read with Indeed, this agrees with the perturbative energy of the Harper-Hofstadter model (3.2) by the replacement (5.4). Note we have adapted the series of E pert F 0 (ν, ) to be consistent with the normalization of H F 0 used in this paper. To formulate the results of the formal "instanton-anti-instanton" correction, we need some terminology from topological string theory on a local Calabi-Yau manifold and its mirror curve.
The coefficient u in the equation of mirror curve (5.1) parametrizes the complex structure moduli space of the curve. The moduli space has several singular points, one of which of particular interest is called the conifold singularity and it is located at u = 4, as it corresponds to the semi-classical limit E F 0 = 0 of the QM model H F 0 . Let us introduce Then the classical periods of the mirror curve are of which t c can serve as a good local coordinate on the moduli space near the conifold singularity. Here K is the complete elliptic integral of the first kind. Furthermore, for the topological string theory on a local Calabi-Yau threefold X, an important quantity is JHEP01(2019)079 the refined free energy F (t, 1 , 2 ). It encodes the numbers of BPS states of the M-theory compactified on X × (R 4 × S 1 ) 1 , 2 , where the parameters 1 , 2 describe how the R 4 is twisted along S 1 . t is a set of coordinates on the moduli space of X, which due to mirror symmetry is mapped to the complex structure moduli space of the associated mirror curve.
In the application to the spectrum of H F 0 , one is in particular interested in the so-called Nekrasov-Shatashvili limit [60] F NS (t, ) = lim 9) and the free energy in the NS limit enjoys a genus expansion Near the conifold singularity, the NS free energies F NS n are functions of t c with at most logarithmic singularity, and we will use the notation They can be computed recursively by the so-called refined holomorphic anomaly equations [10][11][12] in the NS limit, as explained in detail in [15][16][17]. For local F 0 , the first few NS free energies are  We stress that these results are obtained purely in the framework of topological string theory. We do not need any knowledge of the corresponding quantum mechanics. Our goal is to relate these quantities to the eigen-energy in quantum mechanics. It turns out, the formal "instanton-anti-instanton" correction to the eigen-energy of H F 0 , which controls the asymptotic growth of the coefficients of E pert F 0 (ν, ), is given by [15] where C is the Catalan's constant, and f (1) a free constant. The exponential factor is e 16C/ = e 2A/ , and this indeed corresponds to the 2-instanton sector in our terminology. Using the NS free energies of local F 0 , one can write down the terms in the exponential (5.14)

JHEP01(2019)079
Interestingly, the terms independent of can be resummed to log √ 2π16 ν Γ( 1 2 + ν) . (5.15) Furthermore, let us denote the power series in starting from O( ) by Then the "instanton-anti-instanton" correction can be written as where the components after · is a power series starting from constant 1.
We observe here that this result in terms of topological string free energies also reproduces the imaginary part of the instanton-anti-instanton correction (4.44) for the Harper-Hofstadter model after applying the map (5.4). Indeed, if we do so, we find that the factor in front of · in (5.17) agrees with the prefactor before · in (4.44), if we choose Note that this normalization constant can also be fixed through the path integral calculation in section 4.1. Comparing the remaining part with the numerical result (4.47), one conjectures then the A-function should be identified with the opposite of the derivative of the NS free energy for local F 0 , i.e.
We follow the calculation in [15] of the NS free energies for local F 0 by solving the NS holomorphic anomaly equations and push it to a few orders higher than what is explicitly given in [15]. We find  20) and it agrees completely with the A-function (4.48) from the numerical fit.

JHEP01(2019)079
Finally, since the power series in the 1-instanton sector is given by the A-function as shown in (4.46), we claim that the 1-instanton sector can also be expressed in terms of the NS free energy of local F 0 . In fact, by plugging (5.19) into (4.46) for the fluctation and comparing the prefactor (3.43) with the component (5.15), we find an expression similar to (5.13) Note that after mapping → −φ the exponential becomes purely imaginary, and we take its imaginary value in the expression above.

Conclusion and discussion
In this paper, our goal is to understand the peculiar band structure of the energy spectrum of the Harper-Hofstadter model in the semi-classical limit. According to the general philosophy of resurgence, the energy levels should be written as trans-series summing over contributions from all saddle points coupled differently to the Bloch's angles, and in addition, the large order growth of the perturbative sector is controlled by the ambiguity of the energy, which receives leading order contributions from the instanton-anti-instanton sector. We used various techniques to compute the trans-series energy levels. The perturbative series is computed very conveniently using the extended BenderWu package [8,9]. The 1-loop contributions to the 1-, 2-instanton sectors, and the energy ambiguity are obtained by a path integral calculation, albeit restricted to the ground state level. Higher order corrections in the 1-instanton sector and in the ambiguity (imaginary contributions of the instanton-anti-instanton sector) are computed using refined topological string techniques in connection with the local F 0 geometry inspired by a similar work [15]. All these results can be checked against numerical results, which can be computed exactly when the magnetic flux is 2π times a rational number, and they all agree perfectly. This validates all our techniques. In the process, we find that the perturbative-non-perturbative relation 20 relating the perturbative sector and the 1-instanton sector is not satisfied, 21 which is not that surprising since the Schrödinger equation of the Harper-Hofstadter model is a difference rather than a second-order differential equation. On the other hand there still exists a curious relation between the three sectors: perturbative, 1-instanton, instanton-anti-instanton.
Clearly there are many open questions. The most pressing one is how to understand better the nature of instantons, and in particular the [I x I y ] instanton configuration, the treatment of which is rather ad hoc in this work. Namely the inclusion of the saddles, such as instantons, is expected to be dictated by the Picard-Lefschetz theory, and requires a decomposition of the path-integral cycle into the Lefschetz thimbles. We have not rigorously checked whether instanton configurations we analyze are a part of this decomposition, but have argued that they should contribute on physical grounds. The case of [I x I y ] is JHEP01(2019)079 particularly interesting, as it is an object without an analogue in simple 1D quantum mechanics. Naïve application of the BZJ prescription yields a result in perfect agreement with the numerics. On the other hand the BZJ prescription was interpreted as contributions from saddles at infinity [29]. It would be desirable to understand this better in the case at hand. It would also be nice to have a path integral understanding of the higher order corrections computed using topological string techniques. Furthermore, another real world model, one that describes electrons on a triangular lattice, is revealed to be connected to the topological string theory, with the target space being the canonical bundle of the three-point blow-up of P 2 [5]. One can explore whether the similar analysis can be applied in that model as well.

JHEP01(2019)079
energy E(β) reaches the maximum value in the limit β → ∞, and it corresponds to the oscillation between two neighboring highest points of the inverted potential. In (A.3), we have E(∞) = 2 regardless of the sign in the inverted potential, so we simply take + without loss of generality 1 +ẋ 1 2 + cos x 1 = 2 . (A.5) Solving (A.5), we find the following profile of 1-instanton as well as .
Using the conservation law cos y 1 + cos x 1 = 2 , (A. 8) we find that the action is given by In the last line we performed the change of variables t = sin x/2, and used one of the definitions of the Catalan's constant

B The moduli-space metric
We want to find the moduli-space metric of the one instanton. We can do this by adding a factor λδx 2 to the action (3.25) before integration, so as to lift the zero mode. Upon modifying (3.25) by adding such a term, we can do the Gaussian integral and simply get Now let us write the small deviation around the instanton solution as

JHEP01(2019)079
where δx is orthogonal to cos y 1ẋ1 . The first term is a small deviation from the instanton solution in the direction of the zero mode, and t 0 specifies a shift of its position in time.
In fact t 0 is precisely the coordinate we want to isolate, and over which we will integrate exactly, producing a factor of β. Recall that our goal is to find a way to write the path integral in (3.25), as where the prime indicates that the zero-mode has been excluded from the determinant. The µ above is the measure of the zero-mode moduli t 0 (also referred to as moduli space metric), which is what we wish to find. To find it we will add the term λδx 2 into the action as before, and integrate over t 0 . We should get (B.1), up to a constant, which will precisely correspond to µ −1 . To do this, let us plug in the expression (B.2) for δx into the path integral (3.25). It only amounts to adding the term λδx 2 into the action, since the zero mode is annihilated byÕ. Then it is easy to see that the action contains the term where N is given by (3.27). If we now integrate over t 0 and δx ⊥ we produce a term

C The one-instanton determinant
We will compute the determinant of the one-instanton fluctuation operator using the Gel'fand-Yaglom theorem, explained for instance in [67][68][69][70]. Consider an ordinary differential operator O, with a canonical second derivative term O = −∂ 2 t + . . . . We wish to compute the determinant of the operator. For that purpose we consider the space of functions on which the operator acts to be defined on an interval t ∈ [−β/2, β/2] with the The proportionality identity can be made precise by regularizing the operator determinant with that of a simple operator, for instance, the harmonic oscillator where Ψ 0 (t) is the zero mode of the harmonic oscillator O = −∂ 2 t + 1 with the boundary condition (C.4), and it is simply Ψ 0 (t) = sinh(t + β/2) . (C.6) To treat det O with zero mode removed, we can use the relation Therefore we need to compute the zero mode of O λ satisfying the boundary condition (C.4) up to order λ. Now we could take the operator O to simply be the fluctuation operatorÕ given by (3.26). However notice that we have

JHEP01(2019)079
where f (t) is an arbitrary, nonsingular function with no zeros. By taking the derivative with respect to λ and setting λ = 0 we get . (C.10) If we take f (t) = cos y 1 (t) we can define the operator so that we will compute det (O) instead of det (O).
In order to compute it we first have to consider the determinant of det O λ , where O λ = O + λ, at least for small λ. We already know that O has a zero mode given byẋ 1 .
To use Gel'fand-Yaglom theorem we look for a solution where Ψ λ satisfies (C.4). Now assuming λ is small we can write where OΨ (0) = 0 (C.14) and The first of these equations reduces to OΨ (0) = cos y 1 (t) −∂ t 1 cos y 1 (t) ∂ t + cos x 1 (t) Ψ (0) = 0 . (C. 16) This is a second order ODE, and we already know that one solution is ψ 1 (t) =ẋ 1 (t) , (C.17) although it does not satisfy the boundary condition (C.4). In order to find a second independent solution, we notice that the operator O can be factorized in the following way. We introduce operators Q = 1 cos y 1 ∂ t − i sin x 1 sin y 1 , Q † = 1 cos y 1 ∂ t + i sin We want to find the most general homogeneous solution to the equation Oψ = 0. This is the same as finding such a solution for the operator Q † Q. We observe that Q † annihilates JHEP01(2019)079 1/ẋ 1 . If one can find ψ 2 such that Qψ 2 = 1/ẋ 1 , then one concludes immediately from (C. 19) that ψ 2 is another solution to (C.16). Indeed by making an appropriate ansatz we find ψ 2 (t) =ẋ 1 (t) t dt cos y 1 (t ) x 2 1 (t ) . (C.20) Furthermore since the Wronskian is not identically vanishing W 21 (t) := ψ 2 (t)∂ t ψ 1 (t) − ψ 1 (t)∂ t ψ 2 (t) = − cos y 1 (t) , (C. 21) the two solutions are linearly independent. From ψ 2 (t) we can construct the solution to (C.16) satisfying the boundary condition (C.4) .

D Comparison of Fourier analysis with linear formulas
In the section 4.3 we have used (4.33) to compute the 2-instanton θ x , θ y -dependence of the system. These formulas are expected to have an error exponentially suppressed with the coupling. Here we analyze E pert 0 + E IĪ 0 , E 1-inst 0 , E II 0 , E IImix 0 by a direct Fourier transformation in order to check explicitly the validity of the formulas (4.33) and make an estimate on the error. To be explicit, we should have 0 dθ x dθ y cos(θ x ) cos(θ y )E 0 (θ x , θ y ) .

(D.1)
Here we can integrate over θ ∈ [0, π] instead of [0, 2π] because E 0 is an even function of θ x , θ y . We notice that the energy level E 0 is solved from the equation (cf. (2.28)) F Q (E 0 ) = 2(cos θ x + cos θ y − 2) (D.2) where F Q (E 0 ) for φ = 2π/Q is a polynomial of degree Q in E 0 . It is simpler to perform the integration if we can exchange the integration variable from θ x , θ y to E 0 . This can be done in the following way.
To compute E 1-inst 0 , we take θ y = π/2 and E 0 only depends on cos θ x = cos θ. The relation (D.2) is simplified to where E(•) is the complete elliptic integral of the second kind, and conclude in the end (D.14) Equations (D.11), (D.4), (D.5), (D.14) then provide us the means to compute the instanton contributions as Fourier coefficients of the ground state energy. We compare the results of different 2-instanton corrections computed by the two different methods in figure 4. Since we focus on the 2-instanton sector, we only include the corrections of E pert 0 + E IĪ 0 , E II 0 , and E IImix 0 . We find the agreement to be remarkable, with the relative difference to be at 2-instanton or even 4-instanton levels, and is thus negligible.