Creation mechanism of Devil’s Staircase surface and unstable and stable periodic orbits in the anisotropic Kepler problem

In the two-dimensional Anisotropic Kepler Problem (AKP), a longstanding question concerns the uniqueness of an unstable periodic orbit (PO) for a given binary code (modulo symmetry equivalence). In this paper, a finite-level (N) surface defined by the binary coding of the orbit is considered over the initial-value domain D0. A tiling of D0 by base ribbons of the surface steps is shown to be proper, i.e., the surface height increases monotonously when the ribbons are traversed from left to right. The mechanism of creating a level-(N+1) tiling from the level-N tiling is clarified in the course of the proof. There are two possible cases depending on the code and the anisotropy. (A) Every ribbon shrinks to a line as N → ∞ . Here, the uniqueness holds. (B) When future (F) and past (P) ribbons become tangential to each other, they escape from the shrinking. Then, the initial values of a stable PO (S) and an unstable PO (U) sharing the same code co-exist inside the overlap of the F and P non-shrinking ribbons. This corresponds to Broucke’s PO. When the anisotropy is high, only case A is observed; however, as the anisotropy decreases, a bifurcation of the form U ( R ) → S ( R ) + U ′ ( NR ) occurs along with the emergence of a non-shrinking ribbon. (Here, R and NR denote self-retracing and non-self-retracing POs, respectively). We conjecture that, from a classification based on topology and symmetry, case B occurs only for odd-rank, Y-symmetric POs. We report two applications. First, the classification is applied successfully to the successive bifurcation of a high-rank PO (n = 15), where the above bifurcation is followed by S ( R ) → S ′ ( R ) + S ″ ( NR ) . Second, enhancing the sensitivity to the co-existence of S and U POs through ribbon tiling, we examine the high-anisotropy region. A new symmetry-type POs (O-type) are found and, at γ = 0.2, all POs are shown to be unstable and unique. An investigation of 13648 POs at rank 10 verifies that Gutzwiller’s action formula works with amazing accuracy.


Abstract
In the two-dimensional Anisotropic Kepler Problem (AKP), a longstanding question concerns the uniqueness of an unstable periodic orbit (PO) for a given binary code (modulo symmetry equivalence). In this paper, a finite-level (N) surface defined by the binary coding of the orbit is considered over the initial-value domain D 0 . A tiling of D 0 by base ribbons of the surface steps is shown to be proper, i.e., the surface height increases monotonously when the ribbons are traversed from left to right. The mechanism of creating a level-(N+1) tiling from the level-N tiling is clarified in the course of the proof. There are two possible cases depending on the code and the anisotropy. (A) Every ribbon shrinks to a line as  ¥ N . Here, the uniqueness holds. (B) When future (F) and past (P) ribbons become tangential to each other, they escape from the shrinking. Then, the initial values of a stable PO (S) and an unstable PO (U) sharing the same code co-exist inside the overlap of the F and P non-shrinking ribbons. This corresponds to Broucke's PO. When the anisotropy is high, only case A is observed; however, as the anisotropy decreases, a bifurcation of the form  + ¢ ( ) ( ) ( ) U R S R U NR occurs along with the emergence of a non-shrinking ribbon. (Here, R and NR denote self-retracing and non-self-retracing POs, respectively). We conjecture that, from a classification based on topology and symmetry, case B occurs only for oddrank, Y-symmetric POs. We report two applications. First, the classification is applied successfully to the successive bifurcation of a high-rank PO (n=15), where the above bifurcation is followed by  ¢ +  ( ) ( ) ( ) S R

