1. Introduction

Light propagation in arrays of weakly coupled waveguides has been studied extensively for fundamental reasons as well as for applications to optical processing devices. Light injected in a single waveguide couples to more and more waveguides as it propagates under the form of Floquet-Bloch (FB) waves. Specific features of these waves have been studied in the linear (discrete diffraction, beam-like behavior, non-divergent beams…) and nonlinear (self-focusing, discrete solitons…) regimes [1,2]. They have also been exploited to demonstrate diffraction engineering [3], Bloch oscillations [4,5], quantum mechanical decay [6], Anderson localization [7,8], Rabi oscillations [9], and optical routing schemes [10–12]. Arrays have then become an important tool of optics and as such require a robust though simple assessment basis, which can treat efficiently in the linear and nonlinear regimes the complex patterns already built or envisioned in the future [14–16].

Waveguide arrays are series of identical equally-spaced (period = S) confining optical structures. Propagation in a single isolated structure is permitted for discrete optical indices corresponding to the guided modes. Coupling between N structures splits those into sets of N sublevels indexed by K_x∈[-π/S, π/S] corresponding to “supermodes” encompassing all structures. For an infinite array bands are formed and allow for the free propagation of FB waves across the whole space [3]. Most experimental situations involve mainly the upper band generated around the fundamental level of the isolated structure, all the more when this structure is designed to be monomode. The band shape n(K_x) known as the diffraction relation controls all the features of FB waves: direction of rays, divergence of beams [3], refraction and reflection at interfaces [17]... Describing n(K_x) is then the key objective of any model.

This basis is at present provided by a long-established approximation, the coupled mode theory (CMT), which describes “weak” coupling situations by assuming (1) the slowly varying envelope scheme and (2) the orthogonality of the modes of the quasi-isolated confining structures even for finite S. CMT predicts a cosine-shaped relation $n\left(K_x\right)=n_{irw}+2C/{\beta }_{vacuum}\cos \left(K_{x}S\right)$ where n is the sublevel effective index, n_irw is the effective index of the fundamental mode of the isolated waveguide, C is the coupling coefficient between neighboring waveguides, and β_vacuum=2π/λ. In addition simple equations describe the propagation of the electric field reduced to the amplitudes of the FB wave in each confining structure. If refined theories give a more accurate description of systems such as pairs of asymmetric waveguides [18,19] they become too cumbersome for the much more complex case of the waveguide array so that CMT is used in most studies because of its very simplicity and in spite of its drastic assumptions.

Using the slowly-varying envelope scheme assumes that C is much smaller than the propagation constant of the upper guided mode β_irw. This is satisfied in most systems since C<1mm⁻¹ whereas β_irw~10 µm⁻¹. On the other hand the orthogonality assumption and deviations from the pure CMT calculations have not been assessed explicitly in waveguide arrays, even though they have an important impact on the behavior of FB waves through deformations of the diffraction relation. Furthermore, the design of functional devices requires as short lengths of active waveguides as possible, i.e. tends to maximize the coupling strength. This trend pushes CMT to the limits of its range of validity, which calls for the check whether CMT holds within a required accuracy for a given waveguide design and for the determination of the necessary corrections to CMT when this is no longer the case.

Only a few papers have tried to assess the validity of CMT in waveguides arrays. Eyges and Wintersteiner [20] have compared sets of supermode levels for arrays of N≤6 circular waveguides obtained by CMT and numerical calculation and they reported large discrepancies (~30%) at very high couplings. At C~0.1mm⁻¹ Kaplan and Ruschin [21] found only small discrepancies between output intensities calculated by CMT and observed in experiments. In a recent and stimulating paper Cooper and Mookherjea [22] have renewed the field by developing a simple analytical model which takes non-orthogonality and further-neighbor effects into account to find out the supermodes. With a modal approach they obtain accurate values for the full CMT matrix from numerical data through an inverse-problem calculation and call this approach numerically-assisted CMT.

In this paper, using an analytical model similar to theirs but limited to first and second neighbors, we adopt a direct approach based on propagation equations as is usual in CMT rather than a modal approach. All the parameters – C and the correction factors which express non-orthogonality and second-neighbor effects – can be determined knowing only the mode of the isolated waveguide. This directly yields the corrections to the diffraction relation, which straightforwardly reflects the impact of the parameters. Our extension of CMT to the “moderate” coupling regime is successfully tested against supermode calculations in widely used GaAs, InP, and Si-based shallow-ridge or rectangular waveguide arrays. Scaling the various results demonstrates very general rules which give a clear insight into the corrections to CMT and can help to optimize waveguide arrays. The approach also delimits the “weak” coupling regime unambiguously, which makes it possible to design structures for which CMT applies within a tolerated correction level. We finally discuss the impact of the corrections on experiments and potential devices.

2. Analytical modeling: a first-order extension of the coupled mode theory

The derivation of the extended CMT model is given in detail in Appendix, and is based on the following assumptions:

Our model thus yields as desired an extension of the CMT with next-order terms of clear-cut significance in an extended diffraction relation: extra coupling and non-orthogonality between isolated waveguide modes alter respectively lateral diffraction and propagation and arise respectively in the numerator and the denominator. Long range couplings to successive neighbors and the corresponding harmonics of the diffraction relation can be derived exactly using Wannier functions, i.e. a basis of orthogonal functions localized at each site and built from the exact Bloch functions (see e.g [23].). We have checked that, to first order in η and using the diffraction relation given by Eq. (5), Eqs. (3) are equivalent to the set of equations obtained with Wannier functions (see Appendix). Actually the amplitudes a_m are not expected to depend strongly on the shape of the localized base functions. It may be stressed that all the parameters entering the diffraction relation can be obtained here from the mode of the isolated waveguide, owing to integrations over either the whole transverse plane or the ridge perturbation section so that the more involved computation of the supermodes is not necessary, in contrast to the Wannier function approach [23] or to Ref [22]. It may also be noted that such a perturbative extension of CMT remained hidden in the modal approach [22] where the mode propagation constants are input data for an inverse problem, although the connection between the numerically-calculated modes and CMT relies upon the same integrals as in our propagative approach.

We now validate our extended CMT model on selected examples by reference to direct numerical solutions of the Maxwell equations. More precisely we compare the predicted diffraction relation with the set of sublevels of the supermodes obtained for a finite number of waveguides in the array, which is actually a sample of this diffraction relation.

3. Numerical simulation

3.1 The waveguide arrays used for validation

The validation is performed mainly on prototypal systems, shallow-ridge structures built in the GaAs/GaAlAs and InP/InGaAs systems, widely used in published works on waveguide arrays [24,15] and operated at the telecom wavelength λ=1.55µm; Si/SiO₂ buried-strip waveguides will also be evaluated in Section 4.1. We shall vary only the ridge width L_r and the array period S, which tune respectively the individual ridge waveguide and its coupling to its neighbors. All other parameters will be kept at common literature values (Fig. 1 ).

 

Fig. 1 Description of prototypal systems considered.

Download Full Size | PDF

Three reference optical structures are of interest: the lower planar waveguide (lpw) corresponding to the absence of ridges, the upper planar waveguide (upw) corresponding to complete coverage by merged ridges (S=L_r), and the isolated ridge waveguide (irw) corresponding to the presence of a single ridge (S=∞). Their effective indices (n_cladding < n_lpw < n_irw < n_upw < n_core), together with mode maps, are calculated by direct numerical solution of the Maxwell equations with the finite-element method (FEM) implemented by the COMSOL^® package. We checked that results do not depend significantly on boundary conditions. Due to the monomode structure chosen for the planar waveguides, n_lpw and n_upw take single values. n_irw can take several values depending on the lateral confinement i.e. on L_r. (Fig. 2 ).Results for both systems and both polarizations are similar. Indeed empirically scaling both coordinates merges all curves (Fig. 2(c)). The monomode character is then identical in all cases and is essentially governed by the evanescence length (horizontal scale) and confinement factor (vertical scale) of the fundamental mode of the isolated waveguides (see Eqs. (A2) and A3). Figure 2(c) thus validates the one-mode model for the waveguides of Fig. 1. As expected, they remain monomode as long as n_irw remains below the middle of [n_lpw, n_upw].

 

Fig. 2 Variation of the effective indices of the isolated ridge waveguide n_irw with ridge width L_r for the standard structure in InP (a) and GaAs (b) systems: fundamental mode (green symbols), mode 1 (orange), mode 2 (red), with horizontal polarization in full lines and vertical one in dashed lines. Planar waveguide levels n_lpw and n_upw are shown by dark green horizontal lines and their middle by light green ones. The supermode bands for C=1.3mm⁻¹, L_r=3µm (InP) and C=0.3mm⁻¹, L_r=4µm (GaAs) together with the levels of a 7-ridge array are illustrated by violet symbols. In (c) all curves are shown to coincide when drawn in reduced coordinates.

Download Full Size | PDF

We will consider in the following typical cases, i.e. L_r=1.5, 3, and 4µm for InP and L_r=4µm for GaAs. All the results for these structures, as well as results for Si/SiO₂ buried-strip waveguides will be merged in Section 4 within a conclusive scaled plot which will emphasize the broad applicability of extended CMT.

3.2 The isolated ridge waveguide

Results for all shallow-ridge structures are similar and will be exemplified by the case of the InP-based structure with L_r=3µm for the horizontal polarization. Index reference values are n_lpw=3.22974, n_irw=3.231632, n_upw=3.234355. Figure 3 shows the horizontal section of the fundamental isolated-ridge mode. For all cases calculated the mode tails outside the ridge can be quite well fitted by an exponential decay curve with a decay length L_decay very near to the estimation of the lateral evanescence length $L_{evan}=\sqrt{n_{iw}^2-n_{lp}^2}/{\beta }_{vacuum}$: L_decay=2.20µm while L_evan=2.23µm for the horizontal polarization in Fig. 3. Based on the map of the isolated ridge mode all parameters of the model (C, η, ξ, and ζ) can be determined by numerical integration (Fig. 4 ).As expected by considering only the mode tails, their variation for S above 3µm follows exponential laws with a characteristic length of L_decay or 2L_decay, depending on the order of the parameter. More precisely η can be quite well fitted by the law $\eta \left(S\right)=\left(1+S/L_{decay}\right)e^{-S/L_{decay}}$. Deviations from these laws at high S can be attributed to numerical errors and at low S by the erroneous description of the mode map. However they apply in the useful range of S.

 

