Enhanced Confinement of Terahertz Surface Plasmon Polaritons in Bulk Dirac Semimetal-Insulator-Metal Waveguides

A subwavelength terahertz plasmonic waveguide based on bulk Dirac semimetal (BDS)-insulator-metal (BIM) structure is investigated, which indicates that there is an optimized frequency range with the better confinement as well as lower loss. A broadband mode confinement up to λ0/15 with a relatively low loss of 1.0 dB/λ0 can be achieved. We also show that two silicon ribbons introduced into the BIM waveguide can form a dynamically tunable filter tailoring terahertz surface plasmon polaritons in deep-subwavelength scale, which can be further exploited for the design of ultra-compact THz plasmonic devices with dynamical tunability. Our results may also provide potential applications in optical filtering.


Background
Terahertz (THz) wave has been extremely witnessed in the past few decades for its innovative applications, such as THz imaging, bio-chemical sensing, and communications [1][2][3]. To improve sensing sensitivity, imaging resolution, and integration level of THz devices, confining THz wave in deep subwavelength scale is urgently desired [4][5][6]. Surface plasmon polaritons (SPPs), the surface electromagnetic modes stimulated by the interaction between electron in conduction band of noble metal and photons in visible wavelengths, propagate along metal-insulator interface and enable the manipulation of light beyond the classical diffraction limit [7]. Sommerfeld-Zenneck modes, analog of SPPs in visible band, can be supported by metal in THz region. Metamaterials and other artificial structures, such as periodical patches, perforated plates, and brass tubes, have been proposed to tailor this loosely bound surface waves [8][9][10]. Unfortunately, poor confinement, high intrinsic loss, and passive tunability of this mode have severely hindered its practical applications.
Graphene plasmons, with relatively low loss, dynamical tunability, and extreme confinement to THz waves, hold promising applications in high-resolution, ultra-compact, and dynamical tunable devices. Duan et al. propose a broadband gate-tunable graphene heterostructure to coherently generate and control terahertz plasmons with dynamical tunability and higher efficiency. A robust difference frequency signal can be generated due to the tight confinement of graphene plasmon field [11]. Duan et al. firstly investigate discrete Talbot effect in dielectric graphene plasmonic waveguide arrays at THz wavelengths, which provides a novel platform for high-resolution self-image of THz waves in nanoscale [12]. Lin et al. propose an ultra-compact plasmon-induced transparency waveguide, which promises potential applications in slow light of THz waves [13,14]. Li et al. propose a series of functional optical filters and absorber based on 2D material plasmons, which demonstrate high-integration [15], low loss, and dynamical tunability [16][17][18]. From these works, we can convince that it is the extreme confinement of surface plasmons that make it possible to manipulate THz waves at deep subwavelength scale.
Recently, bulk Dirac semimetals (BDS), "3-D graphene," is brought into focus due to its ultrahigh carrier mobility up to 9 × 10 6 cm 2 V −1 s −1 , which is much higher than the best graphene of 2 × 10 5 cm 2 V −1 s −1 [19]. In general, the higher carrier mobility is, the lower intrinsic loss of plasmons would be. Furthermore, the dielectric functions of BDS can be actively tuned by changing its Fermi energy. The good news is that BDS, such as Na 3 Bi [19], Cd 3 As 2 [20], and AlCuFe quasicrystals [21], are easier to process and more stable compared with graphene, which is expected to be a new generation of plasmonic material after graphene. However, the mode confinement of SPPs in BDS-insulator interface is not optimistic. Our recent work has investigated the manipulation of the THz SPPs in double-layer BDS sheet waveguide, which indicates that the symmetric coupling mode has better confinement than the plasmonic waveguide mode in monolayer BDS film [22]. The mode index of the symmetric mode is 1.21 at 1.0 THz with the Fermi energy of BDS E F = 70 meV, which is still inadequate to meet the demand of manipulating THz wave in deep-subwavelength scale.
In this paper, we propose a deep-subwavelength BDS-insulator-metal (BIM) waveguide with enhanced confinement, relatively low loss, and desirable tunability. Dispersion relation, propagation loss, and filtering application of this highly confined mode are investigated. Interestingly, there is an optimized frequency range with an enhanced confinement as well as a reduced loss, which has rarely been reported in the traditional SPP mode in metal structure. A broadband mode confinement up to λ 0 /15 with a relatively low loss of 1.0 dB/λ 0 can be achieved. Different from previously studied BDS-based structure, the mode of this BIM waveguide can be efficiently transmitted through ultra-narrow slit with width smaller than λ 0 /2000. By taking two silicon ribbons as reflection mirrors, a dynamically tunable optical resonator has been achieved. The resonant frequency of the resonator can be dynamically tuned by varying the Fermi energy of BDS, which may find applications in THz switching and filtering.