Introduction
We study the two-dimensional Anisotropic Kepler Problem (AKP) (1), which admits a binary coding of an orbit. The focus of our attention is the coexistence of unstable (U) and stable (S) periodic orbits (POs) of the same code. We clarify the mechanism of a salient bifurcation process whereby, under decreasing system anisotropy, U bifurcates into S and a new type of unstable PO ( ¢ U ). This bifurcation has the form  + ¢ ( ) ( ) ( ) U R S R U NR , where R and NR denote self-retracing and non-self-retracing POs.
Let us first look at AKP, which has been a vital testing ground for quantum chaos ever since its cultivation by Gutzwiller [1][2][3][4][5][6][7] to test his periodic orbit theory (POT) in the high-anisotropy regime (  g n m º 0.2). Because of the Coulomb interaction, this is not a KAM system, and reveals an abrupt transition from the integrable limit if γ is reduced from one [8]. Interestingly, the region 0.7γ0.85 exhibits critical statistics, reminiscent of the Anderson localization in condensed matter physics [8][9][10][11]. Recently, we have successfully extracted information of 2-dim low-rank unstable POs (sequence, period, and even the Lyapunov exponent) from quantized 3-dim AKP levels via the inverse use of POT and suitable symmetrization [12]. We have also shown that the quantum scar in AKP is robust and survives under the successive avoidance of level crossings during large variations in anisotropy [13]. Note that AKP can be realized by a semiconductor with donor impurity, and experimental manufacturing toward this goal has been reported [14,15].
There is a longstanding question regarding POs in AKP. It was shown by Gutzwiller [2], and rigorously probed by Devaney [16], that, for an arbitrary given binary sequence (or code), there exists at least one initial point for a PO that evolves realizing the sequence, provided that γ<8/9 3 . Here, the PO is isolated and unstable.
The uniqueness of a PO for a given code, however, was conjectured by Gutzwiller following extensive numerical investigation, whereas Broucke found two stable POs as a counterexample [17] in the low-anisotropy region 4 . An overview up to 1990 is given in Gutzwiller's book [18].
Recently, Contopoulos and Harsoula found that there are many stable POs, not necessarily in the low anisotropy regime. This implies that AKP is not an Anosov system [19]. Unfortunately, their PO search was limited to the case of perpendicular emission from the heavy axis, and is therefore not capable of correctly detecting ¢( ) U NR in the case of the bifurcation  + ¢ ( ) ( ) ( ) U R S R U NR , as ¢( ) U NR does not have any perpendicular crossings of the heavy axis.
The remainder of this paper is organized as follows. In section 2, three basic points are established.
(1) One-time map. A suitably compactified 2-dim rectangle D is chosen as an initial value domain of the AKP flow and the one-time map   D D : 0 1 is determined. The map  | I restricted to the collision manifold Ì ( ) D I was investigated by Gutzwiller [2], and we integrate his wisdom into a convenient form. This is extended to the interior combining numerical analysis (figure 5).
(2) A level-N devil's staircase surface (DSS) and associated tiling of D by the base ribbons of the steps (figure 3).
(3) PO trapping. An initial point of a PO of a given code must lie in the cross-sectional area of a future (F) and a past (P) ribbon of the respective F and P code 5 .
In section 3, the creation of a level-(N+1) DSS from the level-N DSS is clarified. The mechanism gives a simple proof of the properness of the tiling (monotonicity of the surface)-traversing ribbons from left to right, the height of the step always increases. This mechanism also demonstrates the possibility of a non-shrinking ribbon for large values of N.
In section 4, we explicitly show that there are two fates of ribbons at large N, depending on the code and the anisotropy.
(A) In most cases, both the F and P ribbons shrink at large N. The PO trapped in the cross-sectional area is then singled out at the crossing of the F and P asymptotic curves.
(B) Below a certain anisotropy level, an exceptional ribbon emerges that stops the shrinkage at finite N. 6 We show that this occurs when the F and P ribbons become tangential to each other. The initial points of the U and S POs are then both contained in the cross-sectional area of non-shrinking F and P ribbons. The emergence of such non-shrinking ribbons in turn induces a salient code-preserving bifurcation:  + ¢ ( ) ( ) ( ) U R S R U NR within Y-symmetric POs (figures 9-11). We explain this pattern by considering the symmetry and topology of the PO (figure 14), and conjecture that it should occur for Y-symmetric, rank-odd POs.
Section 5 is devoted to the application of ribbon tiling to the PO search in AKP. To search for POs of all possible symmetry types, brute force shooting would have to be conducted over the vast two-dimensional free parameter space, which is not feasible for high-rank POs. However, a two-step search-to identify the relevant ribbon and then the PO within it-clearly reduces the task to a quasi-one-dimensional problem. Based on this idea, we devise a concrete algorithm for an exhaustive search that is sensitive to the coexistence of S and U POs. First, we take a highrank PO (n = 15) as an example and show that it does bifurcate as  + ¢ ( ) ( ) ( ) U R S R U NR . We pursue the final POs to the very-low-anisotropy region (γ>8/9), and find a second bifurcation figure 16). We find these successive bifurcations are in accordance with the topology and symmetry classification. We then turn to an extensive PO search, extending the rank to n=10 and considering high anisotropy γ=0.2. This is important because there are questions regarding the possible existence of stable POs even in the highanisotropy regime [19]. We find all 19284 isolated unstable POs predicted from symmetry analysis (with a correction for the new symmetric type (O type) PO at ranks seven and nine) and nothing else. This confirms (within the set-up resolution) the uniqueness of a PO for a given code (modulo symmetry equivalence) at γ=0.2. The above unstable POs are used to verify Gutzwiller's approximate action formula. The actions of 13648 POs at rank n=10 are predicted (with the only two parameters of the formula fixed by 44 POs at n = 5) with amazing accuracy (mean squared deviation (MSD)=0.053 6; see figure 23 and table 4).
Part of this paper was reported at the ISQS25, Prague and this is an extended version of the published conference report [20].
2. Devil's staircase surface and the one-time map 2.1. Gutzwiller's rectangle The Hamiltonian of the two-dimensional AKP is given by where º º u p v p , x y and the mass anisotropy is proportional to g -1 with g n m º <1. Because m n = >̈( ) y x y x y x, the orbit tends to cross the heavy x-axis more frequently. Thus, the Poincaré surface of section (PSS) is specified by the condition y=0 in the phase space, and we can encode the orbit using the code a i =±1, the sign of x i at the i-th crossing of the orbit with the x-axis. The future code is (a 0 , a 1 , L) and the past code is (a 0 , a −1 , L) 7 . The constant energy condition determines the kinematically allowed region on the PSS (the physical region for short). As the potential is homogeneous in coordinates, the system has a scaling property. Any value for the energy is equivalent, and we take H=−1/2 by convention [2]. Under y=0 and H=−1/2, the domain of the initial coordinates in the PSS is which is the 'lips'-like region in figure 1(a) 8 . The part of D at which r=0 (the location of the collision) is the core line of D (the u-axis) where = ¥ v and   -¥ ¥ u . Following Gutzwiller, we use an area-preserving map ) in D and D have the same area. (c) Equi-r contours in D approach the backbone at  r 0.
7 To treat the future and past DSS on an equal footing, we include a 0 in both sequences. The choice of the present (i = 0) is immaterial because of the time-translation invariance. 8 In fact, the physical region is a double-cover of figure 1(a), because v can be either positive or negative. Hereafter, we limit our consideration to D with v0, as the case v < is obtained by a simple reflection with respect to the y−axis.
The collision occurs on the I-shaped backbone of D: We call I the collision manifold.

Symbolic coding of orbits and devil's staircase surface
The future and past surfaces over D at level N are respectively described by the height functions å å Figure 2 shows the first three future DSS as examples. These were calculated from a j (X 0 , U 0 ) of each orbit, starting from a site (X 0 , U 0 ) of a fine lattice set on D.
By definition (6), the level-N surface has + 2 N 1 heights. These examples reveal rather a simple structure. Firstly, there is only one step for each height in a surface. We call the base of each step a ribbon (see figure 3) and label it by the code ¼ ( ) a a , , N 0 of the step-height. Then, each ribbon fully extends from U 0 =−B to U 0 =B, and there are no isolated bubbles in the base domain. Therefore, the 2 N+1 ribbons together constitute a tiling of the base domain D. Furthermore, the height of the steps increases monotonously if we traverse the ribbons from left to right. We describe this by saying that the ribbons are properly tiling D.  The reason why the surface is organized in this way at every N is by no means trivial, and we give a proof in section 3. A remark is in order about our approach based on the finite-N (or coarse-grained) DSS. The edges of the ribbons constitute the unstable (stable) manifolds of the collision for the future (past) DSS, because the changes in the code are invoked by the collisions (see figure 6). They can be calculated by Gutzwiller's technique of the collision parameter (section IIIB in [2] ). In a way, our approach is the dual of Gutzwiller's. However, the reason for the monotonic ordering of the manifolds is very hard to see in the latter. More importantly, the proof of proper tiling (by mathematical induction regarding N) in turn clarifies how the higher-level surfaces are created by the repeated application of the one-time map. Through this mechanism, we can then clarify the coexistence of POs having the same code.
Let us now consider how a PO is related to ribbons. The code of a PO is cyclic, that is, a a  a  a  a a  a  , , , ; Here, the first 2n bits are the primary part of the PO and the integer n is its rank 9 . If the initial point of the PO is where the pre-factor takes care of the repetition of the primary cycle. We can then prove the following: A step of the level-N DSS whose height is given by the first N+1 bits of the string (7): is closest to the PO in height. Therefore, its ribbon traps the initial point * * ( ) X U , 0 0 (see figure 3).
To see this, consider the error between z F PO and z N F . This is simply the contribution truncated away from the string (7) under the coarse-graining to level N. Therefore, and we find that  d + | | 1 2 N 1 . From another perspective, the difference in height between neighboring steps is D = 1 2 N N . Therefore, the selected step is the closest 10 . At the large-N limit, the error δ vanishes and the PO asymptotically sits on the selected step, regardless of whether the enclosing ribbon shrinks or not. +

