The well-known Bloch oscillation (BO) was first predicted in 1929 for a periodic potential excited by an external homogeneous static electric field^[1]. As a classical prediction, BO remained an elusive phenomenon for over six decades due to rapid decoherence in electronic systems. This fundamental problem was finally resolved with the groundbreaking experimental observation of BOs in semiconductor superlattices by Feldmann et al.^[2]. Since then, superlattices in waveguide arrays have emerged as an ideal platform for observing both linear and nonlinear optical Bloch oscillations^[3]. In recent years, abundant novel phenomena related to imaging polarimetry have been observed in various systems, including optical analog processing^[4], metalens imaging polarimetry^[5], and angular spectrum differentiation^[6]. To date, BOs have been extensively investigated in a wide range of optical sub-micron and nanoscale optical systems, such as coupled optical waveguides^[7-12], period-doubling topological systems^[13], Floquet-Bloch oscillations^[14], acoustic superlattices^[15-16] and quantum particles in periodic media under the influence of a constant external force^[17]. Efforts are also underway to fabricate a uniform thickness metamaterial structures, such as multilayer metamaterials^[18] and nanoscale optical devices^[19].

Currently, significant resources are being allocated to the development of the BOs in MMWAs due to their potential for carrying optical signals via evanescent electromagnetic waves, which propagate along the interface between metal and dielectric layers in various dispersive systems^[20-24]. For instance, BOs have been realized in three-dimensional waveguide arrays by introducing a refractive index gradient in dielectric layers^[25] and in chirped metal-dielectric structures^[26].

The recognized principle underlying optical BOs is the coupling of evanescent waves in adjacent waveguides, which results in the formation of a collective supermode. Due to the alternating total internal reflection and the transverse Bragg reflection at the boundaries, this collective supermode is confined to a specific region of the waveguide array and propagates along the longitudinal direction. To the best of our knowledge, for the general case, the realization of BOs in waveguide arrays requires either a transverse gradient in the effective refractive index or a monotonic variation in the adjacent layers. This condition is analogous to the constant electric field required for BOs of electrons in a periodic potential. Specifically, the gradient in the effective refractive index induces a variation in the propagation constant, which follows a linear relationship described by $ \beta = {k_0}n $, where $ {k_0} $ represent the wave vector in the volume. Hence, it is evident that both the refractive index and the thickness of the layers play crucial roles in realizing BOs in periodic waveguide arrays.

Quasiperiodic systems represent an intermediate state between periodic and disordered systems, exhibiting long-range order without translational symmetry. Fibonacci is a common and representative quasiperiodic system generated through iterative rules^[27-28]. Unlike the periodic system with a gradient modulation, the optical thickness of quasiperiodic systems is non-monotonic due to the absence of translational symmetry. In this paper, we aim to investigate Bloch-like Oscillations (BLOs) in a planar Fibonacci quasiperiodic MMWA. The fabrication is achieved using an ultrahigh vacuum electron beam evaporation system on a chemically cleaned glass substrate.

In this paper, we propose a planar Fibonacci waveguide array composed of silver and silicon dioxide to theoretically demonstrate spatial BLOs. Periodic oscillations are observed in the ninth- and tenth-generation Fibonacci MMWA using the exact field distributions of the incident Gaussian pulse. In particular, we identify three distinct types of collective oscillation modes when a Gaussian pulse of the same wavelength is incident on different positions of the ninth-generation structure. The oscillation period increases with the redshift of the incident wavelength for both the ninth- and tenth-generation Fibonacci structures. Our work provides a flexible platform for implementing various interesting phenomena based on metamaterials and has potential applications in manipulating the emergence of BLOs.

This study investigates plasmonic BLOs in MMWAs using theoretical analysis and numerical simulation techniques. The simulation models are implemented using the Finite Difference Time Domain (FDTD) method, which provided explicit and dynamic solutions. During the simulation, the Perfect Electric Conductor (PEC) layers are applied to both sides along the X-direction, and the Perfect Magnetic Conductor (PMC) layer are applied to both sides along the Z-direction, serving as numerical absorbing boundary conditions.

The quasiperiodic waveguide arrays are constructed using a Fibonacci structure composed of two types of layers: A and B. Type A consists of dielectric silicon dioxide layers, while Type B consists of metamaterial silver layers. According to the recursion formula, the ninth-generation Fibonacci quasiperiodic MMWA is configured as follows:

$ {{\text{S}}_9} =\text{ABAABABAABAABABAABABAABAABABAABAABABAABABAABAABABAABABA} $ ()

The graphical representation of the ninth-generation MMWA is shown in Fig. 1. The thicknesses of layer A and layer B are $ {d_{\mathrm{A}}} = 28 $ nm and $ {d_{\mathrm{B}}} = 20 $ nm respectively. To simplify the analytical discussion, we ignore the dielectric dispersion and losses in the metal layer. We choose the dielectric permittivity $ {\varepsilon _{{\text{Si}}{{\text{O}}_{\text{2}}}}} = 2.126\;6 $. The real part of the permittivity for the silver in the Drude model is $ {\varepsilon _{\mathrm{S}}} $ ($ {\varepsilon _{\mathrm{S}}} \approx 1 - \varpi _{\mathrm{p}}^2/{\varpi ^2} $) in the high frequency limit. The plasma frequency $ {\varpi _{\mathrm{p}}} $ is that of the reference^[29]. The transport property is effectively defined by the real part, while the imaginary part corresponds to the attenuation. Hence, only the effect of the real part is considered in the FDTD simulation.

The amplitude profile of a Transverse Magnetic (TM) polarized Gaussian pulse, which is normal incident along the Y-direction, has a specific form $ {H_{\textit{z}}} = \exp ( - {(x - {x_0}^2)}/{2{\delta ^2}}) $, where the incident position and Full Width at Half Maximum (FWHM) are $ {x_0} $ and $ \delta {\text{ = }}{200_{}} $ nm, respectively. For convenience, the coordinate system’s origin is set at $ {x_0} = 0.686 $ μm, which is the center of the twenty-eighth slit dielectric layer shown in Fig. 1 (color online). The FDTD program Optiwave OptiFDTD is used to perform the simulation^[30].

The trajectory of the Gaussian pulse is simulated in the quasiperiodic MMWA. It should be mentioned that the wave dispersion relation is calculated by the transfer-matrix method, since the entire Fibonacci structure can be regarded as a unit in the periodic structure^[31], even the fact that the structure is the absence of translational symmetry. Thus, the pseudo-dispersion relation can be written as:

$ \cos[{k_x}({d_{\mathrm{A}}} + {d_{\mathrm{B}}})] = \cos ({k_1}{d_{\mathrm{A}}})\cos ({k_2}{d_{\mathrm{B}}}) - ({k_1^2\varepsilon _{\mathrm{m}}^2 + k_2^2\varepsilon _{\mathrm{a}}^2}/{2{k_1}{k_2}{\varepsilon _{\mathrm{a}}}{\varepsilon _{\mathrm{m}}}})\sin ({k_1}{d_{\mathrm{A}}})\sin ({k_2}{d_{\mathrm{B}}}) \quad, $ (1)

with $ {k_{1,2}} = \sqrt {({\varepsilon _{\mathrm{a},{\mathrm{m}}}}k_0^2 - k_{\textit{z}}^2)} $, where $ {k_0} $ is the free space wave vector, and $ {k_x} $ and $ {k_{\textit{z}}} $ are the Bloch wave vector along the X-direction and the propagation constant along the Z-axis. In general, a simple way to obtain photonic BOs is to introduce a gradient change of the optical thickness between two adjacent waveguide layers^[32]. Therefore, we first introduce a permittivity gradient $ \alpha $ for dielectric layer, where $ \alpha \in (0,0.15) $. In general, there are 34 dielectric layers in the ninth Fibonacci quasiperiodic MMWA, and the permittivity of the Nth layer is $ {\varepsilon _{\mathrm{a}}} = {\varepsilon _0} +\alpha (N - 1) $, where $ N = 1,2 \cdot \cdot \cdot 34 $ and $ {\varepsilon _0}{\text{ = }}1.0 $. Then, the photonic band diagram of the graded MMWA is shown in the Fig. 1(b). An inclined band gap appears for a Gaussian pulse with wavelength ranging from 100 nm to 600 nm when the incident pulse propagates along the X-direction, which is well matched in Bragg mirrors, to some extent^[32]. The lower index layers in the ninth Fibonacci MMWA without gradient permittivity are the equivalent to the barriers, and the lower index layers play the role of quantum wells in conventional semiconductor superlattices. Thus, the adjacent microcavities, which are the layers in the structure, are coupled to each other. And the coupling leads to the photonic Bloch-like oscillations, which are a collected supermode of the localized states with non-commensurate.

