Experimental and theoretical investigation of statistical regimes in random laser emission
Abstract
We present a theoretical and experimental study aimed at characterizing statistical regimes in a random laser. Both the theoretical simulations and the experimental results show the possibility of three region of fluctuations increasing the pumping energy. An initial Gaussian regime is followed by a Lévy statistics and Gaussian statistic is recovered again for high pump pulse energy. These different statistical regimes are possible in a weakly diffusive active medium, while the region of Lévy statistics disappears when the medium is strongly diffusive presenting always a Gaussian regime with smooth emission spectrum. Experiments and theory agree in identification of the key parameters determining the statistical regimes of the random laser.
pacs:
42.55.Zz, 42.65.Sf,05.40.a,42.25DdI Introduction
In a random laser, as theorized for the first time by Letokhov in the 1960s Letokhov (), amplification via stimulated emission is provided by multiple scattering in a disordered gain medium Wiersma1 (); Wiersma2 (); Cao1 (). In contrast to conventional laser modes, which are determined by the characteristics of the optical cavity, in a random laser light propagates as a random walk through a diffusing material. If the gain experienced along the path inside the active medium overcomes the losses, laser emission can occur.
Despite the large number of experimental studies on the emission characteristics of different random laser materials as well as the large amount of theoretical models reported in literature, a clear connection between experiment and theory is still lacking to date. Because the random nature of light propagation, with a large number of modes competing for the available gain, the emission spectrum and the spatial properties of the output are in general hard to predict; the concept of “mode” itself can refer to different definitions or interpretations. Moreover, together with narrowing of the emission spectrum above laser threshold, narrow spikes at random frequencies were also reported in some experimental works. The presence of these narrow spikes in emission spectra led many researchers to suggest different theoretical possible explanations. For instance, Anderson localization was invoked for a material where the scattering mean free path becomes smaller than the radiation wavelength Jiang (); Pradhan (); Cao2 (); Vanneste. Other models were proposed in order to explain the presence of narrow spikes in a weakly diffusive regime (), where the Anderson localization constrain is not satisfied: waveguide structure in a material with random variable refractive index Apalkov (), FabryPerot type cavity determined by the geometry of the active volume within the sample Wu (), single frequency eigenmodes with a long decay time Molen () and lasing mode consisting of traveling waves extended through the whole system Vanneste2007 (). A recent work Ramy () suggests that the crossover from smooth to spiky spectrum, in a weakly scattering medium, can be interpreted in terms of modes that become uncoupled from mode competition. Some experimental and theoretical studies suggest that random fluctuations in these spectra may be due to rare extended modes, identified as possible long lifetime random paths within the active medium; these socalled lucky photons are able to experience a very large gain, giving rise to spikes in total intensity spectrum Mujandar (); Mujandar2 (). It has also been reported that, in the presence of strong fluctuations, the spikes’ intensity distribution follows a Lévylike behavior, i.e. a distribution that exhibits a fat tail associated to infrequent but large events Sharma (). As a results emission can display a large statistical variability even under well controlled experimental conditions.
In a previous theoretical work Lepri (), based on a phaseinsensitive intensity feedback model, it was demonstrated that the fluctuations of the emitted light can change, upon changing of a control parameter, from Gaussian to Lévy statistical regimes, the latter being characterized by a powerlaw tail in the intensity distributions. Actually, the Lévy regime is limited to a specific range of gain length, while Gaussian statistics occur with both low and high pumping energy. Three different statistical regimes are then possible depending on the gain Lepri (). Such regimes have been experimentally identified in weakly scattering media with a shorttime pump pulse duration ( 30 ps) Uppu (). A more recent work presents a theoretical to experimental convergence in determination of the first transition in the statistical regimes as a function of pumping energy for different scatterers concentrations Uppu2 (). A related experimental investigation of Lévy statistical fluctuations for samples with immobile light scatterers can be found in Ref. Zhu2012 ().
A short excitation time has been often considered a requirement for observation of spiky spectra. A long duration may lead to an averaging mechanism among modes, thus reducing the fluctuations.
In this work we show experimental evidence of different statistical regimes of fluctuations and different spectra, presenting smooth or spiky structure, as a function of the pumping energy in a nanosecond excitation. We present also a Monte Carlo simulation of a twodimensional (2D) model, that includes the medium’s population and parallel processing of a large number of random walkers, in either a long or short excitation pulse. In both experiments and numerical simulations, three different statistical regimes by increasing the pumping energy are possible: Gaussian transition to Lévy regime and back again for high pumping energy. We show that the value of the scattering mean free path plays a key role in determining the transition through the different statistical regimes or even in closing the region of Lévy statistics both in theoretical simulations and experimental results. We then present a clear quantitative link between experiment and theory.
The organization of the paper is the following: In Sec. II we describe the theoretical model and the structure of the Monte Carlo simulation and results are discussed in Sec. III. In Sec. IV we describe the experimental set up and the measurement method designed and performed in order to optimize spectra analysis and statistics of fluctuations. The results obtained are presented in Sec. V. Finally, Sec. VI summarizes and comments on the experimental results and theoretical analysis.
Ii Theoretical Model
The numerical model used in this work is based on an intensity feedback that does not include the effects of the field interference and thus treats light propagation through a MonteCarlo dynamics Mujandar (); Lepri (). The sample consists of a twodimensional square region partitioned in equal square cells of linear size given by the distance covered by light in a time step . Each cell is labeled by the vector index , with integers , so that is the total number of cells. The amplifying medium is modeled by a fourlevel system and a population of excited atoms, providing the population inversion . The spontaneous emission is reproduced by a spectrum of 1001 arbitrary channels (labeled from 500 to 500 channel) centered on the resonance of the transition (channel 0) with a linewidth of channels. The light dynamics is described in terms of a large number of random walkers each carrying an energy (or number of photons) with an assigned frequency channel around the central frequency .
The simulation consists in a loop of three distinct steps for each :

