Application of Shannon's entropy to classify emergent behaviors in a simulation of laser dynamics

Application of Shannon's Entropy to Classify Emergent Behaviors in a Simulation of Laser Dynamics J. L. GUISADO Centro Universitario de Mérida Universidad de Extremadura 06800 Mérida (Badajoz), Spain jlguisado©unex.es F. JIMÉNEZ-MORALES Departamento de Fisica de la Materia Condensada Universidad de Sevilla P.O. Box 1065, 41080 Sevilla, Spain jimenez©us.es J. M. GUERRA Departamento de Optica, Facultad de CC. Fisicas Universidad Complutense de Madrid 28040 Madrid, Spain jmguerra©fis.ucm.es Abstract-Laser dynamics simulations have been carried out using a cellular automata model. The Shannon's entropy has been used to study the different emergent behaviors exhibited by the system, mainly the laser spiking and the laser constant operation. lt is also shown that the Shannon's entropy of the distribution of the populations of photons and electrons reproduces the laser stability curve, in agreement with the theoretical predictions from the laser rate equations and with the experimental results. © 2005 Elsevier Ltd. Ali rights reserved. Keywords-Computational physics, Parallel modelling and simulation, Cellular automata, Com plex systems, Laser physics. l. INTRODUCTION Laser dynamics has been traditionally analyzed by solving a set of coupled differential rate equations which describe the interrelationships and transition rates among the electronic states in the laser active medium and the laser photons [1,2]. But rccently, an alternative approach based Partially supported by the Junta de Extremadura (ConsejerÍa de Educación, Ciencia y TecnologÍa) and the European Social Fund, Grant MOV03Bl 70. Typeset by A_MS-TEX 848 J.L. GUISADO et al. on a cellular automata (CA) model has been proposed [3], in which the laser properties arise as emergent phenomena from the collective effect of many locally interacting, simple components. In this paper we focus on the application of Shannon's entropy to recognize and classify the different types of behavior shown by the CA model. Cellular automata are a class of spatially and temporally discrete mathematical systems, characterized by local interaction and synchronous dynamical evolution. They were first introduced in the 1950s by von Neumann [4] to investigate self-reproduction and since then have attracted much interest because of its capability to generate complex behavior from sets of components which interact locally with relatively simple rules [5,6]. In the last two decades, CA have been extensively used to build models of a wide variety of physical systems [7,8], for example reactiondiffusion processes [9], fluid dynamics [10], magnetization in solids [11], growth phenomena [12], molecular excited-state dynamics [13], etc. Recently, CA have become more attractive because the inherent parallelism makes them very suitable to be naturally and efficiently implemented in parallel computers. High performance simulations of physical systems [14,15] can be carried out in this way. A laser cellular automata model can be useful for two kinds of application: first, for cases in which the standard treatment cannot be applied--as happens in lasers governed by stiff differential equations, with convergence problems--second, to take advantage of the cellular automata intrinsic parallel nature in order to efficiently implement three-dimensional simulations of laser devices in parallel computers. In this work, we study the following problem. In the laser cellular automata model, laser field and population inversion may perform different kinds of temporal behavior, depending on the values of the pumping probability, the life time of laser photons, and the upper laser level life time. The problem is how to find a magnitude to characterize all these possible outcomes of the system to be compared with the predictions of the traditional approach for laser dynamics (the laser rate equations) and with the experimental results. In order to address this problem, we choose the Shannon's entropy as such magnitude. In the laser CA model, two main kinds of temporal behavior are found: an overdamped dynamics with almost null Shannon's entropy and an oscillating one with larger Shannon's entropy. With the Shannon's entropy, we locate these different kinds of behavior and check the agreement between the laser CA model results and the predictions of the laser rate equations. This paper is organized as follows. Section 2 describes the cellular automata model that has been used for the simulations. In Section 3, the results of the simulations are presented and analyzed using the Shannon's entropy. Finally, the conclusions are summarized in Section 4. 2. CELLULAR AUTOMATA MODEL The model that has been used is a cellular automaton defined on a two-dimensional square lattice of Arc -- 200 × 200 cells with periodic boundary conditions. Two variables ai(t) and ci(t) are associated to each node of this lattice. The first one, ai(t), represents the state of the electron in node i at a given time t. An electron in the laser ground state takes the value a~(t) = 0, while an electron in the upper laser state takes the value ai(t) -- 1. The second variable ci(t) C {0, 1,2,..., M} represents the number of photons in node i at time t. A large enough upper value of M is taken to avoid saturation of the system. The neighborhood considered is the Moore neighborhood, each cell having nine neighbors: the ceil itself, its four nearest neighbors (situated in the positions north, south, east, and west) and the four next neighbors (in the positions northeast, southeast, northwest and southwest). The time evolution of the CA is given by a set of transition rules which determine the state of any particular cell of the system at time t + 1 depending on the state of the cells included in its neighborhood at time t. These rules represent the different physical processes that work at the microscopic level in a laser system. Application of Shannon's Entropy 849 R1. Pumping: If the electronic state of a cell has a value of ai(t) = 0 in time t, then in time t + 1 that state will have a value of ai(t + 1) = 1 with a probability ,k. R2. Stimulated Emission: If, in time t, the electronic state of a cell has a value of ai(t) = 1 and the sum of the values of the laser photons states in the nine neighbor ceils is greater than a certain threshold (which in our simulations has been taken to be 1), then in time t + 1 a new photon will be created in that cell: ci(t + 1) = ci(t) + 1 and the electron will decay to the ground leveh ai(t + 1) = 0. R3. Photon Decay: A finite life time ~-~ is assigned to each photon when it is created. The photon will be destroyed ~-~ time steps after it is created. R4. Electron Decay: A finite life time ~-a is assigned to each electron that is promoted from the ground level to the upper laser level. That electron will decay to the ground level again ~-~ time steps after it was promoted, if it has not yet decayed by stimulated emission. To simplify the model as much as possible, we consider this decay is entirely nonradiative, i.e., we do not take into account spontaneous emission. Also, as in an ideal four level laser the population of level E1 is negligible, stimulated absorption has not been considered. In addition, a small continuous noise level of random photons in the laser mode is introduced at every time step, in order to represent the experimentally observed noise level, responsible of the initial laser start-up. This is done by making c~(t + 1) = c~(t) + 1 for a small number of cells (< 0.01% of total) with randomly chosen positions. . SIMULATION RESULTS AND SHANNON'S ENTROPY ANALYSIS Three parameters determine the response of the system: the pumping probability (A), the life time of photons (~-c) and the life time of excited electrons (Ta). Initially, ai(O) = O, ci(O) = O, Vi, except a small fraction 0.01% of noise photons. We let the system evolve for 500 time steps. In each time step, the total number of laser photons n(t) N = ~i21 ci(t), and the total number of electrons in the upper laser state (population inversion) N(t) = ~N~ 1 ai(t) are measured. Running test simulations for different values of the three parameters of the system, two main types of behavior are observed: after a transient time, n(t) and N(t) show either a constant value or correlated damped oscillations. We are interested in making a more systematic exploration of the behavior of the laser CA model in the whole parameters space, and comparing it with the response of the laser rate equations. To this end, it should be interesting to find a magnitude to characterize the type of behavior shown by the system for each particular triad of values of the parameters. As such a magnitude, the Shannon's entropy S of the distribution of values taken by n(t) and N(t), after running the simulation for a time interval, is calculated. This Shannon's entropy is computed by dividing the range of values taken by n or N in 10 a equally spaced bins, and computing the frequency (fi) at which this value lies inside every particular nonvoid bin i. Then, S is calculated as s(a, = - logs (1) i=1 where rn is the number of nonvoid bins. The Shannon's entropy S measures the dispersion of the distribution of values taken by the nmnber of laser photons ~ (or the population inversion N). If this number is approximately constant (i.e., all the bins except one are void), S will tend to zero; if oscillations appear, the probability distribution becomes wider and S takes higher values; and finally the maximum value of S would result from an equiprobable distribution. Therefore, S can be used as a good indicator of the presence of oscillations in the system. The dependence of the Shannon's entropy on two of 850 0,10 J. L. C-UISADO et al. =10 C 0,08 0,06 0,04 0,02 0,10 ,,,.,., ,,,,:',... , :, , 20 40 60 80 100 120 140 160 180 200 1; a ,c =10 c S C 2,0 ~ 2,5 'i,nl 3,o 4'0 0,08 S a °~°81i~iii~U~2ii~iiiiiiiiiiiiii~iiiiii~iiiiiiiiiiiiiii~ii~ii~iiiiiiiii~iiii~iiiiiiiiiiiii~iii~i~ii;iiii~i~ii~i~:~ ......... ~ii~i~!ii~!~i~iii~:~i~iiii~i~i~i~i~ii!ii~i~i~iii~iii~i~!i?~;~i~il :!~,,~i~!ii~i~,: ~° iiiiiiiiiiiiiiiiiiiiiiii)!i~iiiiii 3,o ~,!i!i!i!i!i!!!:~4,0 00, i=iliN t 6,0 7,0 0,02 8,0 20 40 60 80 100 120 140 160 180 200 I; a Figure 1. Contour plots of the Shannon's entropy of the distribution of values taken by the number of laser photons Sc (top) and the population inversion Sa (bottom), showing the dependence of S with the pumping probability ,k and the upper laser level life time ~-a for a fixed value of the cavity life time of Tc = 10. ~-a and Tc are measured in time steps. S has been computed for a period of time: t E [250, 1000] time steps. the parameters of the system when the third parameter is fixed is shown in Figure 1. Similar plots can be obtained for other values of the parameter that is fixed (~'c), showing areas in which S tends to zero, indicating a constant behavior (bright zones), and other areas in which S is higher, indicating an oscillatory behavior (dark zones). The results obtained by measuring S can be compared with the predictions of the laser rate equations for a wide range of values of the three parameters of the system. The laser rate equations can be put in their simplest--but still meaningful--form as a system of two coupled Application of Shannon's Entropy 851 20 18 iiiiil 16 iliiil 14 12 10 8 6 4 2 0 -- 0 Figure 2. S 0 °,6667 1,333 ~ ,~,~,~,~,;':-'..!~ ~ 2 ooo ~ 2,667 ~ 3,333 2 4 6 8 10 12 14 16 18 20 R/R t Contour plot of the Shannon's entropy of the distribution of the number of laser photons for a fixed value of % = 10 time steps. Low values of S~ (bright zones) indicate that the response of the system is nonoscillatory, while high values (dark zones) correspond to an oscillatory response. The black line is the theoretical stability curve. Points a and b correspond to the parameters for which the temporal evolution is shown in Figure 3. (a): -R/Rt = 5, Ta/~-c = 18. (b): R/Rt = 16, ~-a/~-c = 3. differential equations [1,2]. The first one gives the variation with time of the number of laser photons n(t) and the other one reflects the temporal variation of the population inversion N(t), dn(t) - K N(t) n(t) - ~(nJ--Z, (2) dt ~'~ dN(t) _ R N(t) KN(t) n(t). (3) dt ~-~ Here, R is the pumping rate and K is a constant called "coupling constant". For the case of small amplitude fluctuations, a linearized small-signal analysis of these equations can be carried out, and two different situations arise. 1. Overdamped System: The solutions for n(t) and N(t) are real exponentials. The response to any perturbation dies out exponentially towards the steady state. 2. Relaxation Oscillations: The solutions have an exponentially damped sinusoidal form. The response to any perturbation shows relaxation oscillations towards the steady state, with a fi'equency (usually called spiking frequency), a~sp= \27aRt] TaTc Rt -1 ' (4) Here, R is the laser pumping rate and Rt is the threshold laser pumping rate, which are linearly related to the pumping probability A and the threshold pumping probability At that appear in the CA model [3], so that R/Rt = A/At. The threshold pumping probability At is calculated in the CA model as the smallest value of the pumping probability A for which after a transient time the number of laser photons is clearly greater than the number of noise photons introduced. 852 J.L. GUISADO et al. o = Q. O.. 100000 90000 80000 70000 60000 50000 40000 30000 20000 10000 0 0 300000 250000 x~ × 200000 x × × X -°°°° I00000 r ,x 200 400 600 800 1000 Time Steps (a) I Population Inversion + I Laser Photons x 50000 0~ I I I I 0 200 400 600 800 Time Steps 1000 (b) Figure 3. Temporal evolution of the system (number of laser photons and population inversion) for the values of the parameters marked as a and b points in Figure 2. The values of the parameters are (a): A = 0.0125, rc = 10, r~ -- 180. (b): A = 0.192, ~-c = 10, ~-a = 30. Here, re and T~ are measured in time steps. Lattice: 400 x 400 cells. The necessary condition for the occurrence of the relaxation oscillations is that Wsp be real, or equivalently, that the quantity inside the square root in equation (4) be positive. This leads to the following condition, __ (R/Rt) 2 ~° > (5) T~ 4 (R/Rt - 1)" Application of Shannon's Entropy 853 Damped oscillations appear for values of the parameters fulfilling condition (5), and a constant behavior appears if this condition is not satisfied. The values of the parameters for which these two situations appear can be visualized by plotting the equation derived from condition (5) with an equal sign, using R/Rt as x-axis and ra/rc as y-axis. The result is the black line in Figure 2. This function is the theoretical stability curve which separates areas of different behavior for the solution of the linearized small-signal laser rate equations: relaxation oscillations appear above and to the right of this curve and an overdamped behavior below and to the left of it. In order to compare the behavior exhibited by the results of the laser cellular automata model with these predictions from the laser rate equations, we have plotted in the same graph in Figure 2, the Shannon's entropy for the number of laser photons, using R/Rt as x-axis and ra/r~ as y-axis. It can be observed that the Shannon's entropy follows the form of the theoretical stability curve. Areas of high S (dark zones)--which indicate an oscillatory behavior in the results of the laser CA mode~appear above and to the right of this curve, and areas of low S (bright zones)--which indicate a constant behavior in the results of the laser CA model--appear below and to the left of it. This is in agreement with the solutions of the laser rate equations. Figure 3 shows the time evolution of n(t) and N(t) corresponding to the points a and b in Figure 2. Relaxation oscillations (also known as laser spiking) are found for point a, included in the zone of high entropy, and a constant behavior for point b, included in the zone of low entropy. The information obtained by the Shannon's entropy not only is in agreement with the results of the laser rate equations, but is richer than that obtained only from the stability analysis of those equations, because it is able to quantify the amplitude of the oscillatory behavior. Therefore, the Shannon's entropy can be used as a parameter to characterize the behavior exhibited by the laser CA model and to compare it with the dynamics predicted by the laser rate equations. 4. CONCLUSIONS Along this paper it was shown the usefulness of the Shannon's entropy concept to test the agreement between the different behaviors observed in the laser cellular automata model with those predicted by the stability analysis performed in the laser rate equations. This agreement reinforces the confidence in the capability of laser cellular automata as a modelling tool alternative to differential equations--to describe laser physics. REFERENCES 1. A.E. Siegman, Lasers, University Science Books, (1986). 2. O. Svelto, Principles of Lasers, Plenum Press, (1989). 3. J.L. Guisado, F. Jim~nez-Morales and J.M. Guerra, A cellular automaton model for the simulation of laser dynamics, Phys. Rev. E 67:066708, (2003). 4. J. von Neumann, Theory of Self-Reproducing Automata, University of Illinois Press, Urbana, IL, (1966). 5. S. Wolfram, Universality and complexity in cellular automata, Physica D 10 (1), (1984). 6. S. Wolfram, Cellular Automata and Complexity, Addison-Wesley, (1994). 7. B. Chopard and M. Droz, Cellular Automata Modelling of Physical Systems, Cambridge University Press, (1998). 8. T. Toffoli and N. Margolus, Cellular Automata Machines: A New Environment for Modelling, The MIT Press, (1987). 9. B. Chopard, P. Luthi and M. Droz, Reaction-diffusion cellular automata model for the formation of Liesegang patterns, Phys. Rev. Lett. 72 (9), 1384 1387, (1994). 10. A.K. Gunstensen, D.H. Rothman, S. Zaleski and G. Zanetti, Lattice Boltzmann model of inmiscible fluids, Phys. Rev. A 43, 4320-4327, (1991). 11. G. Vichniac, Simulating physics with cellular automata, Physica D 10, 96-115, (1984). 12. J. Gravner and D. Griffeath, Cellular automaton growth on Z2: Theorems, examples and problems, Adv. Appl. Math. 21, 241-304, (1998). 13. P.G. Seybold, L.B. Kier and C.-K. Cheng, Aurora borealis: Stochastic cellular automata simulations of the excited-state dynamics of oxygen atoms, Int. Y. Quant. Chem. 75, 751-756, (1999). 854 J.L. GUISADO et al. 14. P.M.A. Sloot, J.A. Kaandorp, A.G. Hoekstra and B.J. Overeinder, Distributed simulation with cellular automata: Architecture and applications, In SOFSEM'99: Theory and Practice of Informatics, Lecture Notes on Computer Science, Volume 1725, (Edited by J. Pavelka, G. Tel and M. Barto~ek), pp. 203-248, (1999). 15. D. Talia, Cellular processing tools for high-performance simulation, IEEE Computer 33 (9), 44-52, (September 2000).