The electromagnetic evolution of the pulse in the waveguide array can be simulated as shown in Fig. 2 (color online).

The results regarding optical BLOs are shown in Fig. 2(a), where a Gaussian pulse with FWHM $ \delta {\text{ = }}200 $ and wavelength $ \lambda {\text{ = }}460 $ nm is incident on the right side of the waveguide arrays at the position $ {x_0} = 0 $. The simulated contour of the magnetic field $ \left| {{H_y}} \right| $ for a Gaussian pulse with $ \lambda = 460 $ nm is shown in Fig. 2(a). It can be observed that an excited oscillating mode spreads over several waveguides and undergoes a discrete diffraction. Finally, the oscillation periods plotted in Fig. 2(b) are about 3.54 μm and 4.19 μm, respectively, while the corresponding oscillation amplitudes have approximately the same value of 0.5 μm. In view of the fact that, quasiperiodic incommensurate MMWA deals with inhomogeneous refractive index, that the maximum value of the dielectric constant is 2.13 in the region below the position $ {x_0} $, while the permittivity of silver is approximately −11.54. We associate the localization of the electromagnetic field with the layers adjacent to the quasiperiodic structure, which have similarities to the electronic field influenced by a non-uniform force in a non-periodic potential.

When the Gaussian pulse impinges on the positions $ x = 0.34 $, $ x = 0 $ and $ x = - 0.34 $, respectively, the distribution of magnetic field intensity is shown in Fig. 3 (color online).

By comparing the Fig. 3(a)−3(c), we can identify some remarkable features. Firstly, there are three different types of collective oscillation modes in which the electromagnetic field is initially inclined to oscillate in the positive and negative directions along the X-axis. The oscillating direction depends on the incident position. For x₀ = 0.314, the excited mode oscillates along the negative direction of the X-axis. Conversely, for $ x = 0 $ and $ x = - 0.34 $, the excited mode oscillates along the positive direction of the X-axis and some of the electromagnetic field also spills out in the reverse direction. Secondly, there are slight distinctions between the periods of each oscillating mode corresponding to 5.92 μm, 6.04 μm and 5.85 μm from Fig. 3(a) to Fig. 3(c). Lastly, the pattern of the magnetic field distribution differs significantly, with a noticeable increase in the field energy leakage around $ x = - 0.5 $ in Fig. 3(c) due to tunneling between the neighboring layers.

Fig. 4 shows the simulated results of the magnetic field evolution within the ninth-generation Fibonacci quasiperiodic waveguide arrays for the Gaussian pulse impinging on the fixed position $ x = 0.21 $. The incident wavelengths of the Gaussian pulse are 488 nm, 514 nm and 532 nm, respectively. According to the above analyses, the coordinates of the trough for the oscillation shift gradually from $ \textit{z} = 3.85 $ to $ \textit{z} = 4.92 $, while the distribution of the magnetic field is mainly confined to the area between $ x = 0.48 $ and $ x = - 0.22 $.

Additionally, in another general scenario, several different Gaussian pulses are selected to impinge on the tenth generation Fibonacci MEWA. The field distributions with the same pattern which are the Bloch-like oscillations are shown in Figs. 5(a)−(c). The input port of the Gaussian pulse is $ x = 0.34 $. From Figs. 5(a)−5(c), it can be observed that all the oscillations initially start with along the positive direction of the X-axis. The periods of oscillation for the input wavelengths of $ \lambda $= 488 nm, 514 nm, 532 nm are 4.89 μm, 5.56 μm, 5.89 μm, respectively. This observation also confirms a redshift, consistent with the one found in the ninth-generation Fibonacci MEWA.

In summary, we have investigated spatial BLOs in the ninth- and tenth-generation Fibonacci planar quasiperiodic waveguide arrays using FDTD simulations. Three distinct collective oscillation modes were observed in the ninth-generation quasiperiodic structure based on the evolution of the electromagnetic field. Additionally, we found that the oscillation period increases with the redshift of the incident wavelength for both Fibonacci structures. The non-monotonic optical thickness also influences BLOs, which are distinctly different from those observed in periodic structures. As a result, the findings presented in this paper provide a flexible platform for implementing a variety of interesting phenomena related to metamaterials and may have potential applications in manipulating the emergence of BLOs.