One-time map
We now investigate the one-time map from a rectangle D 0 on a certain PSS onto D 1 on the next PSS: This is a challenging problem, because the AKP flow involves both hyperbolic and elliptic singularities. Separation and blow-up inevitably occur, and brute-force numerical integration alone is not sufficient to identify the salient features. However, a map restricted to the collision manifold,  | I , was derived by Gutzwiller in [2].
where (ϑ 0 , ψ 0 ) and (ϑ 1 , ψ 1 ) parameterize the initial and final I, respectively (ψ 0 and ψ 1 are either 0 or π, because y=0 on PSS. ψ 0 is taken to be 0 as a choice of the fundamental initial domain). Thus, all we need to obtain  | I in terms of (X, U) is to properly pull back and push forward (15): which is valid on I. The above outline for obtaining  | I is substantiated in appendix A. A careful examination is necessary because, firstly, the initial domain of the map  is divided into three regions (1, 2, 3) under separation by hyperbolic singularities (see the flow in figure A1), and furthermore, the pulling back and pushing forward procedure in (16) introduces additional criticality (when ϑ 0 passes through π/2 or when J 1 passes through 3π/2). Consequently, the collision manifold I is eventually divided into five regions A B C 1, 2 , 2 , 2 , 3, and we have to consider how each is mapped by  | I . This issue is discussed in appendix A and the results are summarized in table 1. Figure 4 is a graphical representation of the table.
We observe that (i) the boundary map  | I acts on the boundaries of ++ ( )and +-( ) separately, sending them to the left and right, respectively, and that (ii)  +-| ( ) rotates the boundary of +-( ) by one sub-region. The virtue of collision manifold analysis is that  | I extends to the interior and helps to grasp the behavior of the internal map  -the point emphasized by Devaney [16].
In figure 5, we combine figure 4 with the numerically calculated interior map. Indeed, the above two characteristics of  | I clearly control the interior map. The middle regions +-( ) and -+ ( ) swap order to become +-¢ ( ) and -+ ¢ ( ) in the image, corresponding to the separation property (i) of  | I . Remarkably, the distortion of the orthogonal lattice in D 0 shows that there exist focus points in the interior map, for instance, ++ v 1 and +-v 1 , which are the images of the sides of the separator curve C:   ++ ++ For the notation of critical objects, see appendix B. Furthermore, the rotation property (ii) implies that the vertical lines connecting 2A and 2B in D 0 must be bent to connect the images ¢ A 2 and ¢ B 2 , which are both in the bottom boundary of D 1 . This produces wing-shaped curves in D 1 that exhibit a folding property in  . The wings contract to a single point + c 0 , which is also a focus point to which the X 0 =0 ends of the horizontal lines in +-( ), namely  + , are mapped: c . 1 We investigate these phenomena in more detail in appendix B. The occurrence of blow-up and contraction may appear ad-hoc, but they are quite systematic. This can be seen if we follow the seven parts on the boundary of +-( ) (and their images) in table 1 counter-clockwise: Here, ℓ ( ) C and  ℓ ( ) denote the length of the separator curve C and of  , , , ). The blow-up and contraction of a part are indicated by ∧and ∨symbols, respectively. Note that +-( ) and +-¢ ( ) are congruent to each other, because +-= +-¢ ( ) ( (( ) )) P T , see (C.5). Indeed, the blow-up and contraction are involved a symmetrically equal number of times, so that the perimeters of +-( ) and +-¢ ( ) intriguingly remain the same. The reason why a curve is contracted into a point and a point is blown up to a curve is as follows. In figure 6, we show the result of our collision parameter analysis, which follows that of Gutzwiller [5,7].
The collision orbits are created as a one-parameter family with parameter A. The first short-time interval is analytically calculated and continued by numerical integration. The positions of the first arrivals on the Poincaré surface form a curve in the ( ) X U , 1 1 -rectangle developed by the parameter A-this is nothing but the separator curve C. This is a typical sample of the blow-up process. Note that one can reversely track back by numerical integration (with slow-down) starting from +-C 1 , but it is difficult to follow the motion after closely approaching the upper end of  .
Continuing the collision trajectory for many crossings of PSS numerically produces the stable and unstable manifolds, as found by Gutzwiller (figure 2 in [2] and figure 26 in [5]). In a way, our approach and collision analysis are complementary to each other. We focus on finite-N ribbons, whereas the collision analysis focuses on ribbon-ribbon curves. That is, the longitudinal line separating two ribbons is nothing but the stable/unstable manifold, because the code change is just induced by the collision. The difference is that, although determining which hierarchy a line is subject to is a difficult problem in collision parameter analysis, our approach has the advantage of simple book-keeping. Hence, it can easily be used to locate the PO directly, and helps to resolve the coexistence of stable and unstable POs of the same code. We discuss these points in detail below.
3. Creation mechanism of proper tiling by ribbons 3.1. Proper tiling by ribbons of the initial-value domain D 0 In this section, we prove the following theorem using mathematical induction with respect to N.
For any N, the ribbons are tiled properly.
The case N=0 is special, and the one-time-map scheme in figure 5 is responsible for N1. The future case is considered here; the past case can be proved similarly. In the proof, we clarify how the level-(N+1) tiling is created from the level-N tiling, and this in turn shows how a non-shrinking ribbon can appear in special cases.

N=1
Now, a 1 is determined by the separator curves + C 0 and - C 0 in the one-time map  . As a result, the height distribution over the ribbons is Ribbon , : : Each ribbon extends from top to bottom of D 0 (U 0 =B to -B). Thus, the tiling at N=1 is also proper.
In the same way, z + ( ) X U , a a , , , , can only be calculated via the full  + ; N 1 it cannot be obtained by simply combining the one-time map  with the given level-N height function z ( ) X U , N 0 0 . One may resort to a numerical calculation, but then the analytic understanding of the system is lost. A resolution can be achieved if one calls up the inverse of  , and maps once backward from D 0 onto -D 1 to obtain One then obtains the necessary , and can examine whether the corresponding level-( 1 gives a proper tiling. This is actually a tiling on -D 1 , but, from time-shift symmetry, it is equally a tiling on D 0 . Let us calculate the tiling explicitly: , , 2 where, in the second line, j j 1 is applied for the sum variable, and the third line is obtained by (21) and (20). Now, to prove that (23) retains the properness of the tiling, we clarify how  -1 acts on the ribbons on D 0 .
The level-N tiling is made up of + 2 N 1 ribbons. Let us label them by k from left to right, so that the left-and right-half of D 0 are tiled, respectively, as By time-reversal symmetry (appendix C), the tilings of + D 0 and -D 0 are symmetric with each other under  -- , . The height z N has a common value for any ( ) X U , 0 0 within a ribbon. Hence, it can be regarded as a class function if we consider a ribbon as an equivalent class of initial points with the same step-height: we write the and, allowing the extended use of a function on a set to a class function, In figure 7, the scheme of  -1 is laid over the level-N tiling of D 0 . To account for the backward time shift ).

DSS creation mechanism
We consider the ribbons , and by the separation property of . Every half-ribbon is elongated to a full height ribbon (from = -U B 0 to +B), and thus the number of ribbons is doubled.
(ii) Because  -1 is orientation-preserving, the order of the ribbons is maintained region by region. Therefore, the properness of the tiling of + ¢ ( ) is inherited by + ( ).
(iii) Because  -1 is area-preserving, the following holds (denoting the area of T k as ( ) S T k ): all have equal height. Therefore, (26) implies that the created ribbons are generally finer than their parent. This is case A from the introduction. However, there is a remarkable exception: This is case B, leading to the Broucke-type stable PO in a non-shrinking ribbon. This point is further examined in section 4.2.

