Floquet-engineered moire quasicrystal patterns of ultracold Bose gases in twisted bilayer optical lattices

TL;DR

This paper demonstrates how Floquet engineering can induce and manipulate moire quasicrystal patterns in ultracold Bose gases within twisted bilayer optical lattices, revealing new symmetry properties and collective modes.

Contribution

It introduces a method to generate and control moire quasicrystal patterns in ultracold gases using Floquet modulation, a novel approach in quantum simulation.

Findings

Emergence of twelve-fold quasicrystal pattern under specific parameters

Pattern symmetry can be tuned by modulation frequencies and amplitudes

High rotational symmetry observed in momentum-space distributions

Abstract

We investigate the formation of moire quasicrystal patterns in Bose gasses confined in twisted bilayer optical lattices via Floquet-engineered intralayer atomic interactions. Dynamical evolutions of the total density wave amplitude exhibit the stage for the emergence of moire quasicrystal patterns, where the pattern formation is closely associated with the momenta of collective modes excited by the weak periodic drive. Through analyzing the radial and angular density wave amplitude, we find that these new collective modes are only coupled radially and cannot be decoupled eventually. The symmetry of quasicrystal patterns can be easily manipulated by the modulation frequencies and amplitudes. Reducing the frequencies and increasing the amplitudes can both facilitate lattice symmetry breaking and the subsequent emergence of rotational symmetry. Notably, a twelve-fold quasicrystal pattern emerges under specific parameters, closely resembling the moire quasicrystal in twisted bilayer graphene. The momentum-space distributions also exhibit high rotational symmetry, which is consistent with the real-space patterns at specific evolution times. Our findings establish a new quantum platform for exploring quasicrystals and their symmetry properties in ultracold bosonic systems.

Equations

$$i\hslash \frac{\partial {\psi }_1}{\partial t}=$$

$${\psi }_1+\hslash \Omega {\psi }_2,$$

$$i\hslash \frac{\partial {\psi }_2}{\partial t}=$$

$${\psi }_2+\hslash \Omega {\psi }_1,$$

$$\text{IPR}=\frac{\int |\psi {|}^{4}d\mathbf{r}}{(\int |\psi {|}^{2}d\mathbf{r}{)}^2},$$

Full text

I Introduction

The recent appearance of twisted two-dimensional (2D) materials has profoundly impacted condensed-matter physics. In particular, the discovery of unconventional superconductivity [7] and correlated insulating behavior [6] in twisted bilayer graphene has revolutionized the field of moiré physics. A moiré lattice refers to a new periodic lattice structure formed when two identical periodic lattices are rotated by a specific angle. Moiré lattices with a small twist angle host rich physical properties, such as flat bands [44, 25, 24], quantum anomalous Hall effect [37, 33, 59, 52, 5, 28], moiré excitons [38, 3, 4, 18], magnetism [39, 36, 26, 8, 35, 47] and strongly correlated insulators [57, 27, 53, 23].

Ultracold Bose gases in optical lattices are becoming an ideal platform for studying twisted bilayer systems, owing to their purity and high tunability. Quantum simulations in twisted bilayer optical lattices are realized by interfering multiple sets of laser beams to create the desired lattice geometries. Different atomic hyperfine spin states are employed to construct spin-dependent lattices of two-component bosonic systems. The schemes utilizing spin-dependent square [17] and hexagonal [29] optical lattices to simulate twisted bilayers have been proposed in the past few years. Notably, Zhang* et al.* have experimentally realized atomic Bose-Einstein condensates (BECs) in twisted bilayer lattices [31]. Based on this experimental foundation, many new research findings concerning twisted bilayer lattices have been reported, including fractal structures [34, 50], interaction-induced interlayer coupling [45, 58], solitons [12, 48], topology [49, 22], dipolar bosons [60] and quasicrystal optical lattices [10].

Previous studies have primarily focused on the ground state of square lattices, in which quasicrystal optical lattices are obtained by adding an external quasiperiodic potential. Interestingly, when the twist angle of twisted bilayer graphene is $30^{\circ}$, its electronic state distributions form a moiré quasicrystal with twelve-fold ($D_{12}$) symmetries [2, 32, 15]. In bosonic systems, the competition between disorder and interaction is the key to the emergence of quasicrystal structures, which typically results in a new quantum phase named Bose glass [16, 13, 20, 56]. Quasicrystal structures are controlled by externally applied quasiperiodic potentials in ultracold bosonic systems, while graphene moiré quasicrystals are primarily governed by the twist angle.

