Level spacing statistics for light in two-dimensional disordered photonic crystals

TL;DR

This paper investigates the statistical distribution of eigenfrequency spacings in 2D disordered photonic crystals, revealing Poisson statistics for localized modes within the band gap and Wigner-Dyson statistics for extended modes outside the gap.

Contribution

It provides the first detailed analysis of level spacing statistics for light in disordered 2D photonic crystals, distinguishing between localized and extended modes.

Findings

Localized TM modes in the band gap follow Poisson distribution.

Extended TM and TE modes follow Wigner-Dyson distribution.

Disorder induces localization and affects spectral statistics.

Abstract

We study the distribution of eigenfrequency spacings (the so-called level spacing statistics) for light in a two-dimensional (2D) disordered photonic crystal composed of circular dielectric (silicon) rods in air. Disorder introduces localized transverse-magnetic (TM) modes into the band gap of the ideal crystal. The level spacing statistics is found to approach the Poisson distribution for these modes. In contrast, for TM modes outside the band gap and for transverse-electric (TE) modes at all frequencies, the level spacing statistics follows the Wigner-Dyson distribution.

Tables

Table 1. Table 1: Ranges of frequencies ω ​ a / 2 ​ π ​ c 𝜔 𝑎 2 𝜋 𝑐 \omega a/2\pi c defined in figure  5 .

Polarization	TM	TE
I	0.1–0.28	0.35–0.55
II	0.28–0.29
III	0.29–0.41	0.55–0.66
IV	0.41–0.43
V	0.43–0.5	0.66–0.83

Equations

$$P_{\text{WD}}(s)=\frac{\pi s}{2}\exp \left(-\frac{\pi s^2}{4}\right).$$

$$P_{\text{P}}(s)=\exp (-s).$$

$$\left[{\nabla }^2+\frac{{\omega }^2}{c^2}ϵ(\mathbf{\rho })\right]E_z(\mathbf{\rho })=0,$$

$$\nabla \cdot \left[\frac{1}{ϵ(\mathbf{\rho })}\nabla H_z(\mathbf{\rho })\right]+\frac{{\omega }^2}{c^2}{H}_z(\mathbf{\rho })=0.$$

$${\psi }_{n,\mathbf{k}}(\mathbf{\rho })=u_{n,\mathbf{k}}(\mathbf{\rho })\exp (i\mathbf{k}\mathbf{\rho }),$$

$${\mathrm{IPR}}_{n,\mathbf{k}}=\frac{{\int }_A{\left|u_{n,\mathbf{k}}(\mathbf{\rho })\right|}^4{d}^2\mathbf{\rho }}{{\left[{\int }_A{\left|u_{n,\mathbf{k}}(\mathbf{\rho })\right|}^2{d}^2\mathbf{\rho }\right]}^2},$$

$${\xi }_{n,\mathbf{k}}=\frac{1}{2}\sqrt{\frac{1}{{\mathrm{IPR}}_{n,\mathbf{k}}}},$$

$$s_{n,\mathbf{k}}=\frac{\Delta {\omega }_n(\mathbf{k})}{⟨\Delta {\omega }_n(\mathbf{k})⟩},$$

$$\mathcal{P}(s)=\underset{0}{\overset{s}{\int }}P(s^{\prime })d{s}^{\prime },$$

Full text

Introduction

