1. Introduction

Photonic crystals are artificial materials consisting of a periodic dielectric structure, designed to engineer the characteristics of the light propagating within [1]. They have received quite some interest. On the other hand, nonlinear optical systems promise many interesting all-optical signal and processing applications [2]. When a nonlinearity is included into the photonic crystal structure, the possibility of truly molding light with light opens up. The most considered nonlinear effect is the Kerr effect. Using this nonlinearity, the propagation of modes can be changed by tuning the band-gap dynamically [3, 4, 5]. Moreover, there also other effects to be studied such as the appearance of localized nonlinear modes and gap solitons [6, 7, 8, 9, 10, 11, 12, 13].

Photonic crystals embedded with nonlinear material are an ideal environment to generate and to observe nonlinear localized modes. In [11], Mingaleev et al. consider waveguides created by an array of dielectric rods introduced into an otherwise perfect two-dimensional photonic crystal. In this work, it is assumed that the dielectric constant of the waveguide rods depends on the field intensity due to the nonlinear Kerr effect. Using a nonlinear lattice model, they are able to demonstrate the existence of 2D spatially localized modes with a frequency in the band-gap.

Recently, in [14], some of us demonstrated numerically the existence of self-localized waveguides in a two-dimensional photonic crystal consisting of a square lattice of Kerr-type rods without linear defects. These self-localized waveguides are Bloch modes with frequency in the band-gap. They are confined in the transverse direction due to the band-gap, but they are able to propagate longitudinally as they induce a waveguide in the material by locally reducing the refractive index. Because a self-localized waveguide overcomes the longitudinal band-gap, it can be seen as a variant of the gap soliton [7, 9]. However, they can also be interpreted as a kind of intrinsic localized mode [6, 8] such as discrete solitons in waveguide arrays [10, 12].

Our numerical method is based on a linear mode expansion method [15], where one chooses a main propagation direction and then divides the structure in sections invariant along this direction. The field in one such invariant section can be described as a superposition of eigenmodes. To combine the different sections the mode-matching technique is used. In this way, the fields throughout the entire structure can be found by matrix manipulations. However, things get a bit more complicated when the Kerr nonlinearity is taken into account and the refractive index is dependent on the local field intensity. The nonlinear material in the structure is then divided into smaller sections using a spatial grid. Using the linear eigenmode expansion, the field intensity can be found. Using this calculated field profile the local refractive index is updated. The new field profile can be calculated and the refractive index brought up to date again. In this way, the nonlinearity is taken into account iteratively [16]. Both on-site nonlinear guided soliton modes, modes with a maximum centered on a nonlinear rod, and inter-site modes, have been reported [14].

The properties of photonic crystals and photonic-crystal waveguides are usually studied by numerical integration of Maxwell’s equation. Such calculations can be quite time-consuming and do not always provide good physical insight. Mingaleev et al. show a novel conceptual approach based on effective discrete equations including long-range interaction from site to site, which has proven to be successful for photonic crystal waveguides [17], for nonlinear localized modes [11] and for nonlinear waveguides [18]. The purpose of this work is to suggest a similar semi-analytical approach to describe the properties of self-localized waveguides based on a nonlinear lattice model. By exploiting the clear one dimensional periodicity of the Bloch modes, it is possible to reduce the infinite two dimensional crystal to a strip of material which is finite in the propagation direction and infinite in the orthogonal one. From the Green’s function of the strip, it is possible to find a set of discrete equations which model the interaction between the field centered at the rods of nonlinear material.

The paper is organized as follows. In Section II, we present the nonlinear photonic crystal studied here and review the numerical results found in [14]. In Section III, we study a strip of photonic crystal. This enables us to restrict our study to one-dimensional longitudinal Bloch waves. In Section IV, we introduce effective discrete equations, which model the site-to-site nonlinear interaction of the scattered fields. In Section V, we find the solutions of the nonlinear lattice model and compare these semi-analytical results with the numerical ones found in [14]. Section VI contains conclusions.

2. The model and self-localized waveguides

The structure under consideration is a two-dimensional photonic crystal consisting of a square lattice of parallel, infinitely long dielectric square rods in a homogeneous dielectric with lower refractive index (in this case air). We assume that the rods are oriented along the z-axis. In this sense, the dielectric permittivity can be written as a two-dimensional function ε(r)=ε(x,y). The rods are arranged in a perfect square lattice with lattice constant a (a_x=a x and a_y=a y). The rods have a refractive index of 3.4 and the ratio of the side of the rods to the crystal period is d/a=0.25. The electric field is polarized parallel to the rods (TM-polarization) and propagates in the xy-plane. The wave equation then simplifies to its scalar form. We only consider monochromatic light, E(r)=E(r|ω)e ^-iωt, which reduces the wave equation to the well known Helmholtz equation

