Generation of terahertz vector beams using dielectric metasurfaces via spin-decoupled phase control

Cylindrical vector beams (CVBs), being a special kind of beams with spatially variant states of polarizations, are promising in photonics applications, including highresolution imaging, plasmon excitation, optical trapping, and laser machining. Recently, generating CVBs using metasurfaces has drawn enormous interest owing to their highly designable, multifunctional, and integratable features. However, related studies remain unexplored in the terahertz regime. Here, a generic method for efficiently generating terahertz CVBs carrying orbital angularmomentums (OAMs) is proposed and experimentally demonstrated using transmission-type spatial-variant dielectric metasurfaces, which is realized by designing the interference between the two circularly polarized transmission components. This method is based on spin-decoupled phase control allowed by simultaneously manipulating the dynamic phase and geometric phase of each structure, endowing more degree of freedom in designing the vector beams. Two types of metasurfaces which respectively generate polarization-dependent terahertz vector vortex beams (VVBs) and vector Bessel beams (VBBs) are experimentally characterized. The proposed method opens a new window to generate versatile vector beams, providing new capabilities in developing novel, compact, andhigh-performancedevicesapplicable to broad electromagnetic spectral regimes.


Introduction
As the vector nature of light, polarization plays an important role in light-matter interactions. Traditional optics usually deal with light beams with homogeneous state of polarizations, revealing many intriguing and important phenomena and leading to various optical devices and systems. Recently, increasing attention has been brought to beams with inhomogeneous polarization features, such as cylindrical vector beams (CVBs) with axially symmetric polarization distributions [1,2]. Similar to vortex beams, CVBs can also carry orbital angular momentums (OAMs) [3]. Such expanded polarization and phase properties over homogeneously polarized beams provide more plentiful interaction manners with matter, thus may enable novel effects and potential applications. One particular and attractive application is using CVBs in high numerical aperture (NA) focusing. It has been demonstrated that radially polarized beams can be focused more tightly than homogeneously polarized beams due to a large longitudinally polarized field component in the focus region [4,5], making those very potential in high-resolution imaging [6], lithography [7], optical trapping [8,9], and sensing [10]. Meanwhile, the polarization distributions of CVBs can match with particular modes in waveguides and spatial polarity distributions of many surface plasmon couplers, thus promising the development of efficient waveguide and plasmonic couplers [11,12]. CVBs can be generated via various approaches involving polarization and phase control, such as segmented spatially variant retardation plates [5,13], radially biased electro-optic crystals [14], spiral phase elements [15,16], subwavelength ring gratings [17,18], and liquid crystal spatial light modulators [19,20], etc. However, these traditional components are very generally bulky and/ or need to work along with other components, which hinders the development of CVBs in integrated applications.
Recently, metasurfaces have attracted immense interest owing to the unprecedented light manipulating ability in a subwavelength resolution [21], which allows pixelated and simultaneous control over the interfacial phase, polarization and amplitude distributions of light at will, serving as a new and promising platform for designing next-generation ultrathin devices, such as anomalous deflectors [22][23][24], lenses [25][26][27], holograms [28][29][30], polarization generators and analyzers [31][32][33][34], surface plasmon couplers [35,36], and special beam generators [22,25,37,38], including CVB generators [39][40][41][42][43]. Among the reported metasurface-based CVB generators, geometric phase method is widely applied owing to its flexible phase control mechanism and convenient polarization reorganization ability in circular polarization basis, such as generating vector vortex beams (VVBs) [44]. Though it allows OAM of controllable value, the inter-related phase responses of the geometric phase method to some extent limit the designing and controlling freedom. For example, it cannot be used to generate VVBs composed by arbitrarily different circularly polarized vortex beams carrying different OAMs. As another example, when radial phase gradient is required in generating vector Bessel beams (VBBs), there will be always one diverging circularly polarized output under certain circularly polarized incidence due to reverse phase sign, which makes the two circularly polarized output components separated with each other and unable to compose the VBBs. Multiplexed method can accomplish these tasks, however, the cost is the reduction of the efficiency. To overcome this problem, it is crucial to introduce dynamic phase control at the same time to decouple the phase responses of the two circularly polarized components [32,42,45].
CVBs are essential in a number of terahertz-related applications. One of which is using radially polarized beams to efficiently excite surface plasmon mode on metal wire for terahertz communications [46][47][48]. So far, besides segmented half-wave plate [49], terahertz CVBs are usually generated actively by employing photoconductive antennas with concentric electrodes [50], illuminating femtosecond CVBs onto nonlinear crystals [51], and applying segmented nonlinear crystals with rotated crystal axis [52]. However, compact and integratable methods to generate terahertz CVBs using metasurfaces are seldom reported [53].
In this article, we propose and experimentally investigate an approach to efficiently generate CVBs using transmissiontype all-dielectric metasurfaces via spin-decoupled phase control method in the terahertz regime [32]. Being free of Ohmic loss, dielectric metasurfaces have been demonstrated to be able to more efficiently manipulate the wavefront compared to the plasmonic counterparts [54]. To illustrate the extra flexibility of the spin-decoupled phase control method, two terahertz metasurface CVB generators which cannot be realized using sole geometric phase method are designed, as schematically illustrated in Figure 1. One metasurface is designed for generating VVBs carrying OAM with topological charge of l = +1 ,a ss h o w ni nF i g u r e1 A .U n d e rl e f t -h a n d e d circularly polarized (LCP) incidence, it generates a righthanded circularly polarized (RCP) vortex beam with OAM of l = +2. Under the RCP incidence, it solely functions as a perfect half-wave plate. Under the linearly polarized incidence, the output wavefront will be a superposition of the RCP vortex wave and the LCP plane wave. By designing the spatial phase difference distributions between these two components, radially and azimuthally polarized VVBs with OAM of l = +1 can be generated under the x-a n dy-polarized incidences, respectively [40,43]. The other metasurface is designed with a similarwayforgeneratingVBBswithOAMof0,asshownin Figure 1B. Under the LCP and RCP incidences, it generates R C Pa n dL C PB e s s e lb e a m sw i t hO A M so fl = +1a n dl = −1, respectively. Under the x-a n dy-polarized incidences, the polarization superposition results in a radially and an azimuthally polarized Bessel beam, respectively. The simulated averaging efficiencies of the sole structured surfaces from the substrate to the air of the two metasurfaces are ∼58 and ∼60%, respectively. The proposed method provides a new way to develop novel CVB-relatedcompactandefficient devices essential in terahertz photonics.
2 Results and discussion

