Metasurface-assisted orbital angular momentum carrying Bessel-Gaussian Laser: proposal and simulation

Bessel-Gaussian beams have distinct properties of suppressed diffraction divergence and self-reconstruction. In this paper, we propose and simulate metasurface-assisted orbital angular momentum (OAM) carrying Bessel-Gaussian laser. The laser can be regarded as a Fabry-Perot cavity formed by one partially transparent output plane mirror and the other metasurface-based reflector mirror. The gain medium of Nd:YVO4 enables the lasing wavelength at 1064 nm with a 808 nm laser serving as the pump. The sub-wavelength structure of metasurface facilitates flexible spatial light manipulation. The compact metasurface-based reflector provides combined phase functions of an axicon and a spherical mirror. By appropriately selecting the size of output mirror and inserting mode-selection element in the laser cavity, different orders of OAM-carrying Bessel-Gaussian lasing modes are achievable. The lasing Bessel-Gaussian0, Bessel-Gaussian01+, Bessel-Gaussian02+ and Bessel-Gaussian03+ modes have high fidelities of ~0.889, ~0.889, ~0.881 and ~0.879, respectively. The metasurface fabrication tolerance and the dependence of threshold power and output lasing power on the length of gain medium, beam radius of pump and transmittance of output mirror are also discussed. The obtained results show successful implementation of metasurface-assisted OAM-carrying Bessel-Gaussian laser with favorable performance. The metasurface-assisted OAM-carrying Bessel-Gaussian laser may find wide OAM-enabled communication and non-communication applications.

the thin array of sub-wavelength structured resonators [32][33][34] . Compared to conventional focusing elements such as optical lens and planar focusing devices demonstrated at optical frequencies [35][36][37][38][39] , metasurface has compact structure and flexible design that provide the capability to tailor the phase of the transmitted light from 0 to 2π to realize full control of optical wavefront. Metasurface facilitates volume reduction of optical elements and flexible spatial structure manipulation of light beams. Most of the metasurfaces can be divided into two categories, i.e. plasmonic metasurfaces based on plasmon resonances [40][41][42] and dielectric metasurfaces based on electric/magnetic dipole (Mie resonances) 43,44 . Dielectric metasurface benefits the alleviation of loss compared to plasmonic metasurface. Beyond all-dielectric elliptical metasurface array 44 , hybrid dielectric/metal metasurface structure adding metal at the bottom as the reflection layer may further improve the efficiency and bring superior performance. In addition to the passive method of Bessel-Gaussian beam generation through axicon-enabled light beam conversion, an active method of Bessel-Gaussian beam laser is of great interest [45][46][47][48][49] . A digital laser with an intra-cavity digitally addressed holographic mirror was proposed and demonstrated with impressive performance for on-demand laser modes, i.e. customizing the output beam shape 49 . Additionally, in order to meet the requirement of compact Bessel-Gaussian beam laser, metasurface-based reflector can be employed to replace relatively complicated volume structure components. In this scenario, one would expect to see a metasurface-assisted OAM-carrying Bessel-Gaussian laser.
In this paper, a Bessel-Gaussian laser structure using metasurface as the resonator reflector is designed. The gain material of Nd:YVO 4 enables the lasing wavelength at 1064 nm. Taking advantages of Fox-Li algorithm 50 , the generation of Bessel-Gaussian modes carrying different OAM values and the laser properties are comprehensively analyzed. The metasurface fabrication tolerance and the laser performance dependence are also discussed. The simulation results show favorable performance of the presented metasurface-assisted OAM-carrying Bessel-Gaussian laser.