Fig. 3 Horizontal section of the fundamental mode of the isolated waveguide around its center, in the InP system with L_r=3µm. Blue (resp. violet) symbols for horizontal (resp. vertical) polarization. O = left linear scale; Δ = right logarithmic scale. Lines are adjustments by exponential decay curves outside the ridge, with L_decay=2.2µm (resp. 1.5µm) for horizontal (resp. vertical) polarization.

Download Full Size | PDF

 

Fig. 4 Model parameters as a function of ridge spacing S for horizontal polarization in InP-based 7-ridge structure with L_r=3µm. Calc, mod, fit in the captions indicate respectively the results of the calculation from the isolated mode, an exponential approximation of this mode, and the values giving the best fit of the FEM levels by the model diffraction relation.

Download Full Size | PDF

3.3 Supermodes of finite arrays

We now introduce the waveguide array by considering a set of N identical ridges with a center-to-center separation S. As stated above this system now exhibits N sublevels around the center level n_irw corresponding to supermodes$a_m{e}^{i{k}_{z}z}$ [14]. Their maps are quite well reproduced by superposition of isolated-ridge modes (see Fig. 5 ) with weights given by the CMT $a_m=\sin \left(m\frac{p+1}{N+1}\pi +\left(p+1\right)\frac{\pi }{2}\right)$ where p is the order of the supermode (p=0 for the fundamental mode). No clear evidence of deviation from mere mode superposition has been observed even for strongly-coupled arrays. Polarization hybridization remains small: for instance the TE mode in dense InP-based arrays with L_r=3µm, S=4µm has only |E_y|/|E_x|<0.05.

 

Fig. 5 X section of the fundamental mode of the 7-ridge array in the InP system (blue circles), with L_r=3µm and S=10µm (a) or S=4µm (b). The blue line is an adjustment by a combination of individual modes (violet lines) merely described by an exponential tail with weights given by the CMT model (red squares).

Download Full Size | PDF

