Universal Exponent for Transport in Mixed Hamiltonian Dynamics
UUniversal Exponent for Transport in Mixed Hamiltonian Dynamics
Or Alus ∗ and Shmuel Fishman † Physics DepartmentTechnion-Israel Institute of TechnologyHaifa 3200, Israel
James D. Meiss ‡ Department of Applied MathematicsUniversity of Colorado, Boulder,Colorado 80309-0526 USA (Dated: November 17, 2018)We compute universal distributions for the transition probabilities of a Markov model for transportin the mixed phase space of area-preserving maps and verify that the survival probability distributionfor trajectories near an infinite island-around-island hierarchy exhibits, on average, a power law decaywith exponent γ = 1 .
57. This exponent agrees with that found from simulations of the H´enon andChirikov-Taylor maps. This provides evidence that the Meiss-Ott Markov tree model describes thetransport for mixed systems.
I. INTRODUCTION
A typical Hamiltonian system with more than one de-gree of freedom has a phase space that consists of regularand chaotic regions intertwined in a fractal structure. Inthis paper we focus on two-dimensional, area-preservingmaps that may arise from Hamiltonian dynamics by Pon-car´e section. For such 2D maps, phase space is parti-tioned by invariant circles that are absolute barriers, aswell as by partial barriers formed from hyperbolic invari-ant sets such as homoclinic trajectories or cantori [1].Invariant circles can enclose elliptic islands of stability,and these are typically embedded in chaotic zones in acomplex, island-around-island structure like that shownin Fig. 1(a).Two paradigmatic models of this dynamics areChirikov-Taylor’s standard map:( θ (cid:48) , J (cid:48) ) = ( θ + J + K sin θ, J + K sin θ ) , (1)and H´enon’s quadratic map:( x (cid:48) , y (cid:48) ) = ( − y + 2( κ − x ) , x ) . (2)The standard map was introduced by Chirikov and Tay-lor as a model for interaction between plasmas andelectromagnetic radiation [2]; it models dynamics nearany rank-one resonance. H´enon proposed his map asa paradigm for local behavior near an elliptic point [3];Karney et al. showed that it is a normal form for a genericsaddle-center bifurcation [4].Invariant structures in chaotic systems can be “sticky”,i.e., nearby trajectories may spend a long time in a neigh-borhood. More precisely, we say a region of phase space is ∗ Electronic address: [email protected] † Electronic address: fi[email protected] ‡ Electronic address: [email protected] sticky if its survival probability distribution—the prob-ability that a randomly chosen initial condition in theregion remains up to time t —has a power-law decay [5]: P sur ( t ) ∼ t − γ . (3)As Karney(and many others, subsequently) showed, theouter boundary of an elliptic island is sticky in this sense,though his numerical experiments showed strong fluctu-ations around what was inferred to be a power-law [6].MacKay, Meiss, and Percival (MMP) proposed thattransport through a connected chaotic componentbounded by elliptic islands could be described by aMarkov model [7]. It was later noted that the “states”in this model should be connected to form a tree (e.g.,Fig. 1), and transport on a Markov tree was studied byMeiss and Ott [8, 9]. Their model assumed that the treewas self-similar, based on renormalization theory both forthe flux through cantori near boundary circles [10] andfor islands-around-islands [11]. These calculations gave apower law (3) with γ = 1 . γ that decay only slowlyas t → ∞ [12]. Nevertheless, Cristadoro and Ketzmer-ick showed that correlated fluctuations in the self-similarscalings of Markov rates will result in a mean decay ex-ponent, (cid:104) γ (cid:105) , that depends upon the ensemble;moreover,if this ensemble is universal for mixed phase spaces, thenthe mean exponent well be universal as well [13]. Theirnumerical simulations of the dynamics of the H´enon map(without using the Markov tree model) give (cid:104) γ (cid:105) ≈ . a r X i v : . [ n li n . C D ] J u l ∅
1 11
S=10
011 111001 1011000 1001 1010 1011 1100 1101 1110 1111 p S → DS p DS → S p S → S0 p S → S1 p S0 → S Boundary Circle p S1 →S FIG. 1:
An elliptic island of (2) and a Markov tree. Each node is labeled by the state S . Several transition probabilities p S → S (cid:48) relatedto (4) and (5) are also indicated. The illustrative transition probabilities correspond to state S = 10. true, dynamical-system-based ensemble, (cid:104) γ (cid:105) ≈ . § II. In particular the resultsfor the survival exponent γ are presented there. In § IIIthe survival exponent of the H´enon map is calculatedfrom the standard map. The results are summarized anddiscussed in § IV II. THE MARKOV TREE MODEL
Consider a phase space with a sticky region formedfrom an island surrounding an elliptic fixed point such asthat depicted in Fig. 1(a). Here we recall the ideas andnotation for the Markov tree model for transport in theconnected chaotic component outside such an island [8].The fixed point is enclosed by a family of “class-zero”invariant circles, the outermost of which is the “bound-ary circle”; this circle is one component of the boundaryof the chaotic region. Typically there will be a family ofbroken circles,“cantori”, that are outside the boundarycircle and that limit upon it [1]. The flux of trajectoriesthrough these cantori limits to zero at the boundary cir-cle. This gives rise to a set of states in the chaotic regionencircling each island, called “levels”, that are boundedby “partial barriers.” In the Meiss-Ott model, these lay-ers correspond to successive rational approximations ofthe rotation number of the boundary circle.For the tree of states depicted in Fig. 1(b), the chaoticregion “far” from the sticky region corresponds to the“root” of the tree, denoted S = ∅ . For calculations ofthe survival probabiity (3) ∅ is viewed as absorbing. Theoutermost layer surrounding the class-zero boundary cir-cle gives rise to the state denoted by S = 1. Succes-sive layers are denoted by adding 1’s to the state, e.g., S = 111 denotes the third layer. Within each chaotic layer there will be a largest islandchain. The Meiss-Ott model assumes there is only onesuch island chain in each layer. Each chain also has aboundary circle, a “class-one” circle. The cantori sur-rounding a class-one circle gives rise to an additional setof chaotic layers. The outermost of these class-one lay-ers is denoted S = 10. Successive layers approachingthe class-one boundary in state 10 again correspond toadding one’s to the state, e.g., S = 10111 . . . .This construction generalizes to each layer: near aclass-one boundary there are encircling periodic orbitsgiving rise to class-two islands, etc. The assumption thatthere is one island chain in each layer implies that the treeis binary.Transport in the connected chaotic region outside all ofthe boundary circles is thus represented by a sequence oftransitions on the tree (levels and classes). If the trans-port is Markov, it is defined by transition probabilities p S → S (cid:48) for each pair of connected states, recall Fig. 1. Theprobability of such a transition is determined by the fluxof trajectories, i.e., the area of the turnstile in the can-torus that separates the states [7]. We denote this fluxby ∆ W S,S (cid:48) = ∆ W S (cid:48) ,S ; it is symmetric because the netflux through any region of phase space must be zero. Theflux through a cantorus can be computed by the MMPaction principle [7].The average transit time through a state bounded bysuch partial barriers is exactly equal to the area of theaccessible portion of phase space in the state S , A S , di-vided by the exiting flux [14]. If these transit times arelong enough that correlations are unimportant, one canassume that the transition probability is p S → S (cid:48) = ∆ W S,S (cid:48) A S and that a Markovian approximation is valid.The only nodes that are connected on the tree areparent-daughter nodes. The daughters of a state S = s s s . . . s j are denoted by concatenation: S S S , obtained by deleting the lastsymbol, is denoted DS . There are two important transi-tion probabilities, p S → DS for moving “up” from state S to its parent, and p S → Si for moving “down” from a state S to its i th daughter. It is convenient to categorize thechange in transition probabilities from state to state bythe two ratios w ( i ) S = p S → Si p S → DS = ∆ W S,Si ∆ W S,DS , (4) a ( i ) S = p S → Si p Si → S = A Si A S . (5)When the tree is self-similar the ratios (4)-(5) are inde-pendent of the state S , though they depend on the choiceof class, i = 0, or level, i = 1 [8].We previously computed the flux ratios (4) for a num-ber of states and a range of parameter values of the H´enonmap in [15, Eq. (22) and Fig. 11]. The ratios were cal-culated using fluxes through periodic orbits of (2) as aproxy for the cantorus fluxes. The distributions of fluxratios are different for class and level transitions, so welabel them as f i ( w ), see Fig. 2(a). Each distributiondoes not depend systematically on the parameter κ of(2). In [15] we compared the f ( i ) ( w ) distributions ob-tained by choosing the parameters uniformly two inter-vals, − . < κ < .
25, and 0 . < κ < .
75, (see Fig. 10and Eq. (22) there). For the current paper, we repeatedthis computation using nonuniform κ distributions (notshown here): the new f ( i ) ( w ) do not differ significantlyfrom those in Fig. 2(a). This gives us more confidencein the assertion that there are a “universal” distributionsfor class and level flux ratios.Here we also extend these results by computing thearea ratios (5) for a number of orbits of the H´enon map.Areas were computed for the periodic orbits giving outerrational approximations to the boundary circle rotationnumbers up to the states in the fourth generation on thetree; ( S = 1001, 1010, etc.); for (2) with κ ∈ [ − . , . f ( i ) ( a ) again differ significantly forclasses and levels, but they still appear to be universal inthe sense that they do not depend systematically on κ .Indeed, Fig. 3 shows separate area distributions for theintervals − . < κ < .
25 and 0 . < κ < . f ( i ) ( a, w ), see Fig. 4. Note that the area and flux ratios exhibit significant correlations, since the probabili-ties are concentrated on irregular regions in ( a, w )-space.Below we use these joint densities to draw values of a and w to give Markov trees with random scalings thatcorrespond, at least according to these first-order statis-tics, to those of the true map. That is, we assume thatthe scaling factors on the different branches of the treeare independent random variables drawn from the em-pirically computed f ( i ) ( a, w ) found from the first threegenerations of islands and levels for (2). This contrastswith [8] where the ratios for each level and each classbranch do not vary with depth on the tree.It is important to note that we did not compute thetrue flux through cantori, nor the true accessible area inany state: we assume that the transition probability ratesthrough the turnstiles of the cantori scale in the sameway as those through the numerically computed periodicorbits. Computing the true cantorus flux is considerablymore difficult since it must be done using a high-periodapproximation to the unstable, quasiperiodic cantorus.In the next subsection, we compute the survival prob-ability exponent γ from Monte Carlo simulations on ran-dom trees. In § II B we compare these results to a masterequation approach.
A. Monte Carlo Simulations
For a Markov tree model, the vector of densities ateach state on the tree can be denoted by an infinite vector (cid:126)ρ = ( ρ ∅ , ρ , ρ , ρ , . . . , ρ S , . . . ), where ρ S is the densityat state S . If the per-step transition probability is small,transport on the tree is governed by the master equation d(cid:126)ρdt = W (cid:126)ρ , W S,S (cid:48) = p S (cid:48) → S − δ S,S (cid:48) (cid:88) S (cid:48)(cid:48) p S → S (cid:48)(cid:48) . (6)The absorbing state, ∅ , is treated by setting p ∅ → S = 0for the daughter states S = 1 or 0. To set an overall time-scale we choose p → ∅ = 0 .
1. The remaining probabilitiesare determined by the ratios (4)-(5), which are drawnfrom the distributions f ( i ) ( a, w ) shown in Fig. 4.Though the tree is infinite, the probabilities for tran-sitions decrease rapidly with level and class, and thus itis reasonable to truncate the tree at a finite number, B ,of branches or generations. The states in the B th gener-ation are connected only to their parents: only p S → DS isnonzero. This gives a finite tree with 2 B states.To perform the Monte Carlo experiment, we chose 10 particles with initial states drawn from a distributionsatisfying detailed balance on the tree [13]. This is anequilibrium of (6) when the absorbing state is removed,and since transient behavior is absent, algebraic decay iseasier to observe. For such a distribution, the survivalprobability exponent is γ − B = 17, and av-eraging over 70 realizations of the tree we find γ ≈ . B = 10 and 50 realizations we find γ ≈ . (a) w -8 -7 -6 -5 -4 -3 -2 -1 f ( i ) ( w ) ClassLevel (b) a -6 -5 -4 -3 -2 -1 f ( i ) ( a ) ClassLevel
FIG. 2: (Color online) Distribution densities of (a) area scalings a and (b) flux scalings w for the H´enon map with − . < κ < . (a) a -4 -3 -2 -1 f ( C l a ss ) ( a ) ! : < < : : < < : (b) a -6 -5 -4 -3 -2 -1 f ( L e v e l ) ( a ) ! : < < : : < < : FIG. 3: (Color online) Distribution densities of (a) class and (b) level area scalings for the H´enon map with parameters chosen in − . < κ < .
25, in gray, and in 0 . < κ < .
75, in red (light gray). (a) a Class -3 -2 -1 w C l a ss -6 -5 -4 -3 -2 -1 (b) a Level -5 -4 -3 -2 -1 w L e v e l -7 -6 -5 -4 -3 -2 -1 FIG. 4: (Color online) Histograms of the joint probability distributions for w and a . (a) f ( Class ) ( a, w ) taken from 2629 class transitions.(b) f ( Level ) ( a, w ) taken from 3608 level transitions. B. The master equation on the tree
Here we will compare the Monte Carlo simulationswith a direct calculation using eigenvalues λ n and eigen- states (cid:126)ρ n of the 2 B × B transition matrix W . The evo-lution of given initial state (cid:126)ρ (0) then becomes (cid:126)ρ ( t ) = B (cid:88) n =1 A n (cid:126)ρ n e − λ n t , A n = (cid:104) (cid:126)ρ n † | (cid:126)ρ (0) (cid:105) where (cid:126)ρ n † is the left eigenvector of W . The survivalprobability is P sur ( t ) = (cid:88) S (cid:54) = ∅ ρ S ( t ) = (cid:88) S (cid:54) = ∅ B (cid:88) n =1 A n ρ n S e − λ n t . (7)where ρ S ( t ) and ρ n S are the S th component of (cid:126)ρ ( t ) andof (cid:126)ρ n , respectively. To compute (7), a reasonable initialcondition is (cid:126)ρ = (0 , , , , ... ) . As before we use the empirical distributions f ( i ) ( a, w )for the ratios (4)-(5) to generate N = 200 realizations ofa Markov tree. Choosing B = 10, the decay of P sur ( t )appears to be a power law up to t = 10 , see Fig. 5.The exponent, computed using a least-squares fit fromthe average (cid:104) log ( P sur ( t )) (cid:105) for 10 . ≤ t ≤ (withequally spaced points on a logarithmic scale) is γ = 1 . ± . γ is estimated from individual realizations:for the upper (lower) bound the product t γ ± P sur ( t ) ex-hibits an increasing (decreasing) behavior on a log-logscale for all realizations but one, see Fig. 6. The com-puted value of γ does not change significantly for larger B . The same result is found if one first computes γ foreach realization, recall Fig. 5, and then average the re-sults.In Appendix A, we demonstrate how a power law canarise from a sum of infinitely many exponential decaysaccumulating on λ = 0. t P s u r -30 -25 -20 -15 -10 -5 FIG. 5:
Plot of the survival probability P sur vs time for 200 re-alizations of the sum (7). The heavy line is the average, decayingasymptotically with the slope γ ≈ . III. STICKINESS OF ACCELERATOR MODES
For large enough values of K , the standard map (1)exhibits special, accelerator orbits for which the momen-tum increases by a multiple of 2 π each period [2, 16, 17]. These are due to the vertical 2 π periodicity of (1). In-deed taking J mod 2 π , accelerator modes are periodicorbits created in saddle-center bifurcations. The simplestof these, at K = 2 πn for integer n , creates two saddle-center pairs; one pair accelerates upward and the otherdownward. Near a saddle-center bifurcation the local dy-namics are modeled by the H´enon map [4]. The ellipticpoints created in these bifurcations remain stable for asmall range of K , and their neighborhoods are thereforeislands like that in Fig. 1(a).In a regime where there are accelerator islands, thevertical transport in the chaotic component outside theislands is dominated by the stickiness of the islands: tra-jectories are trapped near the islands with a survivalprobability (3). This results in super-diffusion of the mo-mentum [6, 18–20].Whenever there is an island with positive acceleration,there is a one with negative acceleration, and the momen-tum transport can be treated as a random walk betweenthese modes; statistically this is a drunkard’s [6, 19] ora L´evy [18] walk. When trajectories are not stuck, theydiffuse in momentum, but the contribution of this givesa negligible contribution to momentum transport.To estimate the exponent γ , we divide each trajectoryinto segments that are trapped either near an upwardor a downward propagating accelerator island. Since theupward propagating mode occurs near θ = π/ θ = − π/
2, if we take − π ≤ θ < π , atransition between the upward and downward motion iscorresponds to a change in sign of θ . When a trajectory isnot trapped near an accelerator island, the probability tostay in the same half of the cylinder decays exponentially;therefore the long-time survival probabilities in each halfof the cylinder will be dominated by the power-law decaydue to the accelerator islands. Computations averagedover 45 parameter values give a survival exponent γ ≈ . IV. RESULTS AND DISCUSSION
We have computed—for the first time—the joint dis-tribution of flux (4) and area (5) ratios for states definedby the island-around-island structure of the H´enon map(2). To do this, we assumed that the ratios for cantoriscale in the same way as those for periodic orbits. Thesedistributions appear to be universal: they do not dependthe parameter κ of the H´enon map (2) in any system-atic way, and this map is the universal local model for anisland of an area-preserving map.Using the Markov tree model, we computed the result-ing power-law decay for the survival probability (3) bothby Monte Carlo simulations and directly from diagonal-ization of the transition matrix. The mean survival ex-ponent γ depends only on the distributions of the scalingratios and not on the particular realization of the tree, inagreement with [13]. Since the scaling distributions are (a) t t : P s u r -3 -2 (b) t t : P s u r -6 -5 -4 -3 -2 -1 FIG. 6: (Color online) Plot of the survival probability t γ ± P sur a) γ + = 1 . γ − = 1 .
4, vs time for 200 realizations of the sum (7). universal, the survival exponent γ is universal as well.Our results are consistent with γ = 1 .
57. This is alsothe value found by direct simulations of the H´enon mapin [13]. Here, we also found this same exponent for thestickiness of accelerator modes of the standard map.Thus it appears that the Markov model successfullypredicts the algebraic decay exponent observed in sim-ulations. The assumptions of the Markov property, thebinary structure of the tree, and the use of periodic or-bits instead of cantori for the ratios do not negativelyimpact the results. Therefore the Markov tree is an ef-fective model for the long-time dynamics of transport inarea-preserving maps with a mixed phase space.It remains an open question whether the fluctuationsin γ are real: namely, do they survive the t → ∞ limit? Acknowledgments : We would like to thank RolandKetzmerick, Arnd B¨acker, Holger Kantz and Oded Agamfor fruitful discussions. OA and SF would like to acknowl-edge support of Israel Science Foundation (ISF) grants1028/12 and 931/16, and the US-Israel Binational Sci-ence Foundation (BSF) grant number 2010132 and bythe Shlomo Kaplansky academic chair. OA acknowl-edge the support of the Guthwirth foundation excellencefellowship. SF thanks the Kavli Inst. for Theor. Phys.for its hospitality, where this research was supported inpart by the US National Science Foundation (NSF) undergrant NSF PHY11-25915. JDM was partially supportedby NSF grant DMS-1211350.
Appendix A: Eigenvalue Asymptotics
A natural question is: how can the master equation(6) give rise to power-law decay? Indeed, for any finitematrix size, the long-time decay of P sur ( t ) will be expo-nential, at the rate of the smallest eigenvalue of W , say λ . Nevertheless, if B is large enough, then the decaydoes look like the power-law (3) for a finite time, as we saw in Fig. 5.For an infinite chain, a power-law decay over infinitetime can occur. For a Markov chain (e.g., keeping onlythe level transitions), the power law can arise from a sumof the form P sur ( t ) ∼ (cid:80) n δ n e − (cid:15) n t , implying that the eigen-values and weights decrease geometrically [21]. Inspiredby this idea, we note that the long-time behavior of thesurvival probability depends upon the density of smalleigenvalues. Approximating the discrete spectrum by acontinuum, then the sum over eigenstates in (7) becomesan integral, P sur ( t ) ∼ (cid:90) λ (cid:88) S (cid:54) = ∅ ρ n ( λ ) S A n ( λ ) e − λt (cid:12)(cid:12)(cid:12)(cid:12) dλdn (cid:12)(cid:12)(cid:12)(cid:12) − dλ. (A1)We now suppose that for large n , instead of the geometricdecay of [21], we have λ n ∼ n − δ , (cid:88) S (cid:54) = ∅ ρ n S A n ∼ n − δ . (A2)To support this hypothesis, we again use the distributionsof Fig. 4 to compute the eigenvectors and eigenvalues of W . The results, shown in Fig. A.1 for one realizationof the tree, show that both of these quantities decreasealgebraically with the estimates δ = 5 . δ = 8 . P sur ( t ) ∼ (cid:90) λ λ η e − λt dλ ∼ t − η − , η = δ − δ − δ . (A3)For the realization in Fig. A.1 this gives η = 0 . N = 200 realizations. This assumes thatthe distributions of δ , δ , and thus of η , are narrow sothat one can use the average curve to estimate δ and δ .Computing the exponents for the averaged curves from (a) n P s = @ ; n s A n -30 -25 -20 -15 -10 -5 y = - 8.8*x + 2.1 (b) n n -16 -14 -12 -10 -8 -6 -4 -2 y = - 5.1*x + 2.1 FIG. A.1:
Empirical verification of the power laws (A2) for a single realization of the matrix W in Fig. 5. (a) Plot of (cid:80) S (cid:54) = ∅ A n ρ n S vs n , and a fit with slope δ = 8 .
8. (b) Plot of the eigenvalues λ n of W vs n , leading to the slope δ = 5 . points uniformly distributed on a log scale for n , we find (cid:104) η (cid:105) ≈ . ± .
16 where the standard deviation is takenas the error. Therefore for the average exponent, (A3)and (3) imply γ = (cid:104) η (cid:105) + 1 = 1 . ± . . Alternatively ifthe fit is done using all values of n (i.e., uniform on thescale of n ) we obtain (cid:104) η (cid:105) ≈ .
490 and γ = (cid:104) η (cid:105) + 1 ≈ . . (A4)Finally, if we instead calculate slopes for each realizationand then compute (cid:104) δ (cid:105) and (cid:104) δ (cid:105) and then use (A3) tofind (cid:104) η (cid:105) , we find γ ≈ .
70 and γ ≈ .
48 for the two fit-ting methods described above, respectively (uniform inlog ( n ), and uniform in n ). Of the two fits, the latterseems to more appropriately weight the long-time behav-ior due to the small eigenvalues. Appendix B: Survival Exponent for the Standardmap
Following Karney in [6], we can compute the survivalprobability from the statistics of the duration of thetrapped segments. For example, for a single trajectoryof length T with N segments, denote the number of seg-ments of duration τ by N τ . Then the probability that asegment has length τ is p τ = N τ /N . However, to cor-rect for the finite time of the simulations—which over-estimates the probability of observing a short trajectory,Karney showed that one should use p τ = N τ N TT + 1 − τ . The cumulative survival probability is then P sur ( t ) = T (cid:88) τ = t +1 p τ . (B1) Following the method discussed in the main text to com-pute N τ , we computed the P sur for 30 values of K chosenfrom equal steps of 0 .
025 in the interval [2 π, .
8] for thosecases which had well established super-diffusion; that is,for which no “singular” islands were present [20, 22]. Sin-gular islands correspond to parameters near the saddle-node, tripling (twistless) and period doubling bifurca-tions. The omitted parameters also correspond to casesin which the calculation of boundary circles in [15] failed.Using a fit with points chosen uniformly in log t in theinterval 10 ≤ t ≤ gives γ values that range over[1 . , .
7] with an average γ = 1 .
604 . Adding 15 morevalues of K in the interval [6 . , . γ ≈ . . The results for all 45 parameter values are shownFig. B.2. t P s u r -12 -10 -8 -6 -4 -2 FIG. B.2:
The survival probability (B1) for simulations of (1).The heavy line is the average over 45 parameter values (see text)resulting in γ ≈ . [1] J. Meiss, Rev. Mod. Phys. , 795 (1992).[2] B. Chirikov, Physics Reports , 263 (1979).[3] M. H´enon, Q. J. Appl. Math. , 291 (1969).[4] C. Karney, A. Rechester, and R. White, Physica D: Non-linear Phenomena , 425 (1982).[5] S. R. Channon and J. L. Lebowitz, Annals of the NewYork Academy of Sciences , 108 (1980).[6] C. Karney, Physica D , 360 (1983).[7] R. MacKay, J. Meiss, and I. Percival, Physica D , 55(1984).[8] J. Meiss and E. Ott, Physica D , 387 (1986).[9] J. Meiss and E. Ott, Phys. Rev. Lett. , 2741 (1985).[10] J. Greene, R. MacKay, and J. Stark, Physica D , 267(1986).[11] J. Meiss, Phys. Rev. A , 2375 (1986).[12] R. Ceder and O. Agam, Phys. Rev. E , 012918 (2013).[13] G. Cristadoro and R. Ketzmerick, Phys. Rev. Lett. ,184101 (2008). [14] J. Meiss, Chaos , 139 (1997).[15] O. Alus, S. Fishman, and J. Meiss, Phys. Rev. E ,062923 (2014).[16] A. Rechester and R. White, Phys. Rev. Lett. , 1586(1980).[17] A. Rechester, M. Rosenbluth, and R. White, Phys. Rev.A , 2664 (1981).[18] G. Zumofen and J. Klafter, Europhys. Lett.) , 565(1994).[19] R. Ishizaki, T. Horita, T. Kobayashi, and H. Mori, Prog.Theor. Phys. , 1013 (1991).[20] G. Zaslavsky, M. Edelman, and B. Niyazov, Chaos ,159 (1997).[21] J. Hanson, J. Cary, and J. Meiss, J. Stat. Phys. , 327(1985).[22] H. R. Dullin, J. D. Meiss, and D. Sterling, Nonlinearity13