Quantum Criticality of Generalized Aubry-André Models with Exact Mobility Edges Using Fidelity Susceptibility

In past decades, the study of phase transitions in disordered systems has gained upsurging interest in condensed matter physics. Anderson’s pioneering theory of localization in noninteracting disordered systems has achieved remarkable success [1]. Such systems may exhibit many-body localization (MBL) in the presence of many-body interactions [2]. Recently, the phenomenon of MBL was investigated extensively from both theoretical and experimental perspectives [3, 4, 5]. A key insight of current phenomenological theory is that MBL systems exist an emergent integrability, which is characterized by a set of quasilocal integrals of motion. The emergent integrability plays a similar role of exact integrability in integral systems, and leads to a breakdown of thermalization in isolated systems. Interestingly, the strong ergodicity breaking arising from integrability and disorder fails to describe the recent discovery of quantum many-body scars, which triggered intensive investigations of weak ergodicity breaking in a variety of quantum systems  [6, 7, 8]. Compelling theoretical and numerical evidence suggests the phenomenon of localization and non-ergodicity can exist in quasiperiodic models. Compared to random disorder, which poses challenges in experimental engineering [9], incommensurate potential systems are more accessible for experimental platforms, particularly in the context of ultra-cold atoms [10], photonic crystals [11] and cavity polaritons [12, 13], providing a fertile ground for investigating incommensurate potential driven localization. The dispensation of disorder averaging has conduced to an in-depth theoretical investigations in quasiperiodic systems. Despite many physical features in common, such as multifractal wave functions and exponentially-diverging localization length, the critical properties of the localization transition in random and incommensurate models are rather different. As the random matrix theory is not applicable to describing incommensurate models [14], there are no rigorous bounds on the critical exponents concerning the localization transition [15]. Considering all eigenstates of Anderson model will be localized for an infinitesimal disorder strength in one and two dimensions (1D, 2D), a distinct feature of incommensurate model is the access to a localization-delocalization transition in 1D. Therefore, quasiperiodic models open the door to effectively simulate some properties of Anderson localization that occur generically in three dimensions (3D) and offer a fascinating platform to explore unconventional phenomena, including fractal bands [16, 17], mobility edge (ME) [18], Bose glass [19], topological phases [20], and non-Hermitian Floquet quasicrystals [21].

Among a diversity of quasiperiodic models, the Aubry-André (AA) model and its generalizations are drawing increasing attention nowadays due to its experimental realizations [22, 23]. The 1D AA model was originally derived from the reduction of a square-lattice Hofstadter Hamiltonian describing conduction electrons in the presence of a uniform magnetic field [24]. A dual transformation between the coordinate and momentum spaces will yield the same Hamiltonian with hopping and potential amplitudes interchanged. As a consequence of the inherent self-duality, all the eigenstates of the system undergo Anderson localization transition from an extended phase (EP) to a localized phase (LP) at a critical value of the incommensurate potential strength [25, 26]. At the transition, the energy spectrum showcases the renowned Hofstadter butterfly pattern [27]. The multifractal gap structure that ultimately evolves into the localization of the particles above the transition point. During the past few years, inspiring progress is made by tailoring the incommensurate potential to introduce a single-particle ME into the AA model [28], which separates localized and delocalized eigenstates at a critical energy. Typical variants encompass slowly varying Aubry-André model [29, 30], mosaic Aubry-André model [31] and generalized Aubry-André model (GAA) [32].

