Broadrange tunable slow and fast light in quantum dot photonic crystal structure
Lotfian Alireza1, †, Yadipour Reza1, Baghban Hamed2
Department of Electrical and Computer Engineering, University of Tabriz, Tabriz 5166614761, Iran
School of Engineering-Emerging Technologies, University of Tabriz, Tabriz 5166614761, Iran

 

† Corresponding author. E-mail: Lotfian@tabrizu.ac.ir

Abstract

Slow and fast light processes, based on both structural and material dispersions, are realized in a wide tuning range in this article. Coherent population oscillations (CPO) in electrically tunable quantum dot semiconductor optical amplifiers lead to a variable group index ranging from the background index (nbgd) to ∼ 30. A photonic crystal waveguide is then dispersion engineered and a group index of 260 with the normalized delay-bandwidth product (NDBP) of 0.65 is achieved in the proposed waveguide. Using comprehensive numerical simulations, we show that a considerable enhancement of slow light effect can be achieved by combining both the material and the structural dispersions in the proposed active QDPCW structure. We compare our developed FDTD results with analytical results and show that there is good agreement between the results, which demonstrates that the proposed electrically-tunable slow light idea is obtainable in the QDPCW structure. We achieve a total group index in a wide tuning range from nbgd to ∼ 1500 at the operation bandwidth, which shows a significant enhancement compared with the schemes based only on material or structural dispersions. The tuning range and also NDBP of the slow light scheme are much larger than those of the electrically tunable CPO process.

1. Introduction

Light velocity control has become one of the most important topics in optical engineering. “Slow light” refers to a significant decrease in light velocity whereas “fast light” refers to its velocity greater than the speed of light in vacuum (c). Today, optoelectronics plays the most important role in modern communication and information technologies, whereas there still exist some challenges such as the ability to store an optical signal in the optical domain without converting to the electrical format or in brief the optical buffering.[1] By realizing optical buffers, there would be no need to optical/electrical/optical conversions[2] and as a result, routing delay will significantly be reduced due to the complexity of protocols, size, and power of the routers. Data packets' collision is one of the still remaining challenges in optical networks. To solve this, one needs a controllable amount of delay or advance. Therefore, we are interested in tunable speed for light signal. Tuning ability is a critical feature of slow light schemes and in general, the range in which light speed can be tuned is of great importance. Light signal velocity is the same as group velocity and has been defined by vg = ∂ω/∂k which can be expressed by [cω∂n(k,ω)/∂k]/[n(k,ω)+ω∂n(k,ω)/∂ω], where n is the real part of the refractive index and k is the propagation constant. Two terms in the group velocity relation are capable of controlling the group velocity. The first term is ∂n(k,ω)/∂k (structural resonances), also called waveguide dispersion which utilizes the frequency-dependent phase variation created by strong transmission or reflection characteristics, which results in a constant group index that cannot be tuned actively by the bias current. The second term, ∂n(k,ω)/∂ω (material resonances), also called material dispersion, utilizes the refractive index change created by a narrow resonance in the absorption or gain spectra of the medium. The later can be controlled and then tuned electrically or optically. On the other hand structurally engineered slow light demonstrates a very larger group index than that from the material dispersion method and then, is very effective in slow light generation. Initial experiments were based on electromagnetically induced transparency (EIT) to achieve a large group index (low group velocities). In an EIT scheme, a narrow transmission window is created due to the destructive interference between a strong pump beam and a weak probe one, therefore, a large group index can be obtained over this transmission window. Using EIT, the light velocity decreased to c/165 in a 10-cm lead vapor cell.[3] Later, by utilizing sodium atoms at low temperatures, group velocity decreased to 17 m/s.[4] Coherent population oscillations (CPO) is also utilized in an alexandrite crystal to realize a low group velocity of 91 m/s at room temperature.[5] When a strong pump field and a weak signal beam co-propagate in a semiconductor structure, the beating between these two fields causes the carrier density to oscillate. These oscillations then lead to the scattering of pump beam into the signal beam which in turn results in a refractive index change for the signal. Several schemes have been proposed to achieve a large delay-bandwidth product (DBP) including stimulated Raman scattering[6] and stimulated Brillouin scattering.[7,8] Also, the generated light is delayed compared with the reference using a four-wave mixing experiment in a hot rubidium atomic vapor.[9] On the other hand, structural resonances in photonic crystals (PCs) have been widely used to achieve the slow light effect. Due to the photonic band gap (PBG) property of PCs, the slow light phenomenon is investigated theoretically and experimentally near the Brillouin zone edge of the photonic crystal waveguide (PCW).[10,11] In addition, the problem of pulse broadening due to group velocity dispersion (GVD), reduces the peak intensity of a pulse[12] and thus makes it less effective for nonlinear applications. In general, a trade-off has to be made between the bandwidth and effective group index of the slow light.[13] Therefore, it is critical to linearize the dispersion relation of PCW over a large range, which corresponds to the regime of slow light with a wide band and low GVD, but more importantly, we also need to enhance the slowdown factor (ng) of slow light as much as possible. Very recently, a number of methods have been proposed to optimize the slow light features of PCWs, such as modifying the waveguide width,[14] varying the hole size,[15] adjusting the hole position,[16] changing the hole shape,[17] and varying the hole refractive index through infiltrating a liquid.[18] For appropriate optimization, the realization of slow light with enhanced group index, enhanced NDBP, and reduced GVD is required.