The optical response of any dielectric medium can be calculated from the eigenmodes $\psi_{n}$ and eigenfrequencies $\omega_{n}$ of Maxwell equations in the medium, with appropriate boundary conditions. In a disordered medium where the dielectric function $\varepsilon$ varies randomly in space, $\omega_{n}$ are random numbers that should be analyzed in statistical terms. Whereas the probability distribution of $\omega_{n}$ can strongly depend on the particular medium under study, the distribution $P(s)$ of normalized nearest-neighbor eigenfrequency spacings $s_{n}=\Delta\omega_{n}/\langle\Delta\omega_{n}\rangle$ often takes one of several universal forms independent of details of the medium. Here $\Delta\omega_{n}=\omega_{n}-\omega_{n-1}$, the eigenfrequencies are assumed to be ordered in the ascending order ($\omega_{n}\geq\omega_{n-1}$), and the angular brackets $\langle\cdots\rangle$ denote ensemble averaging. $P(s)$ first came into light for energy level spacings of complex nuclei [1, 2] and is therefore called the “(nearest-neighbor) level-spacing distribution” in the physics literature whenever one speaks of either quantum-mechanical systems or eigenfrequencies of Maxwell or other wave equations [3, 4, 5, 6]; we adopt this terminology here as well. For light or any other wave in a disordered medium, the shape of $P(s)$ depends on the spatial extent of the modes $\psi_{n}$. For weak disorder, the modes $\psi_{n}$ are extended in space over the whole disordered sample and the system is classically chaotic. The normalized eigenvalue spacings are then expected to obey the Wigner-Dyson (WD) distribution [3, 4, 5, 6]:

[TABLE]

This equation turns out to apply to many types of classically chaotic systems, such as, e.g., chaotic billiards [6]. Its most important feature is the so-called level-repulsion phenomenon: $P_{\textrm{WD}}(s)\to 0$ for $s\to 0$. The level repulsion is directly related to the extended nature of eigenmodes because the latter should be mutually orthogonal and thus cannot have the same frequency given that they overlap in space. Therefore, it is impossible to have $s=0$. In contrast, the frequencies of two spatially localized modes are uncorrelated if the modes are localized in two distant regions of space and have no spatial overlap. This leads to the Poisson distribution of normalized eigenfrequency spacings [3, 4, 5, 6]:

[TABLE]

The level repulsion is absent here and the frequencies of two modes can be arbitrary close to each other or even coincide: $P_{\textrm{P}}(0)>0$.

In this paper, we study the level-spacing distribution $P(s)$ for light in a disordered two-dimensional (2D) photonic crystal. Whereas the eigenmodes of the ideal crystal are extended, disorder induces localized eigenmodes due to the phenomenon of Anderson localization [7, 8]. The latter is widespread in physics and under appropriate circumstances, takes place for electrons in disordered solids [7, 9], atoms in random potentials [10, 11, 12], sound [13] and light [14] in disordered materials. The optical case is of special interest because of numerous existing and emerging technological applications of disordered photonics [15, 16, 17, 18]. Early theoretical work on Anderson localization of light exploited the analogy between Maxwell and Schrödinger equations to adapt to light the results first obtained for electrons in disordered solids [19, 21, 20]. Differences between the two cases were also identified, originating mainly from a different dispersion relation and the resonant nature of scattering for optical waves [22, 23]. Experiments were quite successful in low-dimensional (1D [24], quasi-1D [25] and 2D [26]) systems but inconclusive in 3D [27, 28, 29, 30] where an undisputable evidence of disorder-induced localization is still lacking [31]. Recent work suggests that the vector nature of light (i.e., the presence of the polarization degree of freedom) and its associated strong near-field coupling between nearby scatterers by longitudinal fields (the so-called dipole-dipole interactions) play a crucial role for Anderson localization of light in 1D [32], 2D [33], and 3D [34] systems. In particular, they prevent Anderson localization in a simple 3D model in which light is scattered by randomly distributed resonant point scatterers [34, 35]. Although it is unclear whether the same mechanism counteracts localization of light in other disordered photonic systems (powders of small dielectric or semiconductor particles, porous media, etc.), estimations show that longitudinal near fields in light scattering from a single dielectric sphere are of the same order as for a point-like scatterer [36].