Motivated by the recent experimental progress, specifically the successful simulation of the GAA model by using quantum Cs gas [33] and using synthetic lattices of laser-coupled atomic momentum modes [34], we consider critical properties of the GAA model that has a single-particle ME at a finite filling. Our primary focus is on the scaling behavior of the fidelity susceptibility across the localization transition between the EP and the EP at different filling fractions. In the context of the Widom scaling hypothesis [35], which is foundational to the field of critical phenomena, it is posited that a complete characterization of a universality class of a phase transition can be achieved using just two independent critical exponents. The finite-size scaling of the generalized participation indicated that the correlation-length critical exponent $\nu$ = 1 for the AA model, in which an incommensurate on-site potential is characterized by an irrational number $\alpha$. A second independent critical exponent thus plays a decisive role in determining the universality class. The ground state universality class of quasiperiodic models has been thoroughly investigated using a renormalization group scheme [36, 37]. It has been recognized that the universality class of the ground state depended solely on the continued-fraction expansion of $\alpha$ [38]. Furthermore, the dynamic critical exponent $z$ exhibits a significant dependence on irrational numbers. The value of $z$ at different filling fractions dictate the scaling behavior of the corresponding energy spectrum [39, 40, 41], which proves valuable in comprehending both the multifractal characteristics of the system [42] and the properties associated with multiparticle ground-state transitions. Notably, the relationship between $z$, $\alpha$, and the filling fractions $\rho$ is less thoroughly understood. A strategy using generalized fidelity susceptibility (GFS) to obtain the dynamical exponent $z$ has been developed for the lowest state in the both CP-LP [43] and EP-CP transitions [44] of the $p$-wave-paired AA model.

This paper is organized as follows. In Sec. II, we briefly present the GAA models and the phase diagram. In Sec. III, we show the ME of GAA model by inverse participation ratio and fidelity susceptibility. In Sec. IV, we review the scaling hypothesis of fidelity susceptibility and present the numerical results of the critical behaviour in the GAA model. We also unravel the relationship between the Diophantine equation and the filling fractions. The summary and conclusion are given in Sec. V.

We consider a family of GAA models described by the following Hamiltonian

in which the operator $c_{j}^{\dagger}$ ($c_{j}$) represents the creation (annihilation) operator at site $j$ among total $N$ lattice sites, and the hopping amplitude $t$ is set to be the unit energy scale, i.e., $t=1$. The on-site potential $V_{j}$ is characterized by a modulation strength $V$, a wave number $\alpha$, a deformation parameter $b$ and a random phase $\phi$, given by

When $\alpha$ takes an irrational value, the modulation becomes quasi-periodic. In this work, we employ $\alpha=(\sqrt{5}-1)/{2}$. The on-site potential is a continuous function of the parameter $b$, where $b\in[0,1)$. When $b=0$, the model (1) reduces to the standard AA model. For $V<t$ all the single-particle eigenstates are extended, while they become localized simultaneity for $V>t$. At $V=t$, the system undergoes the Anderson localization transition, leading to all eigenstates becoming critical and exhibiting multifractal behavior. This localization transition at $V=t$ is a consequence of a self-dual symmetry present in the system, which indicates the absence of the ME within the system. For $b\neq 0$, the deformed AA Hamiltonian in Eq. (1) respects an energy-dependent duality symmetry, which ensures the existence of localization transition. The GAA duality warrants the occurrence of the ME described by a closed-form analytical expression:

where $E$ represents the single-particle energy. According to Eq. (3), we depict the phase diagram in Fig. 1, in which the single-particle eigenstates are categorized into three types: completely localized, completely extended and mobility-edge phases. In particular, low-energy states are extended in the mobility-edge phase. In this study, we investigate the phase transition from extended eigenstates to localized eigenstates by varying the filling fractions, denoted as $\rho={n}/N$, in the mobility-edge phase, where $n$ represents the index of the highest occupied eigenstate by fermions.

In order to obtain the eigenenergies of the system and the associated eigenstates on a finite lattice containing $N$ sites, we numerically diagonalize the model Hamiltonian in Eq. (1). We impose antiperiodic boundary conditions to alleviate the boundary effects [45, 46]. In the actual numerical calculations, it is wise to replace the inverse golden ratio with a rational approximant being $F_{k-1}/F_{k}$, which will approach the inverse golden ratio for $k\to\infty$, i.e., $\alpha=\lim_{k\to\infty}F_{k-1}/F_{k}$. Here $F_{k}$ represents the $k$-th Fibonacci number, which is recursively defined by $F_{k+1}=F_{k}+F_{k-1}$, with $F_{0}=F_{1}=1$. The Fibonacci sequence, ${F_{k}}$, demonstrates a three-cycle parity pattern: odd, odd, even. As a result, the sequence is partitioned into three distinct subsequences, $\{F_{3l+1}|$ $\cdots,55,233,987,\cdots\}$, $\{F_{3l+2}|$ $\cdots,89,377,1597,\cdots\}$, $\{F_{3l+3}|$ $\cdots,144,610,2584,\cdots\}$.