Structure design
The basic unit cell of the dielectric metasurface is schematically illustrated in Figure 2A, which is a rectangular-shape silicon pillar on the same silicon substrate. The pillar has a period of P f = P s = 150 μm and a height of h = 200 μm. Such a pillar has been demonstrated to function as an anisotropic effective truncated waveguide, whose phase responses of the two eigen orthogonally polarized waves can be freely controlled to individually cover the whole 2π range by changing the width D f and length D s [55]. This kind of phase caused by the physical dimension of the pillar is called dynamic phase. To further incorporate geometric phase, it is required to introduce rotation operation. The inset of Figure 2A illustrates a schematic of the relation of the pillar and the xy basis. The fs basis represents the local coordinate of the pillar which has a θ rotation angle with respect to the xy basis, where f and s axes are along the two symmetric axes of the pillar. Suppose the eigen transmission matrix of the pillar in the fs basis is

T fs
A f e iφ f 0 0 A s e iφ s , where A i and φ i represent the transmission amplitude and phase shift of the i-polarized wave with i ∈ {f, s}, respectively. To increase the efficiency of polarization control, we consider a special case of A f = A s = A and φ s -φ f = π, which indicates that the pillar functions as a half-wave plate. The corresponding transmission matrix in the xy basis can be calculated by applying rotation operation of coordinate transformation, which can be expressed as: Equation (1) indicates the polarization orientation of the output wave can be controlled by the rotation angle θ, while the phase of this polarized wave can be controlled by the dynamic phase φ f . To show the spin-decoupled phase control ability, transmissions under the circularly polarized incidences are also calculated. For an incident LCP i , the half-wave plate  property will enable the output wave RCP, which can be expressed as: Similarly, for an incident RCP plane wave of −i , the output LCP wave can be expressed as: Here, E L(R) and δ L(R) are the amplitude and initial phase of the incident waves with L and R representing the LCP and RCP waves, respectively; Ф RL and Ф LR are the corresponding additional phase shifts resulting from the interaction with the pillar, the former and latter subscripts represent the output and incident polarization states, respectively (same hereinafter). The ±2θ is the well-known geometric phase. It is seen that phase responses of the two output circularly polarized waves can be decoupled by further controlling the dynamic phase term φ f . From these points, we could re-understand Equation (1) that the geometric phase determines the polarization rotation, while the dynamic phase determines the phase shift of the output polarized wave relative to the incident. In particular, under the linearly polarized incidence, the output wave will also be linearly polarized but with a polarization rotation of twice of the angle formed by the pillar's fast axis and the incident polarization owing to the half-wave-plate property. Equations (2) and (3) indicate that such pillars could be used to simultaneously and individually control the polarization and phase distributions at the interfaces by designing the superposition of the LCP and RCP output components, where the polarization distribution can be further controlled by the incident polarization state, providing a versatile platform in generating various polarization-controlled CVBs.
In order to realize full 2π phase control, eight pillars numbered 1−8 are selected from the structure database developed in our previous work, which is built at 1.0 THz [55]. As shown in the top insets of Figure 2B Figure 2B illustrates the corresponding simulated dynamic phases φ f and φ s . From left to right, both the phases are linearly increased with a step of π/4. Figure 2C illustrates the corresponding simulated amplitude transmission (square root of the power transmittance) of the sole structured surface A f and A s from the silicon substrate to the air. It is seen that the values are all quite near with each other in a relatively high level around 0.85. More importantly, the phase differences Δφ = φ s -φ f are all nearly π and the amplitude differences ΔA = A s -A f are all around 0. These indicate that all the selected pillars function well as half-wave plates at 1.0 THz (see Supplementary section 1), which fulfill the above derivation. As our metasurfaces rely on superposition of the output cross-polarized LCP and RCP waves to generate CVBs, such a feature ensures the co-polarized parts have negligible influence on the metasurface performance. With these eight pillars, two metasurfaces for generating VVBs and VBBs are separately designed.