Disordered photonic crystals have been first proposed as candidate systems for observing Anderson localization by Sajeev John [21] and later studied by many authors [37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50]. (These references are only a few among many works dealing with 2D as well as with 3D disordered photonic crystals.) However, the level spacing statistics of these systems has not been studied until now. Localization of light in a disordered 2D photonic crystal presents some particularities that make it quite different from the “canonical” case of Anderson localization of a scalar wave in a completely disordered system without an underlying periodic structure. First, the band structures of the ideal crystal are different for light polarized perpendicular to the plane of propagation (transverse-magnetic or TM modes) and in the plane (transverse-electric or TE modes). The impact of disorder on the two types of modes turns out to be different as well. Second, localized states inside the bandgap of the ideal crystal and near the band edges are due to two distinct physical mechanisms: the former are defect modes that require disorder to exist but that are destroyed when the disorder is increased, whereas the latter are genuine Anderson localized modes that shrink with increasing disorder. This is a typical situation in disordered photonic-crystal structures [46, 47]. Numerous experimental realizations of 2D disordered photonic crystals exist, allowing for observation of various disorder-related phenomena, including those with a potential for practical applications [37, 38, 39, 40, 41].

To study the level-spacing statistics in disordered 2D photonic crystals and its relation with extended or localized nature of electromagnetic modes, we first estimate the localization length $\xi$ of the latter and perform a finite-size scaling analysis by studying the sensitivity of $\xi$ to the size $L$ of the disordered crystal. The independence of $\xi$ from $L$ signals localized modes. For the disordered photonic crystals studied in this work, we can claim the existence of localized modes with certainty only for TM polarization and frequencies inside the band gap of the crystal without disorder. We then analyze the statistics of spacings between eigenfrequencies of Maxwell equations in our 2D systems and establish a correspondence between the level repulsion phenomenon and the extended character of modes. The level spacing distribution is close to the Wigner-Dyson distribution (1) in this case. In the opposite case of localized modes, the level spacing distribution approaches the Poison distribution (2) expected for wave systems with uncorrelated eigenfrequencies.

The model

In 2D dielectric media, TM- and TE-polarized waves decouple from each other. For a constant magnetic permeability $\mu=1$ and a random dielectric constant $\epsilon(\bm{\rho})$ with $\bm{\rho}=(x,y)$, the electric field $E_{z}$ of a TM-polarized monochromatic wave obeys a scalar wave equation

[TABLE]

where $\bm{\nabla}=\mathbf{e}_{x}\partial/\partial x+\mathbf{e}_{y}\partial/\partial y$. On the other hand, the magnetic field $H_{z}$ of a monochromatic TE-polarized wave obeys

[TABLE]

The crucial difference between the two cases stems from the fact that in the TE case, the electric field $\mathbf{E}(\bm{\rho})$ is a vector lying the $xy$ plane whereas in the TM case, $E_{z}(\bm{\rho})$ is a scalar. This trivial observation can have important consequences as far as Anderson localization is concerned. It was found, for example, that only TM modes become spatially localized in a dense 2D random arrangement of identical two-level atoms whereas TE modes remain extended [33].

We consider a square array of infinitely long dielectric cylinders with a dielectric constant $\epsilon=11.7$ typical for silicon in the infrared [54, 55], embedded in air (see figure 1a). The cylinder radius is set to be $r=0.2a$, where $a$ is the lattice parameter. System sizes $L$ up to $15a$ are considered. The infinite structure ($L\to\infty$) has a complete photonic bandgap for TM modes with frequencies $\omega$ in the range $\omega a/2\pi c\in[0.28,0.43]$ and a partial bandgap for TE modes in the $\Gamma$X direction for $\omega a/2\pi c\in[0.55,0.66]$ [56, 55], see figure 2. These frequency ranges are the most interesting for us because the presence of a spectral gap in the periodic structure indicates strong (destructive) interference effects for light in the dielectric material, creating favorable conditions for the appearance of localized modes upon introducing disorder [21, 22]. Disorder is introduced in the crystal by displacing each cylinder along a randomly chosen diagonal of the unit cell by a distance $d=\alpha d_{\mathrm{max}}$ (see figure 1b), where the maximum displacement is $d_{\mathrm{max}}=\sqrt{2}[a/2-(r+\delta_{\mathrm{min}})]$, ensuring that a cylinder does not approach a cell boundary closer than $\delta=\delta_{\mathrm{min}}=0.02a$. The strength of disorder is quantified by a parameter $\alpha\in[0,1]$, with $\alpha=0$ corresponding to a perfect crystal and $\alpha=1$ to the most random structure. The resulting structure exhibits complicated structural correlations that are different from those considered in the recent work on Anderson localization and bandgap formation in hyperuniform random structures [45].

