Creutz–Hubbard Ladder Insights

• Creutz–Hubbard Ladder is an interacting two-leg lattice model that combines synthetic magnetic flux with local Hubbard interactions to generate flat bands and symmetry-protected edge modes.
• The model exhibits distinctive phenomena including Aharonov–Bohm caging, correlated topological transitions, and pair superfluidity, linking microscopic interactions to macroscopic quantum phases.
• Experimental platforms in ultracold atoms, photonics, and circuit QED enable practical realizations, making the ladder a versatile framework for quantum simulation and advanced research.

The Creutz–Hubbard ladder is an interacting two-leg lattice model obtained by augmenting the Creutz ladder with local Hubbard interactions. In the literature it appears in bosonic, spinless-fermion, spinful-fermion, and photonic forms, and it is used to study the interplay between flat-band dynamics, topology and interactions in a quasi-1D setting with synthetic magnetic flux, compact localized states, and symmetry-protected edge modes (Jünemann et al., 2016, Zurita et al., 2019, Pelegrí et al., 2024).

1. Canonical formulations and parameterizations

A common formulation considers a ladder of length $L$, with sites $\ket{j,\alpha}$, $j=1,\dots,L$, $\alpha=A,B$, and spinless particles which hop along legs, rungs, and diagonals in the presence of a synthetic magnetic flux $\phi$ per plaquette. The Hamiltonian is written as

$H \;=\; H_{0} \;+\; H_{\rm int},$

with

$$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$

For spinless fermions one uses

$H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$

whereas for bosons one uses

$$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$

Here $J$ is the diagonal hopping, $\ket{j,\alpha}$0 the vertical hopping, and $\ket{j,\alpha}$1 is the Peierls phase due to a flux $\ket{j,\alpha}$2 per square plaquette (Zurita et al., 2019).

A frequently studied spinless-fermion variant is the imbalanced Creutz–Hubbard Hamiltonian,

$\ket{j,\alpha}$3

with $\ket{j,\alpha}$4, $\ket{j,\alpha}$5. In this notation $\ket{j,\alpha}$6 sets the $\ket{j,\alpha}$7-flux kinetic scale, $\ket{j,\alpha}$8 is a leg imbalance, and $\ket{j,\alpha}$9 is a rung repulsion (Jünemann et al., 2016).

Bosonic formulations also occur naturally. In the circuit-QED construction of the bosonic Creutz–Hubbard ladder, one introduces bosonic annihilation operators $j=1,\dots,L$0 on the upper and lower legs and writes a Bose–Hubbard Hamiltonian with vertical, diagonal, and leg hoppings $j=1,\dots,L$1, a Peierls phase $j=1,\dots,L$2, and on-site Kerr interaction $j=1,\dots,L$3 inherited from transmon nonlinearity (Alaeian et al., 2018).

These parameterizations are notationally distinct. A plausible implication is that the literature emphasizes different gauges, experimental constraints, and particle statistics rather than different underlying geometries.

2. Flat bands, Aharonov–Bohm caging, and topological structure

In the noninteracting sector the ladder is a paradigmatic flat-band topological system. For the $j=1,\dots,L$4-$j=1,\dots,L$5-$j=1,\dots,L$6 parameterization, the Bloch Hamiltonian can be written as

$j=1,\dots,L$7

with

$j=1,\dots,L$8

and band energies

$j=1,\dots,L$9

At the special point

$\alpha=A,B$0

one obtains

$\alpha=A,B$1

namely two perfectly flat bands (Zurita et al., 2019).

In the $\alpha=A,B$2-parameterization used for the bosonic Creutz ladder,

$\alpha=A,B$3

with

$\alpha=A,B$4

The flat-band condition is

$\alpha=A,B$5

which gives $\alpha=A,B$6 and a lower-band Zak phase

$\alpha=A,B$7

signalling a topologically non-trivial phase with protected edge modes (Pelegrí et al., 2024).

The imbalanced synthetic Creutz–Hubbard model displays the same flat-band limit at $\alpha=A,B$8, where the two bands are $\alpha=A,B$9, the chiral symmetry satisfies $\phi$0, and the winding number is $\phi$1, or equivalently the Zak phase is $\phi$2. For open boundary conditions there are localized Aharonov–Bohm cage bulk orbitals and two exact zero-energy edge modes,

$\phi$3

Turning on $\phi$4 warps the bands, and the Zak phase remains $\phi$5 for $\phi$6 and vanishes for $\phi$7, marking a non-interacting topological phase transition at $\phi$8 (Jünemann et al., 2016).