We have investigated parametric excitations and Faraday pattern formation by periodically modulating the atomic scattering lengths in binary BECs [51]. On the other hand, we note that Chin* et al.* have experimentally realized hexagonal lattice density wave (DW) patterns in BECs through Floquet-engineered atomic interactions [61]. Meanwhile, researches on inducing novel quantum phases through Floquet engineering in optical lattices, achieved by periodically shaking the lattice or modulating atomic interactions, have been reported [55, 62, 40].

This work explores the formation of novel moiré quasicrystal DW patterns in twisted bilayer hexagonal lattices via Floquet-engineered intralayer atomic interactions, in which additional quasiperiodic potentials are not necessary. In quantum systems, DW patterns represent a form of spatial order arising from nonequilibrium dynamics induced by driving fields. Investigating their formation mechanism, dynamical evolution and control is essential for understanding collective effects, symmetry breaking and restoration in complex systems.

II Theoretical Model

The twisted hexagonal bilayer optical lattice of ultracold bosonic system can be realized using synthetic dimension techniques in current experiment [31, 41, 43]. Three laser beams of these two wavelengths intersect at $120^{\circ}$ in the x-y plane to produce two hexagonal lattices $V_{1}$ and $V_{2}$, where each beam is linearly polarized in-plane. A hexagonal lattice potential is given by $V(\boldsymbol{r})=-V_{0}\left|\sum_{j=1}^{3}\epsilon_{j}\exp[i\boldsymbol{k}_{j}\cdot(\boldsymbol{r}-\boldsymbol{r_{0}})]\right|^{2}$, with the lattice depth $V_{0}$, the laser wave vector $k_{j}$ and the polarization $\epsilon_{j}$. These two spin-dependent lattices are twisted by $\theta/2$ each other, where the $\theta$ satisfies $\cos\theta=\frac{n^{2}+m^{2}+4mn}{2(n^{2}+m^{2}+mn)}$ with two integers $(m,n)$ [46]. The system is loaded into a harmonic trap and z-direction is tightly confined by a spin-independent potential, reducing the dynamics to quasi-2D. Meanwhile, the interlayer coupling $\Omega$ of two spin states is controlled by microwaves.

The ground state and dynamics of weakly interacting superfluid BECs in shallow optical lattices can be effectively studied by mean-field approximation and Gross-Pitaevskii (GP) equations. Therefore, the system’s coupled GP equations reads

[TABLE]

where $V_{trap}=\frac{1}{2}m(\omega_{x}^{2}x^{2}+\omega_{y}^{2}y^{2}+\omega_{z}z^{2})$ is a harmonic trap and $\psi_{i}(i=1,2)$ are the wave fuctions normalized as $\int\int\left(|\psi_{i}|^{2}\right)dxdy=N$, with $N$ the total atom number. $g_{ii},g_{ij}$ characterize the intralayer and interlayer atomic interactions strength. Focus on the SU(2) symmetric interaction in calculations, we set $g_{11}=g_{22}=g_{12}=g_{0}$ because of the similarity in scattering lengths $a_{11}$, $a_{22}$ and $a_{12}$ for the 87Rb atoms. In our numerics, the trapping frequency is $\omega_{x,y,z}=2\pi\times\{20,20,1000\}$ Hz, hence the system is seen as quasi-2D because of the strong potential in z-direction. The quasi-2D interaction satisfies $g_{2D}=g_{0}/{\sqrt{2\pi}a_{z}}$, where $g_{0}=4\pi\hbar^{2}a/m$ and $a_{z}=\sqrt{{\hbar}/m\omega_{z}}$ represents the characteristic length along the z axis. Here the scattering length for the 87Rb atoms is $a=100a_{0}$ with $a_{0}$ the Bohr radius.

The GP equations is numerically solved via the imaginary time evolution method to obtain the ground state. We make the scattering length oscillate as $a(t)=100a_{0}+a_{m}\sin{\omega t}$ by modulating the magnetic field near a Feshbach resonance [9]. This paper mainly investigates the oscillating intralayer atomic interactions $g_{11}(t)=g_{22}(t)=4\pi\hbar^{2}a(t)/(\sqrt{2\pi}a_{z}m)$, while keeping the interlayer-atomic interactions $g_{12}=g_{0}$. Then we use the dynamical GP equations to investigate the evolution of driven DW patterns in both real and momentum spaces.