In experiments, 2D photonic crystals are often confined between two parallel reflecting planes perpendicular to the cylinders (or holes) forming the periodic structure [37, 38, 39, 40, 41, 46, 47, 48, 49, 50]. These systems, often called “photonic crystal slabs”, support TM- or TE-like modes that can be excited by an appropriate external excitation [56]. The modes of the array of infinite cylinders shown in figure 1a and those of a photonic crystal slab are equivalent as far as only the in-plane propagation and scattering [wave vector in the $(x,y)$ plane] are concerned. Because some amount of out-of-plane scattering is inevitable in real experimental systems, our results may be considered as an approximation of a realistic structure in which cylinders have finite length and in-plane scattering dominates.

Computational approach

Assuming periodic boundary conditions in the $xy$ plane, we consider our disordered system as a unit cell (a “supercell”) of a 2D unbounded photonic crystal with a lattice constant $L$. The eigenfunctions of equations (3) and (4) in such a crystal can be written in the form [56, 57]

[TABLE]

where $\psi$ denotes either $E_{z}$ or $H_{z}$ and $u_{n,\mathbf{k}}(\bm{\rho})$ are periodic functions of $x$ and $y$. Sufficiently strong disorder is expected to make the statistical properties of eigenfunctions and eigenfrequencies of the disordered crystal that we are going to study in the following independent from the direction of the wave vector $\mathbf{k}$. It is therefore sufficient to consider $\mathbf{k}$ parallel to the $x$ axis and restricted to a half of the first Brillouin zone $0\leq k\leq\pi/L$ (the irreducible zone). We use the standard software package FreeFem++ implementing the finite element method [58] and employ the Galerkin method [59] to calculate the band structure $\omega_{n}(k)$ and the eigenfunctions (modes) $u_{n,\mathbf{k}}(\bm{\rho})$. The localization properties of the modes are characterized by their inverse participation ratio (IPR)

[TABLE]

where $A=L^{2}$ is the area of our disordered system. As can be easily deduced from equation (6), a mode extended over the whole system will have IPR $\sim 1/A$, whereas a localized mode will have a larger IPR. We define the localization length $\xi_{n,\mathbf{k}}$ of a mode as

[TABLE]

where the numerical prefactor is somewhat arbitrary and will be of no importance in the following because we will compare $\xi_{n,\mathbf{k}}$ for different disorder strengths and system sizes but will not base any important conclusion on the magnitude of $\xi_{n,\mathbf{k}}$.

In a disordered system, $\xi_{n,\mathbf{k}}$ are random quantities and only their statistical properties make physical sense. To obtain the average localization length $\langle\xi(\omega)\rangle$ as a function of frequency, we employ a procedure illustrated in figure 3. A large number of $\xi_{n,\mathbf{k}}$ obtained for 50 different realizations of disorder and corresponding to frequencies $\omega_{n}(\mathbf{k})$ within a narrow interval around $\omega$, are averaged to obtain $\langle\xi(\omega)\rangle$ which, in its turn, is then fitted by a polynomial function. We will represent our results by the latter in the following.

Whereas the localization length characterizes the spatial extent of eigenmodes, the statistical properties of the spectrum of a disordered system can be characterized by the level spacing statistics $P(s)$, as we already discussed in the Introduction. To compute $P(s)$, we define a dimensionless spacing between adjacent bands $n$ and $n+1$:

[TABLE]

where $\Delta\omega_{n}(\mathbf{k})=\omega_{n+1}(\mathbf{k})-\omega_{n}(\mathbf{k})$. Examples of $P(s)$ obtained for TM modes are shown in figure 4. Previous studies of $P(s)$ in photonic [51, 52] and electronic [53] crystals without disorder have shown that the degeneracy of the spectrum in symmetry points $\Gamma$, $X$ or $M$ (see figure 2) can have a strong influence on $P(s)$ by suppressing the phenomenon of level repulsion. This suppression can be erroneously interpreted as a signature of localized modes. We observe remnants of this phenomenon at weak disorder $\alpha\leq 0.2$ (see figure 4a). To avoid the confusion between the phenomena due to the high symmetry of the perfect crystal from which we derive our disordered system and the localized nature of eigenmodes, for weak disorder $\alpha\leq 0.2$ we will use only $s_{n,\mathbf{k}}$ for $\mathbf{k}$ near the center of the irreducible Brillouin zone. More precisely, for $\alpha\leq 0.2$, $P(s)$ will be calculated using $s_{n,\mathbf{k}}$ with $k$ within an interval of width $\pi/2L$ around $k=\pi/2L$. In contrast, for $\alpha>0.2$, we will use $s_{n,\mathbf{k}}$ corresponding to all $k\in[0,\pi/L]$ because in this case $P(s)$ turns out to be independent of $k$ within statistical errors of our analysis, as illustrated by figures 4b and 4c (similar results were obtained for TE modes, but we do not present them here). Finally, very low frequencies $\omega_{n}(k)a/2\pi c<0.1$ will be eliminated from the analysis to reduce finite size effects. Such low frequencies correspond to long wavelengths $\lambda>10a$ (wavelength in the air), which are comparable with the size of our disordered system $L$. For TE modes we restrict the analysis to even larger $\omega_{n}(k)a/2\pi c>0.35$ which is made possible by a higher position of the most interesting frequency range—the one corresponding to a partial spectral gap shown in pink in figure 2b—than for TM modes.

Results

Localization length

Figure 5 shows the frequency dependencies of the average localization lengths $\langle\xi(\omega)\rangle$ for TM (figure 5a) and TE (figure 5b) modes, calculated for different disorder strengths $\alpha$. Several important differences are obvious between TM and TE cases. First, $\langle\xi(\omega)\rangle$ reaches smaller values for TM than for TE modes, with the shortest localization lengths achieved inside the bandgap of the initial periodic structure. Second, the localization length of TM modes decays towards the center of the bandgap from both sides, whereas for TE modes, such a decay is clearly visible only from the low-frequency side, the localization length exhibiting a complicated oscillatory behavior on the high-frequency side of the partial bandgap. An interesting feature that we observe in figure 5a is the crossing of lines corresponding to different $\alpha$ at frequencies that are close to band edges but do not coincide with them, lying slightly inside the bandgap of the structure without disorder. We marked the approximate positions of these crossings by cyan vertical lines in figure 5a. Together with the band edges of the ideal crystal without disorder, the latter split the spectrum of TM modes in 5 different regions: the frequencies below (I) and above (V) the band gap of the ideal crystal, the middle of the bandgap (III), and the frequencies near band edges (II and IV). Analysis of TE modes (figure 5b) does not suggest the same fine structure in this case, leaving us with 3 different regions that can be clearly identified: the frequencies below (I) and above (V) the partial bandgap, and the partial bandgap itself (III). We summarize this classification of spectral regions in table 1. Note that the emergence of electromagnetic modes in the spectral regions II–IV leads to progressive narrowing of the gap in the spectrum of TM modes with increasing $\alpha$ and ultimately to its complete closing for $\alpha=1$, as we illustrate in figure 6a. This figure also illustrates the role of periodicity in the formation of the first band gap in the photonic crystal that we consider. It is well known that coupling between scattering (Mie) resonances of dielectric cylinders can give rise to gaps in the electromagnetic spectrum of a structure that they form, even when the structure is not periodic [60, 61, 62, 45]. The TM band gap that we consider is certainly also sensitive to the properties (dielectric constant, radius) of individual cylinders, but closes when sufficient disorder is introduced and thus would not exist in a fully random structure.