The aim of this article is to use both mechanisms simultaneously in the same device to achieve a large group index in a wide tuning range. This device should be compatible and easy to be integrated with other optical subsystems and should be also GVD-engineered and loss-engineered for practical applications. In Section 2 of this paper, a theoretical model is presented for slow light generation based on three different schemes including optimized PCW, QD medium, and the novel QDPCW structures. The PCW is also dispersion and loss-engineered using elliptical air holes and adjusting some structural parameters of the waveguide. In Section 3, the group index, bandwidth, and NDBP of the optimized PCW are calculated using the finite difference time domain (FDTD) method. Slow and fast light generations based on the phenomenon of CPO in the QD semiconductor optical amplifier (QD-SOA) are demonstrated using rate equations coupled with the density matrix formalism. Electrically tunable slow and fast light results are achieved, then the combined waveguide/material dispersion system is discussed and compared with the two individual systems in Section 3. Finally some conclusions are drawn in Section 4.

2. Theoretical model

In the following subsections, a theoretical model is developed to investigate the effect of both waveguide and material dispersions on the total group index.

2.1. Optimized PCW slow light

The 2D PC slab has been widely used for slow light generation due to the easy fabrication process and some excellent characteristics. This PC slab includes two types: the air holes in the dielectric slab and the dielectric rods in the air. The former exhibits a transverse electric (TE) band gap and the latter shows a transverse magnetic (TM) band gap.[19] The TM-like and TE-like polarizations can approximately be distinguished by observing whether an electric field is perpendicular to or parallel with the surface of the wafer, respectively. The waves whose frequencies are within the photonic band gap (PBG) cannot propagate through the structure, and as a result when these waves are incident on the PC, they will be rejected back due to the absence of the corresponding mode. The most important parameter for slow light generation is the group index ng which can be written as[20] where k is the wave vector and vg is the group velocity in the PC waveguide. Assuming that the group index is much larger than the background index of the material, equation (1) can be written in the form of ng = a/2π(dk/df). Therefore, we could expect low dispersion slow light (nearly constant ng) in specific frequency regions where f varies linearly with k.

Ellipse shaped holes are probably one of the easiest air hole shapes for practical fabrication. The scheme showing this optimized PCW can be seen in Fig. 1. The letters A and B denote the lengths of the holes along the x and y-axis, respectively. Due to the deformation of the elliptical unit cell, the periodic symmetry of the structure breaks and waves which satisfy the Bragg criterion do not resonate in various directions. This results in less energy leaking in turn. Therefore, an improvement in power transmission of the PCW consisting of the elliptical unit cell is observed in comparison with the circular unit cell.

Fig. 1. (color online) Optimized PCW with ellipse-shaped air holes and structural optimization parameters.
2.2. CPO slow light

Realization of the CPO phenomenon in the QD medium can be described as follows. In the presence of an intense optical pump field with the photon energy adjusted to the transition energy of a two-level system, saturation of absorption occurs since the population of the lower state depletes. If a weak and slightly detuned signal beam is present as shown in Fig. 2(a), the excited state population will oscillate at the frequency determined by the pump-signal detuning. Significant population beating occurs if the population can closely follow the intensity profile due to the time-dependent interference between the optical fields of the pump and the signal beams. The coupling between the strong beam and the oscillation of the population results in a reduction in the absorption of the probe beam as shown in Fig. 2(b), which then leads to a rapid index variation around the resonant frequency through the Kramers–Kronig relations. If frequency of the probe beam fits in the spectral hole, the probe beam will propagate subluminally along the medium.