At the $k$-th rational approximation, the GAA potential and the Schrödinger equation have a periodicity of $F_{k}$ sites, allowing for experimental realizations with available laser wavelengths. Then, it is convenient to calculate the physical quantities of interest on finite lattices. The eigenspectra and eigenstates can be determined by solving the Schrödinger equation

where $\ket{\Psi_{n}}=\sum_{j}\psi_{n,j}|j\rangle$ expressed in terms of Wannier basis $|j\rangle$, representing the state of a single particle localized at the site $j$ of the lattice. The Hamiltonian (1) can be recast into

where $\hat{\mathbf{c}}=[\hat{c}_{1},\cdots,\hat{c}_{F_{k}}]^{T}$ and $\mathcal{H}$ is a $F_{k}\times F_{k}$ single-particle Hamiltonian matrix, given by

In order to visually characterize the localized and extended nature of the entire energy spectrum, we evaluate the normalized inverse participation ratio (IPR) for $n$-th eigenstate $\ket{\Psi_{n}}$ of the model (1), given by [47, 48, 49]

The IPRs can quantify the extent of distribution over the Fock space. For a set of single-particle states in real space, the IPR scales inversely with the system size $N$ in extended states, while appears to be independent of $N$ in the localized states. We plot the logarithm of normalized IPR for all eigenmodes of Eq. (1) as functions of $V$ with the deformation parameter $b=0.1$ in Fig. 2. Hereafter, unless otherwise specified, we set $b=0.1$. The localization transition from the EP to the LP occurs at mobility edges in Eq. (3). In the following, our primary emphasis will be on two specific cases: one where the fermion number is ${n}=\lfloor N/2\rceil$, corresponding to a filling fraction of $\rho=1/2$, and another where the fermion number is ${n}=1$, corresponding to a filling fraction of $\rho=1/N$. Here $\lfloor\cdot\rceil$ denotes the round function.

As the IPR is bounded, it is not easy to accurately identify the critical point, as is revealed in Fig. 3(a). Here we use quantum fidelity susceptibility (QFS) as a versatile indicator of localization transition points. The fidelity susceptibility quantifies the sensitivity of an eigenstate to variations of the driving parameter, which is defined as the second-order derivative of the quantum fidelity with respect to the driving parameter [50]. Being a purely geometric measure, the fidelity approach offers a notable advantage by not requiring any prior knowledge of the order parameter. A GFS of order $2r+2$ associated with the eigenstate $\ket{\Psi_{n}}$ with the tuning parameter $V$ is given by [51]

The numerator in Eq.(15) denotes the probability of exciting the system away from the state $|\Psi_{n}\rangle$ through a relevant (or marginal) perturbation $\partial_{V}\hat{H}$. The GFS of different orders is encoded by the power of the denominator. For the case of $r=-1/2$, Eq. (15) simplifies to the second derivative of the eigenenergy [52], denoted as $\chi_{1}^{(n)}\equiv\partial^{2}E_{n}/\partial V^{2}$. Similarly, when $r=0$, Eq. (15) reverts to the QFS, i.e., $\chi_{2}^{(n)}\equiv\partial^{2}(1-F)/\partial V^{2}$, where the fidelity $F$ is defined as the absolute value of the inner product between $|\Psi_{n}(V)\rangle$ and $|\Psi_{n}(V+\delta V)\rangle$. This fidelity quantifies the Bures distance between two states as their tuning parameter slightly changes.In the following, we focus on the GFS of the eigenstate $|\Psi_{n}(V)\rangle$ at a finite filling fraction for $n=\lfloor\rho F_{k}\rceil$, where $\lfloor x\rceil$ rounds $x$ to the nearest integer. To this end, the superscript of $\chi_{2r+2}^{(n)}$ is omitted for abbreviation.