This eigenvalue problem can be solved in the case of a perfect linear photonic crystal with permittivity ε(r)=ε_L(r), which is of course a periodic function

with r and s arbitrary integers, and a_x and a_y are the lattice vectors of a primitive unit cell of the photonic crystal.

 

Fig. 1. The band-gap structure of the photonic crystal consisting of a square lattice of square dielectric rods (n=3.4) in an air background. The ratio of the side of the rods to the crystal period is d/a=0.25.

Download Full Size | PDF

 

Fig. 2. A gap soliton with ω=0.38(2πc/a) and k_x=0.7π/a propagating along the x-direction. [Media 1]

Download Full Size | PDF

In the linear case, the band structure of this configuration is shown in Fig. 1. The photonic band-gap extends from the lower frequency ω=0.35(2πc/a) to the higher ω=0.48(2πc/a). Of course, if the frequency of the electric field lies within the band-gap, the field cannot propagate through the structure and will be damped. But, as shown in [14], when the rods contain a Kerr material, this nonlinearity can create guided modes within the band-gap. Such a spatial guided soliton is the linear guided mode of the refractive index profile of the waveguide it induces through the nonlinearity. The Kerr nonlinear refractive index n ₂ has to be chosen negative to induce the waveguide defect and to enable the self-localization to occur. In Fig. 2, an example of such a gap soliton is shown. We have checked that gap solitons propagating in the M-direction (along the structure diagonal) also exist for a negative Kerr index. However, in the following we limit ourselves to propagation in the X-direction.

3. A photonic crystal strip

The self-localized waveguides possess translational symmetry, and the corresponding guided modes can be characterized by the reciprocal wavenumber k_x. Also, such a guided mode has a periodic profile along the waveguide, and it decays exponentially in the transverse direction. Through the nonlinear interaction, the medium loses its periodicity in the y-direction, while retaining it in the x-direction. We will exploit this fact in the analysis.

We will investigate the periodic problem on a strip of the linear photonic crystal. The strip is denoted by the material between $x=[-\frac{a}{2},\frac{a}{2}]$. Let us consider the Bloch wave $E({{\mathbf{r}|\omega }_k}_x,k_x)$ parameterized in k_x. For each k_x∊[0,2π/a], this one-dimensional Bloch wave satisfies:

together with the Bloch condition $E({x+a,y|{\omega }_k}_x,k_x)=E(x,y|{{\omega }_k}_x,k_x)e^{i{k}_{x}a}$. Of course, the original dispersion relation and (two-dimensional) Bloch waves of the entire photonic crystal can be recovered when the Bloch condition in the perpendicular direction is also introduced.

Instead of integrating Eq. (3) numerically, as was done in [14], we will follow a similar approach as in [11, 17, 18] and transform the problem to its integral form using the strip’s Green function. The strip’s Green’s function for a fixed $({{\omega }_k}_x,k_x)$ was proposed in [19] and satisfies:

with ${{\omega }_k}_x$ inside the band gap. The strip Green’s function is related to the whole-space Green’s function G, which is the Green’s function associated to the full two-dimensional photonic crystal, through [19]

By folding the whole space Green’s function onto itself, we reduce our study to the modes which can be considered as one-dimensional Bloch modes with the Bloch scalar k_x. An example of the Green’s function of the entire photonic crystal is given in Fig. 3(a). One can notice the typical evanescent behavior of the Green’s function due to the band-gap. In Fig. 3(b), we show the corresponding strip’s Green’s function. Note that the Green’s function of a perfect linear two-dimensional photonic crystal is symmetric, G(r₁,r₂|ω)=G(r₂,r₁|ω), and also periodic, G(r₁+r a_x+s a_y,r₂|ω)=G(r₁,r₂|ω). It is possible to show that the strip’s Green’s function is periodic, $g\left({\mathbf{r}}_{\mathbf{1}}+r{\mathbf{a}}_{\mathbf{x}}+s{\mathbf{a}}_{\mathbf{y}},{\mathbf{r}}_{\mathbf{2}}+r{\mathbf{a}}_{\mathbf{x}}+s{\mathbf{a}}_{\mathbf{y}}|{{\omega }_k}_x,k_x\right)=g({\mathbf{r}}_{\mathbf{1}},{\mathbf{r}}_{\mathbf{2}}|{{\omega }_k}_x,k_x)$.

4. Effective discrete equations