III Dynamical behavior

To describe the spatial localization of the superfluid DW, we introduce a physical quantity called Inverse participation ratio (IPR) [11, 54], which is given by

[TABLE]

where we set $\int|\psi|^{2}d\boldsymbol{r}=1$. Figure 2(a) illustrates the dependence of IPR on the lattice depth $V_{0}$ and interlayer coupling $\Omega$. The IPR increases with $V_{0}$, indicating a more localized DW distribution. In contrast, it exhibits a non-monotonic dependence on $\Omega$, first decreasing and then increasing. The minimum of the IPR-$\Omega$ curve is around $\Omega=E_{r}$ with energy unit $E_{r}=\hbar^{2}k^{2}/2m$, indicating the best superfluidity. In addition, the GP equations is used to describe the superfluid BECs in shallow optical lattices, but $V_{0}$ must remain sufficiently large to ensure a well-defined lattice distribution. Based on this condition, we select $\Omega=E_{r}$ and $V_{0}=4E_{r}$ as basic parameters for our numerics (red star in Fig. 2(a)). The other parameters used in our numerical computations are $N=2\times 10^{4}$ and $\theta=9.43^{\circ}$ (m=3, n=4).

However, the IPR is always at the order $N^{-1}$ in the whole dynamical evolution, whose variation is too small to reflect the localization of the system. Therefore, we introduce the mean-square displacement (MSD) $\sigma(t)=[\frac{\int r^{2}|\psi|^{2}d\boldsymbol{r}}{\int|\psi|^{2}d\boldsymbol{r}}]^{\frac{1}{2}}$ [45, 1, 63] to quantitatively distinguish the different dynamical behaviors of the system at a long evolution time. $\sigma(t)$ keep a stable value at early evolution time, indicating the system remains unexcited and exhibits localized behavior. Then $\sigma(t)$ tends to linearly increase in a large evolution time, characterizing its delocalized behavior (Fig. 2(b) red curve). Notably, this linear increase comprises different segments with distinct slopes, implying four evolution stages.

IV Moiré quasicrystal patterns

IV.1 Density wave amplitude

We use dynamical DW $n(\boldsymbol{r},t)=|\psi(\boldsymbol{r},t)|^{2}$ and $n(\boldsymbol{k},t)=|\psi(\boldsymbol{k},t)|^{2}$ to investigate the evolution of patterns. To characterize the stages of DW evolution, we compute the total DW amplitude, $A_{total}(\tilde{n})=\int\tilde{n}(\mathbf{k})d\mathbf{k}$ and divide the dynamical evolution into four stages by the variation of $A_{total}(\tilde{n})$. Here, $\tilde{n}(\mathbf{k})=\int e^{-i\mathbf{k}\cdot\mathbf{r}}{n}(\mathbf{r})d\mathbf{r}/2\pi$ is the the density Fourier transformation.

The evolution of driven DW patterns can be divided into four distinct stages (Fig. 3(a)). In the preparation stage, the system remains unexcited and persistently maintains its initial moiré lattice pattern. In the excitation stage, the drive-induced excitation initiates at the central lattice and propagates outward. This process lays the foundation for subsequent lattice symmetry breaking. In the pattern-forming stage, the lattice symmetry is broken and rotational symmetry emerges. The resulting patterns lack strict periodicity but exhibit long-range order with quasicrystal characteristics, and are therefore termed moiré quasicrystal patterns. In the nonlinear stage, continuous energy accumulation eventually exceeds the threshold required to sustain rotational symmetry, leading to disordered pattern configurations. In our numerics, the energy of ground state is $4.8E_{r}$, compared to an energy threshold of approximately $7.8E_{r}$ for pattern stability.

Chin* et al.* experimentally reported the hexagonal lattice pattern formation in driven BECs, characterized by six distinct modes in momentum space [61]. In this letter, we investigate the system using the coupled GP equations, where the DW distribution in momentum space condenses to a point. Therefore, we amplify the momentum-space DW via a logarithmic function $\ln{n(\boldsymbol{k},t)}$. This amplified DW distribution similarly exhibits six distinct modes under lattice symmetry, with the modes equally spaced by $\pi/3$ in directions. In particular, each mode comprises two sub-modes, which is a direct consequence of the bilayer lattices. The angle between these sub-modes precisely matches the twist angle of the bilayer lattices (Fig. 3(b) bottom left).

