{'fraction_non_alphanumeric': 0.04365160498875485, 'fraction_numerical': 0.01809445921079534, 'mean_word_length': 4.458891369047619, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Nonclassical properties and Anderson localization of quantum states in coupled waveguides

We consider the propagation of light beams through disordered lattices of coupled waveguides searching for Anderson localization and investigating the evolution of nonclassical properties of injected quantum states. We assume that the beam is initially in a variety of states, such as the complementary coherent state, the reciprocal binomial state, and the polynomial state. The statistical properties of the evolved states were analyzed numerically as functions of the localization/delocalization parameters averaged over many realizations of disorder. We also numerically reconstruct the Wigner function of the output state. Interestingly, we find that high values of the disorder tend to preserve quantum properties of some input states when we look at the input waveguide despite the coupling between it and the neighboring waveguides.


I. INTRODUCTION
Nonclassical properties of quantized light field and its generation -i.e., the "quantum states engineering"are essential ingredients of quantum optics and quantum information. Indeed, nonclassical states perform a crucial role in many potential applications such as quantum teleportation [1][2][3], quantum cryptography [4], quantum computation and quantum communication [5], quantum internet [6], etc. In this direction, periodic photonic systems has emerged as a platform to engineer new light field structures, presenting numerous significant technological advances [7][8][9][10]. On the other hand, the study of nonperiodic photonic structures, by using small defects in periodic lattices [11][12][13] or disordered and quasi-periodic structures [14][15][16][17][18], has demonstrated a diversity of optical effects in the presence of Anderson Localization, such as engineering of random lasers [19] and structurally colored materials with precisely controllable wavelength and angular dependence of scattering [20]. Recently, it was demonstrated how to tune and freeze disorder in photonic crystals by using percolation lithography [21].
Recently, Anderson localization of nonclassical light was investigated for propagation in an array of waveguides in which neighboring waveguides are evanescently coupled and disorder can be added in a controlled manner [40]. Specifically, that work investigated the consequences of using both sub-Poissonian and super-Poissonian input light on the characteristics of Ander-son localization, verifying the enhancement in fluctuations of localized light and superbunching due to the medium's disorder. Also, an important consequence of sub-Poissonian statistics of the incoming light is to quench the total fluctuations at the output [40]. The system employed in Ref. [40] is similar to that used in Ref. [16] to experimentally investigate the evolution of linear and nonlinear waves in the presence of Anderson localization. Moreover, the disordered one-dimensional waveguide lattice was also used to experimentally investigate an extensive list of phenomena: the signature of a localization phase transition for light by directly measuring wave transport inside the lattice [41], quantum correlations between noninteracting particles evolving simultaneously in a disordered medium [42], Hanbury Brown and Twiss correlations of Anderson localized waves [43], the control of the polarization state of coherent light propagating through an optically thick multiple scattering medium by controlling only the spatial phase of the incoming field with a spatial light modulator [44], the coherent manipulation of two-photon path-entangled states by multimode interference in multimode waveguides [45], the observation of topological phase transitions in photonic quasicrystals [46], the observation of ensembleaveraged quantum correlations between path-entangled photons undergoing Anderson localization [47], the twophoton Anderson localization in a quadratic waveguide array with the emergence of off-diagonal disorder [48], etc.
Here, inspired by the results obtained in Ref. [40], we numerically investigate the propagation of light beams, previously prepared in nonclassical states of the electromagnetic field, propagating in disordered lattices of waveguides and undergoing Anderson localization. Our goal is to verify the influence of Anderson localization on the statistical properties of previously prepared input light field (nonclassical states). To this end, we assume the beam in a variety of states, namely, the complementary coherent state [49], reciprocal binomial state [50], polynomial state [51], thermal state, coherent state, and squeezed state [52].
The rest of the paper is organized as follows. The theoretical model is considered in next section, where we present the dynamical model of the system in Subsec. II A, the characterization of the quantum states under consideration in Subsec. II B and the numerical methods in Subsec. II C. In Sec. III we present the numerical results and our analyzes. We conclude the paper in Sec. IV.

II. LINEAR ARRAY OF WAVEGUIDES
The system we consider is a one-dimensional finite array of monomodal waveguides where prescribed input states of the electromagnetic field can propagate. Classically, fields propagating in the array of waveguides are coupled through evanescent waves passing over their boundary-barriers; on the quantum level, one says that photons can coherently tunnel between neighboring waveguides so that the quantum state represents the overall-overlapping superposition of the modes of the waveguides [16]. The field can be injected into one or a few waveguides and disorder can be implemented in the array either by randomly adjusting the spacing among the parallel waveguides along the x-direction, the propagation being in the z-direction, or by randomly fixing the thicknesses of the waveguides. This kind of system has been constructed on an AlGaAs substrate [53,54] and direct identification and measurements of Andersonlocalization of states have been performed.

A. Theoretical model
The electromagnetic energy density in the array of waveguides, assuming that all media are linear and nonmagnetic, and that the relevant evanescent overlap occurs only between neighboring waveguides, is given by where n is the refraction index and (j, j ′ ) denote nextneighbor pairs. We consider monochromatic fields propagating along the z-direction with velocity c/n. At the quantum level, the field mode in the j-th waveguide is written in terms of photon-annihilation and photon-creation operators, a j and a † j respectively, and the evolution of the system is dictated by a Hamiltonian in the form [42] where β j is the propagation constant associated with the j-th waveguide and C j+1,j and C j,j+1 are coupling coefficients between nearest neighbor waveguides. Notice that, we are considering fields with lower enough intensities to make non-linear effects negligible; also, along this paper unless stated in contrary, we use = 1 and c = 1. The creation and annihilation operators satisfy the commutation relations and we assume the existence of eigenstates |ψ , such that n j |ψ = a † j a j |ψ = n j |ψ , where n j is the photon number of the j-th waveguide.
We will analyse arrays where we can fix a constant coupling (tunneling rate) between neighbor waveguides, C j+1,j = C j,j+1 = C, and introduce disorder by taking the coefficients β j as random variables with zero-mean Gaussian distributions; the model then becomes isomorphic to the one-dimensional quantum tight-binding model used by Anderson [22], with β j being the on-site energy, and the system should then present localization of states. This assertion is experimentally feasible since these coefficients are related with the guide geometry, which can be appropriately adjusted [41,43] to spatially modulate the index of refraction n(x). In this case, the Heisenberg equations can be written as [16,41,43,47,55] where z = ct/n, i.e. measurements of intensity distribution at position z give the time evolution along the array. Now we search for a solution of the Heisenberg equations (4), depending on the initial input state. Since the Heisenberg equations are linear in the annihilation (or for the creation) operators, it can be solved by finding the Green's function in such way that where a l (0) correspond to the input state (at z = 0) into the l-th waveguide. The Green's function correlates fields in the j-th and the l-th waveguides at all positions z. By inserting Eq. (5) into Eq. (4), one gets the following set of first order differential equations for the Green's functions, From our assumptions, G jl (z) depend on the parameters β j , which vary randomly, and C that remains fixed; this is usually referred to as diagonal disorder. These equations can be solved numerically with great precision; taking specific distributions of {β j } and a given value of C, solutions are obtained just depending on the initial input state.
We can work with a great simplification if we consider that the input field |Ψ in is injected into only one waveguide [40], which we label by j 0 , that is, a j (0)|Ψ in = 0 for all j = j 0 . Thus, considering the light injection only into the j 0 -th waveguide, the mean field intensity output by waveguide j as a function of the lattice length z is given by where the mean value of the Green's functions is a standard statistical mean over several realizations for different values of {β j }, while the mean a † j0 a j0 represents the quantum expectation value of the number operator a † j0 (0)a j0 (0), which depends only on the input state. Also, we can calculate the correlation of the intensities between the output of two waveguides j and l, to be given by Another important quantity to qualify the statistics of photons propagating in the array is the second-order correlation function defined, for the j-th waveguide, by The correlation function g (2) (z) indicates weather states evolve in a Poissonian way, or if they present either bunching or antibunching (g (2) (z) > 1 and g (2) (z) < 1, respectively) in their photon distribution. In any case, all these quantities depend on the initial state in the input of the array of waveguides. The expressions for other quantities used to characterize the output field are presented in the results section.

B. Characterization of the input quantum states
In Ref. [40], a theoretical study of Anderson localization of light in an array of waveguides was presented, using coherent, thermal and squeezed states at the input, to investigate the effects of nonclassicality. Here, we extend the work of Ref. [40] by considering different input states that show more general quantum statistics, and by presenting the evolution in the Wigner representation. We shall consider specifically the complementary coherent state (CCS) [49], the reciprocal binomial state (RBS) [50] and the polynomial state (P S) [51] as input states; these states present peculiar statistical properties and it is interesting to investigate how they evolve along the array.
The complementary coherent state [49] is written in the number basis as where the normalization constant ℵ CCS is given by For these states we find the mean number of photons, and the mean a †2 j0 a 2 j0 is given by The reciprocal binomial state [50], written as has normalization constant ℵ RBS given by For the RBS, we find The polynomial state [51] is defined by where H N −k (y) is a Hermite polynomial, with the normalization constant ℵ P S written as For the polynomial state we find For sake of comparison, we also investigate the input states used in Ref. [40]: thermal states (T S), mixed states with density matrix for which a † j0 a j0 T S =n and a †2 j0 a 2 j0 T S = 2n 2 ; coherent states (CS), for which a † j0 a j0 CS = |α| 2 and a †2 j0 a 2 j0 CS = |α| 4 ; and squeezed-vacuum states, with squeezing parameter ζ = |ζ|e iφ , for which a † j0 a j0 SS = sinh 2 |ζ| and a †2 j0 a 2 j0 CS = sinh 2 |ζ|(1 + 3 sinh 2 |ζ|).

C. Computational method
Since we are assuming that the input light is injected only in the j 0 -th waveguide, all relevant quantities depend only on the Green's functions G j,j0 , which satisfy To get solutions numerically, these equations were discretized and solved by using the Crank-Nicholson method [56] with step size ∆z = 0.001. In the numerical analysis, we consider an array with 101 waveguides, which we label from 1 to 101 with light injected only in the waveguide j 0 = 51. Thus the initial condition is that G 51,51 (0) = 1 while G j,51 (0) = 0 for all j = 51; additionally, since the waveguide array is finite containing M waveguides (here 101), we add the boundary conditions G 0,l (z) = G M+1,l (z) = 0 for all l = 1, 2, ...., M . Following Ref. [40], we assume that the random coefficients β j are independent of each other and follow a zero-mean Gaussian probability distribution of the form with the variance ∆ 2 measuring the disorder in the waveguide array. For sake of simplicity but without loss of generality, we fix the interaction parameter C = 1 and take ∆/C quantifying the arrangement disorder.
In order to get the β j coefficients, we have used the Box-Müller method [57] that requires the generation of two random numbers with uniform distribution, which was done using a congruence method. For each set {β j } of disorder parameters, Eqs. (25) are solved numerically (fixing C) to find the relevant Green's functions. The functions G j,j0 (z) depend on the disorder parameters β j and on the coupling parameter C, and completely describe the dynamical evolution of any state injected into the waveguide array, only through the 51 th waveguide. Naturally, to average over the sets of the random coefficients β j one has to consider a reasonable number of realizations of disorder, different sets {β j }; here we take 1000 realizations of disorder in each simulation.
In our study, we take all the input states with the same energy, that is the same mean number of photons, specifically n in (0) = a † 51 (0)a 51 (0) = 10. To present the results obtained by our simulations, we shall choose five states from the families we have shown in Subsec. II B, namely: two CCS states, named CCS 1 (= |Ψ CCS (0.1414, 10) ) and CCS 2 (= |Ψ CCS (1.916, 11) ); two P S states, referred to as P S 1 (= |Ψ P S (0.374, 12) ) and P S 2 (= |Ψ P S (0.9345, 13) ); and, in the case of the RBS state, |Ψ RBS (φ, N ) , we set N = 20 and φ = 0. For completeness and comparison, we also consider as input states a thermal state (T S) withn = 10, a coherent state (CS) with |α| 2 = 10 and a squeezed-vacuum state (SS) with sinh 2 |ζ| = 10.

III. NUMERICAL RESULTS
As explicitly shown in Eq. (7), which is valid whenever light is injected into the j 0 -th waveguide, the average output intensity of the waveguide array depends on the average number of photons of the input state (here fixed as n in = 10) and on the degree of disorder of the system, carried by the random β-coefficients, which is manifested by the Green's functions; the output intensity, that is the mean number of photons at the end of the array, does not depend on any other characteristics of the input state.
The average output intensity, as distributed among the 101 waveguides of the array, is shown in Fig. 1 for some values of ∆/C and two values of propagation distance z. As mentioned before, ∆/C measures the arrangement disorder, i.e., as the value of ∆ increases the coefficients β j become more distinct from each other. In the absence of disorder, (Fig. 1(a)), the light suffers only the standard dispersion. However, as the disorder is increased ( Fig. 1(b-d)), a narrow peak of intensity around the input (j 0 -th, number 51) waveguide emerges, which characterizes the localization of the solution. For ratios ∆/C greater then 1.5, the output-intensity profile does not change anymore with the increasing of the propagation distance, which is shown by the exact overlapping of the solid-line (z = 5) and dashed-line (z = 20) in Fig. 2.
In Fig. 3 we present, using a logarithmic scale, the mean number of photons for some values of disorder (∆/C) and considering two values of propagation distance, z = 5 and z = 20 respectively. The choice of using a logarithmic scale is justified by the fact that an ex- ponential decay of the light intensity, in the waveguides different of the j 0 -th one, would signalize the appearance of Anderson localization. We find that, in all cases, there exist two regions of decreasing exponentials and, except in case ∆/C = 1.0, the pattern holds for both propagation distances. It should be mentioned, however, that plots of the output intensities (mean number of photons) may not be sufficient to state whether the light beam presents Anderson localization, since the system might have yet a small diffusion, which would become more evident by increasing the propagation distance.
In order to further investigate the occurrence of Anderson location, we can also calculate the participation number defined by, and given in our case by, We see that the participation number is not only independent of the mean number of photons of the input state, but it is actually completely independent of the electromagnetic field mode injected in the j 0 -th waveguide of the array, the number 51; it depends only on the waveguide array itself and its disorder encoded in the Green's functions. The participation number indicates in how many waveguides there are photons as a function of propagation distance; thus, if P(z) increases, it means that the beam remains scattering among the waveguides along propagation. It should be emphasized that, in the case of absence of disorder, the dispersion of P(z) is linear, as indicated by Fig. 4(a). Also, for a very small amount of disorder, ∆/C = 0.1, we find that the diffusion occurs faster than in the absence of disorder, while for ∆/C = 0.5 dispersion still occurs, but with a growth rate smaller than that in the case of no disorder. In all other cases reported in Figs. 4(b-d), for larger values of ∆/C, the plots of P(z) stabilise as z increases, characterizing localization of the light in the array. Notice that some fluctuations are still observed, but on a small scale when compared to 1, which is the smallest value that can characterize a single waveguide. Figures 5 and 6 show the second order correlation function g (2) of the field output by waveguide number 51, as a function of the disorder degree for different input states, calculated with Eq. (9), taking propagation distances z = 5 and z = 20, respectively. We first notice that, for all input states analysed, taking measurements for a short propagation length (z = 5, Fig. 5), the second order correlation function increases, as the disorder degree (∆/C) raises, it reaches a maximum for ∆/C = 1 and then decreases to values greater than the value of g (2) without disorder. Similarly, for longer arrays (z = 20, Fig. 6), one finds peaks of g (2) somewhat higher than those of Fig. 5 but occurring at a much lower disorder degree and rapidly decaying for large disorder degree. Yet interestingly, for the propagation distance z = 20, we find that only the CCS 1 state returns to the g (2) < 1 regime for large values of disorder. Note also that the curves for the CCS 2 and the RBS states are nearly coincident due to the fact that, for these states, a †2 51 (0)a 2 51 (0) are very close, 94.78 and 95, respectively.
The variance of the output field intensity, (∆I) 2 = I 2 − I 2 , of the 51 th waveguide, when the input field is injected only in it, is obtained directly from Eqs. (7) and (8) as In order to investigate the influence of the disorder on the intensity variance of the output states, in Figs. 7(a) and 7(b) we show these variances as functions of the disorder degree ∆/C, at z = 5 and z = 20 respectively, for different input states injected into the 51 th waveguide. We see that, for all states discussed, when the degree of disorder increases, the variances increase, some of them tending to stabilize, as for states P S 1 , P S 2 and CS, while others, like the cases of the states CCS 1 , CCS 2 and RBS, reach a maximum and then decrease for large values of the disorder parameter, the decreasing rate being bigger for the state CCS 1 . It is interesting to notice that the patterns of (∆I) 2 practically do not change when one compares the results for the propagation lengths z = 5 and z = 20, Figs. 7(a) and 7(b) respectively. Here, likewise Figs. 5 and 6, the curves for the CCS 2 and the RBS states are almost coincident. We display the states CCS1 in circles (blue), CCS2 in boxes (black), RBS in pentagon (cyan), CS in triangles (red), P S1 in down triangles (grey), and P S2 in diamond (green), respectively.

A. Wigner representation
We can also study the influence of disorder in the propagation of light in the waveguide array using the Wigner representation of the quantum states. The Wigner function (WF) of the output state in the 51 th waveguide is given by with the symmetrically ordered characteristic function given by where a 51 (z) = G 51,51 (z)a 51 (0) and the input state is given by ρ 51 (0) = |ψ 51 ψ 51 |.
We take |ψ 51 = N n=0 c n |n , which is the form of the states described in Subsec. II B. Then, the integral defining W 51 (α, z) can be performed analytically, and writing α = x + iy, we obtain where with p = l + s and q = k − l + j − s and where the summation is over even integers and e(r) is the largest even integer not greatter than r, i.e. e(r) = r if r is even and e(r) = r − 1, for r odd. We can use the WF, given by Eq. (31) to investigate the propagation of states through the array, looking at the output state in the 51 th waveguide, when the input is a truncated state in the number basis, N n=0 c n |n ; Figs. 8 and 9 show examples of this case. In Fig. 8, we present the WF when the input state is the CCS 1 , in Fig. 8(a) for the input state (z = 0), and in Figs. 8(b-d) for the output states in the 51 th waveguide, for the length z = 20, considering the disorder parameters ∆/C = 0, 1, and 7, respectively. We find that, in absence of disorder (∆/C = 0) ( Fig. 8(b)), the WF is similar to the one of the vacuum state; this is also confirmed by the photon number distribution, as shown below. But, for small values of the disorder parameter, Figs. 8(b-c), the WF profiles present very important differences when compared with that for the input state, showing that dispersion prevails. However, by increasing the disorder parameter, the WF of the output state becomes similar to that one for the input state; clearly this behavior is directly linked to a strong localization, but not necessarily of Anderson type.
In Fig. 9(a), we present the WF for the input state P S 1 and, in Figs. 9(b-d), for the output state of the 51 th waveguide at z = 20 with the same values of the disorder parameter as in Fig. 8, ∆/C = 0, 1, and 7, respectively. Here, we observed that, similarly to the results obtained for the CCS 1 (Fig. 8 ( Fig. 9(b)) or a small disorder (here characterized by the value of the disorder parameter ∆/C = 1, Fig. 9(c)) the output state of the 51 th waveguide at z = 20 presents a configuration for the WF similar to that of the vacuum state. Clearly, the weak disorder is not enough to ensure that the WF remains in the same shape, due to the interaction of the main guide (51 th ) with the neighboring guides. On the other hand, as the disorder parameter increases (for example with ∆/C = 7 displayed in Fig. 9(d)) the WF is now preserved, which can be observed when comparing the panels (a) and (d) of Figs. 8 and 9. In other words, by increasing the disorder parameter, the WFs become more robust to changes due to variations of the coupling between neighbor waveguides. Also, it should be pointed out that, in all cases, due to the computational cost the number of realizations to produce the output state was reduced to 100. Finally, we can also use the WF (Eq. (31)) to obtain the photon number distribution (PND) of the output state of the 51 th waveguide, given by where W n (α) is the WF for the number state |n . The integral in Eq. (34) is calculated numerically and, as before, the number of realizations to get the proper averages was reduced to 100. In Fig. 10 we show the PND in the 51 th waveguide when the input state is the CCS 1 and the length of the array is z = 20. In Fig. 10(a), we plot together the PND of the input CCS 1 state and the output PND (P (N ) 51 (n, 20)) for null disorder, for comparison. Interestingly, the PND of the state CCS 1 , which has n CCS 1 = 10, is very close to that of the number state  Fig. 8, but now considering the P S1 as the input state.  Fig. 10 but now considering the input state being the P S1.
|10 , P |10 (n) = δ n, 10 . On the other hand, the output state of the 51 th waveguide (at z = 20), in absence of disorder, is close to that for the vacuum state; this fact reinforces the results for the Wigner functions presented in Figs. 8(b) and 9(b). In Fig. 10(b) we show the PND for the output state in the 51 th waveguide considering three different values of the disorder parameter, ∆/C = 0, 1, and 7. We clearly see that the increasing of disorder tend to favor the output PND to become closer to that of the input state, that is disorder tends to preserve the PND.
In Fig. 11 we present the results of the PND, for the same parameters as in Fig. 10, but considering the input state given by the P S 1 . Now we see, from Fig. 11(a), that the PND of the input P S 1 state differs significantly from the PND for the number state, distinctly to the case of the CCS 1 , although we still have n P S1 = 10. However, also in agreement with the case of the CCS 1 , we see in Fig. 11(b) that the PND of the output state is clearly closer to that of the input state as greater is the value of ∆/C.

IV. CONCLUSION
We have discussed the evolution of quantum states of the electromagnetic field propagating through a disordered plane waveguide array. Specifically, we have analysed the propagation of three truncated states (in the number basis), namely the complementary coherent state, the reciprocal binomial state and the polynomial state; we also considered a thermal state, a coherent state and a squeezed-vacuum state, for comparison.
In our numerical calculations, we considered arrays with a 101 waveguides focusing in the injection and detection of waves in the middle one, the 51 th -waveguide. The in site energy parameters β j were assumed to be independent of each other and randomly taken following a zero-mean Gaussian distribution with variance ∆ 2 , while the coupling between neighbor waveguides were fixed as constant, C; the disorder parameter was defined as ∆/C.
We have investigated some quantities that qualify the nature of the statistics of the state. First, we verified that increasing the degree of impurity of the lattice, the mean photon-number distribution tends to concentrate around the waveguide where the input state is injected, with an exponential decay of the light intensity in the others, for any injected state with a given mean number of photons, here n in = 10; although this is not a definitive indication of localization, it does signalize it.
We also analysed the participation number, indicating in how many waveguides there were photons, which shows the effects of disorder in the array and is totaly independent of the input state. We find, as the disorder parameter is increased, from no-disorder to a high disorder regime, the participation number changes from a linear increase with the propagation distance, characteristic of dispersion, to a flat behavior with very small fluctuations which represents a localization pattern.
We investigated second order quantities as the g (2) function and the output intensity variance at the central waveguide. We observe for all input states that, although localized, the average output state presents a classical behavior relatively to the bunching feature. However, as the disorder increases, the g (2) function decreases and even presents antibunching again for the complementary coherent state and high disorder. While the variance in the output intensity at the central waveguide increases with disorder, and stabilizes to a final value for most of the input states investigated.
In order to investigate how the propagation through the lattice changes the input state, it is not enough to analyse quantities quantum averaged over the input state. It can be noticed, for example, when we look at the results for CCS2 and RBS, two different states with different features, however presenting nearly the same values of g (2) and intensity variance at the output. To tackle this question we reconstructed the output state by means of its Wigner function, which could also be used to obtain the output probability distribution of number of photons. We observed a preservation of the characteristics of the Wigner function for high disorder parameter values, although we also notice a suppression of negative values in average. Moreover, as the output probability distribution of number of photons changes with the features of the array as well as with the input state, a well designed array of coupled waveguides could be used to produce new states of the electromagnetic field.