Different from the smooth transition of the IPR around the mobility edge, the QFS shows a divergent peak in the vicinity of the EP-LP transition point due to the enhanced sensitivity, as is disclosed in Fig. 3(b). The QFS and GFS have been widely used to witness the quantum phase transitions in quantum systems [53] and even holographic models [54, 55]. They have proved to be particularly useful for detecting the critical points of a symmetry-knowledge unknown systems [56, 57, 58]. It has been demonstrated that the QFS not only facilitates the identification of quantum critical points (QCPs) but also adheres to the scaling ansatz, enabling the retrieval of critical exponents. The retrieved critical exponents conform to the hyperscaling relations, signifying the existence of only two independent exponents. The key implication is that the finite-size scaling behavior of the standard fidelity susceptibility $\chi_{2}$ determines the positions of QCP, along with the correlation-length critical exponent $\nu$ associated. Notably, a second critical exponent assumes a pivotal role in establishing the universality class such as the dynamical exponent $z$. In the context of quasiperiodic systems exhibiting spatial complexity, the scaling behavior of the QFS and the universality class of the localization transition have received relatively less exploration. Specifically, we are able to extract two critical exponents and universal information by observing the finite-size scaling of the GFS, which is particularly remarkable considering the constraints posed by the limited number of system sizes available within quasiperiodic systems.

The finite-size scaling theory of GFS is applicable in the vicinity of continuous quantum phase transitions, as the sensitivity is much more pronounced for the GFS around the QCP than that away from criticality [59]. When $V$ is adiabatically crossed over the QCPs $V_{c}$, the GFS of a finite system displays a peak, effectively indicating the positions of pseudocritical points $V_{m}$. With increasing system sizes $N$, the peaks of the GFS become more prominent. Additionally, it is anticipated that the maximal value of the GFS is expected to obey a power-law relation with the system size as [51]

where $\mu=2/\nu+2zr$ corresponds to the critical adiabatic dimension. Here, $\nu$ denotes the correlation-length exponent, and $z$ represents the dynamical critical exponent. According to single-parameter scaling theory, the correlation length $\xi$ becomes the only relevant length scale in the thermodynamic limit that diverges near the phase transition as [60],

and the single-particle spectral gap $\epsilon_{s}$ is predicted to close following the relation

For a relevant operator such as $\partial_{V}\hat{H}$ applied to sufficiently large one-dimensional systems, where $\nu\textless 2$, the pseudocritical points $V_{m}$ tend to approach the critical points $V_{c}$, thereby satisfying

Consequently, the GFS of a finite system with size $N$ in the vicinity of a QCP shall obey a universal scaling form  [61],

Here, $\phi_{r}(x)$ represents a unknown regular universal scaling function for the GFS of order $2r+2$. In this case, estimates for critical parameters can be deduced by plotting the scaled GFS as $[\chi_{2r+2}(V_{m})-\chi_{2r+2}(V)]/\chi_{2r+2}$ against $N^{1/\nu}(V-V_{m})$, with subtle adjustments made to the values of $V_{m}$, $\nu$, and $z$ until achieving a collapse of data.