The set of sublevels n(p) indexed by p are represented in Fig. 6(a) and corresponds to a discrete sample of the diffraction relation (Eq. (5). Data show that an excellent sampling is indeed obtained whatever the number of waveguides (see an example in Fig. 6(b)). We have mostly used arrays with N=7 as a compromise: at higher numbers more modes and hence more data are obtained at the expense of simplicity of numerical calculation, while 7 points correctly describe the diffraction relation. These results are a strong indication that the perturbative treatment CMT relies upon can be used to derive adequate diffraction relations.

 

Fig. 6 FEM data on InP-based structures with L_r=3µm. (a) Fan-out diagram, effective indices versus ridge spacing S for horizontal (top) and vertical (bottom) polarizations. Dashed lines show the indices of planar waveguides. (b) Effective index versus reduced supermode wavenumber k_x for S=6µm and various numbers of ridges N (2, 3, 4 to 13), horizontal polarization; the blue curve is a fit by the extended CMT model. (c) Effective index versus supermode number p, for N=7, horizontal polarization, and S values of 15 (deep red dots), 10, 8, 6, 4µm (violet dots). Dashed lines are predictions of Eq. (5) using calculated values for η, ξ, and ζ, but C values increased by a factor 1.5; full lines use η values decreased by the factor 0.7.

Download Full Size | PDF

4. Comparison of analytical modeling to numerical simulation and discussion

4.1 Test of the extended CMT model

Within CMT n(p) is expected to follow a cosine law:

The fit of effective index data by the extended CMT model is two-fold: C controls the height of the diffraction relation i.e. the width of the effective index band while η, ξ, and ζ control its shape. The FEM-calculated width is very well fitted by the model over 3 orders of magnitude in all cases considered, provided that all the values of C obtained from the model are increased by a common factor <1.5 (Fig. 4). This small discrepancy which is actually smaller in other cases (e.g. 1.2 for L_r=4µm) can be traced back to the limitations of a perturbative treatment which makes use of isolated waveguide modes (see Appendix). The shape of the band and especially its change at smaller S is fairly well reproduced by the model (Fig. 6(c), dashed lines) and a nearly perfect fit can be obtained (Fig. 6(c), solid lines) by slightly decreasing η. These small discrepancies can be visualized in Fig. 4 through the small shift between directly calculated and FEM-derived values. A similar fit is obtained for all InP or GaAs shallow-ridge waveguides, using smaller corrections of C and η.

Si-SiO₂ rectangular waveguides at λ=1.55 µm are also considered here for comparison with the work of Cooper and Mookherjea [22], using 200 nm × 500 nm Si core strips (n=3.48) embedded in SiO₂ cladding medium (n=1.45). In such a system the index contrast is much larger than in shallow ridge waveguides and the coupling can get much stronger, thus making the validity of CMT and even extended CMT questionable. Figure 7(a) shows the fan-out diagrams of our results – quite similar to those of Ref [22]. – and Fig. 7(b) a comparison with the prediction of extended CMT. In this case the fit is quite good and does not require any adjustment. Extended CMT then applies also well to those strongly-coupled systems.

 

Fig. 7 Si-based arrays of 5 waveguides. (a) Fan-out diagrams for vertical (V, top) and horizontal (H, bottom) polarizations. (b) Fit of FEM data for vertical polarization with calculated parameters C, η, ξ, and ζ.

Download Full Size | PDF

Finally the excellent overall fit on many various systems validates our extended CMT model and in addition the pertinence of the parameters which describe the strength of the coupling.

4.2 The limit of “weak” coupling – towards a general criterion

We first discuss the case of shallow ridge waveguides. We note that C/β is always <10⁻³, which confirms the validity of the slowly-varying envelope assumption. On the other hand the correction imposed by the non-orthogonality of individual ridge modes (mostly η, which overrides ξ and ζ) becomes negligible (<1%) only at large values of S i.e. at very low and indeed useless coupling constants (<0.01mm⁻¹). Hence this correction should be considered in most practical cases. In this framework the limit of “weak” coupling below which CMT can be used is defined as the deviation to the CMT diffraction relation that can be tolerated. This maximal tolerated deviation corresponds to a maximum value of η, ξ, and ζ, and in turn for a given system to a maximum coupling strength C determined by curves as those shown in Fig. 4, which are obtained using only the fundamental mode map of the isolated ridge waveguide.

A still simpler though less accurate way may be suggested. Since shallow-ridge structures are all similar, reducing curves to dimensionless coordinates should merge them into a single one. The relevant scale for S is L_evan. With the reduced abscissa s=S/L_evan, η(s) and ζ(s) curves – and to a lesser extent the noisier ξ(s) curve – for various L_r values and for InP as well as GaAs systems do coincide (Fig. 8 ). For the C scale, we note that within the CMT the sublevel band of width Δn=4C/β_vacuum must remain within the 2D guidance band [n_lpw,n_upw] (see Fig. 2a,b). We then propose to use as a reference the coupling constant C_max defined by C_max=(β_vacuum/4).min(n_upw-n_irw,n_irw-n_lpw). Figure 8 shows that C/C_max(s) curves for all systems also coincide very well. These data can be correctly approximated by the formulas derived from our extended CMT model for narrow ridges (Eqs. (A10-)A13): C/C_max ~8e^-s, η ~(1+s)e^-s, and ξ ~ζ ~2e^-s.

 

Fig. 8 Variation of η, ξ, ζ, and C/C_max with reduced ridge period s=S/L_evan for all waveguide structures considered. Only values validated by FEM data are shown. Lines indicate the predictions of Eqs. (A10) and A13.

Download Full Size | PDF

Within this “general” rule, prediction becomes easy. For a given structure one calculates the indices of the isolated ridge and planar waveguides and deduces L_evan and C_max. For a desired maximal correction, the “universal” η(s) yields the minimal s=S/L_evan and hence the minimal S. Finally this minimal s yields the maximal C/C_max and hence the maximal C. For given planar structure (layer stack) and deviation from CMT, if C has to be maximized, the best choice for L_r is to maximize C_max by bringing n_irw near the middle of [n_upw, n_lpw]. For the planar structures considered here, C_max is limited to ~2.3 and 3.9 mm⁻¹ for InP in respectively the horizontal and vertical polarizations, and 0.6 mm⁻¹ for GaAs in the horizontal polarization.

Data obtained for Si-SiO₂ waveguides in vertical (TE type) polarization which involve much shorter evanescence lengths (<0.2 µm) and much stronger confinement factors fit also quite well in this scheme (Fig. 8). Such an overall agreement is somewhat surprising considering the polarization hybridization reported [22]. This effect – which is nevertheless unexpected and has to be taken with caution due to the strong singularities occurring at the corners of the strips – does not seem to prevent satisfactory predictions of interactions between neighboring strips, at least as long as S is not too small. Thus Fig. 8 emphasizes the “universal” character of the above rules which may still be used with rather strongly coupled systems where C reaches 2 µm⁻¹ and C/β about 0.2.

4.3 Distortion of the diffraction relation at strong couplings – impact on experiments

On this basis we now describe the distortions of the diffraction relation with respect to its CMT cosine shape (Fig. 9a ) with decreasing s at the moderately small values of s where the model has been validated – at lower values the extended model is no longer valid as attested by the apparent divergence of its diffraction relation. The most distorted part of the curve lies at high k_x where a very strong change of negative k_z takes place. At low k_x the change appears only as a moderate drop of the whole curve; indeed if curves are normalized to take the CMT value of 2 at k_x=0 distortion is barely visible (Fig. 9b).

 

Fig. 9 (a) Variation of diffraction relation for various values of the reduced inter-ridge spacing s. (b) same with k_z normalized to take the value of 2 at k_x=0.

Download Full Size | PDF

Hence in most experiments involving only low k_x the main effect will be a reduction of the coupling constant C to an effective one C_eff. Many literature studies have been performed in the GaAs system [14]. For instance a value of C=0.75mm⁻¹ is reported for the planar structure we used above with L_r=4µm and S=8µm, and is applied to diffraction management studies [3]. Our assessment of this structure on the lines described above leads to C=0.58mm⁻¹ but s=1.4 which gives C_eff=0.37mm⁻¹. In our own study of “channels” performed on InP-based arrays involving the planar structure we used above with L_r=1.5µm and S=5µm, the published value of 1.5mm⁻¹ was obtained by fits of n(p) but without taking distortion into account. The correct value is C=1.32mm⁻¹ but s=1.24 so C_eff=0.81mm⁻¹. Those examples show the importance of corrections in such moderately coupled arrays. However in both cases since only low k_x were mainly investigated qualitative agreement was obtained.

If on the other hand large k_x are used, the asymmetry of the diffraction curve can lead to important deviations from the CMT behavior. This happens for instance in refraction experiments involving oblique interfaces or following direct injection in high-k_x states [25,26]. The distortion can lead to considerable errors on the direction of FB wave packets (beams) which is controlled by the derivative of the diffraction relation and on their divergence which is controlled by its second derivative, if one considers only the low-k_x value C_eff. In Ref [3]. the authors observe a strong deviation from CMT at k_x above π/2 that they attribute to the contribution of the second band but that could be explained as well by mere high-coupling distortion. Finally in arrays which involve a patterning of the coupling constant the contrast between neighboring high-C and low-C zones is reduced by the distortion. In the case of channels [15] – high-C stripes surrounded by large low-C zones – the C contrast 2.83 reduces to a C_eff contrast of 2.58, still comfortably large enough to insure superguiding.

5. Conclusion and perspectives

In conclusion the present perturbative treatment of arrayed isolated waveguides when they are increasingly coupled to their neighbors yields a first-order extension of CMT. Such a treatment based on a propagative approach gives a simple analytical expression of the propagation constant, which naturally gets rid of the high-order polynomial dependencies involved in the modal approaches. The associated diffraction relation exhibits the distinct consequences of coupling to next-nearest neighbors on lateral diffraction and of non-orthogonality between nearest neighbors on propagation, respectively through two transfer integrals and one overlap integral. The distortion of the shape of the diffraction relation with respect to the CMT shape at increasingly stronger coupling situations is thus nicely accounted for. Quantitative agreement with fully numerical simulation is obtained even nearby the strong coupling regime for various waveguide structures, requiring only slight adjustments of the ab initio parameters. In a broad range of weakly and moderately coupled systems the distortion of the diffraction relation can be assessed using the non-orthogonality overlap integral as an indicator, so that the designer can confidently avoid excessive distortions. Alternatively he can take advantage of fairly simple empirical formulas to improve propagation and diffraction management in the design of arrayed devices.

Appendix – Analytical derivations

In the above discussion numerical resolution of the Maxwell equations for shallow ridge waveguide arrays is compared to more straightforward results derived from approximate models. As shown on the left hand side of Fig. A1 the array is considered as the periodical juxtaposition of planar waveguide portions.

 

Fig. A1 Schematic decomposition of the dielectric constant profile of shallow ridge waveguide arrays in vacuum in (X,Y) plane. Light blue and dark blue dashed vertical lines on the left hand side indicate interfaces where boundary conditions must be written respectively in “etched” and “ridge” planar waveguide regions.

Download Full Size | PDF

Since the dielectric constant is piecewise constant, the Maxwell equations can be reduced to the Helmholtz equation:

(A1)$${\nabla }^{2}E\left(X,Y,Z\right)+\frac{{\omega }^2}{c^2}\epsilon \left(X,Y\right)E\left(X,Y,Z\right)=0$$

The solutions are superpositions of planar modes$g_m^{R,E}\left(Y\right)$, m=1 … (taken as real) in the “ridge” or “etch” (R,E) regions, which have to be connected at the region boundaries using continuity equations. In the following we assume that only the fundamental planar modes with propagation constants ${\beta }_1^{R,E}$ need to be considered, either in TE or in TM polarization (${\beta }_1^R$ and ${\beta }_1^E$ are denoted β_upw and β_lpw in the text). This is quite legitimate for shallow ridge waveguides in which the core layer is separated from the etched mesas by a thick cladding layer, and provided the mesas exhibit good verticality in order to minimize polarization mode mixing. Then the confined mode of the isolated waveguide and the modes of the array – the super-modes – are built using two leftward and rightward propagating waves, with vertical profiles given by the two fundamental planar modes in either region. The horizontal profile of the fundamental (symmetric) mode of an isolated waveguide can thus be cast into the form derived for three-layer planar waveguides [27]:

(A2)$$\left|x\right|\le \frac{w}{2},E\propto \cos pxx>\frac{w}{2},E\propto e^{-q\left(x-\frac{w}{2}\right)}x<-\frac{w}{2},E\propto e^{q\left(x+\frac{w}{2}\right)}$$

where w is the waveguide width, and q=1/L_evan obeys:

(A3)$$r_{R/E}q\frac{w}{2}=p\frac{w}{2}\tan p\frac{w}{2}\text{and}{p}^2+q^2={\beta }_1^R{}^2-{\beta }_1^E{}^2$$

r_R/E being an averaged value of the ratio of the dielectric profiles, either equal or very close to 1 depending on the polarization.

On the right-hand side of Fig. A1, the array is decomposed into the usual perturbation scheme of an isolated waveguide plus a dielectric perturbation ${\displaystyle \sum _n\delta {\epsilon }_n\left(X,Y\right)}$ brought about by the mesas of all the other waveguides. In this approach the Helmholtz equation is solved using a linear superposition of the fundamental modes M₁ of all the isolated waveguides:

(A4)$$E\left(X,Y,Z\right)={\displaystyle \sum _m{a}_m\left(Z\right)M_1\left(X-ms,Y\right)}{e}^{i\beta Z}$$

The Helmholtz equation is not strictly equivalent to the Maxwell equations because of the abrupt changes of the dielectric constant around the perturbation by the mesas, and the approach is valid only for small perturbations. Another way to envision its limitations is to recognize that expansion of a supermode over the isolated waveguide modes Eq. A4 cannot satisfy all the boundary conditions at the mesa edges because the evanescent tails of modes m’≠m do not at mesa m; this would require a deeper discussion of the polarization properties of the modes, which is beyond the scope of the present analysis. Using the definitions of C and of the overlaps by Eqs. (1) and (2) and approximating up to the second neighbor, the perturbation approach [27] yields:

(A5)$$\begin{array}{l}\frac{\partial a_m}{\partial Z}+\frac{\partial a_{m-1}}{\partial Z}\frac{〈01〉}{〈00〉}+\frac{\partial a_{m+1}}{\partial Z}\frac{〈01〉}{〈00〉}+\frac{\partial a_{m-2}}{\partial Z}\frac{〈02〉}{〈00〉}+\frac{\partial a_{m+2}}{\partial Z}\frac{〈02〉}{〈00〉}=\\ C\left\{a_{m-1}+a_{m+1}+a_m\frac{2\cdot 〈101〉}{〈001〉}+a_{m-2}\frac{〈\overline{1}01〉+〈002〉}{〈001〉}+a_{m-2}\frac{〈\overline{1}01〉+〈002〉}{〈001〉}\right\}\end{array}$$

Since the mode has exponential tails, $〈01〉/〈00〉$,$〈101〉/〈001〉$,$〈\overline{1}01〉/〈001〉$, and$〈002〉/〈001〉$ are first order while$〈02〉/〈00〉$ is second order with respect to evanescent coupling. At 0th order one finds the CMT result. At 1st order, denoting z=CZ and using definitions of η, ξ, and ζ by Eqs. (4), one finds:

(A6)$$\frac{\partial a_m}{\partial z}+\eta \left(\frac{\partial a_{m-1}}{\partial z}+\frac{\partial a_{m+1}}{\partial z}\right)=i\left(a_{m-1}+a_{m+1}+\xi a_m+\zeta a_{m-2}+\zeta a_{m+2}\right)$$

The diffraction relation of FB waves ${\text{e}}^{i{k}_{x}x+i{k}_{z}z}$ then becomes:

(A7)$$k_z=\frac{2\cos \left(k_x\right)+\xi +2\zeta \cos \left(2k_x\right)}{1+2\eta \cos \left(k_x\right)}$$

To first order in η, the propagation equation Eq. A6 is identical to what can be derived from the approach based on Wannier functions [23] when the diffraction relation obeys Eq. A7.

Now, using Eq. A3 and the profile of Eq. A2, we calculate the overlap integrals and the extended CMT parameters C, η, ξ, and ζ, and for instance:

(A8)$$\eta =\frac{{\text{e}}^{-qS}}{q}\times \frac{1+q\left(S-w\right){\text{e}}^{qw}+\frac{4q}{p^2+q^2}{\text{e}}^{q\frac{w}{2}}\left[q\text{\hspace{0.17em} sh}q\frac{w}{2}+p\tan p\frac{w}{2}\text{\hspace{0.17em} ch}q\frac{w}{2}\right]}{\frac{w}{2}+\frac{1}{2q}+\frac{\sin pw}{2p}+\frac{\cos qw}{2q}}\times {\cos }^{2}p\frac{w}{2}$$

(A9)$$C={\text{e}}^{-qS}{\cos }^2\frac{pw}{2}\times \frac{\frac{2}{p^2+q^2}{\text{e}}^{\frac{qw}{2}}\left[q\text{sh}\frac{qw}{2}+p\tan \frac{pw}{2}\text{\hspace{0.17em} ch}\frac{qw}{2}\right]}{\frac{w}{2}+\frac{1}{2q}+\frac{\sin pw}{2p}+\frac{\cos qw}{2q}}\times \frac{{\beta }_{vacuum}^2}{2{\beta }_1}\left({\epsilon }_1-1\right){\displaystyle \underset{\text{mesa}\text{edge}}{\int }dY{g}_1^E{}^2}$$

At low ridge widths, formulas reduce to simple expressions:

(A10)$$\eta \underset{qw\to 0}{\to }\left(1+qS\right){\text{e}}^{-qS},\xi \underset{qw\to 0}{\to }2{\text{e}}^{-qS},\zeta \underset{qw\to 0}{\to }2{\text{e}}^{-qS}$$

(A11)$$\frac{C}{C_{max}}\underset{qw\to 0}{\to }{\text{e}}^{-qS}\times 8r_{R/E}\frac{{\epsilon }_1-1}{{\epsilon }_1^R-{\epsilon }_1^E}\times {\displaystyle \underset{\text{mesa}\text{edge}}{\int }dY{g}_1^E{}^2}$$

where ${\epsilon }_1^R={\beta }_1^R{}^2/{\beta }_{vacuum}{}^2$, ${\epsilon }_1^E={\beta }_1^E{}^2/{\beta }_{vacuum}{}^2$ (in the main text n_upw=$\sqrt{{\epsilon }_1^R}$ and n_lpw=$\sqrt{{\epsilon }_1^E}$). Taking the planar waveguide in the etched zone as an unperturbed system, a perturbative treatment of the Helmholtz equation for the planar waveguide in the ridge zone straightforwardly gives to first perturbation order:

(A12)$${\beta }_1^R{}^2={\beta }_1^E{}^2+{\displaystyle \underset{\text{mesa}\text{edge}}{\int }dY{g}_1^E\left({\epsilon }_1-1\right)g_1^E}$$

so that assuming r_R/E=1:

(A13)$$\frac{C}{C_{max}}\underset{qw\to 0}{\to }8{\text{e}}^{-qS}\equiv 8{\text{e}}^{-S/L_{evan}}$$

The above treatment can also be applied to the case of buried strips [22] though in this case a modal method through the periodical slit and strip arrangement might be more effective. Equivalent roles can be played between the above “ridge” and “etched” regions and, respectively, the “strip” and “aperture” regions. A difficulty arises in the aperture regions which are uniform and do not support localized modes. This difficulty can be circumvented by symmetrically adding new strips to the system in the aperture regions, far from the original strip array in the ±Y directions. Such strips give rise to extra confined modes, the highest of which can be both degenerate and weakly coupled with the planar mode for the strip regions, and hence serve as a planar mode for the aperture regions. As a perturbative approach shows, the new system will exhibit similar strip array modes as the original one near the planar mode for the strip region, due to weak evanescent coupling of the latter with the extra modes. Eqs. A10 and A13 are then also valid estimates for the fundamental band in buried strip arrays as long as higher order planar modes can be neglected.

Acknowledgements

These results are within the scope of C’nano IdF; C’nano IdF is a CNRS, CEA, MESR and Région Ile-de-France Nanosciences Competence Center.

References and links

1. D. N. Christodoulides, F. Lederer, and Y. Silberberg, “Discretizing light behaviour in linear and nonlinear waveguide lattices,” Nature 424(6950), 817–823 (2003). [CrossRef]   [PubMed]  

2. J. Fleischer, G. Bartal, O. Cohen, T. Schwartz, O. Manela, B. Freedman, M. Segev, H. Buljan, and N. Efremidis, “Spatial photonics in nonlinear waveguide arrays,” Opt. Express 13(6), 1780–1796 (2005), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-13-06-1780. [CrossRef]   [PubMed]  

3. H. S. Eisenberg, Y. Silberberg, R. Morandotti, and J. S. Aitchison, “Diffraction management,” Phys. Rev. Lett. 85(9), 1863–1866 (2000). [CrossRef]   [PubMed]  

4. T. Pertsch, P. Dannberg, W. Elflein, A. Brauer, and F. Lederer, “Optical Bloch oscillations in temperature tuned waveguide arrays,” Phys. Rev. Lett. 83(23), 4752–4755 (1999). [CrossRef]  

5. R. Morandotti, U. Peschel, J. S. Aitchison, H. S. Eisenberg, and Y. Silberberg, “Experimental observation of linear and nonlinear optics Bloch oscillations,” Phys. Rev. Lett. 83(23), 4756–4759 (1999). [CrossRef]  

6. F. Dreisow, A. Szameit, M. Heinrich, T. Pertsch, S. Nolte, A. Tünnermann, and S. Longhi, “Decay control via discrete-to-continuum coupling modulation in an optical waveguide system,” Phys. Rev. Lett. 101(14), 143602 (2008). [CrossRef]   [PubMed]  

7. T. Schwartz, G. Bartal, S. Fishman, and M. Segev, “Transport and Anderson localization in disordered two-dimensional photonic lattices,” Nature 446(7131), 52–55 (2007). [CrossRef]   [PubMed]  

8. Y. Lahini, A. Avidan, F. Pozzi, M. Sorel, R. Morandotti, D. N. Christodoulides, and Y. Silberberg, “Anderson localization and nonlinearity in one-dimensional disordered photonic lattices,” Phys. Rev. Lett. 100(1), 013906 (2008). [CrossRef]   [PubMed]  

9. K. G. Makris, D. N. Christodoulides, O. Peleg, M. Segev, and D. Kip, “Optical transitions and Rabi oscillations in waveguide arrays,” Opt. Express 16(14), 10309–10314 (2008), http://www.opticsexpress.org/abstract.cfm?URI=oe-16-14-10309. [CrossRef]   [PubMed]  

10. D. N. Christodoulides and E. D. Eugenieva, “Blocking and routing discrete solitons in two-dimensional networks of nonlinear waveguide arrays,” Phys. Rev. Lett. 87(23), 233901 (2001). [CrossRef]   [PubMed]  

11. A. Fratalocchi, G. Assanto, K. A. Brzdakiewicz, and M. A. Karpierz, “All-optical switching and beam steering in tunable waveguide arrays,” Appl. Phys. Lett. 86(5), 051112 (2005). [CrossRef]  

12. R. A. Vicencio, M. I. Molina, and Y. S. Kivshar, “Switching of discrete optical solitons in engineered waveguide arrays,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 70(2), 026602 (2004). [CrossRef]   [PubMed]  

13. A. L. Jones, “Coupling of optical fibers and scattering in fibers,” J. Opt. Soc. Am. 55(3), 261–269 (1965). [CrossRef]  

14. F. Lederer, G. I. Stegeman, D. N. Christodoulides, G. Assanto, M. Segev, and Y. Silberberg, “Discrete solitons in optics,” Phys. Rep. 463(1-3), 1–126 (2008) (and references therein). [CrossRef]  

15. N. Belabas, S. Bouchoule, I. Sagnes, J. A. Levenson, C. Minot, and J. M. Moison, “Confining light flow in weakly coupled waveguide arrays by structuring the coupling constant: towards discrete diffractive optics,” Opt. Express 17(5), 3148–3156 (2009), http://www.opticsexpress.org/abstract.cfm?URI=oe-17-5-3148. [CrossRef]   [PubMed]  

16. J. M. Moison, N. Belabas, C. Minot, and J. A. Levenson, “Discrete photonics in waveguide arrays,” Opt. Lett. 34(16), 2462–2464 (2009), http://www.opticsinfobase.org/ol/abstract.cfm?URI=ol-34-16-2462. [CrossRef]   [PubMed]  

17. A. Szameit, H. Trompeter, M. Heinrich, F. Dreisow, U. Peschel, T. Pertsch, S. Nolte, F. Lederer, and A. Tünnermann, “Fresnel’s laws in discrete optical media,” N. J. Phys. 10(10), 103020 (2008). [CrossRef]  

18. A. Hardy and W. Streifer, “Coupled mode theory of parallel waveguides,” J. Lightwave Technol. 3(5), 1135–1146 (1985). [CrossRef]  

19. W. P. Huang, “Coupled-mode theory for optical waveguides: an overview,” J. Opt. Soc. Am. A 11(3), 963–983 (1994). [CrossRef]  

20. L. Eyges and P. Wintersteiner, “Modes of an array of dielectric waveguides,” J. Opt. Soc. Am. 71, 1351–1360 (1981), http://www.opticsinfobase.org/abstract.cfm?URI=josa-71-11-1351.

21. A. Kaplan and S. Ruschin, “Characterization and performance evaluation of coupled multiwaveguide arrays,” J. Lightwave Technol. 17(10), 1884–1889 (1999), http://jlt.osa.org/abstract.cfm?URI=JLT-17-10-1884. [CrossRef]  

22. M. L. Cooper and S. Mookherjea, “Numerically-assisted coupled-mode theory for silicon waveguide couplers and arrayed waveguides,” Opt. Express 17(3), 1583–1599 (2009), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-17-03-1583. [CrossRef]   [PubMed]  

23. G. L. Alfimov, P. G. Kevrekidis, V. V. Konotop, and M. Salerno, “Wannier functions analysis of the nonlinear Schrödinger equation with a periodic potential,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 66(4), 1046608 (2002). [CrossRef]  

24. H. S. Eisenberg, Y. Silberberg, R. Morandotti, A. Boyd, and J. S. Aitchison, “Discrete spatial optical solitons in waveguide arrays,” Phys. Rev. Lett. 81(16), 3383–3386 (1998). [CrossRef]  

25. T. Pertsch, T. Zentgraf, U. Perchel, A. Brauer, and F. Lederer, “Anomalous refraction diffraction in discrete optical systems,” Phys. Rev. Lett. 88(9), 093901 (2002). [CrossRef]   [PubMed]  

26. D. Mandelik, H. S. Eisenberg, Y. Silberberg, R. Morandotti, and J. S. Aitchison, “Band-gap structure of waveguide arrays and excitation of Floquet-Bloch solitons,” Phys. Rev. Lett. 83, 4752–4754 (1999).

27. B. E. A. Saleh, and M. C. Teich, “Fundamentals of Photonics,” Wiley-Interscience, 2007.