Ultra-Broadband Nonlinearity Enhancement based on a Novel Graphene-Silicon Hybrid Waveguide: Structure Design and Theoretical Analysis

A graphene-silicon hybrid waveguide with a dielectric spacer is proposed to enhance the nonlinear response in ultra-wide wavelength range by applying graphene’s broadband highly nonlinear optical properties. The chemical potential of the graphene layer is tuned to satisfy the resonance condition and hence a low propagation loss is obtained. The dielectric spacer is used for avoiding additional free-carrier-absorption loss due to carrier interchange between the silicon core and the graphene layer. Aiming at the special waveguide structure with ultra-thin graphene layer, a full-vectorial theoretical model is developed to analyze its nonlinear properties. The waveguide dimensions are optimized in terms of the nonlinear parameter. The proposed hybrid waveguide exhibits high nonlinearity enhancement in an ultra-broad wavelength region covering near-infrared and mid-infrared bands. The conversion efficiency for a degenerate four-wave mixing process reaches −18.5 dB only with a pump power of 0.5 W and a waveguide length of tens of microns. In the wavelength range of 1.3–2.3 μm, the conversion efficiency can be kept stable by adopting suitable waveguide geometry and length. The corresponding 3-dB bandwidth can reach 40–110 nm for each fixed pump. The graphene-silicon hybrid waveguide has the potential to support chip-scale nonlinear applications in both near- and mid-infrared bands.