For frequencies outside the range delimited by cyan lines in figure 5a (spectral regions I, II, IV and V), $\langle\xi(\omega)\rangle$ exhibits a behavior that could be expected because stronger disorder (larger $\alpha$) yields shorter $\langle\xi(\omega)\rangle$, see figure 6b. However, inside the frequency band delimited by cyan lines (spectral region III), the dependence of $\langle\xi(\omega)\rangle$ on $\alpha$ is just the opposite: larger $\alpha$ (i.e. stronger disorder) gives longer $\langle\xi(\omega)\rangle$ (blue dashed and magenta solid lines in figure 6b). We believe that such a behavior is due to the different localization mechanisms at work for modes with frequencies in the middle of the gap and near band edges. On the one hand, the modes with frequencies in the middle of the gap (region III in figure 5a) need disorder to appear (because there are strictly no modes in the gap in the absence of disorder), but their localized nature relies on the periodic structure of the system. In other words, they benefit from Bragg reflections from the parts of the photonic crystal surrounding them. Increasing disorder makes these modes more numerous and hence easier to detect, but reduces their spatial localization making their localization lengths larger because of less efficient Bragg reflections. On the other hand, the modes near band edges (regions II and IV in figure 5a) get more and more localized when the randomness $\alpha$ is increased. Stronger disorder (i.e., larger $\alpha$) provides a shorter localization length for this modes, as expected for the genuine Anderson localization. These two different “modes of localization” have been previously predicted to exist in 1D random periodic-on-average systems [63].

Figure 7 illustrates the intensity patterns of typical TM modes in different frequency ranges. No difference can be identified by eye between the modes with frequencies below (figure 7a) or above (figure 7b) the band gap, which both look extended. The modes with frequencies inside the band gap of the periodic structure but near band edges (frequency ranges II and IV, see figures 7c and 7d) also look extended and only the modes with frequencies inside the range III appear as clearly localized in space (figures 7e and 7f). This may seem to be in a contradiction with figure 5a from which we see that the modes shown in figures 7c and 7e, as well as those in figures 7d and 7f, have approximately the same localization length $\xi\approx 2a$ on average. This apparent contradiction illustrates the ambiguity of our definition (7) of the localization length $\xi$ via the IPR. Indeed, the IPR and $\xi$ characterize, respectively, the inverse area and the spatial extent of a region within which the mode intensity is considerable but they do not change if we split the area in $M$ parts of size $\xi/\sqrt{M}$ that we separate in space by arbitrary large distances. Once we realize this, it becomes easier to reconcile the different appearance of the modes shown in figures 7c and 7e (or in figures 7d and 7f) with very similar values of IPR and $\xi$ that they have. Thus, the localization length $\xi$ defined by equation (7) cannot be taken as the only measure of localization and other quantities, such as, e.g., the level spacing statistics that we will analyze below, have to be considered.

Figures 7c and 7d and, especially, 7e and 7f illustrate an important difference in the spatial structure of disorder-induced modes. The modes with frequencies near the lower edge of the band gap tend to be localized inside the dielectric cylinders whereas the modes with frequencies near the upper band edge—in the space (air) between the cylinders. Disorder-induced modes inherit this behavior from the modes of the photonic crystal without disorder, where the corresponding “dielectric” (for frequencies below the band gap) and “air” (for frequencies above the band gap) bands are known to exist [56].