The six modes remain distinct during the preparation stage. Upon entering the excitation stage, these modes progressively develop interconnections. This transformation culminates in the pattern-forming stage where the six modes connect into a unified and discernible hexagram, indicating lattice symmetry breaking and the emergence of a new symmetry in real space. Finally, this configuration becomes largely indistinguishable in the nonlinear stage (Fig. 3(b) bottom row). Although the modulation parameters ($a_{m}$, $\omega$) supporting the moiré quasicrystal patterns vary, all six modes undergo this transformation. Only the duration of each stage differs. The changes in the patterns in both spaces provide clear validation of these four distinct dynamical stages.

In our calculations, the strength of the drive for one layer is $E_{g}=\left\langle|\psi_{1}|^{2}\right|g_{11}(t=0)\left||\psi_{1}|^{2}\right\rangle\approx 0.02$, which is weak compared with other energy scales ($E_{V_{1}}=\left\langle|\psi_{1}|^{2}\right|V_{1}\left||\psi_{1}|^{2}\right\rangle\approx 2.36$, $E_{\Omega}=\left\langle|\psi_{1}|^{2}\right|\Omega\left||\psi_{1}|^{2}\right\rangle\approx 0.50$). Obviously, the energies keep $E_{V_{1}}>E_{\Omega}\gg E_{g}$ in the whole dynamical evolution. Therefore, we think the formation of Moiré quasicrystal patterns is associated with momenta of collective modes that are resonantly excited by the periodic drive. To quantitatively analyze these momenta, we extract the DW amplitude, $A_{\mathbf{k}}=n_{0}^{-1}\int_{\left|\mathbf{k}\right|}\tilde{n}(\mathbf{k})d\mathbf{k}$ [14]. However, different from the Ref. [14], our system is more complex and not easy to get an analytical derivation of the parametric resonance condition from the Mathieu equation. Therefore, we numerically calculate the dependence of $A_{\mathbf{k}}$ on $\mathbf{k}$ during the dynamical evolution.