Aharonov–Bohm caging is one of the defining dynamical signatures. At $\phi$9, destructive interference cages any initially localized state into a small plaquette. One strictly localized eigenbasis is

$H \;=\; H_{0} \;+\; H_{\rm int},$0

and an initial single-site excitation undergoes perfectly periodic breathing with period

$H \;=\; H_{0} \;+\; H_{\rm int},$1

A common source of confusion is that flatness, topology, and caging are related but distinct: flatness is realized at the special interference point, topology persists in a finite parameter region, and imbalance or rung hopping can drive the system into a trivial phase (Zurita et al., 2019).

3. Repulsive interactions, correlated topology, and magnetic phases

Repulsive interactions convert the Creutz ladder into a genuine Creutz–Hubbard model. In the bosonic circuit-QED setting, turning on the transmon Kerr nonlinearity yields a bosonic Creutz–Hubbard model in which all Hubbard parameters $H \;=\; H_{0} \;+\; H_{\rm int},$2 are fully tunable. By tuning $H \;=\; H_{0} \;+\; H_{\rm int},$3 one expects a superfluid–to–Mott-insulator transition, whose critical value depends on the ladder geometry and flux $H \;=\; H_{0} \;+\; H_{\rm int},$4. No full many-body phase diagram is plotted in that work (Alaeian et al., 2018).

For spinless fermions at half filling, the synthetic Creutz–Hubbard model exhibits a much more detailed correlated phase structure. In the weak-coupling regime $H \;=\; H_{0} \;+\; H_{\rm int},$5, a particle–hole transformation plus Jordan–Wigner mapping gives two coupled transverse-field Ising chains, leading to a renormalized critical line between the topological insulator and the orbital paramagnet,

$H \;=\; H_{0} \;+\; H_{\rm int},$6

In the strong-coupling regime $H \;=\; H_{0} \;+\; H_{\rm int},$7, second-order perturbation theory gives an orbital quantum Ising model on the rungs,

$H \;=\; H_{0} \;+\; H_{\rm int},$8

with orbital ferromagnetic order for $H \;=\; H_{0} \;+\; H_{\rm int},$9. The corresponding strong-coupling phase boundary is

$$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$0

At $$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$1, the balanced model has a topological-insulator to orbital-ferromagnet transition at

$$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$2

The full half-filled phase diagram contains a topological insulator, an orbital ferromagnet, and an orbital paramagnet. The topological-insulator to orbital-paramagnet line is a Dirac CFT with $$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$3, whereas the other two lines are Majorana CFTs with $$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$4 (Jünemann et al., 2016).

The numerical characterization of these phases uses matrix-product-state DMRG and tracks the leg-imbalance susceptibility, single- and two-particle gaps, orbital magnetizations, entanglement entropy scaling, and the entanglement spectrum. In the topological-insulator phase all Schmidt levels occur in even degeneracy, while in the orbital-paramagnetic and orbital-ferromagnetic phases there is no systematic degeneracy (Jünemann et al., 2016).

A later spinful mean-field study identified an additional interaction-driven scenario: as the on-site Hubbard interaction grows, the ground state abruptly switches from an anti-ferromagnetic configuration to a ferromagnetic one, and this magnetic transition coincides with a topological transition marked by a quantized jump in the Zak phase from $$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$5 to $$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$6. In that treatment, the Fubini–Study metric diverges at the gap-closing point, and above the critical interaction there coexist a lower-energy ferromagnetic trivial branch and a higher-energy anti-ferromagnetic topological metastable branch (Espinosa-Champo et al., 30 Sep 2025).

4. Attractive interactions and pair superfluidity

The attractive Creutz–Hubbard problem has been studied for spin-$$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$7 fermions with on-site interaction $$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$8. In one standard form the Hamiltonian is

$$\begin{aligned} H_{0} &= -\sum_{j=1}^{L-1} \sum_{\alpha\in\{A,B\}} \Big[ J_\alpha\,c_{j+1,\alpha}^\dagger\,c_{j,\alpha} +J\,c_{j+1,\alpha}^\dagger\,c_{j,\bar\alpha} +\tfrac m2\,c_{j,\alpha}^\dagger\,c_{j,\bar\alpha} +\mathrm{h.c.} \Big], \ J_\alpha&=J\,e^{i\,\sigma_\alpha\,\phi/2},\quad \sigma_A=+1,\;\sigma_B=-1. \end{aligned}$$9

and the noninteracting spectrum consists of two completely dispersionless bands,

$H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$0

The single-particle eigenstates can be chosen as compact four-site plaquette states, so the flatness follows from exact destructive interference in the hopping network (Mondaini et al., 2018).