Fig. 2. (color online) (a) Band diagram in the presence of a resonant pump and a detuned signal, (b) Top: absorption spectrum of the signal in the absence (dash curve) and presence (solid) of a strong pump. The spectral hole (dotted curve) is caused by the CPO produced by the pump and probe beams. Bottom: the corresponding refractive index spectrum in the presence of CPO.

Coupled rate equations have been used for studying and calculating the optical properties and carrier dynamics of QD SOAs.[21,22] In our model, electrical current is injected into the barrier and the pump and probe beams are applied to the QD ground states as shown in Fig. 2. In the self-assembled quantum dots, carriers are injected into the barrier states. Then, they are captured by the wetting layer (WL) followed by relaxing to the excited state (ES) and the ground state (GS). Then, the carrier population dynamics of the barrier states , the wetting layer , the excited state , and the ground state are described by the following rate equations (2)–(5).

Here, J is the injection current density, e is the unit charge, and A is the active region cross section. , , , and are the electron occupation probabilities, respectively, in the barrier, WL, ES, and GS, which are defined by in which is the number of carriers in the barrier, WL, ES, and GS, and is given by . The terms Γ and L are the optical confinement factor and the cavity length, respectively. The optical power and photon energy of the pump and signal beams are Pp(s) and ħωp(s), respectively. The transition time constants of carriers are denoted by where the superscript c(v) represents the electron (hole) band and the subscript xy refers to the transition from state x to state y (e.g., determines the relaxation time from the excited state to the ground state). The τdr and τwr are the recombination time constants in dot and wetting layer, respectively. After the occupation probabilities are obtained from the rate equations, we use the density matrix formalism to calculate the susceptibility, absorption, and refractive index of the medium. The dynamic equations of the density matrix elements associated with the transition between the two energy levels of conduction and valence bands, shown in Fig. 2(a), can be expressed as where ρcc and ρvv are the state populations in the conduction and valence band, respectively; ρcv are the corresponding off-diagonal density matrix elements; ercv = |Mb|· |Menv| is the corresponding dipole moment; E(t) is the optical field containing different frequency components. The Rabi frequency corresponding to the frequency component a and defined as , is the angular frequency difference between the two states. The subscript inc means the incoherent process and Γ2,cv is the dephasing constant. We need the off-diagonal density matrix elements because they are related to the linear signal response which is directly connected to the susceptibility of the signal. The linear polarization density P(s)(ωs) of the signal can be evaluated as where GE and GW are the probability density functions which are assumed to have Gaussian distributions for simplicity. The integration is done to take into account the effect of inhomogeneous broadening. Then we can calculate the susceptibility χ(s)(ωs), through the following expression:

Then, the absorption coefficient, α(ωs), and the corresponding variation of the refractive index, Δn(ωs), of the signal beam can be obtained from the following equations: where Γconf is the confinement factor. The group index can be defined as

2.3. Active QDPCW structure

We consider the optimized PCW slab geometry with QDs embedded in the core region of the structure as illustrated in Fig. 3. The active region of this waveguide consists of a number of self-assembled InAs/GaAs QD stacked layers where QDs are covered with an In0.15Ga0.85 As layer with a thickness of 5 nm which serves as a strain reducing capping layer. The sheet density of QDs is 500 μm−2. The cavity length is assumed to be 500 μm and the photonic crystal slab thickness supports the single mode operation of the waveguide. Increasing the thickness of the capping layer and indium composition can red shift the resonance wavelength of QDs to communication wavelengths. In this structure, the electrical tuning is implemented by passing a current, using p–n contacts, through a laterally doped p+−p−n+ structure which can be formed by ion implantation.[23] In fact, to minimize the voltage drop across the PC cladding regions, these regions are doped heavily with p+-type and n+-type impurities.

Fig. 3. (color online) Schematic diagram of the active QDPCW system. The purple shading illustrates the QDs in the waveguide region.

As mentioned before, our approach involves the tailoring of both the material and structural dispersions. The engineering of the material dispersion curve involves the using of the semiconductor QDs under the CPO mechanism. The material dispersion at the probe wavelength is affected by the optical transition of the carriers of QDs. When the intense optical pump depletes the population of the QD ground state (GS), the absorption of the probe beam with frequency of ωs is saturated (transitions from the valence band GS to detuned conduction band GS). Since the pump and probe beams are slightly detuned, the population of the GS will oscillate. This population beating results in absorption reduction of the probe beam (a hole in the absorption spectra). By the Kramers–Kronig relation, the induced absorption reduction leads to a rapid index variation (the ∂n/∂ω term in the group velocity relation, which is known as the material dispersion).