3.3.
Step-height distribution and end of the proof ; number of ribbons are thus doubled. Meant to illustrate a general case, but data with g = 0.2, N=3 is used.
The first term is simply a 2 0 , because for all ribbons inside ( ) a a , 0 1 , the sign of - In the second term, and via (25), Therefore, (28) gives The distribution of the range of heights is The ribbon tiling inside each region is already proper because of item (ii) of the creation mechanism, and now (31) shows that the range distribution over regions is also proper. Therefore, the level-( + N 1) tiling is proper. + Note that dropping the last term -( ) 1 4, 1 4 in (30) reduces to the basic height distribution of the N=1 tiling in (19).

Stable and unstable periodic orbits in AKP
4.1. Use of ribbons to locate the initial point of a periodic orbit In section 2.2, we showed that, for a given periodic binary code in (7), the initial point * * ( ) X U , 0 0 of the corresponding PO should be enclosed in the level-N ribbon, which has a step height of the first + N 1 bits of the code, and this holds for any N. The same is true for the past case. In general, the future and past tilings are transverse to one another (see appendix C). Therefore, the initial point * * ( ) X U , 0 0 should be inside the junction of responsible future and past ribbons. Furthermore, in section 3, we proved that the tiling is properly ordered by the height. Therefore, it is possible to locate the appropriate ribbon by its height, and we can locate the PO initialpoint * * ( ) X U , 0 0 precisely inside a junction of the corresponding future and past ribbons. Now, however, one must carefully examine whether, at large N, all ribbons shrink to vanishing width (case A), or whether some ribbons escape from the shrinking (case B)-see (iii) in section 3.2. Case A appeals to our naive intuition-some part of a ribbon is definitely chopped away by the separator, and it would then lose some area. As the elongation keeps its height, it would become thinner. At large N, the junction enclosing * * ( ) X U , 0 0 would converge to a point, and the initial point would be identified 12 . Figure 8 corresponds to this case. However, even though case B looks pathological, it really occurs and is the origin of the violation of uniqueness in AKP. Here, each ribbon is indeed chopped by the separator and the resultant two parts are elongated back to their full length. Strangely, one part can have the full area and the other have an area of zero. This can occur when the latter has already shrunk to a line segment! The separator can then only cut out a measure of zero area, and the survivor remains with non-vanishing area. This is very subtle, but the consequence is crucial, as it gives a possibility for a stable PO in AKP. It is amazing that such an exception leads to a physically basic phenomenon.

Advent of stable periodic orbits 4.2.1. Stable periodic orbits and non-shrinking ribbon
We now investigate how the non-shrinking ribbon occurs, taking the case of Broucke's PO as a good example. Broucke reported two stable POs; we consider the shorter one, which has rank n=3 and code +--+--( ) . We call this 'Broucke's stable PO3-6', as its unstable partner is PO3-6 in Gutzwiller's table [4]. Figure 9 shows how the future ribbons evolve as N increases.
The top row shows the tilings T T T , , 6 12 42 . In each, three particular ribbons are shown, namely, the ribbon that encloses the initial point * * ( ) X U , 0 0 of Broucke's PO (we call this Broucke's ribbon and denote it by 'B') and two neighboring ones. Remarkably, only Broucke's ribbon survives, with the other two diminishing rapidly with width~1 2 N . In the first column   T T T 6 7 8 , we clearly observe the process whereby (i) each ribbon is chopped into two parts by the separator curve, (ii) one part is sent to the left and the other to the right half of the rectangle, and (iii) each part is elongated between the upper and lower boundaries. As a result, the full-height ribbons are duplicated, each thinner than its parent by a factor of roughly a half. Now, in the second column starting T 12 , the separator curve is only chopping away a very fine ribbon near the bottom! Further, in the third column, the separator is now inactive, chopping out only a line segment. The non-shrinking ribbon is protecting itself from shrinkage by changing its tail to a line at an early stage. The same occurs for the past Broucke's ribbon. Therefore, the overlap of Broucke's future and past ribbons constitutes a finite-size domain around the initial point * * ( ) X U , 0 0 , and any orbit starting from a point inside the asymptotic junction evolves to produce the same code as Broucke's PO forever in both the future and the past. Broucke's PO in the center is a stable PO in itself.

The bifurcation process  + ¢
U S U ; threshold behavior Let us now examine a decrease in anisotropy by increasing γ and investigate how the Broucke stable PO3-6 appears. As shown in figure 10, there is a threshold g th at which the unstable PO (U) changes into the stable one (S) following the advent of non-shrinking ribbons enclosing S. At the same time, a new unstable PO ( ¢ U ) is born. Thus, the stable PO emerges in a bifurcation process of the form , .
The bifurcation proceeds in the following way. implies p x =0, and the orbit crosses the x-axis perpendicularly at, say, . Thus, the orbit is symmetric under Y-transformation. Now, an infinitesimal shift from * ( ) X , 0 0 leads to a slip-off from the crossing point of the F and P ribbons; this disables the orbit from repeating the code of PO. This is the way in which the PO (U) is unstable in terms of ribbons below the threshold.
(2) At the threshold, the F and P ribbons become tangential to each other at = U 0 0 . (3) Above the threshold, the ribbons stop shrinking and start extending an overlap around = U 0 0 (see figure 11).
Inside the overlap, the orbit repeats the sequence in both the future and the past. The previously unstable orbit (U) becomes stable (S) and remains at = U 0 0 . Therefore, S is also Y-symmetric. The initial point of the newborn PO ( ¢ U ) is located at the corner of the overlap, so that any slight shift (as was the case for U) leads to a slip-off from the overlap. This produces instability while retaining periodicity. It is remarkable that it is also Y-symmetric even though U 0 no longer vanishes. This is not a contradiction, as = U 0 0 is a sufficient condition for the PO to be Y-symmetric, but is not a necessary condition. Indeed, as insets in figure 11 show, ¢ U belongs to a different symmetry class (non-self-retracing) from that of U and S (both self-retracing). This point is further examined below. We add that there is no other PO of the same code within the overlap. For details, see appendix E.

Lyapunov exponents
Let us follow the transition process through g = 0.1-0.8. The change in the maximum Lyapunov exponent of Broucke's PO is shown in the bifurcation diagram in figure 12.
We observe the following.
(1) Above g th , both the stable and unstable orbits of the same code co-exist. As the associated orbit profiles clearly show, the PO in the stable branch is self-retracing (as for the unstable PO below g th ), whereas that in the unstable branches is non-self-retracing. That is, the bifurcation proceeds according to the process where R and NR denote self-retracing and non-self-retracing 14 . See figure 11 for the location of their initial values in the rectangle.
(2) Below and above the threshold, l g g µ -| | max th 1 2 provides a good description. (4) It is interesting to note that the period of the POs is insensitive to the transition. The period of U(R) as a function of γ below g th smoothly continues to that of S(R) above g th . This may be natural, as both POs are  self-retracing and orbit-shape smoothly varies through g th , but the period of the non-self-retracing orbit ¢( ) U NR is also degenerate to that of S(R) in a very good approximation.

