2.1.

Theoretical Model

Nonlinear evolutional dynamics of a 2D ultracold Fermi gas (Bardeen Cooper Schrieffer superfluid of spin-1/2 Fermionic atoms) loaded onto the moiré optical lattices in mean-field theory is governed by the density-functional equation (theory), whose scaled wave function (complex order parameter), $\psi$, yields^16–18

The 2D moiré optical lattices, consisting of two uniform square lattices with a twisted angle $\theta$ to each other, follow

One can know from Eq. (2) that the moiré optical lattice returns to the normal square lattice if we set $V_2=0$ or $\theta =0$. The extraordinary feature of such a moiré pattern is that, when gradually changing the rotating angle $\theta$ (between the two sublattices), the moiré optical lattice will be continuously transformed into the periodic and aperiodic structures, thereby filling a gap crossing of perfect periodic potentials and aperiodic structures and opening a new strategy for light-field manipulation (bandgap engineering) and its nonlinear dynamic control dependent merely on a new rotating degree of freedom of the sublattices even though under shallow lattice modulation. As pointed out in previous publications,^{30–33,44,49,50} such a moiré optical lattice allows the existence of the widest finite gaps separated by neighboring extremely flat Bloch bands provided that it is the perfect periodic pattern which, in fact, is achieved only if we consider the special case of Pythagorean angles $\theta =\arctan (a/b)$, with $\cos \theta =a/c$, $\sin \theta =b/c$, and natural numbers $(a,b,c)$ constitute Pythagorean triples, $a^2+b^2=c^2$.

The stationary solitary matter-wave solution $\varphi$ under chemical potential $\mu$ can be constructed as $\psi =\varphi e^{-i\mu t}$, and after substituting it into the dynamical fundamental equation [Eq. (1)], one can get the stationary equation,

Equation (4) holds the norm (number of atoms) $N=\iint |\varphi {|}^2\mathrm{d}x\mathrm{d}y$.

2.2.

Linear Bandgap Properties

Although the fact that the moiré patterns could be twisted with arbitrary angles $\theta$, allowing ample band diagrams to emerge accordingly, the Pythagorean lattices (the moiré optical lattice under the unique case of Pythagorean angles $\theta$, as stressed above) support the optimal bandgap tuning, which is the spatial tailoring of the largest width of the finite gaps, within which the localized gap modes of diverse types to be excited are of our fundamental interest. Our attention here is only paid to the moiré optical lattices at different Pythagorean twisted angles $\theta$, which, in fact, are equipped with periodic translational symmetry, and therefore, the general linear Bloch bandgap theory widely used in solid-state physics and other periodic potentials, including photonic crystals and optical lattices, can be readily applied here. Bearing this in mind, one can produce the underlying linear Bloch spectrum after a routine operation.

Contour profiles of Pythagorean moiré optical lattices under three different Pythagorean twisted angles $\theta =\arctan (3/4)$, $\theta =\arctan (5/12)$, and $\theta =\arctan (8/15)$ are depicted in Figs. 1(a)–1(c), respectively. The three twisted angles enlarge the lattice period of moiré optical lattices as compared with the original square lattices (${\epsilon }_2=0$ or $\theta =0$). Ignoring the last nonlinear term of Eq. (4), calculating the corresponding eigenvalues within the first Brillouin zone in the reciprocal space [which is depicted in the inset of Fig. 1(d)], one can get the linear Bloch spectrum under the same strength difference $P=0.5$, as displayed separately in Figs. 1(d)–1(f) and in the bottom line of the figure with red dashed vertical lines. Once again, we can see that, as those discovered in other Pythagorean moiré lattices, despite being shallow, there are abundant flat Bloch bands for the Pythagorean moiré optical lattices; strikingly, the former features a very narrow gap (which we name minigap “$Q$”) between the first and second gaps; the latter two cases exhibit more flat bands, and the higher finite gaps can be wider than their first one. Further research revealed that, for the former case, the widths of both the first gap and minigap decrease when increasing $P$, and the minigap ($Q$) disappears at $P>0.53$, while the second gap’s width grows gradually [Fig. 1(g)]; by comparison, the first gap’s width increases, whereas the decrease is only for the fourth and fifth gaps for the latter two cases, according to Figs. 1(h) and 1(i). These rich linear bandgap properties of moiré optical lattices indicate that very flexible and tunable bandgap engineering can be achieved by simply changing the rotating angle $\theta$ and strength difference $P$, providing a new linear (flat-band) localization mechanism (different from that of Anderson localization) of optical and matter waves,³¹ and robust nonlinear localization of various localized gap modes within those numbered finite gaps of the associated linear Bloch spectra, as will be presented below.