On the other hand, the mechanism behind the designed waveguide dispersion assigned to the periodic PC structure can be explained as the fact that the holes form the boundary of the waveguide. The waveguide is narrow where there is a hole and is wider in the absence of holes. The described narrow-wide-narrow periodic geometry forms a Bragg mirror. Provided that the Bragg condition is fulfilled (i.e., λ/2 = na with a being the period of the structure), a standing wave is created. Otherwise, the light coherently scattered by the mirror plane forms an interference pattern along with the incoming light which moves slowly (proportional to ∂n/∂k). This term is also known as the structural dispersion.

Since the structural dispersion depends both on the geometry and on the background index and considering the fact that the material dispersion is proportional to the CPO process in QDs, both the material and structural dispersion affect the total index of the structure. It is difficult to calculate the group index analytically in this combined structure since both material and structural resonance are utilized simultaneously, as in general the waveguide dispersion features will be affected by the material dispersion. The group velocity can be calculated in such a PCW structure by assuming linear constitutive relations for the material. The H, as a harmonic magnetic eigenmode, is the solution to equation ω2/c2 = 〈H,ΘH〉/〈H,H〉 in which ΘH = ∇ × ε(r,ω)−1∇ × H. The group velocity (the derivative of ω with respect to k) can be expressed as[24] where /∂k represents the derivative with respect to k regarding a fixed dielectric function ε. Also d/dk refers to the derivative which involves the dependency of ε with respect to frequency. Since Θ is Hermitian with H being its eigenstate, one can write

Here, the Hellman–Feynman theorem has been considered. Regarding Maxwell equation ∇ × H = −i(ω/c)D and the Hermiticity of the curl operator and considering a low-loss dielectric, one may obtain the group velocity as where is the PC waveguide group velocity without considering the material dispersion. Also, Ed = 〈D,EQD/〈D,E〉 denotes the amount of electric energy in the dielectric region and is calculated from Ed = ∫ dtD,EQD/∫ dtD,E〉. Using ε = n2, the partial derivative of ln(ε) can be written in terms of , and the total group index may be expressed as

It is clear that the total group index depends linearly on the product of and . Also, the QD ensemble group index, , and the filling factor, Ed, affect the total group index.

3. Simulation results and discussion

By simply adjusting the structural parameters in the PCW structure, flat and wideband slow light with a large ng is obtained. The group velocity dispersion also becomes zero since it benefits an S-shaped dispersion curve which prevents the signal from being distorted during propagation. A constant group index changing from 41 to 265 is obtained by optimizing the major and minor axes of the bulk elliptical holes and adjusting their positions and sizes in the first row adjacent to the defect (see Fig. 1). Figure 4 shows the relationship between f and k for the TE even mode whereas parameters A and B are adjusted for optimization. The linear region of each curve corresponds to a flat slow light region with a specific ng determined by its slope.

Fig. 4. (color online) Relationship between f and k of PCW by changing the optimization parameters.

As depicted in Fig. 5, A, B, and y1 can be adjusted properly to yield high NDBP (ngΔλ/λ = 0.577) with the bandwidth (Δλ) of about 17.2 nm and the group index of ng = 52 at λ = 1550 nm. It is noteworthy to mention that the narrower the bandwidth, the larger the achieved group index is.

Fig. 5. (color online) Calculation results of the optimized PCW structure consisting of elliptical holes, where sizes (A1, B1) and positions (y1) of the holes in the first rows adjacent to the core region are changed.