Metasurface for VVB generation
To realize the function illustrated in Figure 1A, the phase distributions Ф RL and Ф LR should satisfy: respectively, where r =(x 2 + y 2 ) 1/2 and α = atan2(y, x) ∈ [0, 2π] represent the radius and azimuth angle of polar coordinate, respectively. According to Equations (2) and (3), the output waves under the LCP and RCP incidences can be expressed as: respectively. Equation (6) indicates that the transmitted beam is RCP vortex beam carrying OAM of l = +2 with initial phase δ L remained from the incident LCP wave. Equation (7) indicates that the transmitted beam is LCP beam carrying OAM of l = 0 with initial phase δ R remained from the incident RCP wave. Under the linearly polarized incidence, the wave can be decomposed to the LCP and RCP components with equal amplitudes where E L = E R = E, δ L = −γ, and δ R = γ, the output wave is given by: where γ =(δ R − δ L )/2 represents the orientation angle of the incident polarization. Equation (8) indicates that the output wave can have inhomogeneous linear polarization and phase distributions. The output polarization rotates linearly as the azimuth angle α with orientation angle of α − γ, resulting in axially symmetric distributions controlled by γ. The output phase at an arbitrary point is equal to the corresponding azimuth angle α, corresponding to a CVB with OAM of l = 1 regardless of γ. Especially, under the x-polarized incidence (γ = 0), the output is radially polarized VVB, while under the y-polarized incidence (γ = π/2), the output is azimuthally polarized VVB. To realize the above design, the metasurface VVB generator is designed based on Equations (4) and (5) using the eight selected pillars. In the metasurface, the pillars are arranged in a square lattice. The coordinate of each lattice point (x, y) can be transformed to the polar coordinate of (r, α), at which a pillar satisfied φ f = α and θ = α/2 should be placed. To simplify the design, the azimuth angle is divided into eight angle regions evenly. Each region contains only one type of the pillars with φ f , whose orientation angles are accurately determined by α. The sample was fabricated using conventional lithography, together with deep reactive ion etching on a high-resistivity silicon wafer. The microscope image of the center part of the fabricated sample is shown in Figure 3B.
A scanning broadband terahertz time-domain microscopy system was used to experimentally evaluate the performance of the proposed metasurface VVB generator [26,30,38,55]. Figure 3A schematically illustrates the experimental setup. The terahertz beam generated from the emitter was firstly collimated by a lens to form a quasi-Gauss beam, which was almost x-polarized. Two metallic terahertz grid polarizers and a quarter-wave plate were then employed to control the incident polarization state. Polarizer 1 was used to further enable the incident terahertz beam to be x-polarized. The LCP and RCP beams could be obtained by adjusting the quarter-wave plate without placing polarizer 2 behind it, while the x-and y-polarized beams could be obtained by further placing and adjusting polarizer 2. After the polarization of the incident terahertz beam was controlled, it normally illuminates the sample from the substrate side. A photoconductive-antenna probe was then placed above the sample to map the output terahertz beam via raster scanning. Here, the probe could only detect the x-polarized output. Under the circularly polarized incidences, the outputs were cross-circularly-polarized, so measuring the x-polarized component (E xL and E xR ) only was enough to exhibit the performances. Under the linearly polarized incidences, however, two orthogonally polarized components should be measured. We did this by rotating the sample and polarizer 2 in combination in the measurements. For example, suppose the initial configuration was used to measure E xx , then E 45°x , E yx , and E −45°x were measured by simultaneously rotating the sample and polarizer 2 by −45°,9 0°, and 45°, respectively. E yy was measured by solely rotating the sample by 90°as compared to that of E xx , then E 45°y , E xy , and E −45°y were measured by simultaneously rotating the sample and polarizer 2 by −45°,90°, and 45°, respectively. It should be noted that this measuring method cannot rebuild the polarization distributions, since the scanning grid relative to the metasurfaces before and after the rotation cannot perfect overlap with each other. However, this method allows us to show the vector and OAM properties of the generated beams in a similar demonstration way as those in the optical range [40,43].
In the measurement, the output beam was mapped in the xy plane at a height of 5 mm above the sample. The scanning area is 6 × 6 mm with a step of 0.2 mm. Figure 4A,C illustrate the measured normalized intensity distributions |E xL | 2 and |E xR | 2 at 1.1 THz under the LCP and RCP incidences, while Figure 4B,D illustrate the corresponding measured phase distributions, respectively. The clear donut-shape intensity distribution and 4π phase loop indicate a good vortex beam output carrying OAM of l = +2 under the LCP incidence. The Gaussian-like intensity distribution and nearly homogeneous phase distribution indicate that no OAM is attached from the sample under the RCP incidence. Figure 4E-L illustrate the measured normalized intensity distributions |E xx | 2 , |E 45°x | 2 , |E yx | 2 , and | E −45°x | 2 under the x-polarized incidence, as well as |E xy | 2 , | E 45°y | 2 , |E yy | 2 , and |E −45°y | 2 under the y-polarized incidence at 1.1 THz, respectively. Clear "s"-shape intensity distributions are observed whose rotation is determined by the incident and detected polarizations. The rotation feature of the "s"-shape intensity distributions with regard to the detection polarizations under the x-and y-polarized incidences, as well as the azimuth angles of the maximum intensities, are consistent with previous reports [40,43], indicating a radially and an azimuthally polarized VVBs with both OAMs of l = +1. To further verify the measured results, numerical simulations were also carried out. The insets in Figure 4 are the corresponding simulated distributions under the same conditions at 1.0 THz, which are consistent with the measured results. The efficiencies of the VVB generator are also evaluated in simulations, which are calculated by normalizing the integrated intensity of the output beam in the air side above the metasurface to that of the input beam inside the substrate. Under the LCP, RCP, xand y-polarized incidences, the efficiencies of the corresponding desired output beams are ∼60, ∼51, ∼57, and ∼58%, respectively.

