Negative refraction in periodic and random photonic crystals

We present a theory of the artificial magnetic activity in a two-dimensional (2D) dielectric photonic crystal. We show, by using a multiple-scale approach, that if the rods constituting the crystal have Mie resonances at large enough wavelengths, the crystal can be characterized by an effective permeability exhibiting anomalous dispersion. The theoretical results are checked numerically for periodic and random photonic crystals.


Introduction
Metamaterials are artificial structures that can be described, in certain ranges of frequencies, by homogeneous parameters ε and µ presenting unusual properties not encountered in natural materials. In particular, it is possible to design devices in such a way that both their effective permittivity and permeability will be negative, resulting in an effective negative index, as studied years ago by Veselago [1]. Several geometries have been proposed, notably the Splint Ring Resonators structure, that can produce artificial magnetism, associated with a wire mesh structure that can produce a negative epsilon [2]- [12]. In [13], O'Brien and Pendry also suggested that it was possible to obtain an artificial magnetic activity in photonic crystals made of dielectric rods with a strong permittivity. In that case, there are internal resonances in the rods at wavelengths which are large compared with the period of the structure. Around the resonances, the magnetic field is strongly localized inside the rods: this produces a loop of displacement current inducing a microscopic magnetic moment. All the microscopic moments add up to a collective macroscopic moment, resulting in artificial magnetism. The numerical and experimental studies concerning these metamaterials have demonstrated many interesting properties, but we feel that a synthetic theory is still lacking, in particular concerning the artificial magnetic activity. The formulas that have been obtained generally result from an asymptotic study of the scattering matrix of one resonator associated with the Clausius-Mossotti formula. Here, we try to give a more global approach to the problem of artificial materials, by analysing the effective properties with a multiple-scale approach. For the sake of simplicity, we study the photonic crystal suggested by Pendry and O'Brien: the metamaterial is a two-dimensional (2D) photonic crystal made of parallel rods. We obtain a very general expression for the effective permeability, which does not rely on a specific geometry. Next, we study numerically the influence of random variations in the parameters of the rods constituting the metamaterial. Our aim is to quantify the influence of disorder on the negative index property. Two kinds of random variations are investigated: on the positions of the rods and on their permittivity.