As soon as $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$1, two up- and down-fermions on the same site can bind into a local pair. In the weak-coupling regime $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$2, projection onto the lower flat band gives a flat-band gap equation implying $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$3 for infinitesimal $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$4, and the superfluid weight is

$H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$5

In the strong-coupling regime $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$6, tightly bound on-site pairs behave as hardcore bosons with effective Hamiltonian

$H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$7

This yields a smooth BCS–BEC crossover in which pairs evolve from large weak-coupling objects to tightly bound strong-coupling bosons (Mondaini et al., 2018).

Mondaini et al. used exact diagonalization on rings up to $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$8 sites and DMRG on open chains up to $H_{\rm int} =U\sum_{j}n_{j,A}\,n_{j,B},$9 sites. For $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$0, the two-particle gap $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$1 vanishes in the thermodynamic limit, while the pair correlation function exhibits algebraic decay $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$2. In the Creutz ladder, $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$3 decreases from $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$4 at weak coupling toward a saturated value $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$5 for $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$6. In a standard two-leg nonflat ladder at the same filling, the corresponding exponent remains above $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$7. This sharper contrast is one reason the Creutz geometry is used as a reference model for flat-band pairing (Mondaini et al., 2018).

5. Doublons, trions, and composite topological quasiparticles

Repulsive Hubbard-type interactions can also generate bound states. In the photonic Creutz–Hubbard analysis of Zurita et al., a repulsive interaction term gives rise to repulsively bound pairs termed doublons. In the strong-coupling regime $$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$8, a Schrieffer–Wolff projection onto the doubly occupied subspace yields an effective doublon Hamiltonian which is exactly a renormalized Creutz ladder with

$$H_{\rm int} =\frac U2\sum_{j,\alpha}n_{j,\alpha}(n_{j,\alpha}-1).$$9

Accordingly, one recovers flat bands, topology, and Aharonov–Bohm caging, but now for doublons. At $J$0, the doublon bands flatten at

$J$1

This produces a striking sector dependence: at $J$2, single-particle Aharonov–Bohm caging coexists with itinerant doublons, whereas at $J$3 one finds itinerant single-particles coexisting with caged doublons (Zurita et al., 2019).

Pelegrí et al. extended the analysis to few interacting bosons in the flat-band Creutz ladder. In the flat-band limit $J$4, the two-boson sector contains on-site-plaquette doublons with dispersive bands and nearest-neighbour-plaquette doublons forming four strictly flat bands. On an open ladder there are two exact edge doublon states formed by bonding a single-particle edge mode at $J$5 with a neighbouring compact localized state, with energies

$J$6

In the strong-$J$7 limit, the effective doublon and trion models are themselves Creutz ladders,

$J$8

so that doublons and trions inherit flat-band and topological features from the original single-particle model. Trion flat bands occur at $J$9 (Pelegrí et al., 2024).

The same work also shows perfect two-body Aharonov–Bohm caging for arbitrary interaction strengths in the flat-band limit. For the initial state

$\ket{j,\alpha}$00

the time evolution remains perfectly confined within a two-plaquette cage for all $\ket{j,\alpha}$01. Matrix-product-state time evolution further provides numerical evidence that the main features of this effect are preserved in an interacting many-body scenario, resulting in many-body Aharonov–Bohm caging (Pelegrí et al., 2024).

6. Disorder, flat-band many-body localization, and re-entrant thermalization

In the interacting flat-band Creutz ladder with disorder, the Hamiltonian is decomposed as

$\ket{j,\alpha}$02

where $\ket{j,\alpha}$03 is the flat-band Creutz-ladder kinetic term, $\ket{j,\alpha}$04, and

$\ket{j,\alpha}$05

At the flat-band point $\ket{j,\alpha}$06, the single-particle spectrum is

$\ket{j,\alpha}$07

and there is an exact basis of compact localized states $\ket{j,\alpha}$08 satisfying

$\ket{j,\alpha}$09

When $\ket{j,\alpha}$10 and $\ket{j,\alpha}$11, the occupations $\ket{j,\alpha}$12 commute with the Hamiltonian, producing an extensive set of exact local integrals of motion and a stable flat-band many-body localized phase (Orito et al., 2021).

The disorder response is non-generic. Weak disorder first destroys the flat-band many-body localized phase and induces a thermal phase; further increase of disorder leads to the conventional many-body-localized phase. The resulting three-regime structure is

$\ket{j,\alpha}$13

with approximate two-particle boundaries

$\ket{j,\alpha}$14