Orbit symmetry consideration 4.3.1. Three classes of Y-symmetric orbits
As all POs involved in the bifurcation process let us now focus on the Y-symmetric POs. Using simple topology and symmetry considerations, we can obtain an overview of this bifurcation process. Firstly, we prepare two keys.
(1) A retracing (R) PO should be distinguished from a non-retracing (NR) PO. A NR-PO is simply a closed curve and is homotopic to S 1 , but in a R-PO, the particle is going back and forth on the same curve connecting two turning points. To account for this specific feature, let us say that the R-PO is homotopic to a squashed S 1 (see figure 13).
(2)n , the number of perpendicular crossings of a Y-symmetric PO with the heavy x-axis, must be either 2 or 0. This is because an orbit with an odd n ⊥ cannot be closed, while =  n 4, 6,  can close but only in disconnected loops.
With the above preparation, we can prove a remarkable fact.
Any Y-symmetric periodic orbit belongs to one of the following three classes:  The key to the proof is to consider how to realize, under a Y-symmetric orbit, an appropriaten for the topology (retracing or non-retracing) (see figures 14(a)-(c)).
(i) For a retracing PO to be Y-symmetric, at least one perpendicular crossing of the x-axis must be included.
However, being a squashed S 1 , just a single perpendicular crossing already saturates = n 2. This crossing has multiplicity 2 by itself. Crossings other than this are X-type junctions, each with multiplicity 4. This is class (a).
(ii) A non-retracing PO can be Y-symmetric even without a perpendicular crossing: either = n 0 (b) or = n 2 (c). In (b), all the crossings are X-type junctions, each with multiplicity 2. In (c), all the crossings are X-type, but for two distinct perpendicular crossings each with multiplicity 1. (i) S must be class (a) (meaning self-retracing):Under the decrease of anisotropy, the orbit is stabilized and S is born. As shown in figure 10, the orbit is born at the threshold g g = th , when the future and past ribbons become tangential to each other; hence, the initial point * * ( ) X U , 0 0 is located in the middle of the maximum ribbon overlap at * = U 0 0 = ( ) p 0 x . This indicates a perpendicular emission from the heavy x-axis. Hence, S must be Y-symmetric. Then, the classification requires that either = n 2 or 0. As at least one perpendicular crossing exists,n must be 2 and the class of S is either (a) or (c). In (a), the crossing is self-retracing and has multiplicity 2, whereas in (c), two separated perpendicular crossings, each with multiplicity 1, are necessary. However, as we see in figure 10, at the threshold, there is no room for the separation and (c) is excluded. Therefore, S should be in class (a), and must be self-retracing.
(ii) Y-symmetry:The transition  + ¢ U S U is induced by the configuration change of the F and P ribbons at the threshold, and, as these ribbons carry the same code before and after the transition, the PO transition must be code-preserving. As shown by Gutzwiller [4], the code of the PO dictates the symmetry of the PO, . so that U and ¢ U must also be Y-symmetric once S is understood to be Y-symmetric. Therefore, we can apply the classification of Y-symmetric POs to all of them.
(iii) U is also class (a) (self-retracing):The Y-symmetry of U was demonstrated in (ii), but can also be deduced from figure 10. The initial point is located at the crossing of the F and P ribbons at * = U 0 0 , and, for the same reason as (i) (now the crossing is a point), U is Y-symmetric and in class (a). In terms of ribbons, the  ( ) ( ) U R S R transition corresponds to the configuration change from a single-point to a finite-width overlap.
(iv) ¢ U is class (b) and non-self-retracing:As U is now understood to be in class (a), let us consider all possible routes via which a PO in class (a) (self-retracing) may change its feature through the transition to the final class [(a), (b), (c)] (see figure 15), and investigate which route is responsible for the birth of ¢ U from U(R).
The PO remains self-retracing, and only the stability is changed at g th . This applies  ( ) ( ) U R S R in (iii) and corresponds to single-point  finite-width overlap. The initial value remains * = U 0 0 and only a small change in * X 0 occurs. This is not the route for  ¢ ( ) U R U , since ¢ U is born with rapidly increasing p x after the g th . Therefore, , the selfretracing property is broken, and the = n 2 becomes = n 0, that is, all crossings are of X-type. Especially, at the crossing in concern, a non-vanishing U 0 is quickly generated after g th . This agrees with the observed threshold behavior of ¢ U ( * g g µ -( ) U 0 th 1 2 ): unstable with an initial value at the edge of the overlap of the F and P ribbons. In contrast, if the route were via [III], the retracing property would be broken in such a way that a rapid DX 0 would be created without DU 0 . This is totally against the threshold behavior of ribbons. Thus, the responsible route must be [II] and ¢ U must be ¢( ) U NR (class (b)).
The above is a post-diction on the Broucke transition, and is summarized in figure 14(d). Two remarks are in order.
(1) For a PO in class (a), the total number of crossings of the x-axis (i.e., the length) is + n 4 2, where n 4 comes from the n cross-junctions with multiplicity 4 and the 2 comes from a single perpendicular crossing with multiplicity 2. The rank of a PO is half its length; thus, the rank of class (a) POs must be odd ( + n 2 1). For this reason, we conjecture that the Broucke-type transition, associated with the non-shrinking ribbon, will occur in Y-symmetric odd-rank POs 15 .

The two-dimensional AKP PO search
To study the POs in this system, especially with regard to the uniqueness issue, we must first search the POs exhaustively. By exhaustively, we mean avoiding any limitation on the PO, irrespective of its symmetry and whether it is unstable or stable. A recent search by Contopoulos et al [19] was of this type, albeit limited to perpendicular emissions from the heavy axis ( ) also fixes the momentum p y 0, via energy conservation; hence, a one-parameter shooting varying only x 0 is sufficient. Their analysis yielded important information, especially on the existence of stable POs in AKP, including relatively high rank POs. However, the limitation is rather severe: all the POs are  ( ) ( ) a b ( =  n 2 0). Non-vanishing U 0 is generated quickly after g th .
[III]  ( ) ( ) a c No change ofn , but multiplicity at crossing changes as  + 2 1 1. X 0 changes rapidly after g th . 15 Our preliminary result, from numerical calculation of z ( ) for slippery steps, is that, from n=3 up to n=23 and g < 8 9, there is only one advent of non-shrinking ribbon at every odd rank (T Shimada and N Nakajima). We are trying to consolidate this issue. obviously limited to be Y-symmetric, and, furthermore, as there is at least one perpendicular crossing (the initial point of the search), they are with = n 2 [class (a) and (c)]. Thus, this approach cannot detect the class (b) POs (NR) that are created as ¢( ) U NR after the bifurcation To ensure an exhaustive search, we proceed as follows.
(1) The basic flow of the analysis is based on specifying the code of a PO, a PO in (7), at a certain anisotropy.
For instance, if the code +--+--( ) is specified, the routine searches out Broucke's PO3-6; the result is U(R) only if g g < th , but is both S(R) and ¢( ) U NR if g g > th . One can be sure that there are no more POs of this code within the setup resolution. Figure 10 shows the outcome of the code request +--+--( ) .
(2) Given the requested code of the PO, we now take advantage of the PO trapping mechanism. That is, the properness of the tiling of D 0 by ribbons guarantees that * * ( ) X U , 0 0 should be inside the base ribbon of the step whose height z ( ) a N F is calculated by (9). We use N=48 for maximum precision in the doubleprecision calculation.
(3) The above asymptotic ribbon with N=48 may, if subject to case A, be extremely narrow (D~( ) X O 1 2 48 ) and the neighboring steps may be almost degenerate in height. However, for case B, it may retain some finite width. Here, we proceed protectively: rather than directly tackling the asymptotic ribbon, we focus on the leveln 2 (48) ribbon, which should embody the asymptotic one by the properness of the tiling 16  can be calculated by a bi-section method to satisfy Now, the search area is reduced from the vast D 0 to a leveln 2 ribbon extending from = -U B to U=B.
(3a) In this step, we check whether the target ribbon shrinks or not by inspecting its full profile ( Î - ).
(4) We now shoot for a PO inside the leveln 2 ribbon-find a point for which the error between the initial and final n 2 -th crossing of the x-axis vanishes 17 . Practically, we have searched for * X 0 , which minimizes c 2 , at every U 0 . This gives a function *( ) , through which we can regard c 2 as a function of U 0 . That is, Thus, the two-parameter search is effectively reduced to a one-parameter search. With steps (3a) and (4a), the search is organized so as not to miss the violation of uniqueness. In the following, we report two applications: a detailed case study of a high-rank PO (g = -0.85 0.93) and an exhaustive PO search regarding the uniqueness issue at high anisotropy (g = 0.2).