Numerical studies of localization phenomena, including the one reported here, are subject to limitations due to the finite size $L$ of disordered regions that can be analyzed ($L\leq 15a$ in our case). It is often difficult to decide whether a mode with a localization length $\xi$ equal to a fraction of $L$ should be considered spatially localized or not. To resolve this difficulty, one usually studies the evolution of $\xi$ with $L$. $\xi$ (roughly) proportional to $L$ indicates that, most probably, the mode is extended and the surface that it will effectively cover in the limit of $L\to\infty$ will be of order $L^{2}$. Localized modes distinguish themselves by a localization length $\xi$ that is (almost) independent of $L$, suggesting that the mode will cover a negligible part of the surface of the disordered region in the limit of $L\to\infty$. A practical realization of this analysis is illustrated in figure 8. This figure witnesses that localized modes appear for sufficiently strong disorder $\alpha\gtrsim 0.4$ in the central part of the band gap where the localization lengths computed for two different sizes of the disordered region $L=9a$ and $L=15a$, roughly coincide for TM polarization. When disorder is weak (e.g., $\alpha=0.2$ in figure 8a) or for frequencies outside the band gap at any $\alpha$, the average localization length $\langle\xi(\omega)\rangle$ behaves as if the modes were extended. The same is true for TE polarization at any $\alpha$ and throughout the considered frequency range. Because for TE modes, $\langle\xi(\omega)\rangle$ increases with $L$ in a roughly uniform way throughout the spectrum, we cannot conclude about the eventual localization of these modes in the limit $L\to\infty$. It may happen that the modes eventually get localized when $L$ is increased and that their apparent extended nature under conditions of figures 8b, 8d and 8f is due to the fact that their localization lengths are larger than $L$. In any case, however, the analysis of figure 8 indicates that TM modes are easier to localize by disorder than TE ones, which under certain experimental conditions may lead to TE modes that can be considered extended for at least some practical purposes. This conclusion, however, is restricted to a structure composed of high-index dielectric cylinders in a low-index material (air). In an inverted structure (cylindrical holes in a high-index bulk material), TE modes may be stronger localized than TM ones [47, 64, 65].

Level-spacing statistics

Now that we have characterized the localization properties of the modes of our disordered photonic crystal, we are ready to study its level spacing statistics $P(s)$ in a meaningful way. Figure 9 shows $P(s)$ computed for TM modes in frequency ranges defined in table 1. For frequencies $\omega_{n}(k)$ below the band gap, $P(s)$ tends to the Wigner-Dyson distribution (1) when the disorder is increased, see figure 9a. The same is true for $\omega_{n}(k)$ above the band gap where a good agreement of numerical results with equation (1) is found for all $\alpha\geq 0.2$ (figure 9b). This reflects the chaotic nature of our disordered photonic crystal in the classical limit and signals that the modes are extended both below and above the band gap, exhibiting a strong phenomenon of level repulsion: $P(0)=0$. A very similar behavior is found for disorder-induced modes inside the band gap but close to its edges (figure 9c). In contrast, $P(s)$ computed for the modes in the middle of the gap, is closer to the Poisson distribution (2), with the phenomenon of level repulsion being suppressed: $P(0)>0$, see figure 9d. This is in agreement with our previous conclusion about the localized character of modes in the spectral region III.

The level spacing distribution $P(s)$ for TE modes is shown in figure 10. $P(s)$ is close to the Wigner-Dyson result (1) in all frequency ranges, including the partial band gap region III (figure 10c), provided that the disorder is strong enough ($\alpha\gtrsim 0.6$). This is in agreement with the extended nature of eigenmodes following from the analysis of their localization lengths. For weak disorder $\alpha\lesssim 0.2$, $P(s)$ differs from the Wigner-Dyson distribution and has its maximum shifted towards smaller $s$. We attribute this to the remnants of the symmetries of the crystal without disorder which are still partially present at small $\alpha$. Even if the shape of $P(s)$ for TE modes at small $\alpha$ differs from the Wigner-Dyson distribution, $P(s)$ exhibits the phenomenon of level repulsion at all $\alpha$. To illustrate this fact, figure 10d shows the cumulative distribution function of $s$ defined as