Results
Concept and principle. The schematic structure of the proposed Bessel-Gaussian laser is illustrated in Fig. 1. The laser can be regarded as a Fabry-Perot (FP) cavity formed by two mirrors (mirror 1, mirror 2). In order to realize a Bessel-Gaussian laser, mirror 1 can be a partially transparent plane mirror, while mirror 2 is usually a combined structure of an axicon and a spherical mirror (e.g. concave mirror) 51 . The combined functions of an axicon and a spherical mirror can be achieved by a metasurface (sub-wavelength structured plane mirror) benefiting from its flexible spatial light manipulation. Hence, here mirror 2 is replaced by a metasurface reflector to construct the resonator configuration. Nd:YVO 4 provides the cavity gain. A 808 nm laser serves as the pump. The top inset of Fig. 1 displays the schematic diagram of the designed dielectric metasurface (rectangular meta-reflection array), which is structured silicon (Si) antenna array on silica (SiO 2 ) with an Ag bottom layer enhancing the efficiency. The right inset of Fig. 1 depicts the desired typical phase distribution (combined functions of an axicon and a spherical mirror) of mirror 2, which will be replaced by the metasurface structure in the designed metasurface-assisted OAM-carrying Bessel-Gaussian laser.
Besides the metasurface-based reflector, the mode-selection element (MSE) 48 is necessary to be placed in the resonator for direct generation of Bessel-Gaussian modes with controllable handedness of the helical phase. In general, Bessel-Gaussian modes with opposite topological charges (right-handed and left-handed helical phase trajectories) have the same spatial intensity distribution. That is, except for the zero-order mode, the Bessel-Gaussian modes with opposite topological charges always coexist in a laser resonator, leading to petal-shaped intensity distributions [51][52][53][54] . The principle of MSE is to break such degeneracy by exploiting the fact that standing-wave intensity distributions for the Bessel-Gaussian modes with opposite topological charges inside the laser resonator are different 48 . Inside the laser resonator, both Bessel-Gaussian 01 + or Bessel-Gaussian 01 − modes have a two-lobe standing-wave intensity distribution rotating in opposite directions by 2π for every wavelength propagation. By appropriately putting structures with nanoscale thickness, a differential loss can be introduced between Bessel-Gaussian 01 + and Bessel-Gaussian 01 − modes, facilitating selective lasing of Bessel-Gaussian 01 + or Bessel-Gaussian 01 − mode. The detailed operation principle of the MSE can be understood in Fig. 2. For selective lasing of Bessel-Gaussian 01 + or Bessel-Gaussian 01 − mode, as illustrated in Fig. 2a, two nanoscale thickness wires (white line) with an angle of 45° are put at the center of the mode. The first wire determines the available standing-wave intensity distribution of Bessel-Gaussian 01 + and Bessel-Gaussian 01 − modes at the position of the wire in the laser cavity. The initial condition is set by the first wire. The two-lobe standing-wave intensity distribution of Bessel-Gaussian 01 + and Bessel-Gaussian 01 − modes rotates anti-clockwise by an angle of −(n + 1/4)π and clockwise by an angle of (n + 1/4)π after propagation of (n/2 + 1/8)λ, where the second wire (45° with respect to the first wire) is positioned. It is noted that the second wire has almost no touch with Bessel-Gaussian 01 − mode but full contact with Bessel-Gaussian 01 + mode, resulting in minimum and maximum losses for Bessel-Gaussian 01 + and Bessel-Gaussian 01 − modes, respectively. Hence, the lasing of the Bessel-Gaussian 01 + mode with relatively large induced loss by the second wire is suppressed. The opposite situation is achievable (Bessel-Gaussian 01 + lasing and Bessel-Gaussian 01 − suppression) simply by changing the distance between the two wires to be (n/2 + 3/8)λ, as shown in Fig. 2b. Alternatively, one can also obtain the opposite situation by keeping the distance between the two wires to be (n/2 + 1/8)λ while rotating the second wire with an angle of −45° with respect to the first wire, as shown in Fig. 2c. Remarkably, selective lasing of Bessel-Gaussian 01 + or Bessel-Gaussian 01 − mode is also achievable using a different group of angle and distance between the wires. Similarly, selective lasing of Bessel-Gaussian 03 + or Bessel-Gaussian 03 − mode can be realized assisted by MSE with proper angle and distance for the two sets of nanoscale thickness wires, as shown in Fig. 2d,e. Different from the situation for Bessel-Gaussian 01 + and Bessel-Gaussian 01 − modes, here in each position multiple wires are employed instead of a single wire, which is determined by the number of lobes of the standing-wave intensity distribution of high-order Bessel-Gaussian modes (e.g. three crossed wires for six-lobe standing-wave intensity distribution of Bessel-Gaussian 03 + and Bessel-Gaussian 03 − modes). Remarkably, the output beam of a laser is determined by resonance geometric parameters such as resonance cavity length, radius curvature and dimensions of reflector. To describe Bessel-Gaussian modes of the designed laser resonator configuration, the transverse electrical field of a Bessel-Gaussian mode at a distance z from the waist in the paraxial approximation can be expressed as 51 where w(z) is the radius of Bessel-Gaussian mode at distance z from its waist. w 0 is the waist radius. The function Φ(z) = k(1 − γ 2 /2)z-arctan(z/z R ) describes the axial phase of the Bessel-Gaussian mode. k = 2π/λ is the wavenumber with λ the wavelength. ρ is the radial distance from the mode center. γ is the conicity angle of the Bessel beam. z R is the Rayleigh range. R g = z + z 2 R /z is the wavefront curvature of Bessel-Gaussian mode. m is the number of topological charge. ϕ is the azimuthal angle. To properly design the metasurface-based reflector, the phase and amplitude distributions of the lasing mode on two mirrors should be considered. The waist of the Bessel-Gaussian mode is set on mirror 1 (a plane mirror) and the radial phase on mirror 1 can be 0 or π. The mirror 2, having a distance of L away from the mirror 1, provides a phase-only reflection function written by It is known that the phase function of a spherical mirror with a radius of R is −kρ 2 /R and that of an axicon is −kαρ. Hence, the phase function in Eq. (2) can be regarded as a combination of a spherical mirror and an axicon with properly selected R and α. Travelling twice through the axicon and replacing z R with kw 2 0 /2, R and α can be expressed as: Eqs (3,4) give the mirror radius R and the axicon parameter α required to enable the lasing of Bessel-Gaussian mode with conicity angle of γ and waist radius of w 0 . Using Eqs (2-4) one can obtain the desired phase function for metasurface design. In addition, Bessel-Gaussian mode has properties of suppressed diffraction divergence and self-reconstruction. In the designed metasurface-assisted OAM-carrying Bessel-Gaussian laser, when maximizing the gain in the cavity volume, the maximum propagation distance of a Bessel-Gaussian mode Z max is equal to the length of the cavity L, i.e. Z max = L = R/(2tanγ), where R is the radius of the metasurface-based reflector. For the sake of manifesting the physical process of mode lasing and propagation in the laser resonator, Fox-Li algorithm 50 (see Methods) can be applied to analyze the evolution of transverse electrical field, as shown in Fig. 3. To evaluate the quality of the generated Bessel-Gaussian modes, the mode purity or fidelity is introduced and defined as follows   Fig. 4a is 700 nm. The incident light is x-polarized, which can be decomposed into two perpendicular components correspond to the long and short axis of the unit rectangular Si resonator, respectively. The reflection amplitudes in both components are close to each other while the relative phase retardation is around π where the linear polarization conversion occurs, leading to a y-polarized reflected light. The dashed line in the inset of Fig. 4a indicates the symmetry axis of the metasurface unit. The angle between the polarization of incident light and symmetry axis is 45°. In addition, the gain medium is Nd:YVO 4 and the lasing wavelength is 1064 nm. For a typical metasurface unit shown in the inset of Fig. 4a, by varying the geometric parameters of the rectangular unit, i.e. length l x and width l y , the amplitude and phase distributions of the cross-polarized scattered light are calculated and plotted in Fig. 4b,c, respectively. One can clearly see from Fig. 4b,c that the normalized amplitude changes from 0 to 1 and the phase shift varies from −π to π. We choose eight kinds of units with different geometric dimensions as shown in Fig. 4a, providing a full phase coverage of 2π and a nearly constant amplitude. The employed metasurface reflector plays combined roles of an axicon and a spherical mirror. Using proper arrangement of eight selected units in Fig. 4a with a full phase coverage of 2π to construct a metasurface reflector can approximately provide the desired phase function in Eq. (2) for the lasing of Bessel-Gaussian modes. When we use the perfect phase distribution in the right inset of Fig. 1 in the simulations regardless of the real metasurface structure, the obtained results in Fig. 5 show almost ideal multi-ring intensity distribution and helical phase distribution of OAM-carrying Bessel-Gaussian modes.
The desired perfect phase distribution of mirror 2 (combined functions of an axicon and a spherical mirror) in the right inset of Fig. 1 shows smooth phase change (i.e. continuous phase pattern). However, such continuous phase pattern is actually not achievable using metasurface structure which can be regarded as multiple pixels with discrete phase modulations. In view of the discrete characteristic of metasurface structure, it is necessary to replace the desired continuous phase pattern with a discrete one. In the simulations, we first discretize the continuous phase pattern. Figure 6b depicts the discrete phase pattern which discretizes the continuous phase pattern shown in Fig. 6a into eight values. Note that the discrete eight values are linked to eight selected metasurface units which provide a full phase coverage of 2π and a nearly constant amplitude. After the discretization, we then replace the discrete phase pattern by corresponding metasurface units, as shown in Fig. 6c-e. As a consequence, the designed metasurface structure is formed to enable the metasurface-assisted OAM carrying Bessel-Gaussian laser.
Following the discretization process, the real metasurface structure and discontinuities in phase and amplitude responses are then taken into consideration in the simulations of laser with the Fox-Li algorithm. The output mirror in the laser influences the oscillation of Bessel-Gaussian modes of different order. It can be briefly explained with the fact that high-order Bessel-Gaussian modes suffer from relatively high loss. High-order Bessel-Gaussian modes might not be excited with limited size of the output mirror. To generate desired Bessel-Gaussian modes with different order, we adjust the radius of the output mirror as well as the MSE. Figure 7 shows simulation  results with real metasurface structure. In the simulations, w 0 of generated Bessel-Gaussian 0 , Bessel-Gaussian +1 , Bessel-Gaussian +2 and Bessel-Gaussian +3 is 0.14 mm. γ of four Bessel-Gaussian modes is 0.6°. To generate Bessel-Gaussian modes with different orders, we appropriately adjust the radius of mirror 1. When the radius of mirror 1 is chosen to be 2.5 × w 0 , 2.6 × w 0 , 2.7 × w 0 , 2.8 × w 0 and the radius of mirror 2 (metasurface-based reflector) is equal to 2Ltanγ, the Bessel-Gaussian 0 to Bessel-Gaussian +3 modes can be excited, respectively. When carefully comparing the real metasurface structure results in Fig. 7 to those ideal results in Fig. 5, one can see slight performance degradation (circular symmetry) of the laser when considering the real metasurface structure and the discontinuities in phase and amplitude responses. However, the obtained multi-ring intensity distribution and helical phase distribution still show the successful implementation of metasurface-assisted OAM-carrying Bessel-Gaussian laser with favorable operation performance (high laser beam quality). The lasing Bessel-Gaussian 0 , Bessel-Gaussian 01 + , Bessel-Gaussian 02 + and Bessel-Gaussian 03 + modes show high fidelities of ~0.889, ~0.889, ~0.881 and ~0.879, respectively. Figure 8 plots the simulated intensity profiles along the radial direction of Bessel-Gaussian modes corresponding to Fig. 7.
To assess the convergence and stability of the Fox-Li algorithm, we analyze the behavior of the relative power loss per round trip, δ k = (I k − I k+1 )/I k , where I k = ∫|u(x, y)| 2 dxdy is the integral of the field intensity at the kth transit. The behavior of the relative power loss as a function of the transits is shown in Fig. 9. We observe that there is a transient state around which the mode oscillates with high loss. Such transient state is unstable. However, after that the system becomes stable with low loss. During the transient process, the field profile evolves to the Bessel-Gaussian mode profile with m = 0, 1, 2 or 3. Once the Bessel-Gaussian mode is stabilized, the low loss remains constant.
We study the fabrication tolerance of the metasurface-assisted OAM-carrying Bessel-Gaussian laser. Figure 10 shows simulated Bessel-Gaussian mode purity as a function of the fabrication errors (Δl x , Δl y ) of metasurface structure, where Δl x and Δl y are the geometric offset of the fabricated eight metasurface units from the designed ones. We choose Bessel-Gaussian 0 mode and Bessel-Gaussian 01 + mode as typical examples. As shown in Fig. 10a, for Bessel-Gaussian 0 mode the maximum mode purity is ~0.889 (Δl x = 0, Δl y = 0). With fabrication errors varying from −10 to 10 nm, the mode purity of Bessel-Gaussian 0 mode is larger than 0.887. For Bessel-Gaussian +1 mode the maximum mode purity is also ~0.889 (Δl x = 0, Δl y = 0), as shown in Fig. 10b. When the fabrication errors vary from −10 to 10 nm, the mode purity of Bessel-Gaussian +1 mode is larger than 0.873. The obtained results indicate that the designed metasurface-assisted OAM-carrying Bessel-Gaussian laser shows favorable tolerance to the fabrication errors. Additionally, high-order Bessel-Gaussian modes show increased sensitivity to the fabrication errors compared to low-order Bessel-Gaussian modes such as Bessel-Gaussian 0 mode.
We further characterize the lasing properties of the metasurface-assisted OAM-carrying Bessel-Gaussian laser. Figure 11a,b show the simulated threshold power (Pth) as functions of the length of gain medium (l) and beam radius of pump (Wr0). One can clearly see that the threshold power decreases with the increase of the length of gain medium and beam radius of pump. Also, higher-order Bessel-Gaussian modes show relatively higher threshold power. Figure 11c describes the relationship between the output lasing power and the transmittance of mirror 1. It is shown that there exists an optimal transmittance of mirror 1 giving a maximum output lasing power. The maximum output power is around 0.91 W under 5 W pump for the zero-order Bessel-Gaussian mode (Bessel-Gaussian 0 ). In addition, the output lasing power decreases for high-order Bessel-Gaussian modes.