Case study of bifurcations of a high-rank PO (n = 15).
This investigation is motivated by the longest PO (length 30) among the stable POs reported by Contopoulos et al [19]. However, the code of their PO is not given, and there is a possibility of misidentification considering the exponentially large number of POs at such high rank. Thus, the description below may be better read in its 16 This corresponds to the relaxed ribbon with boundaries  L and  R in figure E1. 17  x y , with slowdown around the origin when a close encounter with the origin is involved. 18 In this case, a minimum search using a tri-section method on c ( ) U 2 0 is vital. See our earlier report [22].
own right, irrespective of the motivation 19 . The code of our PO15 is We find that it bifurcates twice:  + ¢ U S U and then  ¢ +  S S S , as seen in figure 16. The bifurcation pattern is summarized in figure 17 using the PO class ( figure 14) and bifurcation route ( figure 15).
As the orbit at this high rank is very complicated, let us examine the bifurcation as a challenge to the general theory considered in section 4.
(1) The first bifurcation  + ¢ U S U : All U, S, ¢ U are Y-symmetric and our classification can be applied. Figure 16 shows that the initial value of U, S is * º U 0 0 , whereas the non-vanishing * U 0 of ¢ U is rapidly created after the bifurcation. This leads us to infer that  U S is via route . This means the bifurcation is  + ¢ ( ) ( ) ( ) U R S R U NR , the same as Broucke's PO3-6 and PO5-40. The ribbon structure in figure 18 consistently shows that the bifurcation occurs simultaneously with the respectively. Curves are from a simple fit l q g g g g q g g g g  . 19 We read x 0 by eye from their figures 4 and 6, and confirmed that the reproduced orbit closes, below the bifurcation threshold, at the stated 30th crossing after a slight adjustment. Past the threshold, however, it gradually fails to close even adjusting x 0 . Their analysis is limited to = p 0 x,0 , and there is a possibility that the observed PO is the same as ours ( . We here quote our initial values ( ) X U , . The orbit profiles are given in figure 19. As the orbits are of length 30, they look at first glance like a cloud, but from the density it appears that ¢ U is non-self-retracing and S is self-retracing. One can also pinpoint the turning-point at the kinematical boundary in the latter. A closer inspection is more intriguing: S(R) has a single perpendicular crossing with multiplicity 2, whereas ¢( ) U NR has no perpendicular crossing. This is exactly the multiplicity prediction from the classification (33) and figure 15.
(2) The second bifurcation  ¢ +  S S S :This is a new case, which occurs in the very-low-anisotropy regime (g > 8 9 2 ), where the flow (A.3) no longer has hyperbolic singularities and the existence of unstable POs for an arbitrary code is no longer guaranteed [16]. The rapid initial value variation is in the X direction. Thus, we infer that  S is produced via route [III], which means it is in class (c) and the bifurcation should be  ¢ +  ( ) ( ) ( ) S R S R S NR . Indeed, the initial positions of ¢ S and  S in figure 18 are horizontally alignedhorizontal because the variation is in the X 0 direction, and two initial points for  S because  S is in class (c) ( = n 2). (See the two multiplicity-1 crossings in [III] in figure 15). This can be directly verified in the orbit  profile in figure 20(c). Let us note a further success of the theory. It implies that ¢ U must be NR, but this can hardly be seen from the orbit, which is almost doubled and looks as though it is self-retracing. However, the magnification in figure 20(a) verifies that it is indeed NR. In summary, the first bifurcation of PO15 follows the theoretical classification precisely, and the theory also explains the second bifurcation, which embodies a class (c) PO.

5.3.
Search of all distinct POs (rank n10) at high anisotropy g = 0.2 This section describes our exhaustive search for all POs up to rank n=10 (length = n 2 20) at high anisotropy (g = 0.2). As noted in the introduction, it was previously considered that all POs were unstable and isolated at such high anisotropy, mainly based on the early numerical analysis [4,7]. Recently, however, the possibility that the existence of ample stable orbits in AKP may indicate stable orbits even at high anisotropy has been expressed [19]. Therefore, we revisit high-anisotropy AKP armed with our twoparameter shooting algorithm, which embodies steps (3a) and (4a) above and is therefore sensitive to the possible violation of uniqueness.

Distinct periodic orbits and distinct binary codes
To challenge this problem, let us first consider the counting of distinct POs. The Hamiltonian (1) is symmetric under the discrete symmetry transformation  - , and  t t : . Thus, any partner orbit, generated from one PO by the symmetry transformations, is again a solution of the equation of motion and a respectable PO.
However, they have the same shape, period, and stability exponent, and it is legitimate to place them into equivalent classes. Two POs are distinct only when they belong to different classes. If a PO itself has none of the symmetries, 2 3 orbits belong to the same class and the class has degeneracy s = 2 sym 3 . In contrast, if the PO is self-symmetric under some transformation, then it does not generate a different partner, and the degeneracy is halved for each self-symmetry 20 . There are ten types of PO symmetry, five for non-self-retracing and five for selfretracing POs. The degeneracy factor for each symmetry-type is listed on the fourth row of table 2. The O-type Table 2. The number of distinct POs predicted under uniqueness assumption in each of ten symmetry classes indexed by k. The fourth row lists the degeneracy s ( ) k sym .
Non-self retracing: NR Self-retracing: R    , that is, the total number of distinct codes at rank n. Under uniqueness assumption, this predicts the total number of distinct POs (see table 2). 20 If a PO is non-self-retracing,  produces degenerate pairs orbiting in opposite direction to each other, but if self-retracing,  is immaterial-it amounts to the freedom to choose one among two turning points as the starting point. Therefore, s ( ) k sym for R are respectively half of those for NR. among NR is newly found at ranks 7 and 9, and is discussed below. The number of distinct POs in a given symmetry class listed in the following rows is predicted (under the uniqueness assumption) by counting distinct codes. The prescription was derived by Gutzwiller in detail, and the code table up to rank n=5 is given in table 1 of [4]. Our table 2 is the extension of this table to n=10. (Explicit codes take pages, and only the number of sequences is tabulated). Let us briefly recapitulate the prescription. First, a rank-n PO is represented by a lengthn 2 binary code; if the orbit traverses the heavy axis upwards n times, then it must also traverse this axis downwards n times to return to the initial point. Next, the symmetry of a PO corresponds to the rule of its code. For instance, if a PO is X-symmetric, its code must satisfy a rule = ---X a a : n i i 2 1 [4]. Thus, as all POs are divided into classes by their equivalence under symmetry transformations, to scrutinize distinct codes, the 2 n 2 binary sequences at rank n should be divided into classes using symmetry rules. However, to reach distinct codes, it is not sufficient to divide according to equivalence under X, Y,  rules, but one must also divide by the equivalence of sequences under the code-shift operation t, that is, the freedom to choose the starting bit among the cyclic binary code. At this final step, there is a slight subtlety.
If the rank n is prime, the code-shift symmetry simply amounts to n-fold degeneracy. (It is n rather than n 2 , because, by definition, the starting bit must be chosen from crossings with > p 0 y ). However, note that the set of rank n POs includes m repetitions of a lower rank primary orbit of rank p, where p is a prime divisor p=n/m. In such a case, the code-shift degeneracy reduces to p. Dividing out the code-shift degeneracy to take account of this, one eventually obtains the distinct sequences. Table 3 describes how to organize this counting.
The uniqueness of a PO for a given code means there is one and only one distinct PO for each distinct code. Under this assumption, the bottom row of table 3 gives the number of distinct POs predicted in table 2.

