Interfacing differently oriented biaxial van der Waals crystals for negative refraction

: Negative refraction has a wide range of applications in diverse ﬁelds such as imaging, sensing, and waveguides and typically entails the fabrication of intricate meta-materials endowed with hyperbolic features. In contrast to artiﬁciallyengineeredhyperbolicmaterials,naturalvander Waals (vdW) materials are more accessible owing to their inherent strong in-plane covalent bonding and weak inter-layer interactions. However, most vdW materials manifest uniaxial crystal properties, which restrict their behavior solelytoout-of-planehyperbolicity.Thischaracteristicposes a considerable challenge to their seamless integration via planar fabrication techniques, unless a suitable pattern is employed. Recent advances have identiﬁed natural biaxial 𝛼 -MoO 3 as a promising vdW material capable of exhibiting in-plane hyperbolicity. In this study, we performed numerical simulations demonstrating that negative refraction could be achieved by interfacing differently oriented 𝛼 -MoO 3 slabs coated with tunable graphene on a gold substrate. Our comprehensive analysis yielded three notable outcomes: negative refraction, simultaneous positive and negative refractions, and diffractionless propagation. These outcomes could be operated in a broad range of frequencies and achieved at all angles to offer a superior platform for the ﬂexible manipulation of mid-infrared polaritons. Our ﬁndings provide valuable insights into the potential application of other two-dimensional vdW materials for advances in nanoscale super-resolution imaging, molecular sensing, and on-chip photonic integrated circuits.


Introduction
Negative refraction (NR) is a counterintuitive physical phenomenon where the refracted ray of light is deflected to the same side of the interface normal as the incident ray.NR can be used to realize a superlens capable of achieving a resolution smaller than the wavelength of light, which gives it a wide range of potential applications in super-resolution imaging, molecular sensing, and on-chip nanophotonic waveguides.Over the past two decades, substantial research efforts using diverse structures, materials, and mechanisms have been dedicated to realizing NR in different spectral regimes such as microwaves [1,2], visible light [3,4], and near-infrared (IR) light [5][6][7][8].Many studies [1][2][3][4][5][6][7][8] have employed subwavelength dielectrics and metallic elements to artificially construct hyperbolic media [9], which are characterized by permittivity tensors with one component along a principal axis having an opposite sign compared to the other two components.However, fabricating hyperbolic media using bulk materials is a major challenge [1][2][3][4][5][6][7][8].
A key advantage of vdW materials is their utilization in creating heterostructures that combine disparate 2D atomic layers without lattice mismatch.These vdW heterostructures [27][28][29][30][31][32] exhibit hybrid properties derived from the unique advantages of their constituent elements, which provide an additional degree of freedom for manipulating optical properties and designing novel photonic devices.For instance, monolayer graphene and hexagonal boron nitride (hBN) [33] were combined to obtain hybrid plasmon-phonon polaritons (HPPhPs) that combine the advantages of graphene surface plasmon polaritons (GSPPs) in graphene and hyperbolic phonon polaritons (HPhPs) in hBN to demonstrate lower ohmic losses and higher mode confinements than the constituent polaritons.hBN can be used to modify the optical properties of graphene to realize secondary Dirac points or new plasmonic states [34].Stacking graphene layers with hBN layers overcomes the absence of a bandgap in graphene [35].Graphene has also been combined with transition metal dichalcogenides (TMDs) [36] to modify the Fermi level and Schottky barrier height, which has resulted in unique properties that are not inherent to either graphene or TMDs.
However, uniaxial hBN and certain TMDs exhibit outof-plane hyperbolic and in-plane isotropic dispersions.This behavior arises from the permittivity tensors, which possess two negative components in the in-plane direction and a positive component in the out-of-plane direction [37,38].The out-of-plane hyperbolicity poses a major challenge for the design of planar photonic circuits.hBN and TMDs require nanostructured patterning to achieve inplane hyperbolicity [39][40][41], but this introduces surface roughness and defects that increase optical losses.Specifically, the reconfigurable topological features of hyperbolic polariton vortices [42] induced in pristine hBN flakes are also reported to tailor spin-orbit interactions.
Recently, a natural biaxial vdW semiconductor has been discovered that exhibits in-plane hyperbolic dispersion: -phase molybdenum trioxide (-MoO 3 ) [43][44][45][46].A twisted bilayer [47][48][49][50][51] and interface engineering [52,53] of the structure have been proposed to further manipulate the properties of -MoO 3 and realize exotic phenomena such as topological transitions, wavefront control, directional canalization, and diffractionless propagation.By covering the -MoO 3 slabs with monolayer graphene, their optical properties can be tuned by varying the charge carrier density of the graphene [54][55][56].Two independent groups recently demonstrated NR at the interface of vdW heterostructures.Hu et al. [57] utilized -MoO 3 slabs partially covered with graphene while Sternbach et al. [58] stacked -MoO 3 slabs with type-II hyperbolicity and h 11 BN slabs with type-I hyperbolicity.In this study, we propose a novel mechanism of interfacing differently oriented -MoO 3 slabs covered with a tunable graphene layer on an Au substrate to achieve broadband and all-angle NR.Notably, the use of an Au substrate provides significantly enhanced mode confinement of the polaritons compared to a SiO 2 substrate.Our proposed mechanism presents an innovative platform for effectively manipulating mid-IR polaritons for further advances in the development of polaritonic devices.