Discussion
In summary, we design an OAM-carrying Bessel-Gaussian laser structure consisting of a FP cavity, a pump, a gain medium, and an MSE. The FP cavity is formed by two mirrors, i.e. one partially transparent output plane mirror and the other metasurface-based reflector mirror. The MSE benefits the direct generation of Bessel-Gaussian modes with controllable handedness of the helical phase. The dielectric metasurface with compact structure, providing flexible spatial light manipulation, replaces the traditional complicated volume axicon and spherical mirror. Zero-order and high-order Bessel-Gaussian lasing modes are achievable by properly changing the size of output mirror and putting MSE in the cavity. We also simulate and characterize the designed metasurface-assisted OAM-carrying Bessel-Gaussian laser using the Fox-Li algorithm, manifesting the physical process of mode  Bessel-Gaussian mode with OAM order of 3 is achieved in the simulations. Remarkably, the order of generated Bessel-Gaussian mode is determined by several factors. The structure of laser resonator, the loss and pump efficiency of high-order Bessel-Gaussian modes are the main determining factors. The structure of laser resonator is affected by the resonance cavity length, radius of curvature, dimensions of reflector (mirror 2) and radius of output mirror 1. The loss of high-order modes includes diffraction loss influenced by the conicity angle of the  Bessel beam and structure loss induced by axial nonuniformity and fabrication errors of metasurface-assisted reflector (mirror 2). For pump efficiency, it is affected by the length of the gain medium and beam radius of pump.
In the designed metasurface-assisted OAM-carrying Bessel-Gaussian laser, the conventional axicon together with the spherical mirror is replaced by the metasurface-based reflector. The simulation results show that the real metasurface structure and its induced discontinuities in phase and amplitude responses have slight impact on the performance degradation of the generated OAM-carrying Bessel-Gaussian modes. All these factors could have influence on the order of generated Bessel-Gaussian modes. The presented metasurface-assisted OAM-carrying Bessel-Gaussian laser may further find interesting applications not only in optical communication but also in non-communication fields.

Methods
Fox-Li algorithm is used to manifest the physical process of mode lasing and propagation in the laser resonator. Firstly, the phase and amplitude of laser beams start with random distribution, i.e. initializing the electrical field with random field distribution. Then, taking the diffraction integral equation method for simulations. The beams travel twice through the resonator and are reflected once by the metasurface-based reflector. After that, comparing the calculated results with the previous transverse electrical field to find whether the field distribution profile of beams tends to be stable. If not, then repeating the iterative process until it becomes stable and the desired Bessel-Gaussian mode will be established.