5.3.2.
The results: uniqueness holds at g = 0. 2 We now describe how the crucial steps (3a), (4a) and (5) in the algorithm work at this high anisotropy.
(3a) The ribbon of the code always shrinks at g = 0.2 as~1 2 n 2 . Specifically, for n=10, the maximum width - does not exceed 10 −7 for any of the 13648 tested codes (see table 2).
(4a) The chi-squared (error) curve is always convex with a single minimum with negligible value. It is mostly under 10 −20 (the best value is~-10 25 ) and the worst is~-10 2 for only a few POs. Subsequent stability tests have proved that all POs are unstable. Therefore, we conclude that, within the above resolution, there are only unstable POs, and that the PO is unique for any given code up to n=10 at g = 0.2.

A new type-O symmetry and other PO samples
Let us comment on the new symmetry class (O-type), figure 21, which was not considered by Gutzwiller [4]. We were first surprised when our search result did not satisfy the sum rule, which states that the number of the PO should be 2 n 2 when division by the symmetry equivalence is removed. This may have indicated a violation of uniqueness, even though we were working in the highly chaotic region g = 0.2. Eventually, we were able to locate the reason as being due to this O-type symmetry; there is only one at rank n=7 and seven at rank 9. After counting these cases correctly, the sum rule was indeed satisfied, and the uniqueness held at g = 0.2 (see table 2). The O-type orbit is neither -X nor -Y symmetric, but is symmetric under the transformation  - -O x x y y : , . The observed orbits are all non-self-retracing and s = 2 sym 2 . The symmetry class number 5 is assigned, and 10 is reserved for the possible occurrence of O-type retracing POs at  n 11. Further examples (figure 22) are taken from the 13648 distinct unstable POs at rank 10. The first is asymmetric under both X and Y, and NR; k=1 in the notation of table 2. The second is X-symmetric but Y-asymmetric, and NR; k=2.
The final one is Y-symmetric and our classification scheme can be used to analyze it. It is apparently NR; therefore it can be either class (b) if = n 0 or (c) if = n 2. It looks that it has perpendicular crossing(s) but hard to see how many. The class (c) PO should have two perpendicular crossings, each of multiplicity 1. Close numerical investigation verifies that this is precisely the case.

Verification of Gutzwiller's action formula
In Gutzwiller's periodic orbit theory, the semi-classical description of the quantum density of state is given by This formula has only two parameters, τ and l A , but for n=5, it can provide Φ for any PO in remarkable agreement with the measured value. 22 We have tested this formula using the action data from rank n=10 POs obtained by our high-accuracy measurement. In one of the tests, we first determined the two parameters by fitting to the action data from 44 distinct POs in the rank n=5 family (g = 0. . We now test whether this formula, with the parameter values in (39), can endure the extrapolation to the n=10 family. In figure 23, we show the result as a scatter plot.
At n=10, there are 13648 distinct POs. Each point in the plot represents one of them. It can be seen that the prediction works for the wide range of action constant values extending up to approximately 60. We show in table 4 the MSD, defined as and, for our canonical choice = -E 1 2, T=S holds. For the Lyapunov exponent, the 'crudest approximation' l l = n 0 with l~1.5 0 was adopted in [3] and we also use it. 22 For the family of rank n POs, F max , F min , and áFñ predicted by (38) are easily calculated as F = 0 min , and (at large n) It is clear that the marvelous success of (38) continues up to n=10. Let us briefly discuss the implications of this success. Following Gutzwiller, let us integrate the density of state up to E [3]. Then, we obtain where a is the code of the PO, and a dimensionless variable s is introduced by By this 'Wick rotation,' the quantum formula can be mapped to a statistical formula-the grand partition function of a statistical spin system 23 . Under this transformation, the action formula (38) corresponds to an Ising spin chain, which has only a two-body interaction with exponential decay [23]. Now that the validity of (38) has been confirmed up to n=10, it may be stated that AKP semi-classical quantum theory at the highanisotropy region is dual to the above spin theory.