and typical numbers $\ket{j,\alpha}$15, $\ket{j,\alpha}$16 in units of $\ket{j,\alpha}$17 (Orito et al., 2021).

The dynamical diagnostics track the same structure. In the FMBL phase, the half-chain entanglement entropy remains bounded, the number entropy stays small, the fourth moment $\ket{j,\alpha}$18 is approximately constant, and local densities show persistent oscillations associated with Aharonov–Bohm caging. In the thermal regime, entanglement grows rapidly to a volume law and the return probability decays to zero. In the conventional MBL regime, entanglement shows the familiar logarithmic growth $\ket{j,\alpha}$19 before saturation (Orito et al., 2021).

This phase sequence is often singled out because it reverses the naive expectation that added disorder should only strengthen localization.

7. Quantum-simulation platforms and experimental realizations

The Creutz–Hubbard ladder has been pursued in several synthetic platforms. In ultracold atoms, one proposal synthesizes the ladder in a one-dimensional spin-independent optical lattice by mapping the two legs onto two hyperfine states of fermionic atoms. Static tilting suppresses bare hopping, Raman-assisted inter-leg tunneling produces diagonal spin-flip links, and lattice shaking restores intra-leg hopping with Peierls phases so that the net flux per square is $\ket{j,\alpha}$20. The Raman detuning yields the leg imbalance $\ket{j,\alpha}$21, while the on-site $\ket{j,\alpha}$22-wave interaction gives $\ket{j,\alpha}$23. Suggested probes include Ramsey interferometry with Bloch oscillations for the Zak phase, localized Bragg spectroscopy for edge states, spin-resolved imaging for leg imbalance, and spin-structure-factor measurements for orbital magnetization (Jünemann et al., 2016).

A related experimental step realized a topological Creutz ladder for ultracold fermionic atoms in a resonantly driven one-dimensional optical lattice. There the two legs are the two lowest orbital states $\ket{j,\alpha}$24 and $\ket{j,\alpha}$25, and the cross inter-leg links are generated via two-photon resonant coupling by periodic lattice shaking. The characteristic pseudo-spin winding was demonstrated using momentum-resolved Ramsey-type interferometric measurements. The same framework includes a two-tone driving method to extend the inter-leg link control and a topological charge pumping protocol with quantized pumped charge

$\ket{j,\alpha}$26

when the adiabatic path encircles a band-gap closing point (Kang et al., 2019).

In photonics, the full Creutz and Creutz–Hubbard ladders can be emulated by arrays of coupled optical waveguides or ring resonators, using either real and synthetic dimensions. One route uses a one-dimensional chain with a synthetic leg degree of freedom carried by transverse modes. Another uses a rhombus chain with detuned spine sites, where the effective diagonal coupling is

$\ket{j,\alpha}$27

For the interacting two-particle sector, the standard mapping from two interacting particles in one dimension to a single particle on a two-dimensional lattice leads to a quasi-2D array or a 3D waveguide network with a synthetic particle dimension. Proposed protocols include direct observation of single-particle Aharonov–Bohm caging, edge-state trapping, and doublon caging or doublon edge modes. Typical parameters quoted for waveguide devices are $\ket{j,\alpha}$28–$\ket{j,\alpha}$29, synthetic flux via $\ket{j,\alpha}$30 modulation or $\ket{j,\alpha}$31, and effective $\ket{j,\alpha}$32 in nonlinear ring-resonator networks (Zurita et al., 2019).

Circuit QED provides a bosonic route. There, coupled transmons have their Josephson energies modulated harmonically in time,

$\ket{j,\alpha}$33

which makes the transmon frequencies time dependent. In the rotating-frame Floquet description, the hopping acquires a time-periodic phase factor, and period averaging yields an effective complex hopping controlled by Bessel functions and by the site-dependent phase offsets $\ket{j,\alpha}$34. These phase offsets act as an artificial gauge potential, and their difference controls the Peierls phase on each link. Because of the intrinsic non-linearity of the transmon qubits, such lattices are an ideal platform for simulating Bose–Hubbard Hamiltonians with non-trivial gauge fields, and the construction is directly realizable in a superconducting circuit and easily extendable to a scalable lattice (Alaeian et al., 2018).

Across these platforms, the Creutz–Hubbard ladder functions as a compact model in which flat bands, edge topology, synthetic gauge structure, and strong correlations can be tuned in a controlled way. Its unusually broad phenomenology spans topological insulators, orbital magnetic phases, pair superfluidity, doublon and trion quasiparticles, Aharonov–Bohm caging, flat-band many-body localization, and interaction-driven topological transitions.