Fig. 1

Linear Bloch-wave band structures of 2D moiré square optical lattices. The geometries, reduced Brillouin zone, and bandgap spectra of the 2D moiré square optical lattices under different Pythagorean twisted angles. Platforms of the moiré lattice (shaded red, lattice potential maxima; shaded blue, lattice potential minima) at different rotation angles: (a) $\theta =\arctan (3/4)$, (b) $\theta =\arctan (5/12)$, and (c) $\theta =\arctan (8/15)$, and the linear Bloch-wave spectrum (indicated as chemical potential $\mu$) of the moiré lattice with different twisted angles: (d) $\theta =\arctan (3/4)$, (e) $\theta =\arctan (5/12)$, and (f) $\theta =\arctan (8/15)$; subplot in panel (d) shows the first reduced Brillouin zone in the reciprocal space, where the associated high symmetry points are marked. The width of the finite gaps versus the relative strength of two layers of the lattice ($P$) under different twisted angles: (g) $\theta =\arctan (3/4)$, (h) $\theta =\arctan (5/12)$, and (i) $\theta =\arctan (8/15)$. The second and third lines: the Roman numerals I, II, III, IV, and V are the first, second, third, fourth, and fifth finite gaps, respectively; $\mathit{Q}$ in panels (d) and (g) is the very narrow minigap. $P=0.5$ for panels (a)–(f).

Before going further, let us introduce the numerical calculation recipes for launching the localized gap modes and testing their stability, which yield the following: the stationary solutions of the localized gap modes under consideration are first constructed through Eq. (4) using the widely used method—modified squared-operator iteration method;⁵³ then, the stability and instability properties of the gap modes thus found will be scrutinized in direct perturbed numerical computation of the dynamical equation [Eq. (1)] by taking a converged and highly accurate computational approach called the fourth-order Runge–Kutta method.

3.1.

Fundamental Gap Solitons

We first survey the excitation of the simplest matter-wave localized gap modes with a single peak and isotropic shapes called fundamental gap solitons, which can be upheld by the Pythagorean moiré optical lattices [Eq. (2)], as done in conventional periodic potentials. As far as the aforementioned three Pythagorean twisted angles [$\theta =\arctan (3/4)$, $\theta =\arctan (5/12)$, and $\theta =\arctan (8/15)$] are concerned, the corresponding relationships among the number of atoms, $N$, and chemical potential $\mu$ for the fundamental gap solitons are, respectively, depicted in Figs. 2(a)–2(c), showing a growing tendency through the first gap to the higher ones, and abiding by the well-known principle of the “anti-Vakhitov–Kolokolov” criterion, $\text{d}N/\text{d}\mu >0$, an indispensable condition for the stable gap solitons under the regime of repulsive (defocusing) nonlinearity.^6,8–15,54 Detailed numerical perturbed calculations in the dynamics equation [Eq. (1)] reveal that these fundamental gap solitons are robustly stable localized modes existing in the middle of the finite gaps; they tend to be unstable near the edges of Bloch bands, as underlined in the figures. Moreover, we have examined the dependence between $N$ and the same sublattice depth $V_0=V_1=V_2$, which is summed up in Fig. 2(d), demonstrating a decreasing trend with the increase of $V_0$; it coincides with the significant experimental finding of low threshold (of the laser power in nonlinear optics) for the formation of optical solitons in photonic moiré lattices. It is necessary to note that the typical profiles of the so-found fundamental gap solitons are shown in Figs. 2(e)–2(g), where the former is a stable gap soliton owning negligible side peaks, while the latter two are unstable ones accompanied by multiple side peaks that are induced by strong Bragg scattering. Such side peaks of the latter two solitons may cause instability in the course of time evolution.

Fig. 2