For the CPO-based slow light in OQ-SOA, we assume that the probe signal is a small signal. First, the coupled rate equations are solved to determine the bias condition of the SOA. Then, using density matrix formalism, small signal analysis for the signal beam is done and the CPO effects on the refractive index, absorption, and group index spectra are demonstrated and plotted in Fig. 6. Since we are interested in a large and tunable group index in a large bandwidth, we should control the depth and width of the dip in the absorption spectrum to achieve a tunable group index and then a tunable delay/advance. In a practical system, electrical control of delay/advance is preferable due to easy implementation. For this purpose, the SOA bias current is varied for investigating its effect on the group index and the resulting time delay. If electron–hole pairs are injected into QDs, they will saturate the background absorption. Thus, the depth and width of the absorption dip created by CPO is limited and reduced and the corresponding variation of the refractive index is also reduced which then can be utilized to tune the slowdown factor (ng). As illustrated in Fig. 7, the group index of the signal is varied in the range of around 70 from the background refractive index of 3.4 by tuning the input current source value. Figure 7 also shows switching between fast (negative group index) and slow light by changing the SOA bias condition from the gain to the absorption region. By increasing the bias current, more carriers enter into the dot states which results in reducing the absorption and then, increasing the gain for the probe signal. So, a dip in the gain spectrum yields the fast light effect. A negative group index yields a condition where the group velocity is greater than the speed of light in vacuum.[25] Such a condition is achieved when the bias current of the SOA increases as shown in Fig. 7. In such a regime, the pulse appears to leave the medium at exactly the same time as it enters. It should be noted that fast light propagation is entirely within the bounds of causality and special relativity.[26]

Fig. 6. (color online) Variation of refractive index (top) creation of absorption dip (middle), and resulting group index (bottom) whereas nbgd = 3.4.
Fig. 7. (color online) Controllable group index via CPO, both slow light and fast light (negative group index) are simultaneously achieved in QD-SOA by changing the bias condition of the SOA.

It should be mentioned that the magnitude of the slow light and fast light increase by increasing the α factor of SOA. This can be understood by noting that the population pulsations in SOA create gain and index gratings. The pump and signal traveling through SOA scatter energy through these gratings. When the linewidth enhancement factor is very large, the contribution due to index gratings dominates which results in a large group index variation for the signal.

To investigate and compare the total group index, , of the QDPCW structure according to both material and structural dispersions, three different schemes of the passive PCW structure, the QDs taking part in the CPO configuration, and the active QDPCW structure are discussed. To further investigate the obtained results of the group index calculations, we note the signal pulses whose spectral width coincide with the absorption peak. Figure 8 demonstrates the group index of the pulses which propagate in the PCW structure, QD-SOA CPO system, and the QDPCW structure considering different applied intensities of the pump Ip = cεnbgd|Ep|2. It is clear from Fig. 8 that the PCW group index is constant and does not depend on the power intensity. As discussed before, increasing the number of carriers reduces the group index of the CPO system and saturates and reaches to the background index nbgd at high intensity. The total group index of the active QDPCW system, which tracks both the structural and material resonances, increases dramatically as Ip is lowered and it also presents a dependency just similar to what the mere CPO system presents. It is demonstrated that the group index increases significantly for the QDPCW scheme in comparison with for the pure CPO or PCW subsystem. As shown in Fig. 8, a small variation in the pump intensity causes a large variation in the group index of the joint system in comparison with of the pure CPO system. Thus, by utilizing both the structural and material resonances in the novel structure, it is possible to dramatically increase the group index, leading to a tunable group index active control in a very small waveguide structure. The prediction of Eq. (18) is also plotted for Ed = 0.535 and in Fig. 8, in which an acceptable accordance is achieved between the numerical results and the theory. As predicted before, the novel structure has the advantages of wide range tuning and very high group index in comparison with the single CPO or PCW schemes. Active control of the group index is not possible for the pure PCW scheme.

Fig. 8. (color online) Group index variations versus pump power intensity for the active QDPCW (pink crosses), analytic CPO in QDs (red dotted line), and the PCW structure (black dash line). The blue dash curve presents the results obtained from Eq. (18) with Ed = 0.535.

To give an acceptable justification for the pulse propagation in a diffractive medium such as the PCW structure, a full vectorial study of the electromagnetic fields is required. Therefore, a semi classical approach is used here for light-matter interaction where the FDTD method[27] is combined with the density matrix method to investigate pulse propagation in the QD slow light medium. In this case, propagation of the electromagnetic wave is affected by the microscopic polarization. Furthermore, the Liouville equation is numerically solved on FDTD grids considering the grid points of the QD carrier dynamics. The polarization density then enters into the Ampere’s law through the calculation of the expectation value of the macroscopic polarization. For achieving the FDTD results, at first, the group index caused by the CPO process ( ) is calculated in a specific bias condition through solving the coupled rate equations along with the density matrix formalism. The corresponding is then substituted as the new background index into the FDTD simulation of the PC structure. The boundary conditions required for the FDTD simulation should guarantee the reduction of reflections from the computational window edge. Hence, we use perfectly matched layer (PML) boundary condition surrounding the structure. The utilized parameters in the FDTD simulations include grid sizes of Δx = Δz = a/32 (with a being the period of the structure). The thickness of the PML is assumed to be 500 nm. Then, simulation is repeated for each bias current to obtain the results presented in Fig. 8.

