{'fraction_non_alphanumeric': 0.04541709320616214, 'fraction_numerical': 0.012059623979426231, 'mean_word_length': 4.692438885730529, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Finite-temperature properties of interacting bosons on a two-leg flux ladder

Quasi-one-dimensional lattice systems such as flux ladders with artificial gauge fields host rich quantum-phase diagrams that have attracted great interest. However, so far, most of the work on these systems has concentrated on zero-temperature phases while the corresponding finite-temperature regime remains largely unexplored. The question if and up to which temperature characteristic features of the zero-temperature phases persist is relevant in experimental realizations. We investigate a two-leg ladder lattice in a uniform magnetic field and concentrate our study on chiral edge currents and momentum-distribution functions, which are key observables in ultracold quantum-gas experiments. These quantities are computed for hard-core bosons as well as non-interacting bosons and spinless fermions at zero and finite temperatures. We employ a matrix-product-state based purification approach for the simulation of strongly interacting bosons at finite temperatures and analyze finite-size effects. Our main results concern the vortex-fluid-to-Meissner crossover of strongly interacting bosons. We demonstrate that signatures of the vortex-fluid phase can still be detected at elevated temperatures from characteristic finite-momentum maxima in the momentum-distribution functions, while the vortex-fluid phase leaves weaker fingerprints in the local rung currents and the chiral edge current. In order to determine the range of temperatures over which these signatures can be observed, we introduce a suitable measure for the contrast of these maxima. The results are condensed into a finite-temperature crossover diagram for hard-core bosons.

In this context, quasi-one-dimensional ladder-like systems (see Fig. 1 for a sketch) are particularly timely. The corresponding geometries can be isolated in superlattice potentials [19] and are naturally realized in setups employing a synthetic dimension [20][21][22][23][24][25]. Moreover, ladder-like lattices host rich physics and allow for a direct comparison between theoretical predictions for interacting systems and experimental results. While originally advocated as a suitable system to measure edge states in experiments [19,26], amazingly, it has also been suggested that Chern numbers can be mapped out experimentally in multi-leg flux ladders with as few as up to five legs and open boundary conditions [27], which was studied in a recent experiment [28].
Analytical as well as numerical studies, with a large fraction of the latter work based on the density-matrix renormalization group method [33,34], have provided extensive theoretical results regarding the ground-state Red arrows indicate the direction and by their length, the strength of the local currents. The size of the dots and the background shading represent the local particle density. The Hamiltonian parameters φ, J ⊥ , J and U are introduced in Eq. (1). A comprehensive analysis of the Meissner and vortex-fluid phases can be found in, e.g., Refs. [26,[29][30][31][32].
properties of interacting quantum gases on such ladderlike systems [26,[29][30][31][32]. The presence of gauge fields clearly enriches the corresponding phase diagrams, which by now have been mapped out in large part. For instance, these phase diagrams host superfluid and Mottinsulating vortex-liquid and Meissner phases [29,30,32], vortex lattices [29,44,46], as well as charge-densitywave [46,48] and biased-ladder phases [40]. Typical particle-current patterns found in the Meissner and vortex-fluid phases are shown in Figs. 1(a) and (b), respectively. Furthermore, the possible existence of frac-tional quantum Hall-like states has attracted great interest [39,41,43,47]. Recently, results regarding the steady state of a two-leg ladder system subjected to a dissipative environment have been put forward in Ref. [57], and the ground-state phase diagram of a related but more complicated two-leg Haldane model has been mapped out in Ref. [58]. We note that currents in mesoscopic ringshaped flux ladder systems, which have been experimentally realized in Ref. [59], have been studied in Ref. [60]. However, most of these studies concentrated on ground states, while experiments are naturally at nonzero energy densities.
In this paper, we investigate a two-leg ladder in a uniform Abelian gauge field at finite temperatures T . Especially, we concentrate on the most prominent transition, the vortex-fluid-to-Meissner crossover of hard-core bosons as well as of non-interacting bosons and a related transition for spinless fermions at T > 0 [26,29,30,32,35,36,42]. We compute chiral edge currents and momentum-distribution functions at zero and finite temperatures and determine the range of temperatures in which clear signatures of the vortex-fluid phase can be detected. For this purpose, we employ a matrix-productstate based purification approach [34,61], which we implement for both, canonical setups with fixed numbers of particles as well as grand-canonical setups with a fixed average number of particles. We note that aspects of finite-temperature physics in ladder systems have previously been investigated in Refs. [44,53,56].
We first revisit the case of non-interacting bosons and fermions [26,31,35] and discuss the corresponding vortex(-fluid)-to-Meissner transition. Moreover, in these limits, we can compute results at finite temperatures for very large systems and therefore, we are able to obtain estimates for the typical finite-size effects. Several different quasi-momentum distribution functions are introduced and we argue that these are better suited to detect the vortex fluid at T > 0 than the chiral edge current or local rung currents. We introduce a quantitative measure that detects the presence of finite-momentum peaks in a suitably chosen momentum-distribution function. Moreover, we compare chiral edge currents observed in canonical and grand-canonical setups. Our main results are for hard-core bosons, for which we present a finite-temperature crossover state diagram for model parameters realized in [19]. Finally, we point out other interesting questions for future finite-temperature studies concerning, e.g., the melting of vortex lattices or of the biased-ladder phase at T > 0 [44,46] in the strongly interacting, low-density regime.
In the following, we first define the Hamiltonian of the ladder system and introduce observables of interest in Sec. II. There, we also recap key aspects of the non-interacting model. The purification approach used for the calculation of thermal states of strongly interacting bosons is introduced in Sec. III. Our main results concerning the vortex-fluid-to-Meissner crossover at finite temperatures are presented in Sec. IV. There, we discuss chiral edge currents and momentum-distribution functions of hard-core bosons, non-interacting spinless fermions as well as non-interacting bosons at zero and finite temperatures. Finally, we summarize our work in Sec. V.

II. MODEL AND OBSERVABLES OF INTEREST
A. Two-leg ladder Hamiltonian This paper mainly focuses on a Bose-Hubbard Hamiltonian defined on a two-leg ladder geometry. Using annihilation (and creation) operators b ( †) l,r referring to bosons or spinless fermions located on site (l, r), the Hamiltonian reads, Here, the indices l = 0, 1 and r = 0, . . . , L − 1 label the legs and rungs, respectively, of a two-leg ladder with L rungs. The first term in Eq. (1) represents nearestneighbor particle hopping along the legs of the ladder. We use the corresponding hopping amplitude J as our unit of energy (we seth = 1 and k B = 1). The second term accounts for inter-leg hopping along the rungs of the ladder. It is designed in such a way that the state of a particle encircling a single plaquette of the ladder gains a phase factor e ±iφ , with the sign depending on the direction of the circulation. The last term in Eq. (1) accounts for repulsive interactions between two or more bosons located at the same site (for spinless fermions, this term is absent). In this paper, we concentrate on hard-core bosons with a maximum occupation of one particle per site (U/J → ∞) and non-interacting particles (U = 0). Moreover, in the examples presented throughout this paper, we consider a value of φ = π/2, which corresponds to the experimental setup of Ref. [19]. Thereby, the considered Hamiltonian exhibits a translational symmetry with a unit cell containing four plaquettes of the ladder in the (rung-)gauge chosen in Eq. (1).

B. Non-interacting model
Next, we recap some of the key aspects of the noninteracting (U = 0) flux-ladder model with periodic boundary conditions. The corresponding Hamiltonian is denoted by H F . A comprehensive analysis of this model can be found in, e.g., Refs. [26,31].
To diagonalize H F , we first rewrite it in terms of leggauge operators denoted by an overhead tilde and defined byb 0,r = e ir φ 2 b 0,r andb 1,r = e −ir φ 2 b 1,r , which moves the complex hopping amplitudes from the rungs to the legs of the system. Second, leg-gauge momentum operators b ( †) l,km are introduced by means of Fourier transforming along the legs of the ladder, with m = 0, . . . , L − 1 and momenta k m = 2πm/L. In terms of leg-gauge momentum operators, the noninteracting Hamiltonian with periodic boundary conditions takes the form with Finally, the operators d with ξ = J ⊥ 2J , B k = sin(k) sin(φ/2) and A k = ξ 2 + B 2 k . Thus, H F takes the form where the eigenvalues corresponding to the lower and upper band, + (k m ) and − (k m ), are given by ± (k m ) = −2J (cos(k) cos(φ/2) ± A km ). For a non-vanishing interleg coupling J ⊥ = 0, the two bands are separated by a finite gap. Depending on the parameters φ and J ⊥ , the lower band has either a single global minimum at momentum k = 0 or two degenerate minima at k = ±k with, For k = 0, the system is in the Meissner phase, while for k = ±k, it is in the vortex-fluid phase [26,31]. The corresponding degenerate single-particle ground states are denoted by |±k . For a given value of φ, the phases are separated by J c ⊥ = 2J sin (φ/2) tan (φ/2). The two-band structure of the non-interacting ladder model is shown for different values of J ⊥ and φ = π/2 in Fig. 2(a).  In anticipation of the discussion of non-interacting fermions, Fig. 2(c) illustrates the range [ ± (k)] = [min k ± (k), max k ± (k)] of the lower and upper band for φ = π/2. From this, a transition-point (for φ = π/2) for the fermionic system can be determined from, We use the same terminology -Meissner and vortex phase -as for bosons since the local current patterns show the same phenomenology: a chiral edge current with no rung currents in the Meissner phase and finite rung currents in the vortex phase for open boundary conditions (data not shown, see the discussion of the chiral edge current in Sec. IV B 1). We will use these equations to elucidate the vortex-to-Meissner transition of noninteracting fermions at filling one-half. The fermionic transition is further related to a metal (="vortex phase") to band-insulator (="Meissner") transition [35].
C. Observables of interest

Particle currents
In the context of the ladder model, local particle currents are of special importance. They have been suc-cessfully employed as a key measure for distinguishing between different zero-temperature phases (see [46] for an overview). Moreover, they are experimentally accessible as demonstrated in various cold-atom experiments, see, for instance, Refs. [19][20][21].
The continuity equation for the occupation of lattice sites implies that operators representing local particle currents read Here, j l,r represents the particle current from site (l, r) to site (l, r + 1) and j ⊥ r represents the particle current from site (0, r) to site (1, r). Furthermore, the chiral-edgecurrent operator accounting for global particle transport along the legs of the system is given by 2. Integrated leg-gauge momentum-distribution function (Quasi-)Momentum-distribution functions exhibit characteristics specific to certain zero-temperature phases of the ladder Hamiltonian, Eq. (1) (see the discussion in [26,46]). Moreover, they can in principle be obtained from experimental time-of-flight measurements [62].
The specific momentum-distribution function (MDF) n(k) that we primarily consider in this paper is defined by using the leg-gauge operators of Eq. (2). We will refer to n(k) as integrated leg-gauge MDF. The expectation value • can either refer to the ground state or to a finitetemperature ensemble.

Leg-resolved MDFs
In the context of the ladder model, Eq. (1), other distribution functions that are related to n(k) have been considered [26,53,54]: (i) Leg-resolved leg-gauge MDFs n l (k) are given in terms of leg-gauge momentum operators, Summing over both legs, they add up to the integrated leg-gauge MDF, n(k) = n 0 (k) + n 1 (k).
(ii) Leg-resolved rung-gauge MDFs n l (k) are defined by means of rung-gauge momentum operators b We refer to these operators b l,r , which are also used in Eq. (1), as rung-gauge operators, indicating that the complex hopping amplitudes of the Hamiltonian are located along the rungs of the ladder. The integrated runggauge MDF is given by n(k) = n 1 (k) + n 2 (k).
Leg-resolved leg-gauge and leg-resolved rung-gauge MDFs are related in the following sense, We compare all four types of MDFs (integrated versus leg-resolved, leg-gauge versus rung-gauge) in Figs. 3(a)-(d). Generally, the number of possible maxima is doubled by going from leg to rung gauge. In the experiment in Ref. [19] the rung gauge was realized which leads to the MDF n(k). For convenience, we will present results for the leg gauge, i.e., n(k), yet the former can be obtained from the latter. The analysis of signatures of the vortexfluid phase proposed in Sec. IV A 3 can straightforwardly be adapted to n(k) by simply restricting k to one half of the Brillouin zone.

III. (GRAND-)CANONICAL PURIFICATION
The core concept of purification is that a mixed state in a physical space can be represented as the partial trace of a pure state in an artificially extended space [34]. Thus, by incorporating an auxiliary counterpart for each physical site, thermal states can be represented as matrixproduct states and propagated in imaginary time using well-developed techniques [61,[63][64][65][66][67][68][69][70][71][72][73][74][75][76]. In this context, suitable initial infinite-temperature states need to be chosen. For grand-canonical calculations these states are straightforward to generate [34]: They exhibit maximum entanglement between each physical site and its auxiliary counterpart but are otherwise of product form. On the other hand, the required initial states for canonical simulations exhibit non-trivial long-range correlations. They have been successfully obtained as the ground states of specifically chosen Hamiltonians using the density-matrix renormalization group method [77][78][79].
Here, we generate matrix-product-state representations of canonical infinite-temperature states by employing a bookkeeping mechanism that explicitly accounts for the occupation of physical sites within the matrixproduct state. This approach renders the generation of canonical infinite-temperature states straightforward. We emphasize that the methodology employed here is very similar to the one introduced in Ref. [66]. Comparison of momentum-distribution functions in the single-particle case. The figure shows the position of the maxima in various momentum-distribution functions that are (a),(b) leg-resolved (l = 0, 1) or (c),(d) integrated over l. The values of momentum-distribution functions for each of the single-particle ground states |±k are indicated by the colorbars. Panels (a) and (c) depict the behavior in the runggauge, while (b) and (d) are for the leg-gauge. The integrated leg-gauge momentum-distribution function n(k) of panel (d) will be discussed in the remainder of the paper. See the text in Secs. II C 2 and II C 3.
Our purification approach can be subdivided into three steps: (i) First, a matrix-product-state representation of the (grand-)canonical infinite-temperature state is constructed. (ii) Second, a suitable matrix-productoperator representation of the propagator exp (−βH) is constructed, where β = 1/T is the inverse temperature. (iii) A third step involves the imaginary-time propagation including the control of errors, to obtain a thermal state corresponding to a finite temperature.
In the remainder of this section, we will describe all steps in greater detail.
A. Matrix-product-state structure The matrix-product states employed during (grand-)canonical purification simulations are comprised of tensors representing physical sites of the two-leg ladder system and their auxiliary counterparts, see Fig. 4(a). As shown in Fig. 4(b), these tensors are aligned side-by-side following a zigzag pattern through the ladder and forming a chain of length 4L.
Since canonical calculations are at a fixed number of particles, which is here conserved by the ladder Hamiltonian, we need the matrix-product state to account for the corresponding symmetry. In practice, this means that on the level of each tensor, a bookkeeping mechanism needs Within the purification approach, each physical site is accompanied by an auxiliary counterpart represented by an empty circle. (b) The matrix-product state is comprised of tensors representing physical sites (filled squares) and of tensors representing auxiliary sites (empty squares). These tensors are aligned side-by-side following a zigzag pattern through the ladder as indicated by the blue arrows in (a) and (b). (c) As discussed in the main text, the canonical purification approach requires the bookkeeping of the total occupation of all physical sites. Therefore, the incoming left leg i of each tensor is split into blocks {si} referring to different occupations of the physical sites to the left of the specific tensor. Similarly, the incoming leg p entering each tensor from above is split into blocks {sp} referring to the occupation of the specific site represented by the tensor, and the outgoing leg o to the right of each tensor is split into blocks {so} referring to the occupation of all physical sites to the left of the specific tensor including the tensor itself.
to be implemented. This mechanism keeps track of the total occupation number of the physical sites when going through the matrix-product state from left to right. At this point, we recap that the conventional building blocks of a matrix-product state are three-legged tensors, here denoted by T p,i o . The incoming leg p, which is usually referred to as physical leg, corresponds to local basis states of the lattice site represented by the tensor. The incoming leg i and the outgoing leg o enable a connection between the tensor and its neighbors within the matrix-product state. Against this background, the employed bookkeeping mechanism can be understood as a subdivision of the vector spaces associated with each of the conventional legs into sectors corresponding to different occupation numbers. More precisely, for each tensor, the sectors of its incoming leg i correspond to occupation numbers of the physical sites that are represented by the tensors to its left. Analogously, for each tensor, the sectors of its outgoing leg o correspond to occupation numbers of the physical sites that are represented by the tensor itself and by the tensors to its left. Moreover, the sectors subdividing an incoming leg p refer to the local occupation numbers of the basis states. Hence, the outgoing leg o of the leftmost tensor in the matrix-product state contains information about the total occupation of all physical sites. This bookkeeping mechanism is illustrated in Fig. 4(c).

B. Canonical infinite-temperature state
Matrix-product operators are required for the preparation of the (grand-)canonical infinite-temperature state as well as for the subsequent imaginary-time propagation. Here, they are constructed from single-site operators following the generic approach presented in Ref. [80]. Considering the ladder model Eq. (1) in the hard-core boson limit, a possible grand-canonical infinite-temperature state |β = 0 g.c. can be readily constructed using single- l,r phys. and b ( †) l,r aux. which act on physical sites and their auxiliary counterparts. This state is spread over all possible sectors corresponding to total occupation numbers of the physical sites ranging from zero to 2L. It is given by with |vac being the vacuum state. We emphasize that a suitable canonical infinite-temperature state with a fixed total occupation number can be obtained from |β = 0 g.c. by means of a projection onto the corresponding particle number sector and subsequent normalization. Moreover, the matrix-product-state structure introduced above renders this projection straightforward: All but the symmetry sector of interest are set to zero at the very right end of the matrix-product state. Finally, the canonical infinite-temperature state with a fixed total occupation number of N particles is given by where P N denotes a projection onto the corresponding subspace and subsequent normalization of the matrixproduct state.

C. Imaginary-time propagation
The imaginary-time evolution requires a matrixproduct-operator representation of the propagator e −βH = e −τ H β/τ . For this, we subdivide the ladder Hamiltonian given in Eq. (1) into two parts, H = H 0 + H 1 , and approximate the time-evolution operator by means of a second-order Trotter decomposition, One possible choice for H 0 and H 1 is given by (j = 1, 2) with I n m = (m+n) mod 2. Note that for hard-core bosons (U/J → ∞), the interaction term in Eq. (1) is implemented by means of the cutoff of the lattice site-local Hilbert space at a maximum of one boson. Thus, it has been ignored in Eq. (20) for simplicity. The Hamiltonian parts H 0 and H 1 are composed of commuting contributions which can be associated with every second plaquette of the ladder, meaning that they can be diagonalized numerically. Hence, the matrix-product operators representing the exponentials in Eq. (19) can be constructed automatically from single-site operators using common matrix-product-operator arithmetic [80]. In terms of the employed matrix-product-state structure, the Hamiltonian H exhibits at most fourth-nearest-neighbor interactions and the propagators e −τ H0 and e −τ H1 exhibit at most sixth-nearest-neighbor interactions.
The imaginary-time propagation itself consists of the sequential application of matrix-product operators representing the propagators on the right-hand-side of Eq. (19). Choosing a suitable step width τ , they are applied repeatedly to the initial (grand-)canonical infinitetemperature state |β = 0 c (|β = 0 g.c. ), which decreases the temperature of the so-evolved thermal state. Between these applications, the bond dimensions of the evolved matrix-product state need to be truncated and, due to the fact that the imaginary-time propagators are not unitary, the matrix-product state also needs to be repeatedly normalized. During the imaginary-time propagation, errors arise due to the Trotter decomposition of the propagator and due to the repeated truncation of the evolved state. Regarding the numerical results presented in this paper, these errors have been controlled independently by comparing results obtained for different truncation thresholds at fixed Trotter-step widths and vice versa. We note that this error control can in principle be based on the monitoring of state overlaps. However, here it is sufficient to focus on the values obtained for the actual observables of interest, which are chiral edge currents and momentum-distribution functions.

IV. VORTEX-FLUID-TO-MEISSNER CROSSOVER AT FINITE TEMPERATURES
In Sec. IV A, we address the grand-canonical statistics of the non-interacting bosonic ladder model. Focusing on momentum-distribution functions, the presence of clearly detectable characteristics of the vortexfluid phase at finite temperatures is examined. In Sec. IV B, we concentrate on the grand-canonical statistics of non-interacting spinless fermions and elucidate the

A. Non-interacting bosons
We recap that the non-interacting bosonic ladder model with periodic boundary conditions exhibits two degenerate ground states in the vortex-fluid phase, corresponding to finite characteristic momenta ±k (see Sec. II B). At the vortex-fluid-to-Meissner transition, this degeneracy is lifted and the system exhibits a unique ground state corresponding to a vanishing characteristic momentum, k = 0, in the Meissner phase.

Integrated leg-gauge momentum-distribution functions
Results for the integrated leg-gauge MDF of noninteracting bosons are shown in Fig. 6. In this context, the ground-state vortex-fluid-to-Meissner transition manifests itself in the transition from two distinct maxima at finite momenta to a single maximum at zero momentum, as shown in Fig. 6(a). The integrated leg-gauge MDFs calculated in grand-canonical thermal states at different temperatures J/4, J/2 and J, which are shown in Figs. 6(b), (c) and (d), respectively, reveal that the typical peaks persist at finite temperatures. We note that the results for ground states in Fig. 6(a) and for the grand-canonical thermal states in Figs. 6(b), (c) and (d) are for a fixed filling f = 1/2 and an average filling of one-half, respectively.

Signatures of the vortex-fluid phase at finite temperatures
Next, we elucidate the finite-temperature vortex-fluidto-Meissner crossover. Considering the integrated leggauge MDF n(k) of thermal states, we specify at which temperatures clear characteristics of the vortex-fluid phase can be detected and where the crossover can be observed. For this purpose, we introduce a measure for the contrast of finite-momentum maxima in n(k) via this definition Thus, on the one hand, values of δ larger than zero indicate that two characteristic peaks at finite momenta can be resolved. Hence, δ > 0 is indicative of the vortexfluid phase. On the other hand, values of δ equal to zero mean that n(k) exhibits a single maximum at zero momentum, and, due to the blurring effect of finite temperatures, it cannot be decided whether the underlying values of J, J ⊥ and φ correspond to the vortex-fluid or Meissner phase. The definition of the contrast, Eq. (21), is also illustrated in Fig. 7(a) using the abbreviations M = max k [n(k)]−n(k = 0) and N = max k [n(k)]. There, it is exemplified that the characteristic peaks indicating the vortex-fluid phase are more distinctive at lower temperatures. Considering grand-canonical thermal states at average filling one-half, a contour plot of δ as a function of inter-leg coupling J ⊥ and inverse temperature β is shown in Fig. 7(b). It becomes apparent that the integrated leg-gauge MDFs show clear signatures of the vortex-fluid phase at sufficiently small values of J ⊥ < ∼ J and up to quite high temperatures T ∼ J. Note that the zero-temperature transition occurs at J ⊥ = √ 2J, as indicated in Fig. 7(b) by the dotted line. The role of finite-size effects for the contrast δ is studied in Fig. 7(c) by considering a cut through the graph in Fig. 7(b) at J ⊥ = J/2. We note that finite-size effects of δ are overall small and for a ladder of length L = 128, the differences to larger L (not shown in the figure) are negligible on the scale of the plot. Furthermore, we only observe quantitative differences as L varies, yet the overall dependence of δ on β is the same for all L.
Finally, we comment on another measure for the vortex phases that is often studied in the literature (see, e.g., [29,32,46]), namely the vortex density l v . In Ref. [32], l v was extracted from a Fourier transformation of the local rung-current patterns. At finite temperatures, however, the local rung currents are quickly washed out such that this approach cannot be used here. Rather, one would need to study the decay of rung-current equal-time autocorrelations. These are not easily accessible in exper- iments and we therefore do not further study measures for the vortex density here.

B. Non-interacting spinless fermions
In this section, we discuss chiral edge currents and momentum-distribution functions of non-interacting spinless fermions.

Chiral edge currents
Chiral edge currents corresponding to the ground states of non-interacting spinless fermions behave fundamentally different from those corresponding to the ground states of non-interacting bosons. For bosons, the behavior close to the band minima is important, while for fermions, the behavior close to the Fermi energy is relevant, or for non-universal objects such as the chiral edge current, the occupation of all quasimomenta [35]. One notable consequence is that already for U = 0, in the Meissner phase, non-interacting fermions feature a reversal of the direction of the chiral edge current with increasing filling f [35,36]. Fermionic ground-state chiral edge currents are shown in Fig. 8(a) as a function of inter-leg coupling J ⊥ for different fillings f ranging from 1/(2L) to 1/2. For each filling, a kink of the chiral edge current clearly reveals the critical value of J ⊥ at the vortex-to-Meissner transition. Figure 8(b) concentrates on fermionic chiral edge currents at finite temperatures and average filling one-half. For temperatures below approximately J/10 and values of J ⊥ < ∼ J or J > ∼ 2J, i.e., deep within the Meissner and vortex phases, the chiral edge currents computed in grand-canonical thermal states coincide with the ones obtained in the fermionic ground state. However, even for the smallest temperature considered in Fig. 8(b), which is J/20, the chiral edge current does not exhibit a clear kink at the critical value, J ⊥ = √ 2J, of the vortex-to-Meissner transi-tion. We note that the minima of the chiral edge current systematically overestimate the critical value of J ⊥ with increasing temperature.

Integrated leg-gauge momentum-distribution functions
Within the Meissner phase, the integrated leg-gauge MDFs corresponding to the ground state of noninteracting spinless fermions at filling f = 1/2 take a constant value equal to one, see Fig. 9(a) [note that white corresponds to an occupation of one in the color-coding of that figure]. This can be readily understood from the overlap of the bands referred to in Eq. (9) and shown in Fig. 2(c). In the Meissner phase, the lower band is energetically well separated from the upper band. Hence, at filling one-half, it is fully occupied which leads to a constant MDF. Note that the bands shown in Fig. 2 and referred to in Eq. (9) correspond to the ladder model with periodic boundary conditions while we show data for open boundary conditions in the figure. However, for systems with a large number of rungs L, boundary effects play a minor role and integrated leg-gauge MDFs for open and periodic boundary conditions coincide. On the other hand, in the vortex phase, the groundstate integrated leg-gauge MDF is piecewise constant and takes on the values zero, one and two, corresponding to non-occupied, singly-occupied and doubly-occupied momenta, see Fig. 9(a). With the onset of finite temperatures, the sharp features of the fermionic ground-state integrated leg-gauge MDF blur out and their discreteness is lost. Figures 9(b), (c) and (d) show integrated leg-gauge MDFs computed in grand-canonical thermal states at average filling one-half and at different temperatures J/4, J/2 and J, respectively. In the case of fermions, finite temperature causes the peak structure characteristic for the vortex phase to extend into the Meissner region, opposite from the behavior of bosons.

C. Hard-core bosons
In this section, we give an account of our main results addressing the properties of hard-core bosons on the ladder system introduced in Eq. (1) at finite temperatures. We remind the reader that a value of φ = π/2 is considered.

Chiral edge currents
Let us start with the presentation of our results concerning the chiral edge current j c . In the thermodynamic limit, a kink at the maximum value of j c indicates the critical value of the inter-leg coupling strength J ⊥ of the ground-state vortex-fluid-to-Meissner transition. Furthermore, in the Meissner phase, j c is expected to decay to zero as j c ∝ 1/J ⊥ [32].   Fig. 10(a) show chiral edge currents corresponding to canonical thermal states of different temperatures ranging from J/4 to 5J/2 at filling one-half; note that inverse temperatures β are indicated by the colorbar. For values of J ⊥ approximately greater than 1.5J, i.e., deep within the Meissner phase, and for temperatures below around J/3, the chiral edge currents corresponding to the canonical thermal states overlap with the ones corresponding to the ground state. Further, for sufficiently small temperatures, the finite-temperature currents show clear maxima. However, the positions of these maxima tend to systematically overestimate the critical value of J ⊥ with increasing temperatures. state only). At zero temperature and in the Meissner phase, the chiral edge currents for the L = 64 system are slightly more pronounced than those for the L = 32 system. However, as can be seen in Fig. 10(b), finite-size effects of the chiral edge currents play a minor quantitative role, in particular at T > 0, where there are small differences between the L = 32 and L = 16 data for the selected values of T .
In Fig. 11, we show that for small ladder systems, chiral edge currents corresponding to canonical setups with fixed filling f = 1/2 quantitatively differ from those corresponding to grand-canonical setups with average filling one-half. On the other hand, for large systems, canonical and grand-canonical setups are expected to feature the same characteristics and the respective chiral edge currents should coincide. It can be seen in Fig. 11(a)  lowing, we show only results obtained in the canonical ensemble.

Integrated leg-gauge momentum-distribution functions
Next, we consider integrated leg-gauge MDFs n(k). In this quantity, the ground-state vortex-fluid-to-Meissner transition manifests itself in the transition from two distinct maxima at finite momenta to a single maximum at zero momentum as shown in Fig. 12(a). There, we show results for the integrated leg-gauge MDF computed in the ground state of a hard-core boson ladder with L = 32 at f = 1/2. Figure 12(b) demonstrates that the characteristic ground-state peaks persist in canonical systems with a finite temperature J/4 and that the two peaks in the vortex-fluid phase can be clearly resolved. However, as expected, these sharp ground-state features blur out with increasing temperatures, see Fig. 12(c) and Fig. 12(d).

Signatures of the vortex-fluid phase at finite temperatures
Here, we concentrate on the measure of contrast introduced in Eq. (21) in order to quantify at which temperatures clear signatures of the vortex-fluid phase can be detected. Considering a canonical half-filled hard-core boson ladder with L = 32 rungs, Fig. 13(a) shows a contour plot of δ as a function of inter-leg coupling strength J ⊥ and inverse temperature β. one finds δ > 1/10. This means that two different peaks, and thus a clear signature of the underlying vortex-fluid phase, can be detected. From Fig. 13(a), it also becomes apparent that in the vortex-fluid regime, the contrast δ increases with decreasing inter-leg coupling strength J ⊥ and decreasing temperatures β −1 .
Moreover, the results for δ shown in Fig. 13(a) do not suffer from significant finite-size effects as can be seen in Fig. 13(b). For the vertical cut at J ⊥ = J/2 indicated in Fig. 13(a) [dashed line], Fig. 13(b) shows the contrast δ as a function of inverse temperature β for ladders of length L = 8, L = 16 and L = 32. Note how sharply the curves for δ drop to zero as signatures of the two maxima in n(k) are fully blurred out. Thus the derivative of δ with respect to β would serve as a good indicator for the vortex-fluid-to-Meissner phase crossover. A comparison of the L = 32 data from Fig. 13(b) and Fig. 7(c) reveals that for non-interacting systems, a finite contrast persists at smaller temperatures.

V. SUMMARY
In this paper, we investigated the properties of strongly-interacting bosons as well as non-interacting bosons and spinless fermions which are confined on a two-leg ladder lattice subjected to a uniform magnetic field at zero and finite temperatures. A particular focus was on the vortex-fluid-to-Meissner quantum-phase transition and the corresponding crossover observed at finite temperatures.
The chiral edge currents and momentum-distribution functions were the key observables analyzed in this paper. In order to introduce the physics, the principal features and the analytic solution of a non-interacting ladder model with periodic boundary conditions were recapitulated [26,31]. Our main results concern the vortex-fluidto-Meissner crossover of hard-core bosons. They were obtained by means of a matrix-product-state based purification approach following [66]. This method is applicable to both, canonical setups with a fixed number of particles as well as grand-canonical setups. Typical chiral currents and momentum-distribution functions computed at zero and finite temperatures were discussed in detail.
We presented a comparison between results for chiral edge currents obtained in the canonical ensemble compared to those found in grand-canonical setups. As a result, even for the relatively small particle numbers typical for quantum gases in optical lattices, (the technically simpler) grand-canonical simulations give quantitatively reliable results for the quantities of interest here and at sufficiently high temperature.
Moreover, we showed that clearly detectable signatures of the underlying bosonic vortex-fluid phase persist in suitably chosen momentum-distribution functions at finite temperatures. For this, a measure of contrast was introduced in Eq. (21) that is sensitive to the presence of finite-momentum peaks in the integrated leg-gauge momentum-distribution function.
Considering the notorious difficulty of cooling a quantum gas to low energy densities [81], we expect our results to provide relevant guidance for future experiments which are naturally at nonzero temperature but allow control over the particle numbers. Our results for the contrast can be compared to the resolution of a given experiment and will thus provide an upper limit for the temperature at which the vortex-fluid phase is visible in the strongly interacting, low-density regime. The methods and results presented here will also be relevant in characterizing the final state in quantum-quench dynamics and state-preparation protocols for the phases on flux ladders [82][83][84].
Our work paves the way for future investigations of finite-temperature properties of flux-ladder systems. Open questions concern the finite-temperature properties of the bosonic ladder model in the presence of finite interaction strengths 0 < U/J < ∞. In this regime, the model supports a variety of zero-temperature phases, including vortex-lattice and biased-ladder phases [40,44,46]. It has not yet been investigated if the characteristic signatures of these phases, such as typical current patterns and leg-population imbalances, persist at finite temperatures in the low-density regime. The same applies to predictions for fingerprints of candidates for fractional quantum Hall states in these ladder systems [39,41,43,47].
We thank M. Aidelsburger, C. Schweizer and B. Wang for useful discussions. We are indebted to C. Schweizer for his comments on a previous version. This work was supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under project number 277974659 via Research Unit FOR 2414 and under Germany's Excellence Strategy -EXC-2111 -390814868. We acknowledge funding through the ExQM graduate school.