For ground state, the DW is controlled by six distinct momenta (Fig. 4a red curve), where the two main momenta ($\mathbf{k}_{1}=0.03k_{r}$, $\mathbf{k}_{2}=2k_{r}$) correspond precisely to the central mode and its surrounding six modes in the momentum space (Fig. 3(b) bottom line). A multitude of additional collective modes emerge during the pattern-formation stage ((Fig. 4a blue curve). We track the evolution of these modes throughout the dynamics and uncover several interesting regularities.

In the region of $0<\mathbf{k}<\mathbf{k}_{1}$, $A_{\mathbf{k}}$ keeps decreasing after excitation, indicating these momentum-associated collective modes play a waning influence on pattern formation. In the region of $\mathbf{k}_{1}<\mathbf{k}<\mathbf{k}_{2}$, $A_{\mathbf{k}}$ of original modes exhibits damped-harmonic envelope oscillations around its initial value. In the region of $\mathbf{k}>\mathbf{k}_{2}$ (except for $\mathbf{k}=\mathbf{k}_{3},\mathbf{k}_{4}$), the trend of $A_{\mathbf{k}}$ is similar to that of $A_{\mathbf{total}}$, indicating that collective modes with these momenta are dominant in the excitation (Fig. 4b right). Normally, in the absence of mode coupling, $A_{\mathbf{k}}$ is expected to return to near its initial value or exhibit a pronounced drop after entering the nonlinear stage. Thus we think these new collective modes are coupled each other and eventually cannot be decoupled. In the above region, the weak and peakless oscillations of $A_{\mathbf{k}}$ also reflect radial coupling among the modes. Notably, the three modes ($\mathbf{k}_{1},\mathbf{k}_{3},\mathbf{k}_{4}$) of the ground state retain distinct peaks throughout the whole dynamical evolution. Obviously these three modes have no radial coupling to other additional collective modes, leading to the decline in $A_{\mathbf{k}_{1,3,4}}$ because of the energy transfer from these modes to the remaining uncoupled modes (Fig. 4b left).

From Fig. 3(b), we observe the new modes emerge along specific directions when $\mathbf{k}\leq\mathbf{k}_{3}$. Inspired by the dynamical patterns in momentum-space, we wonder whether the newly appearing collective modes are coupled angularly. To quantify the distribution of these modes along the angular direction, we also evaluate $A_{\phi}=\int_{\phi}\tilde{n}^{\prime}(\mathbf{\mathbf{k}_{\phi}})d\mathbf{\mathbf{k}_{\phi}}$, where $\tilde{n}^{\prime}(\mathbf{\mathbf{\mathbf{k}_{\phi}}})=\int n(\mathbf{r})e^{-i\mathbf{k}_{\phi}\cdot\mathbf{r}}d\mathbf{r}/2\pi$ with magnitude $\left|\mathbf{k}_{\phi}\right|=\left|\mathbf{k}\right|$ and angle $\phi$. Consistent with our expectations, new modes emerge along the directions ($\phi_{2},\phi_{3}$) defined by $\mathbf{k}_{2}$ and $\mathbf{k}_{3}$ (Fig. 4c left), where $\phi_{2}=\frac{(2j+1)\pi}{6}$ and $\phi_{3}=\frac{j\pi}{3}$ with $j=0,1,...,5$. Interestingly, for $\mathbf{k}>\mathbf{k}_{3}$, $A_{\phi}$ consists of a series of distinct peaks (Fig. 4c right), which means the collective modes are angularly uncoupled.

IV.2 The real-space Moiré quasicrystal patterns

The rotational symmetry of the patterns significantly depends on the modulation amplitude $a_{m}$ and frequency $\omega$. We first fix the modulation amplitude and vary the frequency. Increasing the modulation frequency breaks the rotational symmetry of the system. For example, at $a_{m}=50a_{0}$ and low $\omega$, the pattern exhibits $D_{6}$ symmetries and other distinct $D_{6}$ patterns subsequently emerge during the third stage. Each of these non-subharmonic patterns represent a new symmetrical structure. Upon increasing $\omega$ to approximately 252 Hz, the $A_{total}$ decreases significantly, indicating only a minor DW alteration and patterns with lattice characteristics (Fig. 5(b) blue curve). The blue curve in Fig. 5(b) is notably smoother than the other two. We think this frequency matches exactly with the inherent frequency of a collective excitation mode. Under this frequency, system is excited without a preparation stage. In other words, energy is injected into the system directly and efficiently from the drive field, causing a smoother curve. This phenomenon differs from the collective mode momentum of DW yet does not conflict with it. For frequencies above this frequency, the formation of moiré quasicrystal patterns is suppressed, while the duration of the lattice patterns is extended.

We then fix the modulation frequency at $\omega=200$ Hz and vary the modulation amplitude. Increasing $a_{m}$ yields patterns with increasingly complex symmetry structures. As previously demonstrated, a simple $D_{6}$ pattern forms at $a_{m}=50a_{0}$, which means regions other than the center possess only a single layer of $D_{6}$ symmetries. At $a_{m}=75a_{0}$, two additional layers of $D_{6}$ patterns develop. These patterns are not simple six-petal types but rather a more complex cobweb-like structures (Fig. 6(a)).

Remarkably, the system forms a $D_{12}$ pattern at 30 ms with $a_{m}=100a_{0}$ (Fig. 6(b)), closely resembling the moiré quasicrystal realized in twisted bilayer graphene rotated by exactly $30^{\circ}$. This pattern exhibits three distinct layers, whose structure is nearly identical to graphene quasicrystal. In twisted bilayer graphene, the stability of hexagonal atomic structure yields standard Stampfli quasicrystal constructions [42, 19], including triangles, rhombuses, and squares. However, in our driven system, the excited DW presents apparent fluctuations, preventing the formation of patterns with such regular shapes. For a clearer comparison with graphene quasicrystal, we outline the contours of the second layer in Fig. 6(c). The patterns retain $D_{6}$ symmetries at other times.

In general, increasing the modulation amplitude while appropriately decreasing the frequency promotes lattice symmetry breaking and the emergence of rotational symmetry. This process also increases the multiplicity of rotational symmetry, yielding more complex symmetrical structures. For instance, the pattern maintains its lattice symmetry at $a_{m}=50a_{0}$ and $\omega=252$ Hz. Further increasing the modulation amplitude leads to the emergence of quasicrystal patterns. The stable parameter region for moiré quasicrystal patterns agrees well with the regularity of Floquet stable tongues derived from the Mathieu equation (black dotted line in Fig. 8) [30, 21]. The parameters (blue and green lines in Fig. 5(b)) are beyond this stable region, with no moiré quasicrystal patterns emerging.

IV.3 The momentum-space Moiré quasicrystal patterns

We further analyze the symmetry of the DW distribution in momentum space based on Fourier analysis $\psi(\boldsymbol{k})=\int e^{-i\mathbf{k}\cdot\mathbf{r}}{\psi}(\mathbf{r})d\mathbf{r}/2\pi$. The Fourier DW $n(\boldsymbol{k})$ also exhibits high rotational symmetry, typically forming $D_{6}$ patterns. These momentum-space patterns are consistent with real-space patterns at specific times. At $a_{m}=100a_{0}$ and $\omega=200$ Hz, the momentum-space pattern forms a hexagram nearly identical to the real-space pattern at 36.2 ms, but rotated clockwise by $60^{\circ}$. Similarly, at $a_{m}=75a_{0}$ and $\omega=200$ Hz, the patterns in both spaces match well at 42, 44, and 64.5 ms (Fig. 7(b)). At 64.5 ms in particular, both patterns exhibit a symmetrical butterfly structure. The momentum-space patterns also do not exhibit subharmonicity. However, unlike the real-space patterns, some momentum-space patterns exhibit only $D_{4}$ symmetries at the center.

V Conclusions and Discussions

We numerically study the formation of novel moiré quasicrystal patterns of driven ultracold Bose gasses in twisted bilayer hexagonal lattices. Dynamical evolutions of the MSD and total DW amplitude reveal the stage for the emergence of moiré quasicrystal patterns. In our weakly driven system, pattern formation is closely related to the momenta of collective modes. By quantifying the radial and angular DW amplitudes, we obtain the distributions of the emerging collective modes. We conclude that these new collective modes are coupled only radially and eventually cannot be decoupled. The system undergoes a symmetry transition from lattice symmetry breaking to rotational symmetry emergence during the modulation, ultimately forming novel moiré quasicrystal patterns. Each moiré quasicrystal pattern generally persists for less than 1 ms before transforming into the next symmetrical configuration. This symmetry transition heralds a new quantum dynamical phase transition.

The spatial symmetry of patterns evolves temporally when the frequencies or amplitudes are modulated. Increasing the modulation frequency prevents the formation of moiré quasicrystal patterns while preserves the original lattice characteristics. In contrast, a larger modulation amplitude generates more complex symmetrical structures. In particular, a $D_{12}$ pattern emerges at $a_{m}=100a_{0}$ and $\omega=200$ Hz, closely resembling the moiré quasicrystal in twisted bilayer graphene. Furthermore, the momentum-space patterns also exhibit high rotational symmetries, consistent with the real-space patterns at specific evolution times.

In summary, our results demonstrate that Floquet-engineered intralayer atomic interactions can generate novel quasicrystal patterns in bosonic systems without the need for additional quasiperiodic potentials. This phenomenon arises from the competition between the drive-induced disorder and the superfluidity of the system. The dynamical evolution of quasicrystal patterns provides a unique perspective for investigating quantum dynamical phase transitions in quasicrystalline systems. Furthermore, our numerical results reveal a sensitivity of the pattern formation to the modulation frequency, where a variation of only 1 Hz is sufficient to produce completely distinct moiré quasicrystal patterns. This work provides a new quantum platform for exploring quasicrystals and their symmetry characteristics in ultracold bosonic systems. The anyonic behaviors supported by the quasicrystals could offer a promising avenue for topological quantum computation.

Acknowledgments. Our work is supported by the Natural Science Foundation of Shanghai (Grant No. 23ZR1418700), the National Key Research and Development Program of China (Grant No. 2025YFF0515201), the Joint Funds of the National Natural Science Foundation of China (Grant No. U25D8014), the National Natural Science Foundation of China (Grant No. 11774093), the Program of Chongqing Natural Science Foundation (Grant No. CSTB2022NSCQ-MSX0585), and the Innovation Program of Shanghai Municipal Education Commission (Grant No. 202101070008E00099).