Theory and Simulation
The proposed BIM plasmonic waveguide is schematically presented in Fig. 1(a), where the monolayer BDS film with thickness of 0.2 μm is placed at a gap width g away from the silver substrate separated by the dielectric spacer with permittivity ε r . The silver substrate in THz region can be treated as a perfect electric conductor (PEC) boundary. For the TM-polarized incident light, the plasmonic waveguide mode confined in the metal-insulator interface can propagate along the x direction with a wavevector k SPP and exponentially decay along the y direction into the free space. By combining proper boundary conditions, the wavevector k SPP of the BIM waveguide can be obtained from the following dispersion relation: [23].
where k 0 is the wavevector of the incident light. By solving Eq. (1), we can obtain the effective refractive index n eff = k SPP /k 0 = Re(n eff ) + iIm(n eff ) of the proposed plasmonic waveguide. For the highly confined plasmonic waveguide modes, the real part of effective refractive index Re(n eff ) roughly describe the mode confinement, while the imaginary part Im(n eff ) is directly proportion to the mode propagation loss: The larger Re(n eff ) is, the higher the confinement. When g is large enough such that tanh[g(k SPP 2 − ε r k 0 2 /ε 0 )]~1, Eq. (1) would be reduced to the dispersion relation which depicts plasmonic waveguide mode supported by a single layer of BDS alone. The complex conductivity of BDS is presented in methods Eqs (3)-(4).

