{'fraction_non_alphanumeric': 0.04136582117813521, 'fraction_numerical': 0.018768596932936598, 'mean_word_length': 4.3704324505531345, 'pattern_counts': {'":': 0, '<': 28, '<?xml version=': 0, '>': 2, '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}}
Frequency map analysis of spatiotemporal chaos in the nonlinear disordered Klein-Gordon lattice

We study the characteristics of chaos evolution of initially localized energy excitations in the one-dimensional nonlinear disordered Klein-Gordon lattice of anharmonic oscillators, by computing the time variation of the fundamental frequencies of the motion of each oscillator. We focus our attention on the dynamics of the so-called `weak' and `strong chaos' spreading regimes [2010, EPL 91 30001], for which Anderson localization is destroyed. Based on the fact that large variations of the fundamental frequencies denote strong chaotic behavior, we show that in both regimes chaos is more intense at the central regions of the wave packet, where also the energy content is higher, while the oscillators at the wave packet's edges exhibit regular motion up until the time they gain enough energy to become part of the highly excited portion of the wave packet. Eventually, the percentage of chaotic oscillators remains practically constant, despite the fact that the number of excited sites grows as the wave packet spreads, but the portion of highly chaotic sites decreases in time. We show that the extent of the zones of regular motion at the edges of the wave packet in the strong chaos regime is much smaller than in the weak chaos case. Furthermore, we find that in the strong chaos regime the chaotic component of the wave packet is not only more extended than in the weak chaos one, but in addition the fraction of strongly chaotic oscillators is much higher. Another important difference between the weak and strong chaos regimes is that in the latter case a significantly larger number of frequencies is excited, even from the first stages of the evolution. Moreover, our computations confirmed the shifting of fundamental frequencies outside the normal mode frequency band of the linear system in the case of the so-called `selftrapping' regime where a large part of the wave packet remains localized.

We study the characteristics of chaos evolution of initially localized energy excitations in the onedimensional nonlinear disordered Klein-Gordon lattice of anharmonic oscillators, by computing the time variation of the fundamental frequencies of the motion of each oscillator. We focus our attention on the dynamics of the so-called 'weak' and 'strong chaos' spreading regimes [Laptyeva et al., 2010, Europhys. Lett., 91 30001], for which Anderson localization is destroyed, as the initially restricted excitation at the central region of the lattice propagates in time to more lattice sites. Based on the fact that large variations of the fundamental frequencies denote strong chaotic behavior, we show that in both regimes chaos is more intense at the central regions of the wave packet, where also the energy content is higher, while the oscillators at the wave packet's edges, through which the energy propagation happens, exhibit regular motion up until the time they gain enough energy to become part of the highly excited portion of the wave packet. Eventually, the percentage of chaotic oscillators remains practically constant, despite the fact that the number of excited sites grows as the wave packet spreads, but the portion of highly chaotic sites decreases in time. Thus, albeit the number of chaotic oscillators is constantly growing the strength of their chaotic behavior decreases, indicating that although chaos persists it is becoming weaker in time. We show that the extent of the zones of regular motion at the edges of the wave packet in the strong chaos regime is much smaller than in the weak chaos case. Furthermore, we find that in the strong chaos regime the chaotic component of the wave packet is not only more extended than in the weak chaos one, but in addition the fraction of strongly chaotic oscillators is much higher. Another important difference between the weak and strong chaos regimes is that in the latter case a significantly larger number of frequencies is excited, even from the first stages of the evolution. Moreover, our computations confirmed the shifting of fundamental frequencies outside the normal mode frequency band of the linear system in the case of the so-called 'selftrapping' regime where a large part of the wave packet remains localized.
In particular, for the 1D DKG system, which has been studied in more detail than its 2D counterpart whose long-time numerical integration is a very demanding computational task, the existence of two energy spreading dynamical regimes was theoretically predicted and numerically verified: the so-called 'weak' and 'strong chaos' regimes, which are characterized by different dynamical behaviors. More specifically, the second moment m 2 of the normalized energy distribution increases in time t as m 2 (t) ∝ t am with a m = 1/3 (a m = 1/2) for the weak (strong) chaos case, while the corresponding participation number P grows as P (t) ∝ t ap with a p = 1/6 (a p = 1/4). In addition, the appearance of a 'selftrapping' regime, where the bulk of the wave packet remains localized, was also theorized and observed [Kopidakis et al., 2008;Flach et al., 2009;Skokos et al., 2009;Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a]. The generality of these results was substantiated by the fact that the same dynamical behaviors were observed also for the disordered discrete nonlinear Schrödinger (DDNLS) equation [Molina, 1998;Pikovsky & Shepelyansky, 2008;García-Mata & Shepelyansky, 2009;Flach et al., 2009;Skokos et al., 2009;Laptyeva et al., 2010;Bodyfelt et al., 2011a,b;Kati et al., 2020], whose numerical integration is much more demanding from a computational point of view. Further substantiation was produced from high speed iterations of related nonlinear unitary maps using discrete time quantum walks which not only confirmed the generality of the spreading results, but which also allowed to push time horizons of subdiffusion to record 10 12 [Vakulchyk et al., 2019].
All these studies indicated the chaotic nature of energy spreading in nonlinear disordered media. The characteristics of this chaotic behavior were investigated in [Skokos et al., 2013; and  for, respectively, 1D and 2D disordered lattices, through the time evolution of the most commonly used chaos indicator, the maximum Lyapunov exponent (MLE) Λ 1 [Benettin et al., 1980a,b;Skokos, 2010;Pikovsky & Polity, 2016] and the study of the properties of the so-called 'deviation vector distribution' (DVD) of the tangent space vector used for the computation of the MLE. There it was shown that the overall strength of chaos is decreasing as an initially localized energy excitation spreads to more lattice sites, but without any sign of a potentially crossover to regular behavior as was speculated in [Johansson et al., 2010;Aubry, 2011], at least up to the maximum evolution times reached in those studies. This behavior is reflected in the time decay of the finite time MLE (ftMLE) Λ, i.e. Λ(t) ∝ t a Λ , with a Λ ≈ −0.25 (a Λ ≈ −0.3) for the weak (strong) chaos regime of 1D systems [Skokos et al., 2013;, and with a Λ ≈ −0.37 (a Λ ≈ −0.46) for the weak (strong) chaos case of 2D systems . It is also worth noting the existence of a dimension-independent scaling between the wave packet's spreading and chaoticity, quantified by the same behavior of the ratio Λ(t)/m 2 (t) in 1D and 2D systems .
The scope of the current paper is to shed some additional light on the spatiotemporal characteristics of chaos in 1D disordered, nonlinear lattices. The MLE is an average measure of chaoticity providing information about the global behavior of a dynamical system. As such, the value of the MLE by itself is not enough to understand the spatiotemporal evolution of active chaotic regions in a multidimensional system [Tieleman et al., 2014]. Thus, in order to follow the chaotic or regular behavior of each individual lattice site in time we implement the so-called 'frequency map analysis' (FMA) technique [Laskar, 1990;Laskar et al., 1992;Laskar, 1993] and evaluate the fundamental frequency of an appropriate time series produced by the coordinates of each oscillator. Then, we use the variations of these frequencies to identify chaotic behavior. This approach allows us to visualize the evolution of chaos in the propagated energy distribution and to identify differences between the weak and strong chaos spreading regimes.
The paper is organized as follows. In Sect. 2 we present the Hamiltonian function of the 1D DKG model along with the basic quantities we use to characterize the energy spreading and the system's chaoticity. The numerical results of our study are reported in Sect. 3, initially for some weak chaos cases and afterwards for the strong chaos regime, followed by a brief discussion of the selftrapping behavior. Finally in Sect. 4 we summarize our findings and discuss their significance.

Lattice model and numerical techniques
The Hamiltonian function of the 1D DKG lattice system of N anharmonic oscillators is with q i , p i being respectively the generalized position and momentum of site i. The disorder parameters l take uncorrelated random values in the interval 1 2 , 3 2 following a uniform probability distribution, and W > 0 determines the disorder strength. In our study we consider fixed boundary conditions, i.e. q 0 = p 0 = q N +1 = p N +1 = 0. The Hamiltonian (1) is an integral of motion and its value H (usually referred to as the system's energy E) remains constant in time.
We use the ABA864 symplectic integrator (SI) [Blanes et al., 2013] for numerically solving the equations of motionq of Hamiltonian (1), which has been proved to be a very efficient technique for this task Danieli et al., 2019]. In particular, for the implementation of the ABA864 SI we split Hamiltonian (1) in two integrable parts, namely the system's kinetic and potential energy, and use an integration time step τ = 0.5, which, in general, keeps the absolute value of the relative energy error below 10 −5 . In our simulations we typically excite one or few sites at the central region of the lattice and follow the evolution of this energy excitation in time, taking care that the lattice is large enough so that the energy does not reach its boundaries until the end of the integration. In order to analyze the characteristics of the energy spreading we define the normalized energy distribution and compute the distribution's second moment which measures the extent of the distribution, and the participation number which estimates the number of highly excited sites. We note that in (5)ī = N i=1 iξ i is the distribution's mean position.
Here we are interested in the spatiotemporal evolution of chaos in the DKG system (1). Thus, since the MLE provides information about the collective and overall chaotic behavior of the system, we implement the FMA method for identifying the time evolution of the chaotic behavior of individual lattice sites. Neglecting the non-linear term q 4 /4 in the Hamiltonian (1), the system is integrable and the motion of each site is fully determined by a certain linear combination of the eigenfrequencies of this system. Under the influence of the non-linearity these frequencies will become time-dependent and can be determined by using a refined Fourier technique as for example the NAFF (Numerical Analysis of Fundamental Frequencies) algorithm. For a detailed description of NAFF and its application for FMA in various research fields see for example [Laskar, 1990[Laskar, , 1993Robutel & Laskar, 2000;Papaphilippou, 2014] and references therein.
This algorithm uses as input the numerically integrated trajectory q i (t), p i (t) of oscillator i over a time span T and identifies the fundamental frequenciesin the quasiperiodic approximation of this motion. In this work we compute the relative change of the largest of these fundamental frequencies f 1i and f 2i between two successive time windows of length T = 6 · 10 5 time units. 1 In the case of regular, quasiperiodic motion the two frequencies should theoretically be equal (in practice they slightly differ due to their numerical estimation), while for chaotic motion the two computed frequencies are in general different, as there is no reason for them to remain constant. Thus, the relative change of these two frequencies, quantified by can be used to identify the chaotic or regular nature of motion. More specifically, small D i values denote the practical constancy of the fundamental frequency and consequently indicate regular motion, while large D i values signify chaotic behavior characterized by strong variations in the computed frequency values.

Numerical results
Although in our analysis we considered several individual cases belonging to the weak and strong chaos regimes, as well as the selftrapping regime (more information on the definition of these dynamical behaviors can be found in [Laptyeva et al., 2010;Flach, 2010]), we will only report our findings for some representative cases, which have already being studied in the literature. In all cases we consider a disorder realization, i.e. a set of˜ l values in (1), and excite L central sites of a lattice of N oscillators by setting their positions at zero and adjusting their momenta so that each one of them has the same energy density e = E/L. For single site excitations (L = 1) always a positive momentum is used, while for block excitations (L > 1) the sign of the initial momentum of each excited site is chosen randomly. In particular, we consider the following cases: (1) Weak chaos regime

Weak chaos dynamical regime
In Fig. 1 we see that the evolution of the second moment m 2 (t) [ Fig. 1(a)], the participation number P (t) [ Fig. 1 (b)] and the ftMLE Λ(t) [ Fig. 1(c)] exhibits the expected behavior for both the WCI (red curves) and the WCII case (blue curves). More specifically, in accordance to Skokos et al., 2009;Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a] the asymptotic behavior of the wave packet's dynamics is characterized by m 2 (t) ∝ t 1/3 and P (t) ∝ t 1/6 , while the ftMLE shows a decrease Λ(t) ∝ t −0.25 as was reported in [Skokos et al., 2013;. Let us first study in detail the WCI weak chaos case where initially 1 site is excited at the middle of a lattice with N = 999 oscillators. As time increases the system's constant total energy E (1) is spread to more sites as can be seen from Figs. 2(a) and (c) where we depict the energy evolution in two time windows of length 10 8 time units: immediately after the initial excitation [ Fig. 2(a)] and after 9 · 10 8 time units [Fig. 2(c)]. In order to visualize the change in time of the chaoticity of each lattice site we also perform the FMA of the system and compute the evolution of the quantity D i (8), which captures the variation of the sites' fundamental frequencies in two successive time windows of length 6 · 10 5 attributing the computed D i value of each site at the beginning of these two time intervals. The outcome of this process is presented in Figs. 2(b) and (d), respectively for the time intervals of Figs. 2(a) and (c). There each site is colored according to its log 10 D i value. Chaotic behavior corresponding to large frequency changes and high log 10 D i values is indicated by red and yellow colors, while regions of low log 10 D i values denoting small frequency changes attributed to regular, non-chaotic behavior are colored in purple and black. From the results of Figs. 2(b) and (d) we clearly see that chaos appears at the central regions of the wave packet, where also the energy concentration is higher as Figs. 2(a) and (c) indicate. On the other hand, the edges of the wave packet, through which energy spreading takes place, behave more regularly having also significantly smaller energies than the central regions.
It is worth mentioning here a practical issue. In order to compute the change in the fundamental frequencies using the NAFF algorithm for two successive time windows of length 6 · 10 5 time units, we have to ensure that the time series we analyze is reliable and accurate. For this reason we perform the computation of D i only if the absolute values for each the two successive frequencies f 1i and f 2i are above the adopted threshold of 10 −16 . If this condition is not satisfied then we do not regard the computed frequency as accurate and we dismiss it. For the computation of D i (8) we require the obtainment of the Although the color maps of D i values in Figs. 2(b) and (d) accurately capture subtle changes in the chaotic behavior of lattice sites, we also follow below a different approach in order to clearly visualize the spatiotemporal evolution of chaos in the system. In order to do so, let us first note that the frequencies of the normal modes of the linear DKG model, i.e. the system obtained by neglecting the nonlinear terms q 4 i /4 in Hamiltonian (1), belong to a well-defined interval. More specifically, the normal mode eigenvalues ω 2 ν ∈ 1 2 , 3 4 + 4 W , ν = 1, 2, . . . , N (see e.g.  and references therein), so that The presence of nonlinearity induces a frequency shift to higher frequency values . This shift increases with the growth of the system's energy E, which actually plays the role of nonlinearity's strength. Thus, the fundamental frequencies of all lattice oscillators are expected to belong to a frequency band similar to (9), having somehow raised lower and upper limits. This is actually true as our numerical results reported below clearly verify. In our study we consider set-ups for only two values of W , namely W = 3 (case SCII) and W = 4 (all other cases). For each W value we consider a range of possible frequencies, which is slightly larger than the one reported in (9) in order to be certain that it will contain the actually computed frequencies of the nonlinear system (1). In particular, for W = 3 we consider the frequency band [0.1, 0.28] and for W = 4 the interval [0.1, 0.26]. These ranges are then divided into n = 500 bins of equal length, setting the width of each frequency bin to be of the order of 10 −4 . During the wave packet's evolution the fundamental frequencies are computed and for each site the bin which contains the found frequency is registered. This is called one visit to that bin. In the course of the integration, the number of visits to the bin that received the most visits up to this point for each site is divided by the number of total registered bin visits (i.e. total number of frequency computations) for this site. This gives the relative number of visits of the most visited frequency bin for each site and is used to measure what we call 'frequency locking' of site i (F L i ). This number is always 0 ≤ F L i ≤ 1. The value F L i = 1 for a particular site i indicates that during the evolution this site had always the same frequency (actually this means that although the frequency might change a bit it remained in the same bin, i.e. it does not change much as the bins' width is ∼ 10 −4 ), signifying effectively regular behavior. On the other hand, small F L i values denote many changes of the fundamental frequency as the placement of the computed values in one frequency bin is low, and consequently indicate chaotic behavior, which becomes stronger as F L i decreases.
The computation of F L i is done based on data from a time window of 1.2 · 10 6 time units and repeated each 1.2 · 10 5 time units. To avoid saturation effects the computation is re-initialized after a fixed time interval, which is set to be t F L = 2.4 · 10 7 in order to better follow the time evolution of this quantity in different evolution stages of the wave packet. By doing that we would for example be able to identify potential regularly behaving epochs in the dynamics of a lattice site, i.e. time intervals for which its frequency remains practically constant, although from time to time the value of the practically constant frequency might change.
The results of the computation of F L i for the WCI case are shown in Figs. 3(a) and (b) respectively for the time intervals 0 ≤ t ≤ 10 8 and 9 · 10 8 ≤ t ≤ 10 9 . Points for which the fundamental frequencies of the lattice sites were computed are colored according to their F L i values, with red and yellow colors denoting small F L i values indicating strong chaotic behavior. On the other hand, regions with very high F L i values close to 1 colored in black correspond to regular motion. The re-initialization of the F L i 's computation every t F L = 2.4 · 10 7 time units results in the creation of the equidistant vertical 'lines' in both panels of Fig. 3.
In accordance to the results presented in Figs. 2(b) and (d), Fig. 3 shows that the edges of the wave packet behave regularly, while chaos appears at the central region of the lattice's excited part. The fact that in the F L i plots of Fig. 3 we do not restrict ourselves in computing the exact frequencies, but we allow some frequency changes of the order of the width of the bin, i.e. ∼ 10 −4 , results in the creation of pictures depicting the system's dynamics more clearly than the D i color plots [Figs. 2(b) and (d)], emphasizing especially the more active chaotic regions of the lattice, which exhibit large changes in their frequencies.
The D i plots are more detailed and of course more accurate, but the many details they provide make the picture somewhat more obscure, hiding the main dynamical features of the system. For example, in Fig. 3(b) the footprint of a horizontal 'stripe' of regular motion around i ≈ 480, lasting for almost the whole  (d)] 9 · 10 8 ≤ t ≤ 10 9 (WCI case). The considered intervals are: 0 < F L i < 0.2 (yellow curves), 0.2 ≤ F L i < 0.4 (orange curves), 0.4 ≤ F L i < 0.6 (red curves), 0.6 ≤ F L i < 0.8 (blue curves), 0.8 ≤ F L i < 1 (purple curves). Insets in panels (a) and (c): The time evolution of the percentage of sites with F L i = 1, i.e. points plotted in black in Fig. 3, which correspond to regions of regular motion. All presented results are obtained by the analysis of the outcomes presented in Fig. 3. The re-initialization of the F L i 's computation every t F L = 2.4 · 10 7 time units leads to the vertical 'stripe'-looking feature appearing in all panels.
duration of the presented time interval, is easily identified as points there are mainly colored in black. This region is also present in Fig. 2(d) where points are mostly colored in purple/dark red, indicating values of log 10 D i ≈ 10 −5 , but is not as easily visible as in Fig. 3(b).
The comparison of Figs. 3(a) and (b) shows that as time grows the portion of the wave packet exhibiting chaotic behavior (i.e. points not colored in black) decreases. This is mainly due to the fact that the regularly behaving sites at the black colored edges of the wave packet grow. In order to quantify this observation we plot in Figs. 4(a) and (c) the evolution of the percentage p F L of sites, respectively shown in Figs. 3(a) and (b), having F L i values in various intervals. In particular, we consider the intervals 0 < F L i < 0.2 (yellow curves), 0.2 ≤ F L i < 0.4 (orange curves), 0.4 ≤ F L i < 0.6 (red curves), 0.6 ≤ F L i < 0.8 (blue curves), 0.8 ≤ F L i < 1 (purple curves), while in the insets of Figs. 4(a) and (c) we plot the evolution of p F L for sites having F L i = 1, i.e. sites whose fundamental frequency does not change much as it always remains in one bin indicating regular motion. The presented results are obtained by the analysis of the data shown in Fig. 3, so we see again here a vertical 'stripe' formation of length t F L = 2.4 · 10 7 related to the re-initialization of the F L i computations. In Figs. 4(b) and (d) we depict the accumulated percentages From the results of Figs. 4(a) and (b) we see that at the initial stages of the evolution the percentage of chaotic sites (i.e. the ones having F L i = 1) decreases, and obviously the number of sites showing consistent regular behavior (F L i = 1) increases [see inset of Fig. 4(a)]. It is worth noting that the percentage of strongly chaotic regions, i.e. the ones having F L i < 0.4, remains relative low, e.g. p F L ≈ 3.4% for sites with 0.2 ≤ F L i < 0.4 [orange curve in Fig. 4(a)], while the p F L values for F L i < 0.2 are too small (≈ 0.02%) so that the corresponding yellow curve is practically not visible in Fig. 4(a). Furthermore, the fraction of sites of moderate (0.4 ≤ F L i < 0.6 and 0.6 ≤ F L i < 0.8 -red and blue curves respectively) and of weak chaotic behavior (0.8 ≤ F L i < 1 -purple curve) show signs of saturation to values p F L ≈ 10% for each one of these three F L i ranges. This saturation to some asymptotic p F L and P F L values becomes evident in Figs. 4(c) and (d). More specifically, we see that the red, blue and purple curves in Fig. 4(c) oscillate around values p F L ≈ 10% with red (purple) curves being slightly below (above) this value. Furthermore, the percentages of well-established regular behavior remains practically constant around p F L ≈ 65% [inset of Fig. 4(c)], while the strong chaotic behavior (F L i < 0.4), corresponding to orange and yellow curves [which are practically not visible in Fig. 4 Fig. 3 and data presented in the main panels of Fig. 4], are mainly concentrated at the central regions of the wave packet. These results suggest that chaoticity is not extended, as the large portion of sites at the edges of the wave packet, whose frequency does not practically change, exhibit regular behavior. Thus, in general, points which are slightly excited at the wave packet's boundaries are initially not chaotic (having relative small D i and large F L i values) until they become part of the highly excited component of the wave packet, when more energy reaches them.
Furthermore, during the above-mentioned F L i computations additional information is gathered. More specifically, for each frequency bin the number of sites having their fundamental frequency in the bin is registered. The intention here is to identify frequencies (or more accurately frequency bins, i.e. small frequency ranges) which dominate the dynamics of the system by having many sites oscillating with these particular frequencies. In this way we quantify the importance or the influence of the particular frequency range on the dynamical evolution of the lattice. The outcome of these computations are shown in Fig. 5 for the two time intervals 0 ≤ t ≤ 10 8 [ Fig. 5(a)] and 9 · 10 8 ≤ t ≤ 10 9 [ Fig. 5(b)] we have considered so far. If at a given time the computed fundamental frequencies of some sites are located inside a particular frequency bin then a point marks this occurrence in Fig. 5. This point is colored according to the number of sites belonging in the frequency bin at each time instant (we will refer to this measure as the 'frequency strength', F S, of the particular frequency bin) following the color code depicted at the right side of Fig. 5. Thus, red and yellow colors indicate many sites in a particular frequency range, while the lack of computed frequencies correspond to white regions in Fig. 5.
From the results of Fig. 5 we see that as time grows more frequencies are excited since the extent of the white regions in Fig. 5(b) is smaller that the one observed in Fig. 5(a). Furthermore, in both time windows only a few frequency bins are highly populated (the ones colored in red and yellow). The corresponding frequencies are mainly located near the middle of the frequency band of Fig. 5 (f i ≈ 0.19), although some frequencies with higher F S values are also present in the interval f i = 0.2-0.22 [especially in Fig. 5(b)].
The results presented in Figs. 2-5 for the weak chaos case WCI are quite general and characteristic of the weak chaos spreading regime, as similar behaviors have been observed for all other studied weak chaos cases. As a testimony to that we present in Fig. 6 results obtained for the WCII case in which we initially excite L = 21 neighboring lattice sites, in contrast to the single site excitation (L = 1) of case WCI. For the WCII case we see again that as the wave packet spreads the region of intense chaoticity [i.e. points colored in red and yellow in the D i plots in Figs. 6(a) (h) clearly show that, as in the WCI case (Fig. 5), initially a few frequencies are excited, but later on, as the energy spreads to more sites, the number of excited frequencies increases. Again a rather small number of frequencies, mainly located near the middle of the frequency band [points colored in red and yellow in Figs. 6(g) and (h)] strongly influence the dynamics of the lattice as they are highly populated for most of the integration time.

Strong chaos dynamical regime
Let us now focus on the strong chaos spreading dynamical regime, which is characterized by a potentially long-lasting faster expansion of the wave packet, with respect to what is observed in the weak chaos case, eventually followed by a crossover to a slower subdiffusive weak chaos spreading [Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a;. In Fig. 7 we present results for the two representative strong chaos cases we consider in our study, namely cases SCI (red curves) and SCII (blue curves). All quantities are eventually behaving according to what is expected in the strong chaos regime. More specifically the second moment (5) and the participation number P (t) (6) follow respectively power law growths m 2 (t) ∝ t 1/2 and P (t) ∝ t 1/4 [Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a], while the decrease of the ftMLE (7) is well described by Λ(t) ∝ t −0.3 .
In order to investigate in depth the characteristics of chaos in the strong chaos regime, along with their differences to the weak chaos behaviors presented in Sect. 3.1, we analyze in detail the SCI case. In our comparisons between the strong and weak chaos behaviors we will refer to the WCI case, which, as we already mentioned in Sect. 3.1, is a good representative of the weak chaos spreading regime, since all other tested cases (like for example the WCII case of Fig. 6), exhibited similar behaviors. In the strong chaos regime spreading is faster than in the weak chaos case (something which is quantified by the larger exponents in the power law time evolutions of m 2 and P ), so for its numerical investigation we need larger lattices in order to avoid the appearance of boundary effects. For this reason the lattice size in the SCI and SCII cases is at least doubled with respect to the weak chaos cases of Sect. 3.1. This faster spreading is also evident by the comparison of the energy distributions ξ i (4) in Figs. 8(a) and (c) with the ones in Figs. 2(a) and (c), as in the former case significantly more sites have been excited as same times.
What is more interesting though, is the picture emerging by contrasting the chaotic behavior of the SCI and WCI cases through the comparison of their D i (8) distributions respectively presented in Figs. 8(b), (d) and Figs. 2(b), (d). In all these figures the edges of the wave packet are colored in purple and black, indicating low D i values and consequently non-chaotic behavior. On the other hand, the regions of active chaos colored in red and yellow are significantly more extended in the SCI case, covering a larger portion of the excited part of the lattice than in the WCI case. This clear difference between the strong and weak chaos dynamical behaviors is also vividly seen by the comparison of the frequency locking F L i values presented respectively in Figs. 9 and 3. From the results of Fig. 9(a) it is obvious that from the first stages of the evolution chaos (i.e. regions with F L i = 1) is quite extended, covering almost the whole wave packet  apart from a rather narrow area at its edges where regular behavior (i.e. points corresponding to F L i = 1, which are colored in black) appears. This arrangement persists also at later stages of the evolution as can be seen in Fig. 9(b).
It is worth noting though that the area of strongly chaotic regions in Fig. 9(b), corresponding to very low F L i values (i.e. points colored in orange and yellow) is decreased with respect to Fig. 9(a). In order   Fig. 4. The analysis of these results, as well as their comparison with the outputs presented in Fig. 4 for the WCI case, reveal some interesting characteristics about the chaotic behavior of the strong chaos regime, along with some important differences with respect to the weak chaos case. More specifically, the percentage of sites showing consistent regular behavior, i.e. F L i = 1 after some initial fluctuations [see inset of Fig. 10(a)] more or less stabilize its value at p F L ≈ 20%, something which remains true also at the later stages of the evolution, as is seen in the inset of Fig. 10(c). This behavior is also revealed by the practical constancy of the height of the colored area in Figs. 10(b) and (d), which indicate the P F L values of sites with F L i = 1. A similar behavior was also observed for the weak chaos regime, i.e. the percentage p F L of sites with F L i = 1 saturates to an almost constant value [see insets of Figs. 4(a) and (c)], but that value is significantly larger, p F L ≈ 60%, than the one obtained for the SCI case. This difference in p F L values actually quantifies the observation that the width of the regularly behaving zone at the edges of the wave packet (i.e. points colored in purple and black in Figs. 8(b) and (d) and in black in Fig. 9) is smaller than in the weak chaos case [see Figs. 2(b) and (d), and Fig. 3].
The differences between the weak and strong chaos cases go beyond the decrease of the extent of the regular component in the latter case (or equivalently the increase of the chaotic portion), as also the chaotic parts of the wave packet exhibit disparate characteristics. In particular, from Figs. 10(a) and (b) we see that at the initial stages of the evolution the percentage of strongly chaotic regions, i.e. sites with F L i < 0.4, shows a tendency to slightly decrease. This tendency becomes more evident in Fig. 10(b) where the height reached by the orange colored area moderately lessens attaining values P F L ≈ 30% towards the end of the depicted time interval. These levels are significantly larger than the P F L ≈ 3% values observed in the weak chaos case WCI for the same time interval in Fig. 4(b). At later stages of the evolution [Figs. 10(c) and (d)] the fractions of sites with F L i values at different intervals tend to oscillate around some well defined levels. In contrast to the weak chaos case of Fig. 4(d) where P F L 3% for sites with F L i < 0.4, in the strong chaos case we see that the percentages of these highly chaotic sites oscillates around much higher levels as P F L ≈ 20%.
The higher degree of chaos in the strong chaos case, in comparison to the weak chaos regime, is also evident by the elevated p F L values of sites exhibiting moderate chaotic behavior. More specifically, we see that p F L ≈ 30% and p F L ≈ 20% respectively for the 0.4 ≤ F L i < 0.6 [red curve in Fig. 10(c)] and 0.6 ≤ F L i < 0.8 [blue curve in Fig. 10(c)] cases of the strong chaos regime, while we had p F L ≈ 8% [red curve in Fig. 4(c)] and p F L ≈ 10% [blue curve in Fig. 4(c)] for these two cases in the weak chaos regime. Furthermore, although the percentages of low strength chaos, i.e. sites having 0.8 ≤ F L i < 1, are similar in both the strong and weak chaos cases, with p F L ≈ 17% [purple curves in respectively Figs. 10(c) and 4(c)], it is worth noting that this fraction is lower (higher) than the ones corresponding to the 0.4 ≤ F L i < 0.6 and 0.6 ≤ F L i < 0.8 cases in the strong (weak) chaos regime. This behavior indicates again the higher degree of chaos appearing in the strong chaos case.
In Fig. 11 we see the time evolution of the computed fundamental frequencies f i for the SCI case. More accurately for the creation of this figure we register the frequency bins in which site frequencies reside in the time intervals 0 ≤ t ≤ 10 8 [ Fig. 11(a)] and 9 · 10 8 ≤ t ≤ 10 9 [ Fig. 11(b)], with each occupied bin colored according to its F S value (i.e. the number of sites having frequencies in that bin), in analogy to Fig. 5. Similarly to what was observed in the weak chaos case in Fig. 5, the number of excited frequencies increases in time, although this increase is not as profound as in the WCI case. The striking difference between Figs. 5 and 11 is that in the strong chaos case significantly more frequencies are excited than in the weak chaos regime, for both depicted time intervals. This is in perfect agreement with the theoretical prediction that the strong chaos spreading behavior is related to a widespread excitation of normal mode frequencies Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a]. Another important difference is that in the strong chaos case many more frequencies influence the dynamics than in the weak chaos case. This is vividly seen in Fig. 11(b) where a lot of frequency bins, throughout the whole frequency spectrum, are highly populated (i.e. having large F S values and consequently colored in red and yellow). Nevertheless, again a concentration of bins with high F S values at the central region of the frequency band is observed.
The results presented in Figs. 8-11 for the SCI case are typical for the strong chaos spreading dynamical regime, as similar behaviors have been obtained for all other tested strong chaos cases. The explicit results for one such case, the one we named SCII in Sect. 3, are presented in Fig. 12. It is worth noting that the width of the frequency band (9) of the system's normal modes, denoted by the horizontal dashed green lines in Figs. 12(g) and (h), is larger than the ones related to all other studied cases having W = 4, as W = 3 for the SCII case. Nevertheless, the qualitative features observed in Figs. 12(g) and (h) are similar to the ones seen in Fig. 11.

Selftrapping dynamical regime
Although the main objective of our work is to understand the spatiotemporal behavior of chaos for energy propagation in disordered lattices, as well as to identify the differences between the weak and strong chaos spreading regimes, we also consider in our study, mainly for completeness' sake, a case belonging to the selftrapping regime. This is the STI case mentioned in Sect. 3, whose dynamics has been studied in Skokos et al., 2009].
The existence of the selftrapping dynamical regime was theoretically predicted for the 1D DDNLS system in [Kopidakis et al., 2008] and was further discussed in . According to those studies, in the case of single site excitations, and for strong enough nonlinearities, the induced nonlinear frequency shift moves the frequencies of some of the excited oscillators outside the model's linear spectrum, resulting in the localization of a part of the wave packet, with the remainder spreading subdiffusively. Although the theoretical prediction of this behavior was done only in the case of the 1D DDNLS model, the selftrapping regime was numerically observed both for the 1D DDNLS and DKG systems in several publications, e.g. Skokos et al., 2009;Laptyeva et al., 2010;Bodyfelt et al., 2011a], mainly through the monitoring of wave packet profiles and the computation of m 2 and P , with the latter remaining practically constant. Thus, these numerical studies denoted the generality  Fig. 6 but for the strong chaos case SCII. It is worth noting the different frequency band of the linear DKG system (9), which is indicated by the horizontal thick dashed green lines in panels (g) and (h), as in this case W = 3 in contrast to W = 4 in all other studied cases.
of the selftrapping behavior.
The time evolution of the STI case's fundamental frequencies, along with their F S, given in Fig. 13, clearly show, to the best of our knowledge for the first time, the existence of lattice sites having frequencies outside the frequency band of the linear model (identified by the horizontal, dashed green lines in Fig. 13), something which is a basic ingredient for the theoretical description of the selftrapping regime. Thus, our numerical results provide strong additional evidences of the correctness of the developed theoretical framework.
Furthermore, it is worth noting some similarities between Fig. 13 and the behaviors observed for the weak chaos case [Fig. 5 and Figs. 6(g) and (h)], as in both cases fewer frequencies are excited with respect to the strong chaos case [Fig. 11 and Figs. 12(g) and (h)]. In addition, in both the weak chaos and selftrapping cases a rather small number of frequencies have large F S values, corresponding to red and yellow colored points in Figs. 5, 6(g) and (h) (weak chaos) and in Fig. 13 (selftrapping). Moreover, the number of excited frequencies increase as time grows, as the size of the white region in Fig. 13(b) for 9 · 10 8 ≤ t ≤ 10 9 is smaller than the one of the initial phase of the evolution in Fig. 13(a).

Summary and discussion
Based on results obtained by the FMA of the evolution of initially localized energy excitations we studied the features of spatiotemporal chaos in a prototypical model of 1D disordered nonlinear lattices, namely the DKG system. We focused our investigation on well-established dynamical behaviors leading to the destruction of Anderson localization, which have been investigated in the past by different approaches Skokos et al., 2009;Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a,b;Skokos et al., 2013;Antonopoulos et al., 2014Antonopoulos et al., , 2017: the so-called weak and strong chaos spreading regimes. In particular, we studied in detail the characteristics and the evolution of chaos in these two cases, as well as explored the differences encountered between these regimes.
Using the relative change D i (8) of the computed fundamental frequencies in two successive time windows as an indicator of chaos for each lattice site [Figs. 2(b), (d); 6(a), (b); 8(b), (d) and 12(a), (b)] we showed that, in both the weak and strong chaos cases, chaotic behavior mainly appears at the central regions of the wave packet, where also the energy density ξ i (4) [Figs. 2(a), (c) and 8(a), (c)] is relatively large, while sites at the edges of the wave packet forms a zone of regular motion. The width of this zone is larger in the weak chaos case, which consequently means that the chaotic component of the wave packet is more extended in the strong chaos regime.
In order to quantify further these findings, and to also investigate in more detail the features of the chaotic part of the dynamics, we followed the progression of a quantity we named frequency locking (F L i ), which is measuring the degree of practical frequency constancy (denoting non-chaotic behavior), of each oscillator in time windows of fixed length [Figs. 3;4; 6(c), (d), (e), (f); 9; 10; 12(c), (d), (e), (f)]. The analysis of our F L i computations showed that the percentage of non-chaotic oscillators in the wave packet is about 3 times larger in the weak chaos regime with respect to the strong chaos case. Obviously this means that the fraction of oscillators behaving chaotically is much larger in the strong chaos regime. Apart from this important difference, our study showed that the percentage of strongly chaotic sites, i.e. oscillators continuously exhibiting many and large changes in their frequencies and consequently having low F L i values (typically F L i < 0.4), is much higher in the strong chaos case, usually around 5 times higher than in the weak chaos one. These observations, legitimize and further support, in some sense, the naming of the weak and strong chaos regimes [Flach, 2010;Laptyeva et al., 2010], as in the latter case chaotic behavior is indeed more pronounced and extended.
An additional interesting feature we noticed is that, in both spreading regimes, although the percentage of chaotic oscillators eventually remains almost fixed (because the percentage of regularly behaving sites with F L i = 1 practically saturates to a more or less constant value [insets of Figs. 4(a), (c); 6(e), (f); 10(a), (c) and 12(e), (f)]) the fraction of highly chaotic oscillators (sites having e.g. F L i < 0.4) decreases in time. Thus, chaos persists although it becomes weaker in time. This observation is in agreement with the picture emerging from the computation of the ftMLE Λ(t) [Figs. 1(c) and 7(c)], that as time passes, the strength of chaos [quantified by the value of Λ(t)] decreases, but nevertheless without showing any tendency to crossover to regular dynamics (which corresponds to Λ(t) ∝ t −1 ), as was also shown in various studies through the computation of other quantities  [Skokos et al., 2013;.
All these results allow us to better understand the mechanisms of energy spreading and the process of chaotic destruction of Anderson localization for the case of initially localized excitations in disordered lattices. The constant energy spreads in time to more sites as oscillators at the edges of the wave packet start acquiring energy and begin their motion, initially in a regular, non-chaotic fashion. As time passes these oscillators gain more energy and eventually become chaotic, and get incorporated in the significantly excited part of the wave packet, not belonging any more to the regularly-behaving edges. As additional oscillators get excited, the energy is distributed among more sites reducing in this way the amount of energy per particle in the excited part of the lattice. Thus, although the number of chaotic oscillators is continuously increasing the strength of their chaotic behavior is lessened.
A significant outcome of our study is the finding of the excited frequencies in the dynamics and, more importantly, the identification of the ones which influence the most the system's behavior [Figs. 5; 6(g), (h); 11 and 12(g), (h)], through the computation of the number of oscillators moving according to these frequencies (a quantity we refer to as the related frequency strength, F S). A remarkable difference between the weak and strong chaos regimes is that in the latter case a significantly larger number of frequencies is excited, even from the first stages of the evolution. This contrast remains valid for the whole duration of our numerical simulations, with more frequencies been activated as time passes, although a higher increase rate is observed in the weak chaos regime. The excitation of more frequencies in the strong chaos case is in accordance to the extended chaotic interactions of the system's normal modes, which were theorized to take place in this case Laptyeva et al., 2010;Flach, 2010;Bodyfelt et al., 2011a].
Furthermore, it is worth mentioning that in our study we also provided some significant numerical results for the selftrapping dynamical regime, where a large part of the wave packet remains localized. Our findings (Fig. 13) confirm the shifting of frequencies outside the normal mode frequency band, a phenomenon which is at the core of the theoretical treatment of this dynamical behavior.
Apart from understanding in more depth the chaotic spreading processes in disordered nonlinear lattices, the numerical approaches implementing in our study (FMA, variations of fundamental frequencies, F L and F S) can be used to follow possible subtle changes in the local dynamical evolution of multidimensional systems. The overall characterization of chaos through the computation of some chaos indicators, like the MLE, depends on the collective behavior of the system and fails to identify localized differences in the dynamics, like for example the appearance of localized chaotic hot spots [Tieleman et al., 2014]. The study of the properties of the deviation vector distribution (DVD) related to the computation of the MLE has already been implemented as a mean to tackle this problem [Skokos et al., 2013;Ngapasare et al., 2019;Miranda Filho et al., 2019]. The numerical methods presented here can be implemented as alternative and complementary approaches for the same purpose.