In this Section, we reintroduce the nonlinearity into the dielectric permittivity as an addition to the linear dielectric constant:

Let us reconsider the Helmholtz problem defined in the strip, where the eigenvalue ${{\omega }_k}_x$ lies in the band gap. The mode profile $E(\mathbf{r}|{{\omega }_k}_x,k_x)$ satisfies

Using the theory of the Green’s function, it is possible to transform this partial differential equation into an integral equation

The nonlinearity is a natural part of the rods which make up the photonic crystal. The nonlinear permittivity will therefore have a similar periodicity as the linear dielectric constant. We model the nonlinearity as

where δ_rod is 1 inside the rods containing nonlinear material and zero outside. The negative sign ensures that the local intensity reduces the refractive index of the rod, and as such, creates a waveguide defect. The electric field is renormalized such that n ₂=1. The size of the nonlinear rods in the photonic crystal is assumed to be sufficiently small so that the electric field can be considered constant inside the rods. We number a center rod as 0. The other rods in the strip are then numbered with positive integers in the positive y-direction and negative in the negative y-direction. An approximate discrete nonlinear equation for the electric field in rod n can be written, when we insert Eq. (10) in Eq. (9) and by averaging over the rods

where J_n,m, the coupling coefficient, describes the interaction between rod n and m. Due to the symmetry and periodicity of the Green’s function and considering the symmetry of the photonic

crystal strip, the coupling constant only depends on the distance between two rods, reducing J_n,m=J|_n-m|:

with r₁ denoting the center of rod l.

 

Fig. 3. (a) Photonic crystal geometry with Green’s function at ω=0.38(2πc/a) superimposed. (b) Strip Green’s function at ${{\omega }_k}_x=0.38(2\pi c⁄a)$ and k_x=0.7π/a calculated using Eq. (6). In both cases r ₂=0.

Download Full Size | PDF