[TABLE]

for a small value of $s=0.1$. ${\cal P}(s)$ is equal to the probability of having an eigenfrequency spacing inferior to $s$ and can serve as a measure of level repulsion if we take $s$ sufficiently small. We see from figure 10d that for TE modes, the level repulsion remains weak and roughly the same ${\cal P}(0.1)\simeq 0.05$ independent of the strength of disorder $\alpha$. In contrast, the disorder considerably weakens the level repulsion for TM modes and increases ${\cal P}(0.1)$ by an order of magnitude when $\alpha$ grows from 0.2 to 0.6. Once again, this observation is consistent with our previous conclusions based on the analysis of the localization lengths of TM and TE modes. Note that the suppression of level repulsion weakens for strong disorder $\alpha>0.6$. We relate this phenomenon with the reduction of Bragg scattering and the weakening of the associated Bragg-scattering-assisted localization mechanism at strong disorder.

Conclusions

We performed a detailed statistical analysis of eigenmodes of a 2D disordered photonic crystal composed of parallel, circular dielectric (silicon) rods in air. Dependencies of the average localization lengths of eigenmodes on their frequencies and on the size of the disordered system were analyzed. In the TM case (electric field parallel to the rods), the modes arising in the middle of the band gap of the ideal crystal are spatially localized due to Bragg reflections on the almost periodic structure surrounding them. The localization length of these modes grows with increasing disorder. In contrast, the modes appearing near band edges have a tendency to become more and more localized when the disorder is increased (genuine Anderson localization), although we cannot claim with certainty that these modes are spatially localized in our relatively small samples (up to $15\times 15$ rods). In the TE case (magnetic field parallel to the rods), the band gap of the ideal photonic crystal in the considered frequency range is only partial and no clear sign of localized modes were identified. This indicates that TM modes are easier to get localized by disorder than TE ones, although the conclusion may be different in an inverted structure (cylindrical air holes in a high-index dielectric material).

For TM modes outside the band gap and TE modes at all frequencies, the level spacing statistics is found to be close to the Wigner-Dyson distribution, typical for classically chaotic systems and for disordered wave systems with extended eigenmodes. Significant deviations from the Wigner-Dyson distribution are observed at weak disorder and are attributed to the residual symmetries due to the periodicity of the initial photonic crystal. The latter symmetries partially suppress the level repulsion phenomenon and modify the overall shape of the distribution. For TM modes inside the band gap of the ideal crystal, the phenomenon of level repulsion is significantly suppressed and the level spacing statistics approaches the Poisson distribution, as expected for statistically independent eigenfrequencies corresponding to spatially localized eigenmodes. The level repulsion weakens upon increasing disorder, reaches a minimum at some optimal disorder for which the interplay between the Bragg-scattering-assisted and disorder-induced localization is optimal, and then starts to strengthen again due to the suppression of Bragg scattering by disorder.

An interesting extension of this work may concern a comparison of results obtained here for a system of dielectric rods in air with those for an inverted structure in which circular holes are drilled in a dielectric material. The latter structure provides stronger scattering for TE modes and such a comparison can therefore give an insight into the role of polarization in the phenomenon of Anderson localization of light.

Acknowledgements

This work was funded by the Agence Nationale de la Recherche (project ANR-14-CE26-0032 LOVE). All the computations presented in this paper were performed using the Froggy platform of the CIMENT infrastructure (https://ciment.ujf-grenoble.fr), which is supported by the Rhone-Alpes region (grant CPER07_13 CIRA) and the Equip@Meso project (reference ANR-10-EQPX-29-01) of the programme Investissements d’Avenir supervised by the Agence Nationale de la Recherche.

Author contributions statement

S.E.S. formulated the problem. J.M.E. performed the calculations. Both authors contributed to the discussion of results and writing the paper.

Additional information

Competing interests. The authors declare no competing interests.