Conclusion
In this paper, we have fully utilized the ability of symbolic coding in AKP and considered the level-N DSS and the tiling of the initial-value domain by ribbons introduced by this surface. We have proved the properness of the tiling by ribbons from the mechanism of creating the level-( + N 1) tiling from the level-N case, which clarifies how non-shrinking ribbons can emerge.
Our key points are as follows: (i) Non-shrinking ribbons emerge when the future and past asymptotic curves become tangential to each other at = U 0 0 at the threshold anisotropy value.
(ii) The bifurcation is of the form  + ¢ ( ) ( ) ( ) U R S R U NR . The initial point * * ( ) X U , 0 0 of the stable PO S is located in the maximum overlap of the future and past ribbons at * = U 0 0 , whereas the unstable PO ¢ U is located at the edge of the overlap.
(iii) From topology and symmetry considerations, we have explained the above bifurcation scheme, and conjectured that the stable PO occurs in the Y-symmetric, odd-rank self-retracing orbit. A case study of high-rank PO15 supports this statement.
(iv) An exhaustive PO search has verified that the uniqueness conjecture holds at high anisotropy (g = 0.2), if the newly found type-O orbit is accounted for. Gutzwiller's action approximation formula works amazingly well for all POs up to rank 10 at this anisotropy level.
The topology and symmetry approaches developed in this paper provide a framework for constraining the bifurcation, but we must admit that a direct analytical understanding of the advent of non-shrinking ribbons is urgently required. On a related note, it is tempting to look at this bifurcation from the quantum side-separating S and ¢ U from quantum data, using inverse chaology as a quantum prism. Previously, we could successfully extract low-rank unstable PO data (including Lyapunov exponents) from the AKP quantum spectrum [12]. Tackling this with a bifurcation is a serious challenge, requiring higher-order corrections in ÿ as well as taking satellite POs into account [24][25][26]. Work in this direction is in progress. to remove the singularity. Then, the first two equations take the autonomous form of (14) and the χ-equation decouples (in fact, χ remains ¥ for any finite ¢ t ), and (14) describes the evolution of ϑ and ψ in the limit c  ¥ [2]. The one-time map  in (11) describes how D 0 is mapped onto D 1 . Because the limit c  ¥ restricts D to the collision manifold I With c  ¥, we observe U e 2 cos cos , sign cos 2 1 cos . 2 Therefore, except for the critical case q p = 2, we can use the following relation for the necessary conversion: where for X we have used the fact that, on the Poincaré section y=0, either y = 0 (the case > X 0) or y p = (the case < X 0). This is (17) in the text. To extract  | I from , it is best to follow [2]. Fix y = 0 0 (that is, choose > X 0 0 ), and let J 0 increase from 0 to π. (This is sufficient because of the symmetry of the system). This makes ( ) X U , 0 0 circulates around the boundary of the half-Gutzwiller rectangle (see figure 4). Then, follow the streamline in figure A1 until it reaches the next PSS, that is, until it crosses either y = 0 1 or y p = 1 . In this way, one can read off the final (J y , 1 1 ) for  and, via (A.4), one gets Î ( ) X U I , 1 1 for  | I . Let us now study the flow in figure A1. There are two types of singularities for g < 8 9 : and f p y as the consequence of the time-reversal symmetry of the Hamiltonian (1).
The flow divides the initial domain into three regions, depending on whether the next crossing of they axis occurs with > X 0 ). The division is determined by two critical angles J v and J h [2]. The three regions are as follows.  Figure A1. The flow given by the autonomous ODE (14). This is the same with figure 1 in Gutzwiller [2] and equivalent to figure 2 in Devaney [16]. The anisotropy is g = 0.   Figure B1 shows how sets of vertical and horizontal line segments (sets b and c) meeting at the separator  ++ 0 are mapped by  . We observe that every set turns into a 'beak' with its tip at ++ v 1 , while the full vertical line (a) that does not touch  ++ 0 is simply distorted as line (a'). Thus, we find a rule of contraction By moving a vertical line segment (and its image in ++ ¢ ( ) ) until it reaches ++ v 0 , we find a rule of blow-up: Equations (B.1) and (B.2) constitute a time-reversal pair.
Let us now proceed to +-( ) C 0 . This is adjacent to +-( ) and  | I rotates the boundary 2, as seen in figure 4. Figure B2 clarifies how the rotation affects the interior map.
(i)  maps each set of vertical and horizontal line segments (  p q , , ) meeting at +-C 0 into a beak with its tip at a point +-v 1 . Therefore, This just corresponds to the contraction rule (B.1). The body of the beak reflects the rotation as follows.
(ii) The left-side of the beak (dashed line) always connects +-v 1 with another critical point + c 1 . It is also the image of the line segment horizontally connecting  +-0 to  + 0 (the > X 0 side of the center line of the collision manifold I in D 0 ). Therefore, we find a contraction rule (iii) The right side of the beak, however, is the image of the vertical line segment connecting  +-0 to either B 2 or C 2 . In the former, it is simply a long vertical curve. In the latter, it is a short hook, connecting +-v 1 and a point in ¢ C 2 , now in the upper boundary following the rotation.  The vertical line segment a is critical. Body of a makes a curve ¢ a , its top end +-v 0 makes +-C 1 , and they together make a beak (enveloping wings), whose top +-v 1 is then the image of the +-C 0 . for  is equivalent to (a 0 , a 1 interchanged) -+ ¢ = ( )B for  -1 . Therefore, it follows from time-reversal symmetry that . The direction of increasing DSS height is depicted by arrows, which are  - where P acts on D 0 as well as D 1 as  - -P X X U U : , . Finally, (C.2) and (C.3) together describe that the eight sub-regions ¼ A D , , and 1Ā D , , are grouped into two classes: The direction of increasing height of the future N=1 surface is shown by an arrow in figure C1. This is nothing but the N=1 fact in the proof of the properness of ribbon tiling by mathematical induction. Now, by the above symmetry relation, the arrow for the past N=1 surface is the  -U U mirror of the future arrow. Therefore, the future ribbons and past ribbons are transverse to each other. This property is inherited by higher-level ribbons, as the level  + N N 1 evolution maintains the properness of ribbon tiling. The sole exception is where future and past ribbons become tangential to each other, as seen in section 3.

Appendix D. Transverse chopping and longitudinal splitting
If one is content with just comparing the location of the new and previous ribbons, this is simple: each of the previous ribbons longitudinally splits into two finer ones. The height z ( ) X U , However, as we proved, the tiling at + N 1 is proper. Therefore, the level-N ribbon longitudinally splits into two finer ribbons in such a way that the left of the split-line gets -D + N 1 and the right gets D + N 1 . See figure D1. Note that, in the case of non-shrinking ribbons, the longitudinal split-line occurs on the boundary of a level-N ribbon.
The level-( + N 1) ribbons are created by chopping the N ribbons transversely into two parts by the (fixed) separator curve and mapping each into different halves of a rectangle with elongation from -B to B, as shown in section 3. Hence, picking some ribbon at large N and tracing back to find which part of the initial domain it comes from is a tremendously difficult task (see figure D2).

Appendix E. Uniqueness of PO within the overlap
As for the unstable PO3-6, the initial position is remarkably close to the edge of the overlap. This is just as it should be-the unstable PO must, to be periodic, lie on the union of future and past ribbons of its code, and yet, to be unstable, should not be far inside the junction. This is a subtle point, as the exact corner is a homoclinic point and cannot be a periodic point. The fact is that the unstable PO3-6 turns out to be extremely close to the corner, as shown in figure E1. We thank Tanikawa and Shibayama for pointing out this issue of homocriticality and mentioning that this kind of close proximity often occurs. The unstable PO3-6 is non-self-retracing ('NR'). There is no other PO of the same code on the union of the future and past ribbons, as is also clear in figure E1. We should add that we have observed a stable satellite in Broucke's island. The details are under investigation. Figure D2. From the level N=1 2 2 ribbons, higher level ribbons (2 3 at N = 2, 2 4 at N = 3 are created by repeated application of  -1 and lines to track descendants are shown. It is extremely hard to track-back to the parent already through  -2 (twice the baker-map). Figure E1. c 2 contours showing that no period 6 PO other than PO3-6 (rank 3 id 6) S and ¢ U (figures 11(a) and 12(a)(b).) exists inside Broucke's non-shrinking ribbon. g = 0.6. Union of future and past ribbons has asymptotic boundaries b L and b R (N = 48), which is well inside the N=6 union ( L to  R ). Inside the latter, a fine mesh of initial points (10 10 3 3 ) is set and c 2 ((35), = n 2 6) is calculated for each point. This is a fail-safe procedure as well as a useful devise for magnifying near the narrow edge region. The equicontours focus the S. Inset: The same c 2 contours in pqrs (after a stretch in the X direction). They focus the unstable ¢ U at the edge.