Results and Discussion
First, we demonstrate the dependence of the mode confinement and propagation loss of the BIM waveguide on the BDS-metal gap width g and Fermi energy E F . By taking E F = 70 meV, we calculate the effective refractive indices of the SPP waveguide mode n eff for different values of g, Fig. 1 Schematic illustration of the BIM plasmonic waveguide: a monolayer BDS film is placed at a gap width g away from the silver substrate separated by a dielectric spacer with permittivity ε r . The TM-polarized plasmonic waveguide mode propagates along the x direction and decays along the y direction. Schematic depiction of the E x distribution is depicted by the red line where its real and imaginary parts, Re(n eff ) and Im(n eff ), are plotted in Fig. 2a, b, respectively. As depicted in Fig. 2a, the curves for g = 10 and 100 μm merge at frequencies higher than 0.05 THz, which suggests that the plasmonic waveguide modes are so tightly confined in the BDS-insulator interface that most of the SPP fields are distributed within the scale of 10 μm and the silver would not work at such a large gap width. While the mode confinement is dramatically enhanced after the gap width g is gradually reduced from 1 μm, the smaller g researches, and the stronger confinement can be obtained. The similar tendency can be observed in the dependence of propagation loss on the gap width g, as depicted in Fig. 2b. On the other hand, for a fixed gap width smaller than 1 μm, Re(n eff ) each initially shows a pronounced reduction to a minimum and then exhibits a gradually increasing behavior, while Im(n eff ) each decreases monotonically as the frequency increases. Thus, there is an optimized frequency region where the mode confinement is strongly enhanced while the propagation loss is gradually reduced. This characteristic has rarely been observed in traditional plasmonic waveguide modes at metal-insulator interface. Figure 2c, d depicts the dependence of the mode confinement and propagation loss on the Fermi energy E F of the BDS film, where the gap width g = 1 μm. Similar with the case of a monolayer and double-layer waveguide, the mode confinement and propagation loss continuously decrease with the increase of Fermi energy, which can be attributed to the enhanced metallicity and prolonged carrier relaxation time of BDS. For example, the confinement factor of plasmonic waveguide mode at 2.5 THz can be up to λ 0 /15, where is λ 0 the incident wavelength, with a relatively low loss of 1.0 dB/λ 0 when BDS-metal gap width is 10 nm and the Fermi energy is 70 meV. Therefore, relying on frameworks that have already been discussed above would increase mode confinement with a relatively low loss, which can be utilized to the design of integrated optical filters, buffers, and Mach-Zehnder interferometer.
To examine aforementioned analysis, we perform numerical calculations about the transmission intensity and field distribution of the proposed waveguide structure. The simulation setting is described in methods. Comparing with the monolayer BDS waveguide with same Fermi energy E F = 70 meV, the transmission intensity of the BIM waveguide at the frequency of 1.56 THz is 0.97 which is higher than that of the former, as shown in Fig. 3a, which suggests that the plasmonic waveguide mode in BIM structure suffers lower propagation loss. On the other hand, as indicated in Fig. 2a, the real part of effective refractive index of BIM at 1.56 THz Re(n eff ) = 2.45, which is much higher than that of the monolayer case of 1.002. To visualize this statement, the magnetic field Hz distributions of these modes are presented in Fig. 3b, c. It can be clearly found that the highly confined plasmonic mode in BIM waveguide shows shorter oscillation period than that of the monolayer BDS case. In addition, most of the plasmonic field are localized in such narrow slit~λ 0 /2000, which holds promising applications in near-field enhancement for nonlinear physics.
Among all the applications mentioned above, optical resonator is the essential element for tailoring THz plasmonic waveguide mode. As illustrated in Fig. 4a, two silicon (n Si = 3.4) [24] ribbons are embedded into the dielectric spacer to form the reflective mirrors, where the propagating plasmonic wave can be reflected back and forth at the silicon-air interface forming localized standing wave resonance in the BIM region between the two silicon ribbons. Only the incident frequency satisfies the resonance condition of the standing wave, the plasmonic waves can transmit to the output of the waveguide via coupling with the designed optical resonator. Figure 4a presents the transmission spectrum of the BIM waveguide with two silicon ribbons, where two transmission peaks with FWHM (full width at half maximum) values of 0.12 and 0.09 THz can be obviously found at the frequency of 1.56 and 2.22 THz, which demonstrates novel band-pass filtering effect at terahertz region. The magnetic field distributions (|H z | 2 ) of the transmission peaks are depicted in Fig. 4c, e, which implies that the BIM region sandwiched by two silicon ribbons can be regarded as a Fabry-Perot (FP) cavity. The first-and second-order resonance can be clearly found in the FP cavity. The incident plasmonic wave near the resonance frequency can be coupled into the FP cavity and then transmits through the BIM waveguide, which generates the transmission peak in the spectrum. While, for the non-resonant frequency region, the standing wave cannot be formed and thus the incident waves are prohibited in the left port of BIM waveguide, as shown in Fig. 4d. Moreover, combined with the dispersion relation of BIM waveguide, the transmission intensity can be analytically calculated by coupled mode theory (CMT) [17]: where ω 0 is the resonant frequency of the FP cavity, respectively. Here, κ w = ω 0 /(2Q w ) and κ i = ω 0 /(2Q i ) are decay rates related to the waveguide coupling loss and intrinsic loss of the FP cavity, respectively. The total and intrinsic loss quality factor can be estimated by Q t = ω 0 / FWHM, and Q oi = − Re(n eff )/(2Im(n eff )), respectively. Then, the waveguide coupling loss quality factor can be obtained by subtracting the intrinsic loss from the total loss, namely, Q ei = Q oi Q ti /(Q oi − Q ti ) [17]. The analytical results based on CMT hold good agreement with the numerical simulations, as depicted in Fig. 4b. Figure 5 shows the dependence of resonant frequency on the cavity length L, where g = 1 μm, d = 5 μm, and E F = 70 meV. The transmission peak tends to red-shift with the increase of L, as presented in Fig. 5a, which can be further described by the standing wave resonance condition 2k SPP (ω r )L + θ = 2mπ (m = 1, 2, 3, ...), where θ is the reflective phase shift from the silicon-air interface and k SPP (ω r ) is the wavevector of BIM waveguide at resonant frequency. As shown in Fig. 5b, the resonant frequencies of the first and second modes indeed exhibit a red shift with the increase of L. According to Eq. (1), the mode confinement is affected by gap width g which therefore have impact on the resonant frequency. Figure 6a presents the transmission spectra for different g, where L = 120 μm and E F = 70 meV. With the increase of g, the resonant peak in same order exhibits a blue shift. This phenomenon can be attributed to the dramatic decrease of Re(n eff ) as shown in Fig. 6c. The tuning of Fermi energy of BDS can be realized by alkaline surface doping in experiment. Figure 6b presents the transmission spectra for different Fermi energy, where the other parameters are same as Fig. 4b. As the Fermi energy increases, the transmission peak presents a blue shift, which can also be involved in the standing wave resonance picture. For a fixed length L, the FP cavity supports the resonance with defined SPP wavelength λ SPP = λ 0 /Re(n eff ), where λ 0 is the incident wavelength. As shown in Fig. 6c, Re(n eff ) is reduced with the increase of Fermi energy. As a result, the incident wavelength λ 0 should be decreased as well to keep λ SPP as a constant. That is the reason why the transmission peak tends to blue-shift with the increase of Fermi energy. Meanwhile, the bandwidth of the transmission peak is narrowed, which can be attributed to the decrease of Im(n eff ), i.e., the propagation loss of the plasmonic waveguide mode in BIM waveguide.