Figure 4 shows the logarithm of the GFS, including $\ln\chi_{2}$ and $\ln\chi_{4}$ of the GAA model, as a function of $V$ for various system sizes $N$ with $\phi=\pi$. We observe that $\chi_{4}$ exhibits a much more pronounced divergence compared to $\chi_{2}$ at a critical point, yet there exist distinctions between the cases of $\rho=1/N$ and $\rho=1/2$. The logarithmically scaled maximum values of the GFS, $\chi_{\rm 2,max}$ and $\chi_{\rm 4,max}$, in proximity to the QCP are plotted against various system sizes $N$, which confirm the linear scaling relations by the numerical fittings. In-depth analysis of the insets in Figs. 4(a) and (b) reveals interesting behavior regarding the filling fraction $\rho=1/N$. Specifically, at this fraction, only a single subsequence is apparent. However, an interesting deviation occurs for even subsequences compared to odd ones when the phase $\phi$ is set to 0. To mitigate this odd-even effect and ensure consistent data analysis, we appropriately adjust the value of $\phi$ to $\pi$, thereby aligning the subsequences more uniformly. According to Eq. (20), the linear fitting for various system sizes $N$ gives rise to $\nu=0.998\pm 0.0015$ and $z=2.3705\pm 0.0025$, while the numerical fitting for even number of system sizes $N$ yields $\nu=0.998\pm 0.009$ and $z=2.371\pm 0.018$. On the contrary, when considering $\rho=1/2$, the numerical result reveals a division of the sequence into three distinct subsequences for numerous values of $\phi$. As shown in the inset of Fig. 4(c), the linear fittings of the log-log plot give rise to $\nu=1.004\pm 0.0515$ for $N=F_{3l+1}$, $\nu=0.9985\pm 0.0030$ for $N=F_{3l+2}$ and $\nu=1.0005\pm 0.0015$ for $N=F_{3l+3}$. Likewise, from the inset of Fig. 4(d), the dynamical exponent can be fitted as $z=1.8265\pm 0.0082$ for $N=F_{3l+1}$, $z=1.8305\pm 0.0138$ for $N=F_{3l+2}$ and $z=1.8355\pm 0.078$ for $N=F_{3l+3}$. By calculating their averages, we determine the correlation-length critical exponent of $\nu=1.001\pm 0.0186$ and the dynamical critical exponent of $z=1.83083\pm 0.03332$ for the filling fraction $\rho=1/2$.

Recent studies have highlighted an intimate relation between finite-size effects and Diophantine equation in the AA model [46]. The solutions of the Diophantine equation, obtained through continued fraction expansions, elucidate a mathematical framework that captures the essence of the system’s quasi-periodicity, which is fundamental to understanding its renormalization group dynamics [36]. The Diophantine equation establishes a link between the filling fraction and the disorder’s incommensurate frequency, refining the resonance criterion $Q\alpha-P=\rho$ of commensurate cases [62], where $Q$ and $P$ are integers that do not exceed the system size $N$. This association is based on the continued-fraction expansion of the frequency $\alpha$. When $Q$ and $P$ remain fixed irrespective of system size, and the resonance equation is solved precisely, the system is deemed commensurate. In this state, $Q$ determines the perturbation order where the resonance between the Fermi momentum and the potential frequency $\alpha$ occurs. Conversely, when $Q$ and $P$ scale with the system size, the system is recognized as incommensurate. Under these circumstances, the resonance criterion transcends into a Diophantine equation requiring a solution that accommodates the scaling variables, which then forms the mathematical foundation for understanding the resonance phenomena within the incommensurate framework. The incommensurate filling $\rho$, marked by the growth of $Q$ and $P$ with system size $N$, implies a dynamic interplay between the system’s physical dimensions and its spectral properties. This Diophantine equation not only serves as a tool for solving for these properties but also encapsulates the essence of the incommensurate state within the GAA model. Our investigation focuses on exploring this relationship further within the GAA model, particularly examining the sequence effects.

We note that $\alpha$ can be approximated by $\alpha_{k}=F_{k-1}/F_{k}$ and ${n}=\lfloor\rho F_{k}\rceil$. Consequently, the Diophantine equation is reformulated as

To seek integer solutions ($Q_{k}$, $P_{k}$) for given integer values of $F_{k}$ and ${n}$. For each value of $k$, an infinite spectrum of solutions ($Q_{k}$, $P_{k}$) hat satisfy the resonance condition in Eq. (21). The solutions are constrained to $|Q_{k}|\leq F_{k}/2$,and among these, the smallest $|Q_{k}|$ is selected. When solving the Diophantine equation for different values of $F_{k}$ and $n$, the ratio $|Q_{k}|/F_{k}$ for large $k$ displays a periodicity of $p$ distinct values, thus dividing the sequence of denominators into $p$ subsequences. As illustrated in Table 1, $p$ also corresponds to the number of scaling functions, such as $p=1$ for $\rho=1/N$ and $p=3$ for $\rho=1/2$ at half filling. Further details are available in Appendix A.

