Modal properties of 2D disordered optical structures have been numerically analyzed, in the Mid-IR region, varying the amount of scattering and the disorder level. The modal properties study has been carried out through the use of Finite Element Method, highlighting the localized regime transition and investigating the quality factor. The results have been interpreted in a statistical fashion, investigating light diffusion in these structures with the help of Monte Carlo Method. An alternative measure of randomness weight has been proposed to support the numerical results.
© 2015 Optical Society of America
Light transport in disordered structures is a wide research topic that covers an ample range of applications. It includes lasers technology , spectrometry  and photovoltaic field where the disorder enhancement is promising to improve the absorption efficiency of thin-film solar cells . Multiple scattering plays a key role in the properties of the disordered optical media, where the dielectric constant is a random function of the spatial coordinates. The multiple scattering phenomenon is quite common in nature and appears in all of the opaque materials, for example clouds, milk or sugar powder. The light, traveling through such materials, is scattered thousands of times in a fashion that can be compared to a random walk. Nevertheless, because of the interference between the light paths inside disordered media, the multiple scattering is a complete coherent phenomenon. This fact leads to the formation of electromagnetic modes [1, 4], which are eigen-solutions of Maxwell’s equations in the disordered structure. The nature of these modes and consequently, the light transport, are strongly related to the scattering strength. Varying the scattering contribution, the light propagation shows a transition from a ballistic transport to a classical diffusion and eventually to a localized regime. Diffusion regime is characterized by a linear transmission coefficient with respect to the thickness of the structure, while localized regime is characterized by a transmission coefficient which decays exponentially [5, 6]. As a consequence, the spacial extension of the modes varies significantly, permitting a distinction between extended and localized modes [7, 8]. It is worth noting that such a distinction is not very sharp and it is influenced by the physical mechanisms used to explain the localization phenomena, which can be ascribed both to the backscattering, which brings to a weak localization, and to the generation of closed scattering channels, which are related to a strong localization . An appropriate measure for the scattering can be obtained evaluating the product between the mean free path l and the wavevector k. The onset of strong localization is expected when kl ≤ 1, which is known as the modified Ioffe-Regel criterion .
As theoretically discussed by Letokhov , when a disordered structure is built including an active optical material, it is possible to create an atypical class of lasers called random lasers . In these objects the lasing principle is similar to a conventional laser with an active material that provides the amplification through the stimulated emission. The difference is related to the feedback mechanism, where the light is confined by the effect of multiple scattering and not by the presence of a resonant cavity. This fact has some consequences on the modes configuration driving random lasers to show a wide emission spectrum characterized by a number of sharp peaks, apparently not correlated with each other. Under this point of view, the study of the leaky modes and the mechanisms of light diffusion in the random structure are fundamental to understand the behavior of random lasers [12, 13]. Experiments show that random lasing phenomenon is a complex scenario, where modes can interact with each others , and the emission spectra can be characterized by the co-existence between localized and extended modes [15, 16].
In this contribution a 2D disordered structure, similar to that designed by Liang in order to demonstrate random lasing in Mid-IR region , has been investigated. The modal properties of the structure have been studied with different numerical strategies, varying the amount of the disorder and the level of scattering in the structure. The first numerical approach, based on Finite Element Method (FEM), has been used to evaluate the resonant modes, their spectral distribution, and their confinement degree, considering a lossy structure without any gain effect. The second approach, based on Monte Carlo Method, has been exploited to obtain information about the statistical parameters, such as scattering distance and travel distance, that influence the light diffusion inside the material.
2. Planar structure
The structure under investigation has been modeled on the device, developed in  in order to operate as a random laser source in the Mid-IR region with a wavelength around 10 μm which is approximately equivalent to a frequency of 30 THz. The main structure is a square slab with side length l equal to 150 μm and thickness t equal to 2.2 μm, obtained in a Quantum Cascade Laser (QCL) heterostructure of In0.53Ga0.47As/In0.52Al0.48As with a refractive index ns equal to 3.35. Considering the ratio between the slab surface width and the thickness, and considering that the QCL heterostructure guides primely TM modes, it is possible to model the system as a 2D structure, where x and y axes are in plane. Random structures in 2D approximation play an important role in this research field, being the simplest systems to achieve a complete light transport transition from diffusive to localized regime, and being enough versatile to be engineered [18, 19]. A pattern of circular air holes, the scatterers, with a diameter d equal to 3 μm and a refractive index na equal to 1, has been obtained inside this slab, with different dispositions according to the strength of the scattering and to the amount of disorder. The scattering has been varied modifying the number of holes, which is equivalent to modify the filling factor FF , a parameter defined as the ratio between the holes total area and the square slab area. The disorder has been varied generating the holes pattern exploiting two different algorithms. The first random pattern has been obtained choosing x and y components of the holes centers, uniformly distributed between 0 and l, avoiding superpositions between holes and discarding the holes that partially fall outside the slab. To obtain a more realistic structure, a minimum distance dmin, with a value of 150 nm, has been assigned between adjacent holes. An example of full random pattern structure with FF equal to 0.3 is shown in Fig. 1(a), while some details of full random patterns with FF equal to 0.5 and to 0.1 are shown respectively in Figs. 1(b) and 1(c). The second algorithm is intended to gradually move the holes from a periodic triangular lattice to a random configuration, using the parameter Ks which defines a random variable S as the shift of the holes centers with respect to the holes position in the periodic lattice:Figure 1(d) depicts the triangular lattice while Fig. 1(e) shows an example of pattern shifted with a weight equal to 0.6, both of these structure are generated considering a filling factor equal to 0.3.
3. Modes spectra analysis
This analysis has been performed using the FEM in order to solve for Maxwell’s equations in the frequency domain. The slab, modeled as shown in Fig. 1(a), has been surrounded by an empty region of 30 μm which terminates with a Perfectly Matched Layer (PML) framework. The boundary conditions have been chosen according to the scattering theory. These strategies pledge an irrelevant reflection due to the boundaries, fundamental in this case of study. In order to obtain solutions totally independent from the discretization of the problem, a particular attention has been given to the mesh parameters, setting the maximum elements size less than one seventh of the wavelength. In the setup phase of this FEM simulation work, different levels of mesh fineness have been tested in order to improve the accuracy of the solutions keeping the calculation time reasonably small. The maximum elements size has been chosen equivalent to the reference wavelength inside the media, which is similar to the expected fluctuation of solutions due to the light interference. With a heuristic procedure, the elements size has been decreased evaluating the numerical stability of a reference group of eigenfrequencies near 30 THz region. It has been observed that for elements size smaller than one seventh of the wavelength the solutions remained constant with a large number of digits, enough to obtain a confident values of quality factor, as further defined.
The first study has been focused to investigate the spectral behavior of the structure, performing a sweep in a range between 29 THz and 31 THz, using a step of 5 GHz. Due to the nature of the slab material, a QCL heterostructure, the modes of interest, i.e. TMs with Hz = 0, have been stimulated imposing 1 A of line current, perpendicular to the xy plane. As emerged in , the slab has been uniformly pumped to stimulate TM modes over the whole structure. To obtain a similar behavior the source injection current has been distributed over a large number of points, precisely 5000, located on the active substrate. With this number, the average distance between two points is roughly the wavelength in the medium. The location of these points is carried out by the same algorithm that generates the full random pattern. At the output the electromagnetic energy has been recorded on the area of the square slab. Following this operative method, the impact of the scattering weight on the system has been investigated. Numerous realizations of full random patterns, with filling factor varying in a range from 0.1 to 0.5, have been simulated. Figure 2(a) shows the comparison between five energy spectral distributions starting from a low scattering case with FF equal to 0.1 and ending to strong scattering case with FF equal to 0.5. It is possible to notice the evolution of energy inside the system. With a filling factor smaller than 0.3, the spectrum appears a sequence of smooth peaks and valley, with rare well defined peaks. At filling factor equal to 0.3 the first sharp peaks arise sharing the band with wider gather of energy. At values of filling factor larger than 0.3 peaks become narrow and well separated between each other. To give an interpretation of these spectra, a second study has been performed using FEM solver searching for eigevalues and eigevectors of the electromagnetic problem in frequency domain. Investigating for the modes located near the frequencies corresponding to a maximum of energy it is possible to see that a sharp and narrow peak corresponds to a well definite mode, while a wide and smooth maximum of energy is associated to the superposition of different modes which share the same band and the same area in the slab. Defining the quality factor as:Figs. 2(b) – 2(d). Figure 2(b) shows the principal mode located at 29.60 THz in the structure with FF equal to 0.1. It is clearly a lossy extended mode, occupying all the area of the slab, its calculated Qf is equal to 1247 and its nature can be associated to a diffusive transport. Figure 2(c) shows the mode located at 29.91 THz in the structure with FF equal to 0.3. Unlike the previous one, this mode corresponds to a narrow peak and its quality factor is substantially higher with a value of 2819. It is a mode with a transient behavior toward a localized regime, which is reached in the structure with filling factor equal to 0.5, as depicted in Fig. 2(d), where it is possible to see the mode corresponding to peak at 29.98 THz. This mode exhibits a Qf larger than 4500.
Following the same procedure described previously, the 2D slab has been analyzed shifting from a condition of perfect periodicity to a disordered situation. Triangular lattice has been chosen as the basic periodic pattern, while randomness has been added varying the shift coefficient Ks form 0.1 to 0.9. Figure 3(a) shows the change in the spectrum at different level of Ks, the filling factor of these random patterns is set at 0.3. The spectrum of triangular pattern is almost flat exhibiting some isolate energy peak in correspondence of the frequencies where the solutions of Maxwell’s equations assume a form with an envelope that varies periodically with the pitch of the lattice. These solutions show a quality factor very high, for instance the peak, located at 29.12 THz, has a calculated Qf greater than 6100 as shown in Fig. 3(c), this value is two times larger with respect to the maximum Qf found in the full random pattern. Disturbing the positions of the holes from the periodic lattice it is possible to see evolution of the energy inside the slab area. With Ks equal to 0.1, corresponding to a very small shift, some wide and distinct energy stacks starts to grow in correspondence of the periodic lattice peaks which are absorbed inside the stacks. This trend is more clear at Ks equal to 0.3, while at Ks equal to 0.6 the energy spectrum is uncorrelated with respect to the periodic spectrum. With values of Ks larger than 0.6, the spectrum is actually similar to the one obtained using the full random algorithm. This evolution can be also highlighted examining the evolution of the average quality factor. Figure 3(b) shows the mean Qf calculated at different FFs, referring to the periodic lattice, to the full random pattern, and to different shift values. The mean has been calculated by averaging 400 modes, roughly all the modes in the band under examination, using 10 different realizations of disorder. The Qf curve of periodic lattice gradually gets close to the Qf curve of the full random pattern, increasing the value of Ks. At Ks = 0.6 the average quality factor is practically equal to the full random pattern. This trend is more evident when strong resonant frequencies in the periodic pattern fall inside the band between 29 THz and 31 THz, this is evident for the patterns with filling factors equal to 0.2 and to 0.4.
It is interesting to notice that the random shift, generated by the algorithm previously explained, and described by the coefficient Ks = 0.6, is particular small, remaining highly correlated with the periodic pattern, as it is evident in the comparison between Figs. 1(d) and 1(e). Nevertheless such disturbance with respect the periodic lattice is enough to lead the system to show a spectrum similar to a full random generated pattern. This consideration suggests to give an alternative interpretation of disorder. As previously underlined, the localization, one of the key phenomenon that leads to random lasing, can be explained as the development of closed scattering paths in the case of strong localization, or as backscattering enhanced paths in the case of weak localization. In both cases these paths shall be a multiple of the resonant wavelength, which can be written as:11], m is a positive integer and F(•) is a function that weights all the possible jumps between two scatterers. λ and ls are random variables. Although F(•) can be written explicitly only for special and easy cases, for instance a multilayer random structure [21, 22], Eq. (3) suggests that the resonant modes statistic depends on the statistic of the scattering distance partial sum. Under this point of view a small disturbance in the periodic lattice, i.e. a random variable ls with a small variance, can generate resonant frequencies with a large variance, if the sum is large enough. Effectively, to understand the modes statistic, it is helpful to investigate the statistic of the light travel distance inside the structure, which is equivalent to ∑n ls.
4. Statistical analysis
In order to collect informations about the statistic of light diffusion in this system, it has been chosen to simulate the structure through a Monte Carlo Method, which is particularly suitable to analyze complex disordered structures [15, 23]. A custom software has been implemented assuming that photons are traveling grouped in packets identified by their energy. The scattering has been calculated according with the Mie theory, approximating the phase function with the Henyey-Greenstein function. To better simulate the structure, the program has been provided with a pumping system and with a complete model which implements absorption, spontaneous and stimulated emission. To make the software numerically efficient and reliable, random number generation for Monte Carlo operations, has been implemented with a single instance pseudo-random sequence, using singleton design pattern, and using the Marsenne Twister-64 bit algorithm provided by C++11 standard library.
The simulations has been performed supposing to electronically pump the system at negative time, to stop the pump at t = 0 and to start spontaneous emission at positive time. The system has been set and with a negligible absorption and without gain, to be coherent with the simulations described in the previous section. With these settings, the software has been run in order to diffuse and collect statistics from a large number of photons packets, roughly 100 × 106. The structure has been simulated varying the filling factor and the shift coefficient, in ×agreement with the FEM analysis performed in the previous section, searching for the statistic of the photons scattering distance, ls, and for the statistic of travel distance, ∑n ls. For every disorder configuration and for every scattering configuration, results have been averaged over 10 different realizations. Figure 4(a) shows the product between the the wavenumber k and the scattering mean free path calculated with the structure generated with the help of full random algorithm and varying the FF parameter. It is possible to notice that there is no circumstance where the Ioffe-Regel criterion has been fulfilled, as the minimum calculated value of is equal to 9.2 when the FF is equal to 0.5. This result implies that it is not possible to reach ·a strong localization in these structures, therefore the evolution of localization shown in Figs. 2(b)–2(d) can be explained using the weak localization mechanism. Similar conclusion has been suggested also by Liang .
To better understand how the degree of disorder affects the structure, two important parameters have been defined:Fig. 4(b) the Kl trend vs the filling factor has been depicted for various values of disorder. Considering that the scattering mean free path is weakly dependent on the disorder parameter, it is possible to notice that low disorder structures have a high variance and a Kl roughly flat profile varying the filling factor, while higher disorder structures tent to decrease the variance with the increase of filling factor. This fact can explain the reason why periodic structures support a few number of localized modes with respect to random structures. Encountering a larger scattering step variance, modes tent to spread through the structures and only few frequencies, enveloped with the lattice pitch, can resonate in the system. It is worth noting that for a low filling factor all the disorder configurations have an high Kl, therefore all the modes are extended and the transport regime is diffusive. A similar conclusion can be infered analyzing the coefficient Kd, Fig. 4(c). This parameter is quite useful because it gives informations about localization. When the system begins to support localized modes, the value of Ks increases and the curves at different values of disorder start to separate with each others, showing that high disorder configurations support more resonant frequencies with respect to periodic lattice and low disorder structures. This is evident when the filling factor is high. When filling factor is low all the modes are extended and all the structures have a similar behavior. Furthermore Fig. 4(c) shows that a little increase of the shift from periodic lattice pushes the curves to approach the full random behavior. Kd can measure how fast a disturbance from periodic lattice can behave like a random structure, this means that Kd is a good randomness weight for this planar disorder structures.
A disordered optical planar structure has been investigated through different numerical methods in order to show the modal properties in the Mid-IR region, around 10 μm. It has been also shown the nature of light transport, varying the amount of scattering, through the control of filling factor, and the impact of disorder from a periodic lattice to a full random pattern of scatterers. With the help of a Finite Element Method solver, the modes energy distributions has been investigated in a band between 29 THz and 31 THz, highlighting the relation between peaks of energy and resonant frequencies. A clear pictures of quality factor evolution has be given and the transition between diffusive transport regime and localized regime has been identified when the filling factor overcomes 0.3. A particular behavior of the structure, where a small disturbance of the periodic lattice induces the system to acts like a full random structure, has been exhibited, and therefore an explanation has been proposed, underlining the dependence of modes spectral properties on the statistic of light travel distance in the structure. In order to verify this hypothesis, a custom and versatile Monte Carlo Method based software has been developed. With this software the diffusion of photon inside the system has been simulated and information about the statistic of important parameters, like scattering distance and travel distance, has been collected. Two new parameters, the weighted scattering distance Kl and the weighted travel distance Kd have been defined. It has been shown that these parameters can explain the evolution of the modal properties varying the filling factor and the disorder weight. In particular Kd has been identified as a key parameter to measure randomness in this 2D structure.
This work is supported by the Agency for Science Technology and Research through the Advanced Optics in Engineering Programme (Grant No. 122 360 0011) and A*STAR Research Attachment Program.
References and links
1. D. S. Wiersma, “The physics and applications of random laser,” Nature Phys. 4, 359–367 (2008). [CrossRef]
2. B. Redding, S. F. Liew, R. Sarma, and H. Cao, “Compact spectrometer based on a disordered photonic chip,” Nature Photon. 7, 746–751 (2013). [CrossRef]
4. K. Vynck, M. Burresi, F. Riboli, and D. S. Wiersma, “Photon management in two-dimensional disordered media,” Nature Mater. 11, 1017–1022 (2012).
5. M. Segev, Y. Silberberg, and D. N. Christodoulides, “Anderson lacalization of light,” Nature Photon. 7, 197204 (2013). [CrossRef]
6. D. S. Wiersma, P. Bartolini, A. Lagendijk, and R. Righini, “Localization of light in a disordered medium,” Nature 390, 671–673 (1997). [CrossRef]
8. C. Vanneste and P. Sebbah, “Complexity of 2D random laser modes at the transition from weak scattering to Anderson localization,” Phys. Rev. A 79, 041802 (2009). [CrossRef]
10. V. S. Letokhov, “Light generation by a scattering medium with a negative resonant absorption,” Sov. Phys. JETP 16, 835–840 (1968).
11. D. S. Wiersma and A. Lagendijk, “Light diffusion with gain and random lasers,” Phys. Rev. E 54, 4256 (1996). [CrossRef]
12. H. E. Tureci, A. D. Stone, and B. Collier, “Self-consistent multimode lasing theory for complex or random lasing media,” Phys. Rev. A 74, 043822 (2006). [CrossRef]
14. H. Cao, X. Jiang, Y. Ling, J. Y. Xu, and C. M. Soukoulis, “Mode repulsion and mode coupling in random lasers,” Phys. Rev. B 67, 161101 (2003). [CrossRef]
16. J. Fallert, R. J. B. Dietz, J. Sartor, D. Schneider, C. Klingshirn, and H. Kalt, “Co-existence of strongly and weakly localized random laser modes,” Nature Photon. 3, 279–282 (2009). [CrossRef]
19. F. Riboli, N. Caselli, S. Vignolini, F. Intonti, K. Vynck, P. Barthelemy, A. Gerardino, L. Balet, L. H. Li, A. Fiore, M. Gurioli, and D. S. Wiersma, “Engineering of light confinement in strongly scattering disordered media,” Nature Mater. 13, 720–725 (2014). [CrossRef]
20. P. Sebbah and C. Vanneste, “Random laser in the localized regime,” Phys. Rev. B 66, 144202 (2002). [CrossRef]
21. X. Jiang and C. M. Soukoulis, “Transmission and reflection studies of periodic and random systems with gain,” Phys. Rev. B 59, 6159 (1999). [CrossRef]
23. S. Mujumdar, R. Torre, H. Ramachandran, and D. Wiersma, “Monte Carlo calculations of spectral features in random lasing,” J. Nanophoton. 4(1), 041550 (2010). [CrossRef]