Proposed structure
Figure 1A and B illustrate the proposed structure and its cross-sectional view (x-z plane), respectively, which comprises air/graphene/-MoO 3 on an Au substrate.Here, t is the thickness of the -MoO 3 slab, [100] and [001] are the principal axes of the -MoO 3 crystal, and  is the angle between the x-axis and [100] direction.The employment of a metal substrate helps the proposed structure realize significantly higher mode confinement than the conventional HPhP, which is referred to as a hyperbolic image plasmon-phonon polariton (HIPPhP) [23].The HIPPhP originates from the hybridization of a hyperbolic polariton in the -MoO 3 slab and their mirror image with an inverted charge distribution in the perfect electric conductor (PEC) [21], which behaves as a mirror.At the mid-IR band, Au can be considered an approximate PEC.Therefore, most of the electric field is confined in the -MoO 3 slab, which significantly enhances the mode confinement due to the opposing charge distributions of the original and mirror polaritons.

Numerical simulations
We specifically investigated the reststrahlen band II (RB II) of -MoO 3 , which encompasses the frequency range of 816-972 cm −1 .Within this frequency range, the permittivity () components along the principal axes [100], [001], and [010] exhibit the characteristics of  x < 0,  y > 0, and  z > 0, respectively.The RB is a range of wavelengths where a material exhibits strong absorption due to the resonance between the incident radiation and lattice vibration.Furthermore, an HPhP occurs when the incident radiation couples strongly with the oscillations of the lattice vibrations in the -MoO 3 slab.Therefore, the RB can give rise to hyperbolic dispersion and the formation of hyperbolic polaritons in the -MoO 3 slab.Notably, the negative sign of Re ( x ) × Re ( y ) < 0 implies that the -MoO 3 has in-plane hyperbolicity in the x-y plane.
The complex permittivity of Au can be obtained from the literature [59].Here, the graphene sheet was modeled as an infinitely thin layer with a surface current density.The surface conductivity of graphene can be calculated by using the Kubo formula [60] (see Supplementary Material, §1).A numerical model was constructed in COMSOL Multiphysics with the mesh resolution set sufficiently fine and the computational domain sufficiently large to ensure simulation accuracy.The permittivity of -MoO 3 was obtained by employing the conventional Lorentz model [38,43,44,58] (see Supplementary Material, §2.1).The crystallographic directions [100], [001], and [010] of -MoO 3 coincided with the coordinates x, y, and z, respectively, under the condition  = 0 • .The real and imaginary components of the permittivity of -MoO 3 were plotted (see Supplementary Material, Figure S1).We placed an electric dipole source with polarization along the z-direction at a position 50 nm above the top graphene layer.The electric fields, represented as Re(E z ), were calculated at a probing height of 20 nm above the graphene layer.To obtain the corresponding IFCs in wavevector space (k x , k y ), we performed a Fourier transform of the Re(E z ) fields using spatial Fast Fourier Transform (FFT) with spatial sampling resolutions set to N x = 1000 and N y = 1000.To minimize undesired reflections, the computational window was enclosed by a perfect-matched-layer boundary condition with a thickness of 1 μm, effectively absorbing the transmitted waves.Simulations were performed to analyze the effects of modifying the -MoO 3 /Au structure with graphene, the effects of interfacing α-MoO 3 slabs at different orientations, and the effects of multiple interfaces.