For the CPO based slow light the tuning range is very small (from nbgd to ~ 30). It can be seen (from Fig. 8), as predicted by Eq. (18), that the slow light realized in the active QDPCW can be tuned in a wide range from nbgd to ∼ 1500. The background index of the material which is assumed to be 3.4 increases due to both the structural and material dispersions (as discussed before) and according to relation vg = c/ng, tunable group index ng of 3.4 ∼ 1500 results in slow propagation of light with the corresponding group velocity. It should also be mentioned that the tunability feature of the combined system arises from tunability of the CPO scheme and the NDBP factor is significantly enhanced with respect to the pure tunable CPO system.

4. Conclusion and perspectives

A novel scheme is proposed for light speed control utilizing both material and waveguide dispersions simultaneously in the electrically pumped active QDPCW. We introduce a loss-engineered slow light waveguide which reduces the propagation loss and GVD by an appropriate waveguide design. A flat band slow light PCW is designed for the 1.55-μm operating wavelength, and a theoretical model is developed for the active QDPCW. It is found that the length and dimension of the active QDPCW can be highly integrated due to the significantly enhanced figures of merit of this novel structure. In brief, the slow-down factor ng, NDBP, GVD, and other slow light parameters may be significantly improved compared with those of the systems utilizing only material or structural resonances. This could be very important for realizing the efficient and tunable control of propagation of the pulses in compact semiconductor devices.

Reference
[1] Tucker R Ku P Chang-Hasnain C 2005 IEEE J. Lightwave. Tech. 23 4046
[2] Zhang B Yan L S Yang J Y et al. 2007 IEEE Photon. Tech. Lett. 19 1081
[3] Kapasi A Jain M Yin G Y et al. 1995 Phys. Rev. Lett. 74 2447
[4] HauL V Harris S E Dutton Z et al. 1999 Nature 397 594
[5] Bigelow M Lepeshkin N Boyd R W 2003 Science 301 200
[6] Okawachi Y Foster A Sharping E et al. 2006 Opt. Express 14 2317
[7] Yu C Yan L Willner S P et al. 2007 IEEE Photon. Tech. Lett. 19 861
[8] Okawachi Y Bigelow S Sharping E et al. 2005 Phys. Rev. Lett. 94 153902
[9] Dong D S Zhou Z Y Shi B S 2013 Chin. Phys. 22 114203
[10] Minkov M Savona V 2015 Optica OSA 2 631
[11] Xing M X Zheng W H Zhou W J et al. 2010 Chin. Phys. Lett. 27 024213
[12] Sharma A Kumar M 2015 IET Optoelectron. 40 2666
[13] Casas-Bedoya A Husko C Monat C et al. 2012 Opt. Lett. 37 4215
[14] Schulz S A O’Faolain L Beggs D M et al. 2010 J. Opt. 12 104004
[15] Kurt H Üstin K Ayas L 2010 Opt. Express 18 26965
[16] Liang J Ren L Yun M et al. 2011 Appl. Opt. 50 98
[17] Wan Y Fu K Li C H et al. 2013 Opt. Commun. 286 192
[18] S Y Zhou J L Zhang D 2010 Chin. Phys. Lett. 27 034205
[19] Joannopoulos J D Johnson S G Winn J N 2008 Photonic crystals: molding the flow of light Princeton
[20] Ma J Jiang C 2008 IEEE J. Quantum Electron. 44 763
[21] Rostami A Lotfian A Yadipour R et al. 2013 J. Opt. Commun. 34 1
[22] Puris D Schmidt L C Ludge K et al. 2012 Opt. Express 20 27265
[23] Ellis B Mayer M Shambat G et al. 2011 Nat. Photon. 5 297
[24] Lægsgaard J Bjarklev A Libori S E 2003 J. Opt. Soc. Am. 20 443
[25] Chang-Hasnain C Chuang S L 2006 J. Lightwave. Tech. 24 4642
[26] Garrison C Mitchell M Chiao R et al. 1998 Phys. Lett. 245 19
[27] Taflove A Hagness S C 2005 Computational Electrodynamics Boston Artech House