Summarizing the procedure: first of all one chooses the frequency of the nonlinear mode to obtain the Green’s function. This frequency must lie in the range of the band-gap, because in that case the Green’s function has an evanescent form (see Fig. 3(a)) and it is possible to perform the folding procedure from Eq. (6). By choosing the Bloch parameter k_x and reducing the Green’s function to within the strip, we find g, the strip’s Green’s function (see Fig. 3(b)). Averaging the strip’s Green’s function over a rod, delivers us the coupling parameters from Eq. (12), that model the long-range nonlinear interaction between rods. The last step is of course finding the solution of the nonlinear lattice model in Eq. (11). This system is a closed set of nonlinear algebraic equations in the fields at the center of the rods carrying nonlinear material. Techniques for calculating the Green’s function can be found in the literature. We have calculated the Green’s functions using CAMFR (CAMFR simulation software is freely available from http://camfr.sourceforge.net). The nonlinear system can be solved using an iterative Newton-Raphson method.

The coupling parameters defined in Eq. (12) and a similar lattice model have also been presented in [11] for the case of nonlinear localized modes induced by an array of dielectric nonlinear rods embedded in an otherwise perfect linear two-dimensional photonic crystal. Two important differences between that approach and the procedure presented in this work exist. Firstly, in our case every rod of the photonic crystal is nonlinear, while in [11] the nonlinear rods were embedded into the photonic crystal. Secondly, the solutions of the lattice model presented here will always be nonlinear guided modes which are only localized in the transversal direction, while in the examples in [11] the longitudinal band-gap was not overcome. A related approach to finding guided modes in linear waveguide defects have been presented in [17] and for linear waveguides with embedded nonlinear defects in [18]. We have checked that the theory based on the strip’s Green’s function can reproduce the results found on a linear waveguide defect as presented in [17]. The major goal of the work presented in this paper is, like in [18]: to further develop the nonlinear lattice model allowing the use of fast numerical techniques. In [18] a linear lattice model is used to approximate the linear modes (either propagating or evanescent) of a linear waveguide and then construct the field traveling through embedded nonlinear defects in this waveguide as a linear combination of these linear modes. In the next Section, we will study several examples of nonlinear guided Bloch modes using the procedure presented here.

5. Examples of nonlinear guided Bloch modes

5.1. Self-localized waveguides

By solving the set of nonlinear effective discrete equations Eqs. (11) for a chosen set $({{\omega }_k}_x,k_x)$ by means of a Newton-Raphson algorithm, we can find the field amplitudes at the nonlinear rods of the guided spatial soliton Bloch wave. In Fig. 4, we present the calculated mode profile and the corresponding coupling parameters J_n. These modes can be characterized by their modal energy:

The energy, as defined above, is not the total energy of the electric field, but is merely related to it. It can, however, be used as a measure for the total energy. This energy is sufficient to create a waveguide wherein a Bloch wave with frequency ${{\omega }_k}_x=0.38(2\pi c⁄a)$ and Bloch number k_x can propagate. In Fig. 5, we have plotted Q as a function of frequency for a constant Bloch number. As in [14], we are able to obtain both symmetric (on-site) and anti-symmetric (inter-site) modes. These results are compared with the more rigorous numerical simulations of self-induced waveguides of [14] in Fig. 5. The overall qualitative agreement is good, although there is an underestimation of the field on the central rod with the Green’s function approach. However, the semianalytical approach is able to explore a larger k_x domain. The quantitative error can be reduced by replacing the central nonlinear rod, by several smaller ones.

 

Fig. 4. Coupling parameters J_n and electric field at the center of rod n for a self-localized waveguide with ${{\omega }_k}_x$ and k_x=0.7(π/a).

Download Full Size | PDF

5.2. Diatomic photonic crystal

Now, we consider an alternative photonic crystal geometry, a diatomic crystal depicted in Fig. 6. It consists of two offset square lattices with rods of different size. The structure is able to exhibit larger absolute band-gaps (overlapping gap for both polarizations), compared to the normal square lattice [20]. This is because the reduction of symmetry can lead to the canceling of degeneracies at band edges. We will consider the larger rods to be linear, so they impart most of the linear scattering. On the other hand, we assume that the small rods have a positive nonlinearity, so they may be able to induce nonlinear localization.

Because the nonlinear sections are smaller than in the previous section, we expect a better agreement between the strip Green method and the rigorous calculations. All rods in the photonic crystal have a refractive index of √12 (the background is air), and the diameters are, respectively, 0.2a and 0.1a. This crystal has a TM bandgap between ω=0.389(2πc/a) and 0.436(2πc/a). An example of the Green’s and strip Green’s function is shown in Fig. 6. With both the strip Green method and the numerical simulations we find self-localized waveguides in the diatomic photonic crystal. A comparison of the modal energy is shown in Fig. 7. We see a good agreement between the two methods, as expected. However, with the semi-analytic method it is easier to cover the whole k_x-range. An example of a propagating self-localized waveguide is depicted in Fig. 8. Notice how in this case the staggered Green’s function gives rise to an unstaggered mode, which is also shown by the unstaggered strip Green’s function in Fig. 6b. Indeed, if the rods aside from the central column have a larger field than the center rod, and if they are multiplied by the appropriate phase factor (see Eq. (6)), the sign can change.

 

Fig. 5. Modal energy Q and modal electric field amplitudes of the self-localized waveguide at ${{\omega }_k}_x=0.38(2\pi c⁄a)$, calculated with strip’s Green theory (thin line) and exact simulations (thick line), respectively. The field in the center of rod 0 (black), 1 (red), 2 (green) and 3 (blue) are shown.

Download Full Size | PDF

 

Fig. 6. (a) Diatomic photonic crystal geometry with Green’s function at ω=0.4(2πc/a) superimposed. (b) Strip Green’s function at ${{\omega }_k}_x=0.4(2\pi c⁄a)$ and k_x=0.85π/a.

Download Full Size | PDF

 

Fig. 7. Modal energy Q for the diatomic photonic crystal at ${{\omega }_k}_x=0.4(2\pi c⁄a)$, calculated with strip’s Green theory (thin line) and exact simulations (thick line), respectively.

Download Full Size | PDF

 

Fig. 8. A propagating diatomic gap soliton with ${{\omega }_k}_x=0.4(2\pi c⁄a)$ and k_x=0.908π/a. [Media 2]

Download Full Size | PDF

Stationary 2D modes also appear in the diatomic photonic crystal [21]. However, because these solitons are zero-dimensional, such as point-defect modes, their description needs one degree of freedom less. There, one frequency leads to one mode. This is in contrast with our self-localized waveguides, where one frequency leads to a one-parameter family of modes, characterized by the Bloch constant k_x.

6. Conclusion

We have developed a consistent theory of nonlinear solitonlike Bloch waves, inducing their own waveguide into a monoatomic or a diatomic photonic crystal with a Kerr-type material in the rods. We have shown that it is possible to reduce this effectively two-dimensional problem to one dimension by only studying a strip of the photonic crystal. We have demonstrated that these Bloch modes can be adequately described by a nonlinear lattice model that includes the long-range interaction, necessary to describe diffraction and interference, and an effectively nonlocal nonlinear response. This semi-analytical approach was compared with a more rigorous numerical technique [14]. Overall qualitative agreement was found to be good.