Homogenization of the system
Throughout the paper, we assume a time dependence of exp(−iωt), and the field is p-polarized (the magnetic field reads H = u(x, y)e z ). The structure under study is basically a 2D photonic crystal depicted in figure 1, with basic cell Y , the cross-section of a rod is denoted by D (it is circular on the figure, but we develop the theory for an arbitrary, connected, cross-section). The rods have a permittivity of ε r and they are embedded in a host medium of permittivity ε e . The magnetic field satisfies the following equation: div(ε −1 ∇u) + k 2 u = 0, where k = 2π/λ is the wavenumber in vacuum and ε is equal to ε r inside the rods, to ε e in the host medium and to 1 outside the structure. In order to derive the effective parameters of the photonic crystal, we assume that the wavelength λ is large with respect to the period d and we denote η = d/λ. The permittivity ε r of the rods is supposed to be large and the topological radius r is small compared to the wavelength. The resonances depend mainly on the optical diameter r × √ ε r , therefore, it can be assumed that the existence of a resonance is invariant with respect to a scale transformation. More precisely, if the radius is reduced and the permittivity increased by keeping the optical diameter constant, we expect the physics to be the same. Let us therefore scale the radius by η: r −→ ηr; then if we scale the permittivity by ε r −→ ε r /η 2 , the resonances should not move.  Let us check this hypothesis numerically. To do so, we consider only one rod and we plot its transmission spectrum (i.e. the flux of the Poynting vector through a segment situated below the rod) when we rescale the radius and the permittivity. We start with ε = 200 + 5i and R = 1/4 and we compute the transmission for η = 1, 3/4, 1/2 and 1/4, cf figure 2. It is clearly seen that the first Mie resonance remains unchanged.
We are now in a position to perform a scaling analysis of the medium: the permittivity of the rods is now ε η = ε r /η 2 and their radius r η = ηr. The field satisfies div(ε −1 η ∇u η ) + k 2 u η = 0. We denote χ η = ε −1 ∇u η . The Helmholtz equation is then equivalent to the following system: The field now depends upon the small parameter η, and our point is to find the limit of u η when η tends to 0. In order to do so, we make the following multiple-scale ansatz [14,15]: where the fields (u j , χ j ) depend on two variables: the macroscopic variable x = (x 1 , x 2 ) and the fast variable x/η. These fields are periodic with respect to the second variable: The limit field u 0 depends on the microscopic variable y which is an internal degree of freedom and the physical macroscopic field u eff is obtained by averaging over Y : Let us plug the expansions (3) and (4) into Helmholtz equation. The derivation operators become The identification of the terms in 1/η and 1/η 2 leads to the following relations: div y χ 0 = 0; On Y \D: On D: We obtain directly that u 0 does not depend on y on Y \D. More information can be gained from the transmission condition. Indeed, we have where the superscript ± denotes the inner and outer traces of the field on ∂D and n is the outer normal of ∂D. The multiple-scale expansions leads to the following relations on ∂D: By taking the divergence of equation (12), we get the microscopic equation satisfied by u 0 on D: The microscopic behaviour of the magnetic field u 0 is ruled by the following problem: Indeed, denoting we have the relation Let us now derive the macroscopic effective equation. By taking the divergence of equation (10), we get, with equation (9): to which we have to add the boundary values on ∂D obtained in (15). This shows that ∇ y u 1 depends linearly on ∇ x u 0 . More precisely, let us introduce the following problems: where j = 1, 2 and n = (n 1 , n 2 ). From these functions, we can form the matrix: and we have the relation ∇ y u 1 = A(y)∇ x u 0 . The averaging of equation ( and by averaging equation (10) on Y \D, we get 6 Institute of Physics ⌽ DEUTSCHE PHYSIKALISCHE GESELLSCHAFT the effective permittivity turns out to be Inserting (25) into (24), we get finally the macroscopic equation: The equation clearly shows that the effective medium has both an effective permittivity ε eff and an effective permeability µ eff . In order to find a close form for the effective permeability (19), we have to solve problem (18). To do so, we introduce the following spectral problem: This leads to a family of eigenvalues E p and associated eigenfunctions (possibly degenerate) ψ pl (the index l accounts for the degeneracy). Then, we can solve for m(y) by expanding it on the Hilbertian basis {ψ pl } by noting that m − 1 = 0 on Y \D: By plugging this expression into the Helmholtz equation in (18), we get and An explicit expression for µ eff will be derived in the next section for the case of circular rods. Let us however remark that, near a resonance, the effective permeability is proportional to −k 2 ε r /(k 2 ε r − E p ). Provided that the imaginary part of ε r is not too important, the sign of the real part of µ eff changes when crossing the resonance, therefore there exists a band of wavelengths over which the permeability is negative. Let us now turn to some numerical illustration of this theory.

Numerical results
The numerical study of metamaterials can be done by means of temporal domain or Fourier domain methods. The temporal domain methods (i.e. FDTD) are sometimes a bit unstable when strong resonances occur. We turn to modal methods in the Fourier domain. A 2D photonic crystal can be modelized by a stack of diffraction gratings, which is satisfactory for the study of periodic structures, but is not very suited when some disorder is introduced. We use both a multiple scattering method for parallel rods [18] and a rigorous grating code.

Test of the homogenization result
In this section, we test the theoretical results with a structure made of a stack of three lamellar gratings (period d) made of square rods with permittivity ε = 200 + 5i and side a = d/2. The eigenfunctions are nm (y) = 2 sin(nπy 1 ) sin(mπy 2 ) and the corresponding eigenvalues are k 2 nm = π 2 (n 2 + m 2 ). The expansion of m on this basis leads to the following effective permeability: wherek 2 nm = k 2 nm /a 2 ε i . Besides, we compute the effective permittivity to be ε eff = 1.7. The homogenized structure is thus a slab of height h = 3d with electromagnetic parameters ε eff and µ eff . In figure 3, we plot the real and imaginary parts of the transmission for both the real structure and the homogenized one. A very good fit is observed, except for the resonance at λ/d = 6.5. This is due to the fact that this resonance has a null mean, and therefore vanishes in our approach.
Let us now turn to the study of the negative refraction in photonic crystals. We start with a periodic medium.

Periodic resonators
We use a closely packed set of rods with permittivity −1, in which the resonant rods are embedded periodically. A basic cell is made of nine rods; the rod at the centre has a permittivity of 200 + 5i (cf [13]) and the boundary rods have a permittivity of −1. These rods simulate a homogeneous medium with negative permittivity. All the rods are embedded in air. For circular resonators, the above theory leads to the following spectral problem: − ψ = Eψ for r < R and ψ = 0 for r = R. (32) Let us denote by χ nm the mth zero of the nth Bessel function (i.e., J n (χ nm ) = 0), then the eigenvalues of the above operator are E nm = (χ nm /R) 2 . These eigenvalues are twice degenerate (except for n = 0) and the associated normalized eigenfunctions are ψ ± nm = |πR 2 J n (χ nm )| −1 J n (χ nm (r/R))e ±inθ . It is obvious that 1|ψ = 0 for n = 0, therefore the only remaining resonance is the fundamental one: where k 2 r = χ 2 00 /ε r R 2 and χ 00 ∼ 2.405. The effective permittivity is equal to ε eff = −1.1 (see [16]). We have plotted µ eff in figure 4. Near the resonances, µ eff shows anomalous dispersion and becomes negative. Therefore, the medium can be described by homogeneous parameters that are both negative around the resonances: our theory predicts that it is a left-handed medium near resonances, where the product ε eff × µ eff is positive. We plot the transmission spectrum in figure 5. Near λ/d = 9, the transmission increases, indicating that the waves can propagate (the sharp peak near λ/d = 12 is linked to a leaky mode). Here, there is a mixing of both properties of negative permittivity and negative permeability unlike the wire mesh case [17]. This can be understood by the fact that only the tunnel effect (i.e., the coupling by evanescent waves) is necessary in order to excite resonances inside the rods. The effective magnetic activity exists without any propagating modes, contrary to the case of a negative permittivity, which is the   result of destructive interferences. In a left-handed medium, negative refraction occurs, i.e., an incident beam is refracted on the same side as the incident field. This is due to the fact that the group velocity and the phase velocity are opposed [19]. In order to test the left-handedness of the medium, we use an incident Gaussian beam, with mean angle of incidence 45 • at λ/d = 9. We plot the field amplitude on two lines, one above the structure and one below (cf figure 6) as well as the map of the Poynting vector in figure 7. It is clearly seen that the beam propagates negatively inside the structure, confirming our prediction.

Random resonators
From the experimental point of view, it might be quite difficult to realize a perfect structure, with identical scatterers and a true periodicity. It is therefore interesting to know to what extent the above results still hold with a certain amount of disorder. We investigate the consequences of two kinds of disorder: first, the resonators are all identical but are settled randomly inside the host medium; second, the resonators are disposed periodically but their permittivity is random.
For the random disposition case, we use 26 resonators that are settled randomly among the closely packed rods, with a uniform law inside the photonic crystal. The transmission is averaged over 50 realizations, ensuring a sufficient convergence. We have plotted the averaged transmission in figure 5. It is seen that there still are propagative waves near the resonance. The random disposition of the scatterers has only a small influence on the physical properties of the medium. This can be rather easily foreseen in the theoretical framework developed above. Indeed, the disorder can be taken into account by using the supercell trick: we take a piece of the random medium and use it as a supercell to construct a new periodic medium, to which periodic homogenization applies. The supercell is schematized in figure 8, it contains N resonators D j embedded in a host medium. By periodic homogenization, we get that the auxiliary problems are now multiple: for each resonator D j we have to solve: This leads to N spectral problems like (28) but the resonators being identical, all the eigenvalues are identical (they do not depend on the position of the resonators). Consequently, the average permittivity is exactly the same as obtained in the purely periodic case. This result can be understood in terms of the quality factor of the resonators. It is very high, therefore the coupling between resonators is quite small and consequently the position of the resonators is not important. However, the randomness of the positions of the resonators has an important influence on the propagation of a beam. The structure is illuminated by a Gaussian beam, with the average angle of 40 • and the wavelength is λ/d = 9. The map of the amplitude of the Poynting vector is given in figure 9. For clarity, we have not included the map inside the photonic crystal. Although there are strong scattering events inside the crystal, there still is an indication of a negative refraction. A more detailed study would be necessary to draw firm conclusions, but it seems that the presence of disorder does not prevent the medium from being left-handed.
We turn now to a geometrically periodic structure, in which we have introduced some disorder in the permittivity of the rods. We choose a uniform distribution in the interval [0.8ε r , 1.2ε r ]. The transmission averaged over 50 realizations is given in figure 5. It is seen that the statistical dispersion in the radius makes the left-handed property disappear (the dip near λ/d = 11 is due to the presence of a leaky mode). Indeed, as the resonators all have different resonance frequencies, the local polarizations do not add to build up a global polarization. This effect is all the more important as the statistical dispersion is large. If the distribution law has a sharp maximum near one value, then the medium remains left-handed. However, this indicates that the experimental realizations of such a medium would require a good control over the statistical dispersion of the resonators' size.

Conclusion
We have proposed a general approach to the problem of artificial magnetism in a 2D photonic crystal. In order to produce magnetic activity, we have used rods with a very high permittivity, because they support Mie resonances at wavelength large with respect to the period of the crystal. For the sake of simplicity, we have limited our study to 2D structures but rods of finite length can be studied by the same approach. However, the method can be extended to more general structures (Split Ring Resonators, polaritonic or plasmonic rods). Practically, the realization of a device cannot be perfect, hence the expected properties should present some stability with respect to random deviations from a periodic structure. This is why we have studied the influence of disorder. Firstly, on the position of the resonant rods and secondly on their optical radius. We have shown that the disorder on the position of the rods was not really relevant, but that, in contrast, the statistical dispersion on the radius of the resonant rods could provoke the disappearance of the negative index effect, as noted already in [20]. This shows that the realization of a metamaterial working in the infrared or terahertz range would require a very stable process to make the resonators. It is therefore highly desirable that the resonators be of the most simple geometry and cast doubt on the true usefulness of Omega or Split Ring Resonators at small scales.