. In the linear optical response, the equivalent complex relative permittivity of the Δ-thick graphene can be expressed as , where the linear conductivity σ g (1) is calculated from the Kubo formula 33 . At optical frequencies, graphene's complex refractive index κ = +  n n i g g g is defined as ε =  n g g , where κ g is called the extinction coefficient and determines the linear optical absorption. As shown in Fig. 1a, for ω µ  2 c  , the attenuation constant is very large and vary slightly at high optical frequency. When the photon energy drops below the transition threshold (ω µ < 2 c ), no interband transition is allowed and the absorption sharply decreases. However, in the long wavelength range, the intraband transition becomes dominant and the absorption becomes obvious once more. The location of this low absorption window depends on graphene's chemical potential and can be tuned in a wide wavelength range.
As for the nonlinear case, we only consider the third-order conductivity σ ω ω ω ( , 3 and neglect higher order nonlinearities due to the inversion symmetry of the graphene lattice. According to the analytic results of the third-order conductivity (such as Kerr nonlinearity in Fig. 9a in ref. 27 ), for a given chemical potential there are certain frequency ranges where the nonlinear response of graphene are exceptionally large due to the resonance enhancement of photon absorption. At high (infrared) frequencies, the imaginary part of third-order conductivity σ ω ω ω − Im( ( , , )) g (3) of Kerr nonlinearity achieves its maximum at the resonant condition: ω µ ≈ 2 c . As the phenomenological relaxation parameter Γ is a small value ([0, 60] meV) 24 , the position of high-frequency resonance mainly depends on the Fermi level, i.e. µ c . The equivalent complex third-order susceptibility χ ω ω ω − ( , , ) g of Kerr nonlinearity relates to the third-order conductivity as χ ω ω ω σ ω ω ω εω ( 3) 0 . Based on the analytical theory in ref. 27 , in a broad band (0.9 ~ 3 μm), Fig. 1b shows the inverted hump-shaped response of χ g (3) at different chemical potentials. However, assuming a low chemical potential µ = . 0 1 c eV, the calculated values of χ Re( ) g (3) are smaller than the results reported by Hendry et al. 4 by two orders of magnitude and have weak dependence on the relaxation parameter in the range of [0, 60] meV. In practice, we adopt the nonlinear refractive index n 2 rather than the third-order susceptibility to characterize the Kerr nonlinearity of a particular kind Scientific RepoRts | 7: 12290 | DOI:10.1038/s41598-017-12554-6 of material, which relates to the real part of χ (3) as ε χ = c n Re 3/(4 ) ( ) 2 (3) (ε denotes the material's permittivity). As both graphene's linear and nonlinear conductivities have considerable real parts and imaginary parts, the derivation of nonlinear refractive index n g 2, should follow the results of del Corso and Soles 34 . For a given chemical potential, the complex n g 2, of grapheme depicted in Fig. 1c has different shape from the real part of χ g (3) in Fig. 1b, and both the real and imaginary parts achieve their maxima near the resonant frequency ( ω µ ≈ 2 c  ). Usually, the positive imaginary part of n 2 indicates the TPA process. While in graphene, the imaginary part is negative for ω µ > 2 c  and can be understood as a correction to the simple linear prediction of the absorption as the one-photon absorption is always positive in this frequency range. We also find that the calculated values of n Re( ) g 2, at specific wavelength is in agreement with the experimental measurements 21 . For a monolayer of pristine graphene with µ = .
0 02 c eV (taking into the slight substrate-induced doping), the nonlinear refractive index (red curve in Fig. 1c) is in the order of 10 −13 ∼ 10 −12 and essentially dispersionless over the broad wavelength range. By tuning graphene's Fermi level, the peak response of n Re( ) can be shifted closely to the central frequency ω 0 of the interested frequency band according to the resonance condition ω µ . In the vicinity of the most important feature at ω µ ≈ 2 c 0 , the amplitude of n Re( ) Fig. 1c shows no obvious increase compared with the pristine grapheme (µ = . 0 02 c eV), but the linear absorption is greatly reduced as shown in Fig. 1a. Therefore, shifting the Fermi level to match the resonance condition at ω 0 shows the potential in improving graphene's nonlinear performance by increasing the figure of merit (FOM) 10 . Figure 2a shows the structure of the proposed graphene-silicon hybrid waveguide. Monolayer graphene film is on the top of the SOI nanowire and between them there is a 7-nm thick Al 2 O 3 spacer 6 . This insulating layer blocks the carrier interchange between the silicon core and the graphene layer as the photoexcited electrons may transfer from graphene to silicon and cause additional free carrier absorption 35 . Then a gate voltage can be applied to this hybrid waveguide for actively tuning the Fermi level of the cladding graphene film. Further considering the effective interaction between the light field and the graphene in the particular band centered at ω 0 , the dimensions of the waveguide cross-section is designed and optimized to support the electric field maximized at its top and bottom surfaces (fundamental TM mode) when graphene's Fermi level is tuned to match the resonant condition . Using the finite-element-frequency-domain (FEFD) method, the electric field distribution at 1.55 μm for a waveguide cross-section of 360 nm × 266 nm is simulated in the inset of Fig. 2b. Also, Fig. 2b shows the normalized electric field distribution E norm along the y axis and sufficient light concentration can be expected near the graphene-silicon boundary. As the fact that graphene only interacts with the tangential (in-plane) electric field of electromagnetic waves, we only consider the tangential components of the electric field and magnetic field in the graphene layer in the following simulation.
In order to derive the coupled equations for parametric conversion process in graphene-silicon hybrid waveguides, we start from the conjugated form of the Lorentz reciprocity theorem in frequency domain: (see e.g. ref. 36 The arbitrary integration area s is in the waveguide cross-section plane and ẑ is the unit vector oriented along the longitudinal axis. At a certain frequency ω, we introduce two sets of electromagnetic fields: the vectors ( ω E ( ) 0 , ω H ( ) 0 ) represent the total unperturbed fields of the ideal lossless waveguide and the vectors ( ω E( ), ω H( )) correspond to the situation under which the optical properties are perturbed both by the nonlinear optical effects, such as SPM, XPM, and FWM, and the linear and nonlinear losses. They are related to the vector F c , which is defined as = × + × is the perturbed polarization with δ ω P ( ) lin and ω P ( ) nl representing its linear and nonlinear components, respectively. The change of the permittivity δε ω ( ), induced by both linear effects (such as optical losses and free-carrier generation) and nonlinear effects, is the origination of this perturbed polarization 37 .
Supposing the practical field for graphene-silicon hybrid waveguides consists of a superposition of quasi-monochromatic waves, with ω′ denoting any one of the frequencies involved, we can write the total electric field t E( ) and its Fourier transform ω E( ) as where the monochromatic-wave amplitude ω′ E is a vector phasor and = The vectorial solution of Maxwell's equations for the unperturbed fields in such ideal lossless waveguide results in a complete orthonormal set of forward, backward, and radiation propagation modes with certain boundary conditions applied 36 . In our proposed single-mode hybrid waveguides, the input light is coupled into one particular polarization of the fundamental mode, like the TM mode, and one can see no coupling between these two polarization modes during the light propagation by simply applying the reciprocity theorem above. Thus for the propagating TM mode, the forward term of the unperturbed fields notated in vector phasor form (evaluated at ω′) can be expressed as , , ) ; is a normalization constant that defines the mode power and β is the propagation constant at frequency ω′. Modal field distributions, the vector phasors ω′ x y e( , , ) and ω′ x y h( , , ), of propagating modes in the graphene-silicon hybrid waveguide can be obtained through a variety of numerical methods including the finite element method. Assuming that the linear and nonlinear perturbations are small enough not to affect the modal field distribution, the corresponding perturbed fields ω′ E and ω′ H can be expanded according to the unperturbed ones: i z i z It should be noted that the frequency dependence of perturbed fields is totally contained within the amplitude coefficient ω′ a z ( , ), and the mode optical power is ω ω We only consider the unidirectional wave propagation for which the back scattering of a forward propagating laser beam is neglected, and hence ignore the coupling between the unperturbed field with the backward and radiation modes of the perturbed field 30 . Substituting equations (2)-(5) into the simplified form of equation (1), one can obtain: In frequency domain, the induced n th -order polarization can be written in a general form: Equation (7) shows a complete description that requires the knowledge of the tensorial and dispersive properties of the nonlinear susceptibility χ n ( ) . Since there is no second-order nonlinearity in graphene or silicon, here we only consider the third-order nonlinear processes and neglect higher order nonlinearities, i.e.
. As silicon belongs to the m3m point-group symmetry 28 , in Cartesian coordinates, among 81 elements of χ ijkl (3) (i, j, k, l = x, y, z) only 21 elements are nonzero, which depend on only four independent quantities (χ ( 3) ). In practice, the most relevant nonlinearity is the Kerr nonlinearity involving only , and in the wavelength range considered, photon energy is well below silicon's energy gap E g Si , (≈1.12 eV). Therefore, the last three terms are equal. As a result, χ ijkl (3) can be written as: ( 3) represents the nonlinearity anisotropy and has a nearly constant real value for ω < E g Si ,  . Within the rest of this paper we ignore the last term of equation (8) which affects the polarization dependence of nonlinear effects in silicon waveguides and thus ρ = 1 30 . The real and imaginary parts of χ xxxx (3) can be derived from the Kerr coefficient and the TPA coefficient 28 , so here its wavelength dependence is determined by the experimental measurements 38 through the polynomial fit method. In the following numerical calculation, silicon's third-order susceptibilities for the corresponding SPM, XPM and FWM processes are the same and derived from the experimental Kerr nonlinearity. Due to the two-dimensional nature of graphene 22 , among all eight nonzero components of the third-order susceptibility there are three are independent: ( 3) . Recent research reports have established an initial general theory of third-order nonlinear electrodynamic effects of graphene under some conditions, such as linear-dispersion band structure approximation, two-band tight-binding, single-particle effect and so on 23,27 . Beneath the present theory framework, it is reasonable to treat the three independent components of graphene's nonlinearity as equal.
For the continuous-wave DFWM process, the total propagation field is the superposition of the monochromatic waves with frequencies of ω p , ω s , and Here the labels p, s, and i stands for the pump light, the signal light, and the idler light, respectively. In this case, the nonlinear term for the converted idler lightwave ω i on the right side of equation (6) can be written as 30 : The electric/magnetic-field-related terms in equation (9), different from that of the scalar form 28 , are the direct results of full vectorial calculation. After integrating them over the waveguide cross-section plane with the spatial distribution of the third-order nonlinearity as its weighting factor, one can get the nonlinear parameters: γ ipps describes the four-wave mixing process, γ ippi and γ issi are responsible for XPM and γ iiii relates to SPM. Their real parts characterize the capability of attaining nonlinear phase-shift during the light propagation, while the imaginary parts play an important role in the phenomena of nonlinear loss (TPA and light absorption saturation).
The linear part of the perturbed polarization δP lin in equation (6) is determined by the linear change in the dielectric constant induced by optical losses in silicon and graphene and free carrier generation in silicon core. In the wavelength range considered here, the photon energy ω is below the band gap E g Si , of silicon, so the interband one-photon absorption (1PA) in the silicon core is forbidden. While for the case of Fermi-level-tunable graphene, the linear absorption coefficient can be calculated from its complex surface conductivity as stated above. Using the finite element method, the intrinsic linear absorption coefficient of this graphene-silicon hybrid waveguide can be derived from the imaginary part of its effective mode index. Besides the intrinsic linear absorption coefficient in the hybrid materials, the total linear loss coefficient α lin also includes an empirically additional loss due to the waveguide fabrication such as the sidewall roughness, which is assumed to be 3 dB/cm. When ω > E /2 g Si , , TPA occurs in silicon and the interband transition will generate free carriers, which will cause free carrier absorption α F and affect silicon's refractive index. In practice, it's common to employ the empirical formulas (18)- (20) in ref. 28 to calculate these factors 39 . We find that, compared with the absorption coefficient, the variation of the effective mode index is negligible. As stated above, it's reasonable to ignore higher-order multiphoton absorption under propagation of weak light field. Thus, in graphene-silicon hybrid waveguide, the linear term of equation (6) can be expressed as: With all necessary parameters settled, the coupled equations of the amplitude evolution for the interacting frequencies (ω ω ω , , i s p ) are obtained from equation (6): The above system of coupled equations are general and can be applied to an arbitrary waveguide. In the following, we will numerically investigate the wave dynamics upon propagation in the proposed graphene-silicon hybrid waveguide and discuss its nonlinear performance (nonlinear parameter, conversion efficiency and conversion bandwidth) in the near-and mid-infrared wavelength range. From the derivation above, all the effective parameters of the coupled equations can be considered as the optical properties of the corresponding material averaged over an inhomogeneous cross-section weighted by the field distribution. In order to achieve the best nonlinear performance of a DFWM process in the graphene-silicon hybrid waveguide, the cross section dimensions should be designed to maximize the nonlinear parameter at the pump wavelength λ p , and then the waveguide length and pump power should be optimized for the highest conversion efficiency. According to equation (9), the effective nonlinear parameter γ for Kerr nonlinearity at frequency ω can be written as During the optimization, we take the telecom band as an example and set the pump wavelength λ p at 1.55 μm. Though graphene and silicon have shown their broadband nature of both linear and nonlinear properties in nearand mid-infrared spectrum, for a specified sectional dimension of the graphene-silicon hybrid waveguide, simulation results show that γ Re( ) has a smooth inverted hump-shaped distribution with respect to the light wavelength. The Fermi level of graphene satisfy the resonant condition ω µ ≈ 2 p c  for the pump at 1.55 μm, i.e. µ = .
0 4 c eV, and as depictedin Fig. 3a and b, the position of the peak γ Re( ) depends on the dimensions of waveguide's cross-section and will shift to longer/shorter wavelength as the width W (height H) increases/decreases. So the bandwidth of nonlinear response is limited by the wavelength-selected performance of the waveguide structure itself. Therefore, the best nonlinear parameter in the wavelength band interested can be achieved through: firstly, tuning the Fermi level of graphene to match the resonant condition at the pump wavelength λ p , and secondly, sweeping the waveguide's width and height simultaneously to get the maximum γ at λ p . Under this determined chemical potential for the pump at 1.55 μm, Fig. 3c simulates the nonlinear parameter as a function of the waveguide width and height. The hybrid waveguide is optimized to be 360 nm × 266 nm at 1.55 μm to acquire the maximum absolute nonlinear parameter of 6246.9 W −1 m −1 , an order of magnitude larger than that in graphene-silicon hybrid cavity 11 . The waveguide's linear loss α lin is given in Fig. 3d, where one can find that the linear loss is relatively low for the optimized waveguide and the value is about 0.07 dB/μm. As the wavelength range for DFWM varies, the pump wavelength should be tuned correspondingly. As a result, the chemical potential and the cross-sectional dimensions of the hybrid waveguide should also be optimized once more.
Conversion efficiency of the DFWM process, which is defined as the ratio of the output idler power with respect to the incident signal power, is an important parameter which needs to be analyzed under various constrains. According to the coupled equations (11)-(13), the wavelength conversion ability of a specified waveguide structure mainly depends on two factors: the optical loss and the phase-matching condition. With the optimal cross-sectional dimensions (360 nm × 266 nm) at the pump wavelength of 1.55 μm, one can numerically solve these differential equations and achieve the maximum conversion efficiency by optimizing the input pump power and the waveguide length. Figure 4a shows the simulated conversion efficiency η versus the pump power (P p ) and the waveguide length (L), where the signal wavelength λ s = 1.549 μm is set very closely to the pump wavelength. In this case, the linear phase mismatch is reasonable to be neglected. One can see that, for a given input pump power, the smooth conversion efficiency curve in Fig. 4b is hump-shaped along the propagating length and its maximum increases with the input pump power. The conversion efficiency will reach its maximum of −18.43 dB when the hybrid waveguide is 82 μm long and the pump power is 0.5 W. As the Fermi level of graphene ( µ c = 0.4 eV) is tuned to match the resonant condition at the pump wavelength, graphene exhibits the saturable absorption because the imaginary part of its nonlinear refractive index is negative from Fig. 1c. The estimated saturation field strength is about 3 × 10 7 V/m 24 which corresponds to about 1.2 W input light power in our case. For the validity of perturbation assumption of the whole theoretical model and avoid higher order nonlinearity in much tenser light intensity regime, we limit the pump power below 0.5 W in our numerical calculation.
Above the nonlinearity enhancement is demonstrated in near-infrared (telecommunication) band centering at 1.55 μm. Actually, according to the tunable low-loss bands of graphene and silicon, the designed graphene-silicon hybrid waveguide structure is expected to enhance the nonlinearity in a wide wavelength range, covering both near-and mid-infrared bands. For the applications in different bands, the pump wavelengths should be different, and the chemical potential of the graphene layer needs to be tuned to match the resonant condition (ω µ ≈ 2 p c ) in Fig. 5a. Correspondingly, the cross-sectional dimensions and the waveguide length of the graphene-silicon hybrid waveguide also need to be optimized for various pump wavelengths. By using the  Fig. 5b. One can find that the required waveguide width and height almost increase linearly with the working wavelength. In comparison, we also perform the cross-sectional-dimension optimization in the case of pristine graphene ( µ c = 0.02 eV) (hollow marks), and find that, for the same pump wavelength, the optimal dimensions are in agreement with the case of Fermi-level-tuned graphene (solid marks) as the marks in Fig. 5b of both cases overlap each other completely. Based on the DFWM process, the nonlinearity enhancement is evaluated with the optimal waveguide cross-section for each pump wavelength λ p by setting the wavelength spacing between the  pump and the signal as 1 nm. Figure 5c shows the required waveguide length L 0 and the maximum conversion efficiency that supported for the Fermi-level-tuned graphene (solid marks) and pristine graphene (hollow marks) cases when the pump and signal powers are 0.5 W and 1 mW, respectively. When the pump wavelength λ p varying from 1.31 μm to 2.25 μm, in both cases the optimal waveguide length L 0 increases linearly and the corresponding maximum conversion efficiency η max keeps in a high level, about −18.4 dB for the Fermi-level-tuned hybrid waveguide and −20.75 dB for the pristine-graphene one, and both fluctuates slightly in the range of 0.6 dB and 0.16 dB, respectively. Moreover, the conversion efficiency can be obviously improved by tuning the chemical potential of graphene. The maximum conversion efficiency of the Fermi-level-tuned hybrid waveguide is averagely 2-dB higher than that of the pristine-graphene hybrid waveguide since the propagation loss is reduced when the resonant condition of graphene is satisfied.
Bandwidth performance is an important indication of the parametric conversion process in the graphene-silicon hybrid waveguide. According to the graphene's Fermi level and the optimized waveguide cross-section dimensions for each pump wavelength, the conversion response for the DFWM process is simulated, as shown in Fig. 6a, where the conversion efficiency is calculated as the signal wavelength varies around the pump wavelength and the pump and signal powers are kept as 0.5 W and 1 mW. Here, a series of pump wavelengths in different bands are considered and the waveguide lengths are selected as their optimal values to support the maximum conversion efficiencies. The corresponding 3-dB bandwidth is shown in Fig. 6b as a function of the pump wavelength λ p . One can find that the DFWM in the graphene-silicon hybrid waveguide shows a broadband performance and the conversion bandwidth increases with the pump wavelength. In telecommunication band near 1.55 μm, a bandwidth of 62.5 nm can be obtained in an 82-μm-long 360 nm × 266 nm hybrid waveguide with the graphene's chemical potential of 0.4 eV. While, in the mid-infrared 2 μm band, the conversion bandwidth reaches 97.7 nm, where a 105-μm-long 470 nm × 346 nm waveguide is needed and the graphene's chemical potential need to be 0.31 eV. The results show that the designed graphene-silicon hybrid waveguide can support broadband nonlinearity enhancement in an ultra-broad wavelength range covering both near-and mid-infrared band.
In summary, a graphene-silicon hybrid waveguide with a dielectric spacer is proposed to enhance the nonlinearity using the high third-order susceptibility of graphene in an ultra-broad wavelength range covering near-and mid-infrared bands. The propagation loss of the hybrid waveguide is effectively reduced by tuning the graphene's chemical potential to meet the resonance condition, and the dielectric spacer blocks the free-carrier exchange to avoid the additional loss. A full-vectorial theoretical model is established to analyze the nonlinear DFWM effect in the proposed hybrid waveguide. The waveguide's cross-sectional dimensions are optimized in terms of the nonlinear parameter and the working wavelength. With a 0.5 W input pump power, the maximum conversion efficiency in the Fermi-level-tuned hybrid waveguide with optimal geometry reaches about −18.5 dB and the required waveguide length is only tens of microns. The conversion bandwidth is as broad as 62.5 nm in 1.55 μm band, and it reaches 97.7 nm in 2 μm band. The proposed hybrid waveguide structure has a stable nonlinear enhancement from 1.3-2.3 μm, which exhibits that the micron-scale graphene-silicon hybrid waveguide has the potential to support chip-scale nonlinear applications in both near-and mid-infrared bands.

Methods
By considering the graphene layer as a thin film with thickness of 0.5 nm 33 , its linear conductivity σ g (1) is obtained from the Kubo formula 32 and the third-order conductivity σ ω ω ω ( , , ) 3 is calculated based on Mikhailov's theory 27 . The simulation results of the effective mode index and the field distribution are found by using the two-dimensional finite-element-frequency-domain (FEFD) method. The coupled equations (11)-(13) is obtained using the reciprocal theorem 36 and solved by a fourth order Runge-Kutta method with a length step Δz = 0.3 μm. The simulation parameters have been given in our paper. Data availability. All data generated or analyzed during this study are included in this article.