The universality class, central to the study of critical phenomena, posits that associated critical exponents and scaling functions are universal, subject only to symmetries and space dimensionality. We validate the obtained critical exponents using the data collapse technique, as demonstrated in Fig. 5. The scaling hypothesis suggests that appropriately scaled curves can converge into a single unified curve. To facilitate this data collapse, we employ the rescaled fidelity susceptibility as defined in Eq.(20):

where $\delta=V-V_{m}$ is the deviation from the pseudocritical point, and $\mathcal{F}_{q}$ ($q=1,2,\cdots p$) is scaling functions. Based on Eq. (22), we plot $\left[\chi_{2}(V_{m})-\chi_{2}(V)\right]/\chi_{2}(V)$ versus $N^{1/\nu}(V-V_{m})$ for different values of $N$ as depicted in Fig.5. For $\rho=1/N$ all curves of distinct chain sizes in the vicinity of $V_{m}$ collapse onto a single curve in Fig. 5 (a), while for $\rho=1/2$ the scaling of $\chi_{2}$ collapses all the curves onto three scaling functions. A very good quality data collapse confirms that the scaling hypothesis is in place and the estimated value of $V_{c}$, $\nu$, $z$ obtained through the prescriptions in Eq.(16) and (19) are self-consistent.

Similar strategies are applicable across various filling fractions. Our main findings in this work are the numerical results of the dynamical exponent $z$ for different filling $\rho$ and various values of $b$, as summarized in Tab. 1. We observe that the value of $z$ and the number of scaling functions $p$ are influenced by the filling fraction $\rho$, yet they remain independent of the deformation parameter $b$. As indicated in Tab.1, for a filling fraction of $\rho=1/N$, the dynamical critical exponent $z$ is approximately 2.734. At $\rho=1/3$, $z$ is around 2.00, and for $\rho=1/2$, $z$ is approximately 1.825. Intriguingly, at $\rho=1$, despite the localization transition not occurring at the same critical points as in the case of $\rho=1/N$, a symmetry persists, and the dynamical critical exponent is approximately 2.734 in both scenarios. This consistency of the dynamical critical exponents across similar filling states under varying mobility edges demonstrates that the Diophantine equation is applicable to the AA model, even when self-duality is broken.

In this study, we employ quantum fidelity susceptibility (QFS) to pinpoint the mobility edge in the Anderson localization transition. Recognized as a potent tool for probing quantum criticalities, the fidelity susceptibility has demonstrated its efficacy in one-dimensional disordered lattice models, as well as in the p-wave-paired Aubry-André-Harper model. Here our research thus focuses on the generalized Aubry-André (GAA) model, which is characterized by an exact mobility edge (ME) leading to a transition between the localized phase and the extended phase within its energy spectrum. A significant aspect of our approach involves utilizing generalized fidelity susceptibility (GFS) to deduce both the correlation-length critical exponent $\nu$ and the dynamical critical exponent $z$. A notable observation in our study is the manifestation of finite-size effects on the critical behavior of GFS. Furthermore, the number of scaling functions, denoted as $p$, conforms to the Diophantine equation, taking into account the varying filling fractions and system sizes. This led us to conduct extensive tests across different filling fractions and under diverse ME conditions, culminating in the derivation of the correlation-length and dynamical critical exponents, as shown in Table 1. The recent experimental verification of the existence of MEs [33] underscores the timeliness and relevance of our theoretical approach. By elucidating the localization transformations and critical exponents in the GAA model via QFS, our work lays a robust theoretical foundation for understanding the underlying mechanisms of ME phenomena. In conclusion, this research makes a substantial contribution to the field of quantum physics, particularly in understanding phenomena like the localization transition, mobility edge, and critical behavior in quantum systems. It introduces novel analytical methodologies and provides theoretical insights that align closely with recent experimental findings, thus enriching the current scientific discourse in this domain.