Polariton modes with and without graphene
To gain insights into the mode characteristics, the dispersion properties of the -MoO 3 /Au structure at t = 150 nm were analyzed.The principal crystallographic direction [100] was aligned with the x-direction (i.e.,  = 0 • ).The results indicated an electromagnetic field with unprecedented confinement (see Supplementary Material, Figure S2A).The relationship between the propagation length (L p ) and t (see Figure S2B) was found to follow the conventional tradeoff between mode confinement and propagation loss, where an increase in t leads to looser mode confinement and longer L p .The figure of merit [FOM = Re(k x )/Im(k x )] was obtained to assess the benefit-cost ratio between mode confinement and propagation loss [22,46].High FOM values were obtained at the wavenumbers  = 880-920 cm −1 (see Figure S2C), which indicates a favorable tradeoff in this frequency range.The group velocity (V g ) (see Figure S2D) exhibited a similar trend to L p .The field distributions of the HIPPhP (see Supplementary Material, Figure S3) clearly demonstrated tighter confinement and a shorter wavelength as the thickness decreased.
Next, the dispersion properties of the graphene/ -MoO 3 /Au structure were obtained.Increasing the Fermi energy (E f ) was found to decrease Re(k x ) with a more pronounced dependence at higher frequencies (see Figure 2A).In contrast, L p , FOM, and V g all increased (see Figure 2B-D).
A higher FOM indicates a stronger mode confinement or lower mode loss.Interestingly, when E f is less than 0.1 eV, the case without graphene exhibits the strongest mode confinement but also the largest energy loss.Conversely, when E f exceeds 0.2 eV, the situation reverses, resulting in a smaller FOM (see Figure 2C).This phenomenon can be attributed to plasmons reaching the interband transition threshold at very low E f for specific wavevectors, leading to overdamped outcomes [61,62].Notably, the value of E f reaching the interband transition threshold decreases as the wavevector reduces.In contrast, plasmons transition to the intraband transition as E f increases.As a result, the HIPPhP modes have larger FOM when E f > 0.2 eV and smaller when E f < 0.1 eV while comparing to the case without graphene for the condition given in the proposed system.The Re(E z ) distributions at different values of E f (see Figure 2E-H) are also shown that both polariton wavelength and L p increased with increasing E f .These results offer valuable insights into the influence of graphene doping and its effects on the field propagation characteristics of the HIPPhP.

Effect of different orientations
We first analyzed a single -MoO 3 slab rotated counterclockwise around the z-axis by the angle , as shown in Figure 1A.Here,  corresponds to the orientation of the crystallographic direction [100] with respect to the x-axis.
We used this coordinate system to obtain the anisotropic permittivity tensor of the -MoO 3 slab (see Supplementary Material, §2.1) under the conditions of  = 910 cm −1 , t = 150 nm, and E f = 0.15 eV.The Re(E z ) field distributions varied with  and demonstrated directional canalization along the corresponding rotation angles (see Figure S4).These results provide compelling evidence of the in-plane hyperbolicity of the HIPPhP.To analyze the wavevector properties of the Re(E z ) fields, we performed a Fourier transform and plotted the resulting isofrequency contours (IFCs) in wavevector space (k x , k y ) (see Figure S5).The amplitude of the IFCs reflected the energy distribution in the corresponding wavevector space.Notably, a stronger amplitude indicates a higher concentration of energy.At  = 0 • , a significant portion of the energy was concentrated around k y /k 0 = −100 to 100, while k x /k 0 = ±20.The symmetric axis of the IFCs (i.e., k x /k 0 = 0) exhibited a tilt as  was varied, which implies that the directional propagation of energy can be flexibly tuned by adjusting .To validate the accuracy of the numerically calculated IFCs, we superimposed them with analytically calculated IFCs obtained by solving a four-layer waveguide system (air/graphene/MoO3/Au), as described in reference [58], for the first three cases (refer to Figure S5A-C).The inner elliptical and outer hyperbolic IFCs, indicated by green dashed lines, represent the fundamental and first-order HIPPhP modes.It is evident that there is excellent agreement between the results obtained from the numerical and analytical IFCs.To explore the tunability of the system, we obtained the Re(E z ) field and corresponding IFC for a bare -MoO 3 slab and graphene/-MoO 3 slab with different values of E f at  = 20 • (see Figure S6).We observed canalization of the field distributions at different doping levels.Notably, the opening angle of the hyperbola asymptote of the IFC increased with E f , and the bare -MoO 3 slab exhibited the smallest opening angle.This behavior aligns consistently with the field distribution patterns.The Re(E z ) field profiles showed the emergence of a hyperboloidshaped profile along the [100] direction and ellipse-shaped profile along the [001] direction.As E f was increased to 0.3 eV, the elliptical distribution became dominant.The variations in the field profiles were also reflected in the IFCs, which indicates a doping-driven topological transition from the manipulation of E f [58].Notably, increasing E f reduced the size of the elliptical profile in the IFCs.
Next, we investigated NR when two -MoO 3 slabs, both covered by monolayer graphene, were interfaced with different orientations.Figure 3A  ) are also represented by cyan dashed lines that is the same color used for its inner elliptical IFCs.The NR angle increased moderately with increasing t.This behavior can be explained by considering the formation of the resultant HIPPhP, which arises from the coupling of a GSPP with an elliptical IFC and HIPhP with a hyperbolic IFC.The contribution of the GSPP to the resultant mode of the -MoO 3 slab increases with decreasing t.However, decreasing t also lowers the threshold E f required for the topological transition to occur [56], which weakens the hyperbolicity of the resultant HIPPhP at a fixed E f .Thus, the largest NR angle was observed at t = 300 nm, as shown in Figure S7L.In addition to E f , t provides an additional parameter for controlling the occurrence of the topological transition.Therefore, while mode canalization for any given thickness of -MoO 3 can be adjusted through graphene doping, this adjustment is accompanied by a notable reduction in field strength along the canalization direction, while simultaneously resulting in a stronger field propagating in all directions.To study the dependence of NR on frequency, the Re(E z ) fields for different operating frequencies within the RB II of -MoO 3 along are shown in Figure S8A-F with the corresponding IFCs (see Figure S8G-L).Increasing  from 830 to 930 cm −1 cause the IFC patterns to transform gradually from a closedform ellipse to an open-form hyperbola.The broadband and the all-angle NRs highlight the potential of this system for realizing NR over a broad range of frequencies and incident angles.and [100] directions of the left, upper-right, and lower-right slabs are denoted as  1 ,  2 , and  3 , respectively.This configuration allowed the crystal orientations of each -MoO 3 slab to be manipulated and introduced additional degrees of freedom for controlling the optical properties of the system.