Spontaneous emission
For each lattice cell, a spontaneous emission event can randomly occur with a probability given by . If this event occurs, the local population of the cell is reduced by a unity and a new walker is created with an initial energy , a random initial direction of motion and a frequency . The latter is drawn at random from a Cauchy distribution of width centered at the resonance of the transition [see Eq.(3 below]. 
Diffusion
Each walker present in the slab can undergo a scattering event with a probability , that randomly changes its direction. The position of each walker is thus updated according the trajectory of its motion. If a walker reaches a slab boundary it is destroyed and its number of photons and frequency are recorded for the statistical analysis. 
Stimulated emission
The population of each lattice cell and the number of photons carried by each walker are deterministically updated according to the following rules:(1) (2) where is population of the lattice cell where the th walker is localized, whereas the stimulated emission coefficient depends on frequency:
(3)
It should be emphasized that, although each walker evolves independently from all the others, they all interact with the same population distribution, which, in turn, determines the photon number distributions. It is thus mandatory to evolve in parallel the trajectory of all the walkers that are present in the simulation box. This requirement is computationally demanding but is necessary to take into account the effect of gain saturation and spatiotemporal coupling among different walkers in a selfconsistent way. Note also that the gain length (i.e. the typical length needed to increase the walker energy by a factor ) is not constant being a function of both time and space.
In the following, the initial population distribution at is assumed to have a Gaussian spatial shape
(4) 
The normalization factor is assigned to have a prescribed initial energy :
(5) 
This assignment corresponds to considering an infinitely short excitation during which the sample absorbs and energy from the pump beam. In the following, we also investigate the case in which pumping occurs on a finite time ; in such a case the updating of the population follow this rule:
(6) 
where the pumping rate is given by
(7) 
Thus, the important parameters to be set for each Monte Carlo are: the
spontaneous emission rate , the linewidth of the transition ,
and the scattering
probability per unit time . The latter determines
the actual mean free path
. This that can be compared with the linear size
of the extension of the initially excited population,
as defined by Eq.(4).
It should be noted that gives a measure of the effective
size of the
gain volume, i.e. the actual region where amplification occurs.
In numerical simulations spatial lengths, time intervals and
frequencies are expressed in units of , and
respectively, while the energy is reported as the
number of initial excited molecules in the lattice and in units
of , where is the central
frequency of the transition.
In the following we fix = 40 and = 150. We will study the statistic regimes as a function of other parameters and in particular of .
Iii Numerical Statistical Analysis
The central limit theorem states that the sum of a large number of independent, uniformly distributed random variables tends to a normal distribution; nonetheless there are many systems and observations for which a nonGaussian stable model is to be expected. In these cases, a distribution with infinite variance can be described, according to the generalized central limit theorem, by a Lévystable distribution, that takes into account extreme changes in a variables, associated with strong fluctuations with rare and very large values having an infinite variance alpha0 (); alpha00 (). Such distributions are parametrized by the exponent (the index of stability or tail index) that describes the rate at which the tail of the distribution tapers off: At the tail exhibits a asymptotically powerlaw behavior. Upon approaching , the distribution tends to a Gaussian one.
In the following we characterize the measured data by fitting them to a Lévystable distribution. Before entering in the presentation of the results, we emphasize that this procedure deserves a word of caution. Due to the finite number of measurements it is in general rather difficult to reliably resolve the tail of the distributions. In particular, for marginal cases where is close to 2 there is a certain degree of arbitrariness in discriminating between a Gaussian distribution and a Lévy one. We thus decided to fix an arbitrary threshold value, , above which the characteristics of the intensity distribution and emission spectra converge toward a Gaussian behavior leading to an attribution of Gaussian regime.
Figure 1 shows the simulation results for three different pump energies. The chosen parameter values are , a emission linewidth of channels and a , whereas the pumping temporal type is impulsive, with the whole amount of energy supplied at the time . Some typical emission spectra are reported in the first column for three values of the pumping energies, while in the second column are reported the histograms of the intensity distribution of the spectrum peak as obtained by 10000 simulations with identical initial conditions. The values of the pumping energy are chosen to be sub threshold, near threshold and above the random laser threshold. The threshold is the energy of the change of slope of the peak intensity behavior vs pump pulse energy averaged on more than 1000 simulations. The fluctuations of the spectral intensity at different frequencies near the maximum were completely uncorrelated. The characteristic parameters of the Lévystable are estimated via a regressiontype fitting alpha1 (); alpha2 (); alpha3 ().
It is important to note that in correspondence of the intermediate pump energy value, the emission spectrum is characterized by the presence of spikes at random frequencies with strong fluctuations of the peak intensity. The fitting of the histogram yields indeed , demonstrating that for this parameters choice the system is in the regime of Lévy fluctuations. Conversely, for the lower and the higher pumping levels, fluctuations are modest and is approaching the Gaussian value of 2. For this simulation , i.e. comparable to .
In Fig. 2 typical spectra and two extreme values of the are reported: The data on the left are calculated with a , that is larger than pumped region dimension ( = 40), and the data on the right with which is much smaller. For the same energies are also reported the values of index. It is evident that the spectra are drastically different in the two cases: In the weakly diffusive case () the intensity fluctuations are Lévydistributed even at larger energy while in the strongly diffusive case () the distribution is always Gaussian.
For a more quantitative analysis, in Fig. 3, we report the index calculated as a function of pump energy for four different values of . Assuming, as discussed above, the threshold value of for the passage from Gaussian to Lévy distribution, we note that a strongly diffusive medium displays a Gaussian intensity distribution at any value of the pump pulse. In the case of weaker diffusive media three different statistical regimes are possible increasing the pumping energy: The intensity distribution is Gaussian at low pump energies, becomes a Lévytype distribution at increasing pumping energy and then becomes Gaussian again for higher pumping energies. A key parameter for the statistical behavior is the ratio of over . When the random paths in the excited part of the medium have a low probability of being spatially overlapped for low or intermediate pumping energies. In this condition a first transition of the statistical distribution from Gaussian to Lévy is possible, corresponding to a transition from a spectra dominated by spontaneous emission to one where it is possible that rare long random paths within the active medium, with a very large gain, give rise to spikes in the spectrum. By increasing the pump pulse energy, the number of random walkers increases, as does the probability of having spatial overlapping and temporal coincidence. The random walks are coupled by the population inversion, with a general amplification with limited fluctuations,leading to a second transition to a Gaussian distribution. The line is reported in figure for emphasize the crossovers between statistical regimes (Gaussian for or Lévy for ). For ( = 4), the random paths are strongly coupled and the intensity fluctuations are always limited presenting a Gaussian distribution both in the regime of spontaneous or stimulated emission ( always above 1.8).
From a temporal point of view, a critical issue for the emergence of Lévy statistics is the comparison of two characteristic times: the mean scattering time , which is 40 in the case of Fig. 1 and 100 and 4 in the cases reported on Fig. 2 (in units of ), and the ballistic time span required to cross the gain volume , which can be estimated as , which is 40 in units of for all cases. The presence of Lévy statistics in the whole range of pump energy is possible when becomes of the same order or larger than , as it is showed in cases with a mean scattering time of 40 and 100, while a Gaussian statistics is always observed when is much smaller than (=4).
Another set of numerical simulations was performed in order to clarify the role of the pumping time duration on the spectrum characteristics. One may expect that a larger could lead to an averaging mechanism inhibiting the Lévy regime of fluctuations. The results reported in Fig. 4 show two different trends of values as a function of pumping energy supplied to the medium in a time equal to the spontaneous emission lifetime. In the less diffusive medium case, a large zone of Lévy distribution is present while, in a more diffusive medium, the Lévy region disappears and the is approaching 2 for any value of the energy. This result indicates that the time duration of the pumping mechanism does not determine the possibility of emerging of strong fluctuations while it is again crucial the condition of the weak or strong diffusive medium.
One can wonder how this model is connected to the concept of mode and whether our results can be interpreted in that theoretical frame. The theoretical approach followed in this work suggests an interpretation of “open mode” as possible random path with a certain spatial profile within the active medium; under this picture, the competition of modes in the gain has a critical role, determining the condition of different spectral profiles and the statistical regime of fluctuations. The condition of modes spatially uncoupled, possible in a weakly scattering active medium, can lead to strong fluctuations and isolated peaks in the spectrum (Lévy regime) that disappear in the high pump energy condition when many modes are populated simultaneously. In the case of strong scattering active medium, the modes are always coupled, leading to weak intensity fluctuations and a smooth spectrum (Gaussian regime): In this condition the Lévy regime disappears for any pump energy.
Iv Experimental Set Up and Samples
In our experiments we utilized samples made up a solution of 1 mM concentration of Rodhamine 6G dye dissolved in methanol, with different concentrations of a TiO nanopowder, with a mean radius of 300 nm, suspended in the solution to create disorder. The nanopowder concentration was varied from 0.1 to 3 mg/cm respectively corresponding to scatterer concentrations from 2.2 10 to 6.7 10 cc,with an independently measured that ranged from 4 to 0.13 mm. The experimental values are estimated by extrapolation from the values obtained by a measure of intensity attenuation of the ballistic beam at very low scatterer concentrations (from 0.003 to 0.03 mg/ml), in samples where the singlescattering condition is well satisfied.
In order to prevent inhomogeneity and clustering of TiO nanopowder, the samples are subjected to a ultrasonic bath before and shuffling by magnetic agitator during measurements. The sample was excited by the second harmonic of a monomode Nd:Yag laser ( = 532.8 nm) with a 10 Hz repetition rate and a time duration of 4 ns, a time comparable to the spontaneous emission decay time of Rodhamine 6G. The pump beam is focused on a spot with diameter of 30 m and with an energy up to 360 J. The shottoshot emission from the sample was collected by an Ocean Optics Spectrometer with spectral resolution of 0.3 nm. An optical fiber can collect the radiation emitted from the sample at different tilting angles with respect to the pump beam direction. The spectra reported in this paper are relative to an angle of 15.
To perform the statistical analysis of the emitted intensity, for each laser shot we recorded six values of the spectra, spaced by 0.5 nm, in a small spectral region around the average peak wavelength for a large number of consecutive spectra, up to 2000 for each pumping energy interval. Great care has been used to monitor the shotbyshot pump pulse energy to measure the emission characteristics as a function of pump energy limited to a small energy interval with a posteriori analysis in order to eliminate the fluctuations of the laser pump pulse as spurious cause of fluctuations in the random laser emission intensity. The energy intervals is variable from 1 to 15 J to follow the critical behavior of the statistical fluctuation characteristics versus pump pulse energy.
As a further experimental condition during the measurements we check also the condition of monomode laser pulse excitation which assures a smooth temporal shape. This issue is also important in order to discriminate intrinsic emission statistics of the random laser emission from strong temporal laser pulse intensity fluctuations that are present in a multimode pump laser.
V Experimental Results
In Fig. 5 typical spectra at the concentration of 0.3 mg/cm (corresponding to a = 1.3 mm) for three different energies of the pumping pulse are reported. These energies are chosen to be sub threshold, near threshold and above threshold. The threshold is determined as the energy value at which a large change of slope of the average peak intensity occurs as a function of the pump pulse energy.
The corresponding intensity histograms are also reported with the fit of the Lévystable distribution and corresponding parameter as done for the theoretical computed spectra. The spectra show that the three different regimes for the intensity fluctuation statistics: Gaussian for the lowest energy ( = 1.96), Lévy with spiky spectrum ( = 1.64) at intermediate energy, and Gaussian again with a smooth spectrum for the highest energy as it was in the numerical simulations. In Fig. 6 are reported spectra for the same three energy levels in the cases of 0.1 mg/cm and 1 mg/cm, corresponding to = 3.9 mm and = 0.4 mm, respectively. In the case = 3.9 mm the system display a transition from a lowfluctuations regime at low energy to a highfluctuations one for both the other two energy values. Is important to note that the value of 360 is the maximum that we could use to prevent damage to the sample, and then a possible return to Gaussian statistics for larger pulse energy. In the case = 0.4 mm the spectra are always smoothed with a corresponding value of approaching 2. The Lévy region is then absent in this case as in the numerical simulations for small . As done in the numerical simulation we quantify this behavior in the experimental results reporting in Fig. 7 the value versus pump pulse energy for three value of the . The curve corresponding to = 1.3 mm shows the three different statistical regimes, from Gaussian for low pump pulse energy to Lévy statistics () and again Gaussian statistics for high pump pulse energy. The more diffusive case, = 0.4 mm, presents a closure of the Lévy regime region with a Gaussian behavior throughout the whole energy range. The case of weak diffusion, = 3.9 mm, presents a rapid transition to Lévy regime that remains for all the values of the pump pulse energy, and is limited by the constraint of preventing sample damage for greater pumping intensities. The trend of in this last case seems indicate that a transition to a Gaussian region could be reached again if a larger energy of the pump pulse was achievable without sample damaging. The behavior of the spectra and value versus pump pulse energy for different values of the show common characteristics and the same evolution throughout statistical regimes that emerged from the results of our theoretical model.
Vi Conclusion
We have studied the statistical regimes of a random laser performing both simulations and experiment. Both the numerical and experimental results consistently demonstrate the presence of three fluctuation regimes upon increasing the pumping energy in a weakly diffusive medium. An initial Gaussian regime is followed by a Lévy one, and Gaussian statistics is recovered again for higher pump pulse energy. This behavior is always present in the condition of weakly diffusive medium (small ).
The situation is different for the strongly diffusive case (large ): Here the region of Lévy statistics is not detectable and the medium always yield a Gaussian regime with smooth emission spectra.
In our study the linkage between theory and experiment emerges in the framework of spectral analysis of random lasers, in particular regarding the origin of random spikes and statistical regimes crossovers, finding that the role of the and the dimension of the active medium are crucial. From this study it emerges also that the excitation time is not crucial for the statistical regime characteristics The experimental and numerical data in Figs. 2 and 6 are theoretically explained by the simple model which has been recently been discussed in Ref. Lepri13 (). Such model arises by a number of simplifying assumptions and amounts to a stochastic partial differential equation for the energy field with contains a multiplicative randomadvection term, yielding intermittency and powerlaw distributions of the field itself Lepri13 (). The analysis of such equation indicate that Lévytype of fluctuations are more likely to be observed close to threshold and for small samples when the mean free path of photons is not too short with respect to the sample size (and more generally to the size of the gain volume). In this case the random advection term dominates the diffusive and gain terms and is responsible for unconventional fluctuations. Its relevance is gauged by the effective noise strength which in the present case should be roughly of the order of . It is clear that the difference observed between samples with small and large reported in Figs. 2 and 6 nicely fit with this prediction. A more quantitative comparison will be addressed in the future. Such an fact is, however, a promising indication that connections between the physics of random lasers and the theory of nonequilibrium phenomena in spatially extended systems are likely to be established.
References
 (1) V.S. Letokhov, Éksp. Teor. Fiz. 53, 1442 (1967) [Sov. Phys. JETP 26, 835 (1968)].
 (2) D.S. Wiersma, Nat. Phys.4, 359 (2008)
 (3) D.S. Wiersma and A. Lagendijk, Phys. Rev. E 54, 4256 (1996)
 (4) H. Cao, J. Phys.38, 10497 (2005)
 (5) X.Y. Jiang and C.M. Soukoulis, Phys. Rev. Lett. 85, 70 (2000)
 (6) P. Pradhan and N. Kumar, Phys. Rev. B 50, 9644 (2001)
 (7) H. Cao, J.Y. Xu, D.Z. Zhang, S.H: Chang, S.T. Ho, E.W. Seelig, X. Liu, R.P.H. Chang Phys. Rev. Lett. 84, 5584 (2000)
 (8) V.M. Apalkov, M.E.Raikh, B. Shapiro Phys. Rev. Lett. 89, 016802 (2002)
 (9) X. Wu, W. Fang, A.Yamilov, A.A. Chabanov, A.A. Asatryan, L.C. Botten and H.Cao Phys. Rev. A 74, 053812 (2006)
 (10) K.L. van der Molen, A.P.Mosk, A. Lagendijk Phys. Rev. A 74, 053808 (2006)
 (11) C. Vanneste, P. Sebbah and H. Cao Phys. Rev. Lett 98, 143902 (2007)
 (12) R. ElDardiry, R. Mooiweer and A. Lagendijk, New J. Phys.14, 113031 (2012)
 (13) S. Mujumdar, M. Ricci, R. Torre and D.S. Wiersma, Phys. Rev. Lett. 93, 053903 (2004)
 (14) S. Mujumdar, V. Türck, R. Torre and D.S. Wiersma, Phys. Rev. A 76, 033807 (2007)
 (15) D. Sharma, H. Ramachandran and N.Kumar, Fluc. Noise Lett. 6, L95 (2006)
 (16) S. Lepri, S. Cavalieri, G. Oppo, D.S. Wiersma, Phys. Rev. A 75, 063820 (2007)
 (17) R. Uppu, A.K. Tiwari and S. Mujumdar, Opt. Lett. 37, 662 (2012)
 (18) R. Uppu and S. Mujumdar., Phys. Rev. A 87, 013822 (2013)
 (19) G. Zhu, L. Gu and M.A. Noginov, Phys. Rev. A 85, 043801 (2012).
 (20) J.P. Bouchaud and A. Georges, Phys. Rep, 195, Issues 45, 127 (1990)
 (21) G. Samorodnitsky and M.S. Taqqu, Stable nonGaussian random processes: stochastic models with infinite variance, Chapman and Hall/CRC, New York (1994)
 (22) I. A. Koutrouvelis , JASA 75, 918 (1980)
 (23) I. A. Koutrouvelis, Commun. Stat.  Simul. Comput. 10, 17 (1981)
 (24) J. P. Nolan, Commun. Stat.  Stochastic Models. Comput. 13, 759 (1997)
 (25) S. Lepri, Phys. Rev. Lett., 110, 230603 (2013)