Metasurface for VBB generation
According to the above designing strategy, this method should be able to generate Bessel beams under the circularly polarized incidences. Compared to the vortex beam, the generation of Bessel beams further requires a phase gradient along the radial direction. For an l-th order Bessel beam, the interfacial phase distribution in the polar basis can be expressed as [37]: Ф(r, α)=−2πrNA/λ + lα, where NA = λ/P is the numerical aperture, λ is the working wavelength, P is the overall period of 2π phase coverage along the radial direction, and l is also the topological charge which represents the value of the carried OAM. To realize the function illustrated in Figure 1B, the phase distributions Ф RL and Ф LR should satisfy: respectively. According to Equations (2) and (3), the output waves under the LCP and RCP incidences can be expressed as: respectively. Equation (11) indicates that the transmitted beam is RCP Bessel beam carrying OAM of l = +1 with initial phase δ L remained from the incident LCP wave. Equation (12) indicates that the transmitted beam is LCP beam carrying OAM of l = −1 with initial phase δ R remained from the incident RCP wave. Under the linearly polarized incidence, the output wave is given by: It can be seen that the output beam still has Bessel phase distribution but no longer carries any OAM (l = 0). The polarization response is the same as that of the VVB generator as indicated by Equations (8) and (13). Under the x-polarized incidence (γ = 0), the output is radially polarized VBB, while under y-polarized incidence (γ = π/2), the output is azimuthally polarized VBB.
To realize the above design, the azimuth angle is also divided into eight angle regions evenly, meanwhile, we also divide the metasurface into several regions along the radial direction with a phase interval of π/4. In this case, according to Equations (9) and (10), each newly formed intersecting region contains one type of the pillars with φ f = −2πr/P and P = 1240 μm, whose orientation angles are also accurately determined by α with θ = α/2. Figure 3C illustrates the microscope image of the center part of the fabricated sample.
To experimentally characterize the proposed metasurface VBB generator, the same method was utilized as for the VVB generator. In the measurement, the output beam was mapped in the xy plane at a height of 7 mm above the sample. The scanning area is 6 × 6 mm with a step of 0.2 mm. Figure 5A,D illustrate the measured normalized intensity distributions |E xL | 2 and |E xR | 2 at 1.1 THz under the LCP and RCP incidences, while Figure 5B,E illustrate the corresponding phase distributions, respectively. The Bessel distributions with donut-shape intensity profiles at the center and the spiral phase distributions with ±2π phase loops indicate that the output beams are +1-th order and −1th order Bessel beams, respectively. To illustrate the diffraction-free feature of the Bessel beams, the intensity distributions in the yz plane under the LCP and RCP incidences were also measured, as illustrated in Figure 5C,F, respectively. The distributions are measured from −2.5 to 2.5 mm with a step of 0.1 mm along the y direction, and from 2 to 12 mm with a step of 0.5 mm along the z direction. It is seen the diffraction-free distances are longer than 20λ. Figures 5G to 5N illustrate the measured normalized intensity distributions |E xx | 2 , |E 45°x | 2 , |E yx | 2 ,a n d|E −45°x | 2 under the x-polarized incidence, as well as |E xy | 2 , |E 45°y | 2 , |E yy | 2 ,and |E −45°y | 2 under the y-polarized incidence at 1.1 THz, respectively. Clear split-lobe patterns (including center main lobe and outer sidelobes) are observed, indicating vector beams carrying OAM of l = 0. Under the x-polarized incidence, the output lobes are the strongest along the detection polarization directions, indicating a radially polarized VBB. Under the y-polarized incidence, the output lobes are the strongest along the directions perpendicular to the detection polarization directions, indicating an azimuthally polarized VBB. The insets in Figure 5 illustrate the corresponding simulated distributions at 1.0 THz, which are consistent well with the measured results. The simulated polarization distributions of the metasurface VBB generator can be seen in Supplementary section 2. Under the LCP, RCP, x-a n dy-polarized incidences, the simulated efficiencies of the corresponding desired output beams are ∼59, ∼59, ∼60, and ∼61%, respectively.
It is seen from Figures 4 and 5 that the measured phase distributions are not as good as those in the simulations in the areas where the intensities tend to be zero. This can be attributed to the limited signal-to-noise ratio of our system, which cannot measure the accurate phase values in these areas. Besides, the applied 0.2 mm scanning step in the measurements also reduced the phase resolution to some extent. Nevertheless, our measured phase distributions can well demonstrate the OAM feature of the generated beams, which is an important parameter of CVBs. The working frequency of the fabricated metasurfaces (1.1 THz) here is slightly different from that in the design (1.0 THz). This can be attributed to the fabrication errors and the refractive index deviation from the design owing to the dispersive phase responses of the eight selected silicon pillars by the dynamic phase contributions (see Supplementary sections 3 and 4). Taking the overall half-waveplate performances of the eight selected silicon pillars into account, the effective working bandwidth is around 0.04 THz. Besides, it is more interesting to note that there are two phase singularities in both the measured and simulated phase distributions of the generated vortex beam carrying OAM of l = +2 by the metasurface VVB generator under the LCP incidence, as illustrated in Figure 4B. This is caused by the phase discontinuous feature of the metasurface, where the ideal single phase singularity will be divided into two. The distance between the two phase singularities increases as the pixel size of the structure increases, and also as the number of the angle regions divided in the azimuthal direction decreases. Another reason should be the incomplete circular-polarization conversion efficiency of the composed structures, the interference effect between the output RCP and residual LCP components also contributes to the emergence of the two phase singularities since the x-polarized component was measured here. Detail information can be found in the calculated results in the Supplementary section 5. Owing to the fabrication error induced working frequency shift, the actual performance of the structures' half-wave plate property reduces, so the two phase singularities are more obvious in the measured results.

Conclusions
We propose a general method for generating CVBs using all dielectric metasurfaces in the terahertz regime. This is achieved by controlling the superposition of the two crosscircularly-polarized output beams generated under the circularly polarized incidences. To realize more degree of freedom to manipulate the two output beams, a spindecoupled phase control method is applied by simultaneously manipulating the dynamic phase and geometric phase of the unit cells. Two types of metasurface CVB generators are experimentally realized, one for generating VVBs while the other for generating VBBs under the linearly polarized incidences. Both of them are composed of structures with spatially variant dimensions and orientations, which are eight basic rectangular-shape silicon pillars possessing half-wave plate property and 2π dynamic phase coverage. The averaging efficiencies of the two VVB and VBB generators are ∼58 and ∼60%, respectively. Such efficiencies can be further enhanced by using a substrates with lower refractive index to decrease the reflection loss. Comparing with the transmission efficiency of a bare silicon interface, the relative efficiencies of the two generators achieve ∼83 and ∼86%, respectively. The ultrathin and high-efficiency properties of the proposed CVB generators may find broad applications in integrated systems. Though the metasurfaces are demonstrated in the terahertz regime, the designing method can be extended to other frequency ranges of the electromagnetic spectrum.