Effect of multiple interfaces
The angles  1 ,  2 , and  3 can be tuned to engineer the dispersion, wave propagation, and NR characteristics to realize various intriguing optical phenomena and functionalities.slabs.Figure 4D shows a top view of the Re(E z ) field distribution and highlights the power flows represented by the Poynting vectors P i , P pt , and P nt .In this structure, positive refraction and NR coexist, which allows it to function as a nanoscale beam splitter in photonic integrated circuits.The proposed design can serve as a functional Mach-Zehnder modulator [63] while interfacing seven -MoO 3 slabs with various orientations.This demonstrates the versatility and applicability of the proposed structure for different types of photonic devices.
To further control the polaritons, we interfaced two -MoO 3 slabs partially covered with graphene, as shown in Figure 5A.The left -MoO 3 slab had a crystal orientation of  1 = 0 • , as shown in Figure 5B, and was completely covered with monolayer graphene, which was electrically gated with E f 1 .The right slab had the crystal orientation  2 , and it was partially covered with graphene, which was electrically gated with E f 2 .The separation between the graphene layers was denoted as S. A z-polarized electric dipole was positioned at a distance of d = 1 μm from the left interface on the x-y plane.The other conditions were set to  = 910 cm −1 , t = 150 nm, S = 1 μm, E f 1 = 0.3 eV,  1 = 0 • , and E f 2 = 0.25 eV.
Figure 5D and E display the Re(E z ) field distributions for  2 = 0 • and −20 • , respectively.As a comparison, Figure 5C illustrates the results when only the left -MoO 3 slab was covered with graphene with both  1 and  2 set to 0 • [58].
The |E| fields and IFCs for the three cases are also shown.Figure 5C corresponds to the structure proposed by Hu et al. [58], which exhibited NR (Figure 5I) and planar focusing (Figure 5F).However, after passing through the focal point, the field experienced conventional diffraction.In contrast, partially covering the graphene on the right -MoO 3 slab with a suitable E f 2 = 0.25 eV overcame the diffraction of the light beam, and diffractionless propagation was achieved (Figure 5D and G).In addition, the limitation of a fixed focal point [58] was eliminated.This achievement opens up the possibility for diffractionless propagation in imaging systems and steering the diffractionless beam in wide angles to any desired direction by presetting the crystal orientation (Figure 5E and H).The IFCs for Figure 5D and E are also shown in Figure 5J and K, respectively.The parameters d, S, E f 1 , and E f 2 can be adjusted to realize different imaging capabilities.