Conclusions
In summary, we have demonstrated the highly confined terahertz plasmonic mode supported by a BIM waveguide. The mode confinement and loss characteristics have been discussed with the variations of BDS-metal separation and Fermi energy, which indicates that there is an optimized frequency range with enhanced mode confinement as well as reduced propagation loss, which has rarely been reported in traditional SPP mode in metal structure. Different from previously studied BDS-based structure, the mode of this BIM waveguide can be efficiently Here, g = 1 μm, d = 5 μm, and E F = 70 meV supported in very narrow slit with width smaller than λ 0 / 2000. By taking two silicon ribbons as reflective mirrors, a dynamically tunable band-pass filter has been achieved, where the resonant frequency can be actively controlled by adjusting the Fermi energy of BDS film without re-optimization of its structural parameters.

Methods
Numerical results are obtained by using 2D finite-difference time-domain (FDTD) method, where the perfectly matched layers are set to absorb the scattering light in the x and y directions. The mesh size of the BDS film are set as dx × dy = 1 μm × 0.02 μm to achieve good convergence.
The frequency-dependent conductivity of BDS is described by Kubo formula with random phase approximation [12,25]. Imσ where G(E) = n(−E) − n(E) and n(E) is the Fermi-Dirac distribution function, E F is the Fermi energy of BDS, k F = E F /ћv F is its Fermi momentum, and v F = 10 6 m/s is the Fermi velocity. ε = E/E F , Ω = ћω/E F + iћτ −1 /E F , where ћτ −1 = v F /(k F μ) is the electron scattering rate which shows strong dependence on the carrier mobility μ. ε c = E c /E F (E c is the cutoff energy beyond which the Dirac spectrum is no longer linear), and t is the quantum degeneracy factor. Taking AlCuFe as an example, the fitting parameters in our calculations are set as follows: t = 40, ε c = 3, μ = 3 × 10 4 cm 2 V −1 s −1 and E F = 70 meV. No human participants, data, or tissue or animals are involved in this research.