Properties of fundamental gap solitons of ultracold Fermi gases controlled by moiré square optical lattices. Number of atoms, $N$, versus chemical potential $\mu$ (a)–(c), respectively, and the same sublattice depth ($P=1$) under $V_0=V_1=V_2$ (d) for the fundamental gap solitons, whose profiles for marked points ($A$, $B$, and $C$) in panel (a) are shown in the third line, in 2D moiré square optical lattices at twisted angles: (a), (d) $\theta =\arctan (3/4)$, (b) $\theta =\arctan (5/12)$, and (c) $\theta =\arctan (8/15)$. Other parameters: (d) $\mu =4$; (e) $\mu =4$, $N=14.7$; (f) $\mu =4.7$, $N=36.8$; (g) $\mu =5.18$, $N=66.4$. Here and below, red dashed lines in panels (a)–(c) represent unstable regions.

In the perturbed time evolution, as shown in Fig. 3, the stable fundamental gap soliton sustains its shape [Fig. 3(a)], while there are oscillating or growing side peaks (lobes) for the unstable fundamental modes [Figs. 3(b) and 3(c)]. Further insights into their evolutions confirm that the unstable fundamental gap solitons are only in the weakly unstable mode; this once again verifies the excellent platform for soliton formation and control in physical settings with moiré optical lattices.

Fig. 3

Perturbed evolutions of fundamental gap solitons. Top-view profiles (a)–(c) and the perturbed evolution in real time (d)–(f) of typical fundamental gap solitons for the marked points ($A$, $B$, and $C$) in Fig. 2(a): stable (a), weakly unstable (b), and strongly unstable (c) gap solitons.

3.2.

Vortex Gap Solitons with Topological Charge $S=1$

A more interesting thing is that whether such Pythagorean moiré optical lattices can give rise to robust vortex gap solitons and how to stabilize them at will. It is thus the target in this section to reveal these issues for gap vortices, the vortex gap solitons with topological charge $S$; for the sake of discussion, we only consider the case of $S=1$. As reported elsewhere, such gap vortices are usually featured by a null value in the middle and consist of four bright (fundamental gap solitons) modes/peaks that are entangled with $\pi /2$ phase shifts between them, enabling the whole $2\pi$ phase circulation ($S=1$). The stability regions of the gap vortices reduce a little bit compared with their fundamental gap solitons counterpart because of the existence of destructive interaction among the four solitons, as can be seen from the curve $N(\mu )$ and the associated stability regions portrayed in Figs. 4(a)–4(c) for the three Pythagorean angles. Also, a decreasing relation $N(V_0)$ proceeds for the vortex gap solitons, according to Fig. 4(d). The schematic shapes of the vortex gap solitons are displayed in Figs. 4(h)–4(j), the cross-sectional cut at $y=1.5$ is displayed in Figs. 4(e)–4(g), and the associated spatial phase structures are depicted in the bottom line of Fig. 4 [see Figs. 4(k)–4(m)], showcasing the twisted modes unique for moiré lattices.

Fig. 4

Properties of gap vortices of ultracold Fermi gases loaded upon moiré square optical lattices. Number of atoms, $N$, versus chemical potential $\mu$ (a)–(c) and the same sublattice depth ($P=1$) under $V_0=V_1=V_2$ (d) for the vortex gap solitons at topological charge $S=1$ in 2D moiré square optical lattices at twisted angles: (a), (d) $\theta =\arctan (3/4)$, (b) $\theta =\arctan (5/12)$, and (c) $\theta =\arctan (8/15)$. Other parameters: (d) $\mu =4$; (e) $\mu =3.5$, $N=9.4$; (f) $\mu =4.6$, $N=120.5$; (g) $\mu =5.25$, $N=256.1$. Red dashed lines in panels (a)–(c) represent unstable regions. The cross section (taken at $y=1.5$), the associated profiles, and phase structures are depicted in the third, fourth, and fifth lines, respectively.

In direct perturbed simulations by solving the dynamical equation [Eq. (1)], we can see that the stable vortex gap soliton remains unchanged in long time evolution [Fig. 5(a)]; and for the unstable ones, as those displayed in Figs. 5(b) and 5(c), only their side peaks undergo regular oscillations and small strength growth, demonstrating weak instability.

Fig. 5

Perturbed evolutions of gap vortices. Top-view profiles (a)–(c) and the perturbed evolution in real time (d)–(f) of typical vortex gap solitons for the marked points ($D$, $E$, and $F$) in Fig. 4(a): stable (a), weakly unstable (b), and strongly unstable (c) gap vortices.