Conclusions
In this study, we demonstrated the feasibility of interfacing differently oriented -MoO 3 slabs covered with tunable graphene on an Au substrate to realize NR.We provide compelling evidence for the presence of NR across various structures.By tailoring the orientations of the -MoO 3 slabs and modulating the Fermi energy of the graphene layer, three remarkable outcomes were achieved: broadband and all-angle NR, simultaneous positive refraction and NR, and diffractionless propagation for flexible manipulation of mid-infrared polaritons.The integration of multiple -MoO 3 slabs opens up exciting approaches for designing and implementing photonic devices such as beam splitters and Mach-Zehnder modulators.The significance of these findings extends beyond the materials and structures considered in this study as they are potentially applicable to other types of vdW materials.
have accepted responsibility for the entire content of this manuscript and approved its submission.Conflict of interest: Authors state no conflict of interest.

Figure 1 :
Figure 1: The structure with single -MoO 3 slab.(A) Illustration of the air/graphene/-MoO 3 /Au structure where [100], [001], and [010] are the principal axes of the biaxial -MoO 3 slab and  is the angle between the x-axis and [100] direction.(B) The x-z plane of the structure, denoted as (A), where t is the thickness of the -MoO 3 slab.

Figure 2 :
Figure 2: Wavenumber  of an -MoO 3 slab at t = 150 nm and covered with monolayer graphene at several values of E f .against (A) Re(k x ), (B) L p , (C) FOM, and (D) V g .E z field profiles of an -MoO 3 slab with t = 150 nm: (E) no graphene and covered by monolayer graphene at (F) E f = 0.1 eV, (G) E f = 0.3 eV, and (H) E f = 0.5 eV.

Figure 3 :
Figure 3D, E, and I-K.In the case of the left -MoO 3 slab with a fixed orientation angle of  1 = 20 • , we have indicated the inner elliptical and outer hyperbolic IFCs using green dashed lines.For the right -MoO 3 slab, we have indicated the elliptical IFCs with different angles of  2 using cyan dashed lines.Note that the hyperbolic IFCs associated with the first-order HIPPhP modes of the right -MoO 3 slab are not visible due to their significantly weaker amplitudes.The calculated IFCs for the right slab exhibited inclined orientations that varied with  2 .The trajectories of the incident (k i ) and refracted (k nt ) wavevectors are depicted with dashed green and red lines, respectively, while the incident (P i )

Figure
Figure 4A shows the extension of the proposed structure to multiple interfaces comprising three -MoO 3 slabs with different crystal orientations.The angles between the x-axis

Figure 4 :
Figure 4: The structure with three differently oriented -MoO 3 slabs.(A) Interfacing three -MoO 3 slabs with crystal orientations  1 ,  2 , and  3 between the x-axis and [100] directions of the individual slabs.(B) 3D view of the Re(E z ) field distribution.(C) IFC with the incident wavevector (k i ), two transmitted wavevectors (k pt and k nt ), and the incident (P i ), positive refractive (P pt ), and negative refractive (P nt ) poynting vectors, where the green, red, and cyan dashed lines indicate the analytical calculated IFCs of the left, upper-right, and lower-right -MoO 3 slabs, respectively.(D) Top view of the Re(E z ) field distribution.

Figure 5 :
Figure 5: The structure with two differently oriented -MoO 3 slabs partially covered with graphene sheet.(A) Interfacing two -MoO 3 slabs, in which the left slab is entirely covered with graphene (E f 1 ) but the right slab is partially covered with graphene (E f 2 ).(B) The same as (A), except peeling off the graphene layer to clearly show the orientations of -MoO 3 slabs.The two graphene layers have a separation distance of S = 1 μm, where  1 = 0 • (left slab) and  2 = 0 • or −20 • (right slab) indicate the angles between the x-axis and [100] directions of the individual slabs.Re(E z ) field distributions for various cases: (C)  1 = 0 • and  2 = 0 • [56] without the graphene layer (E f 2 ) covering the right -MoO 3 slab, (D)  1 = 0 • and  2 = 0 • , and (E)  1 = 0 • and  2 = −20 • .(F)-(H) Corresponding |E| fields.(I)-(K) Corresponding IFCs, where the green and cyan dashed lines indicate the analytical calculated IFCs of the left and right -MoO 3 slabs, respectively.

Figure
Figure 4C visualizes k i and the two transmitted wavevectors k pt and k nt , where k pt is the positive refraction in the IFCs derived from the three differently oriented -MoO 3 slabs.The analytical IFCs representing the fundamental modes of the left, upper-right, and lower-right α-MoO3 slabs are denoted by green, red, and cyan, respectively.Additionally, only the hyperbolic IFCs associated with the first-order mode of the left -MoO3 slab are observable.The corresponding P i , P pt , and P nt are also indicated, where P pt is the positive refractive Poynting vector.The three overlapping IFCs resulted from the three interfaces between the -MoO 3