Rapidly converging chaos indicator for studying dynamic aperture in a storage ring with space charge

The determination of dynamic aperture in storage rings and colliders is a numerically intensive procedure. When realistic space-charge forces come into consideration, the numerical load becomes even heavier. Furthermore, dynamic aperture estimation using chaos indicators like frequency map analysis (FMA) raises reliability issues when the dynamical system has a time-dependent perturbation like the space-charge force. In this article, we apply a rapidly converging chaos indicator called reversibility error method (REM) to study the space-charge contribution to the dynamic aperture of the integrable optics test accelerator (IOTA) storage ring at a small value of the space charge tune shift. The strength of REM is addressed through examples, including a particle-core model of halo formation. We also develop a toy model of the IOTA lattice to further reduce the computing time required to estimate the dynamic aperture, and we compare this model with a realistic space-charge simulation for verification.


I. INTRODUCTION
The integrable optics test accelerator (IOTA) is a storage ring at Fermi National Accelerator Laboratory dedicated to testing innovative beam-physics ideas. One of the key operation modes of the IOTA ring is to test the novel concept of nonlinear integrable optics. In this concept, the single-particle transverse optics is highly nonlinear but integrable, so that coherent collective instabilities are expected to be Landau damped [1][2][3][4][5] due to the large frequency spread, while maximizing the dynamic aperture. Since integrability can be broken by various factors including perturbative nonlinearities, fringe fields [6][7][8][9][10], and space-charge forces [11], it follows that the dynamic aperture becomes finite.
We consider the following Hamiltonian system, Hðξ; tÞ ¼ H 0 ðξÞ þ εVðξ; tÞ ð 1Þ where ξ ≡ ðx; p x ; y; p y Þ is the vector of transverse phasespace coordinates, H 0 is an autonomous integrable Hamiltonian, and ε ≪ 1 represents the size of the perturbation V, which can be time-dependent. The independent time variable t often represents longitudinal location or betatron phase advance over the ring [12]. Since an equilibrium beam distribution is difficult to obtain (or may not exist), the presence of space charge typically adds a nonperiodic time-dependent perturbation. Furthermore, the finite number of simulation particles and the numerical discretization of the space-charge model add sources of time-dependent random numerical noise. The brute force determination of the dynamic aperture involves a numerically heavy simulation of many particles sampled in phase-space for a tracking time of interest that can be very long. In practice, chaos indicators like frequency map analysis (FMA) [13] are widely used for dynamic aperture estimation using a relatively short term simulation [12,14,15]. In the core of FMA, the numerical analysis of the fundamental frequency (NAFF) [13] algorithm attempts to obtain the following harmonic decomposition of an orbit zðtÞ by assuming that the orbit is quasiperiodic: where d is the number of degrees-of-freedom of the Hamiltonian system and ν ≡ ðν 1 ; …; ν d Þ is the vector of fundamental tunes that depends on oscillation amplitude [13]. If the measured tunes of an orbit diffuse over time, the attempt to construct a harmonic decomposition fails and FMA indicates the corresponding orbit to be chaotic. However, when the perturbation V is time-dependent, the reliability and convergence of FMA for dynamic aperture estimation is questionable. This is because the measured tunes can diffuse not only on chaotic orbits but also on regular orbits, due to the time dependence of the perturbation. It is unclear how to define the dynamic aperture of the system Eq. (1) when the time dependence is non-periodic. Thus, in the presence of space-charge, it makes sense to define the dynamic aperture only after the beam has relaxed to a matched periodic equilibrium. However, in many cases, convergence to an equilibrium is not guaranteed or may require a very long time (compared to the time of interest or the numerically manageable time). Therefore, we consider dynamic aperture to be the stable phase-space volume from a specified time to the time of interest. If the time dependence is weak (or becomes weak due to very fast beam filamentation), such that the evolution of the dynamic aperture is slow compared to the time of interest, the chaos indicators may be able to predict long term dynamic aperture using data sampled over a shorter time.
In order to apply chaos detection methods to the IOTA ring with space charge, we address these difficulties using two approaches in combination. One approach is to construct a simplified toy-model to represent particle dynamics in the ring using a symplectic map. In this model, the space-charge effect is parametrized by the space-charge induced tune depression. The second approach is to use a rapidly-converging and reliable chaos indicator called reversibility error method (REM) [16][17][18][19] which will be introduced in the following section. Throughout this paper, we often compare REM against FMA, which is one of the most popular tools for dynamic aperture study in the accelerator community. Comparison against other chaos indicators can be found elsewhere [17].
The goals of this paper include (1) introducing REM along with its strengths: rapid convergence, numerical efficiency, and reliability for systems with time-dependent perturbations and (2) studying dynamic aperture of the IOTA ring with weak space-charge near the nominal design choice for the nonlinear magnet parameters, using a simplified toy-model for comparison with a realistic space-charge simulation result. This paper is organized as follows. In Sec. II, we present an overview of REM and illustrating examples. In Sec. III, the IOTA toy-model is defined. In Sec. IV, REM is applied to study the dynamic aperture of the IOTA toy-model, and the results are compared to those obtained using FMA. This is followed by a discussion comparing the two indicators. In Sec. V, REM is applied to a realistic IOTA lattice with space-charge, and the results are compared with FMA. In Sec. VI, we investigate how the dynamic aperture responds to a change of the nonlinear magnet parameter. Finally, a conclusion follows in Sec. VII. Additional details are contained in the three Appendixes.

II. REVERSIBILITY ERROR METHOD
When the particle dynamics is symmetric under timereversal, REM uses numerical integration forward and backward in time to indicate how sensitive the orbit is to a small initial condition perturbation. Consider a numerical integration forward in time. During each integration step, a numerical round-off error associated with the finite digits of precision occurs. If the round-off errors from different integration steps can be modeled by a random walk, the accumulated error in the orbit will be order of ffiffiffi n p after n steps [20,21]. However, there exist correlations of round-off errors between different steps, and the accumulated error grows as a power law in time [17] for regular orbits. For a chaotic orbit, the accumulated error can grow exponentially over time due to the presence of nonvanishing Lyapunov exponents. The accumulated error can be measured by comparing the difference between the initial point of the orbit and the return point obtained after first integrating forward over time t, and then integrating backward over time −t, as illustrated in Fig. 1. Throughout this paper, we use double precision for all the numerical operations.
In this section, we benchmark REM in comparison with FMA using two examples. An additional example illustrating the sensitivity of FMA on resonance is presented in Appendix B.
A. Example: Henon-Heilies potential As a first example, we applied REM to the Henon-Heilies problem [22] whose Hamiltonian is given by: Through this example, we would like to illustrate how the reversibility error grows over time for chaotic and regular orbits and define the value of the chaos indicator using the reversibility error. The black curve is the exact orbit, the red curve is the orbit obtained using numerical integration forward in time, and the blue curve is the orbit obtained using numerical integration backward in time. Figure 2 shows a Poincaré section illustrating two orbits with energies (i.e., Hamiltonian values) E ¼ 1=12 (regular orbit) and E ¼ 1=9 (chaotic orbit). For each orbit, Fig. 3 illustrates the phase space distance between the initial point and the return point after integrating forward and backward over time t, shown as a function of t. For the chaotic orbit, note that this distance grows exponentially with t, exceeding the corresponding result for the regular orbit near t ¼ 20, and it becomes 10 10 times larger at t ¼ 500. This motivates us to define the following measure of chaos using forward-backward integration. For a reference, we also define the measure of the frequency diffusion we will use throughout this paper: Here, Δν x denotes the change in the horizontal tune when it is computed using two successive time intervals of the same length. (Throughout, "sampling period" refers to the time between successive data points in the time series, and "time interval" T refers to the length of the time series. That is, to calculate FMA we track time T forward and another time T forward to evaluate Δ FMA .) A similar definition applies for Δν y . We define ν x and ν y to be the tunes corresponding to largest spectral peak obtained using the horizontal and vertical time series data, respectively. (It is worth mentioning that these are not necessarily equal to the fundamental tunes, as each time series contains many harmonics of the fundamental tunes, as in Eq. (2). It is also worth mentioning that the order of the first largest and the second largest peak can change over the two successive time intervals when two or more spectral peaks have comparable amplitude. Therefore, we carefully choose the peak during the second time interval so that the chosen peak corresponds to the largest peak during the first time interval.) The quantity Δx is the difference between the initial and return coordinate x after forward and backward integration, and Δp x , Δy, and Δp y are similarly defined. When computing Eqs. (4)-(5), we typically work in a canonically transformed system of coordinates in which all phase space variables are dimensionless and of comparable size. Now we apply the chaos indicators with the defined measures to the Henon-Heilies example. A set of initial particles was prepared in the surface defined by x ¼ 0 and Hðx; p x ; y; p y Þ ¼ E. Particles are numerically tracked forward a time T ¼ 2048, then another time T ¼ 2048 to compute Δ FMA (with sampling period t ¼ 1)  other in the overall structure of the chaotic regions, except for some fine bright lines of the FMA result in Fig. 4. Such delicate differences can be partly understood from the fact that FMA is sensitive to the presence of resonant orbits, as illustrated in Appendix B. It appears that all initial conditions with E ¼ 1=12 shown in Fig. 4 correspond to regular orbits, while the set of initial conditions with E ¼ 1=9 shown Fig. 5 contains chaotic bands. Throughout this paper, we often refer to similar figures as dynamic aperture plots.
B. Example: Particle-core model As an example of a time dependent Hamiltonian system, we use the 1-dimensional particle-core model [23][24][25][26][27] in a constant focusing channel. The particle and the envelope (i.e., beam core radius) can be described by the following Hamiltonian [24], HðX where ðx; pÞ and ðX ; PÞ are phase space variables describing the particles and the beam envelope, respectively, and η ¼ ν=ν 0 is the tune depression ratio. (Here ν 0 is the zero current tune and ν is the depressed tune of particles that remain within the beam core.) The normalization is taken such that t ¼ 2π is the betatron period in the limit of zero current, and the matched beam envelope size is X ¼ 1.
When the 1-dimensional envelope Hamiltonian is time dependent there is a rich content of envelope dynamics [26]. Since the envelope Hamiltonian Eq. (7) is autonomous, the envelope oscillation is integrable. In Eq. (6), note that there can be a numerical error associated with the finite integration step size when particles cross the hardedge envelope. Luckily, such an error can be made time reversible under the split-composition method [28,29] of symplectic map construction. Figures 6 and 7 show stroboscopic plots of particle orbits moving in the field of an emittance dominated beam η ¼ 0.8 with X ð0Þ ¼ 0.651 and a space-charge dominated beam η ¼ 0.1 with X ð0Þ ¼ 0.577, respectively. They show particles of an initial array of 32 particle coordinates uniformly distributed along the positive horizontal and vertical axes (as in [24]) with sampling at every minimum of the mismatched envelope. Note that orbits are regular in Fig. 6 while there are chaotic bands in Fig. 7.
In order to illustrate the effect of time dependent potential on the tune diffusion measurement Δ FMA , we sample the orbit data with two different sampling periods. One sampling period is chosen to be equal to the stroboscopic period (equivalently, the envelope oscillation period) as we did for Figs. 6 and 7 so that the envelope motion between each data time appears to be constant to the observer. With such a sampling period, Fig. 8 shows FMA and REM dynamic aperture plot when tune depression ratio and mismatch are same as the Fig. 6 case. Both methods well indicated that the separatrix line to be chaotic. Although FMA reveals much rich structure, the meaning of such complex structure does not seems to be intuitively interpretable. Similarly, Fig. 9 shows dynamic aperture plot when tune depression ratio and mismatch are same as the Fig 7. Both methods well indicated the chaotic bands as well as resonance islands that appeared in Fig 7. Another sampling period is chosen to be slight larger than stroboscopic period (by 0.1%) so that slow envelope motion can be observed as shown in Fig 10 while REM result looks consistent. Under this sampling period, FMA indicates the regions bounded by separatices or chaotic band to be chaotic. It can be understood from the fact that there is spurious tune diffusion due to envelope oscillation.

III. IOTA TOY-MODEL
The single particle dynamics of the IOTA lattice can be described by the following on-momentum autonomous Hamiltonian [1,30,31] whose canonical variables are the (dimensionless) normal coordinates x n , y n , p x;n , and p y;n and the independent time variable is the phase advance ψ: Here, τ is the (dimensionless) strength of the nonlinear potential and, with z ≡ x n þ iy n , is the nonlinear integrable potential, which can be achieved by a specially designed nonlinear magnet [32]. Note that the Hamiltonian is dimensionless. Then a second invariant of motion functionally independent of H is given by: Note that the two invariants are symmetric under the reflections ðx n ; p x;n Þ → ð−x n ; −p x;n Þ and ðy n ; p y;n Þ → ð−y n ; −p y;n Þ. Note also that the potential U has singular points at ðx n ; y n Þ ¼ ðAE1; 0Þ.
Since the Hamiltonian is autonomous, the Lie map can be written as where 2πν is the phase-advance over the nonlinear potential element, and the colons are used to represent the Poisson bracket operator following Alex Dragt's notation [33], so that ∶f∶g ≡ ½f; g. Due to the reflection symmetry of the two invariants, the integrability will be preserved when the map N is concatenated by a linear rotation map R with a half-integer tune advance in each plane: Here 2μ x ; 2μ y ∈ Z and The concatenated map M ¼ N R approximates the oneturn map of the actual IOTA ring, which consists of the nonlinear magnet insert section, represented by N , and the linear arc section, represented by R, as shown in Fig. 13. We introduce a crude approximation of weak spacecharge effects into this model by introducing small perturbations δμ x and δμ y to the horizontal and vertical tune advance over the arc section, denoted here by μ x and μ y , respectively. The map of the nonlinear magnet section is not affected. This is justified by the fact that the length of the linear arc section is an order of magnitude longer than the nonlinear magnet insert section, so that the spacecharge perturbation to R is expected to be dominant over the space charge perturbation to N . (See Appendix C for a more detailed analysis.) We benchmark how well this simplified model can predict dynamic aperture against the prediction from a realistic space-charge simulation in Sec. V.

IV. COMPARISON OF CHAOS INDICATORS USING THE IOTA TOY-MODEL
This section has three goals. The first goal is to address the numerical reliability of REM. The second goal is to illustrate the relatively rapid convergence of REM in comparison with FMA. The third goal is to present and understand the dynamic aperture of the IOTA toy-model at a nominal design choice of nonlinear magnet tune advance ν ¼ 0.3, magnet strength τ ¼ −0.4, and various tune advance errors δμ x , δμ y , to compare with the realistic space-charge model in later sections. The incoherent tune depression can be modeled if we know δμ x and δμ x for each particles. However, for simplicity, we use same tune error δμ ¼ δμ x ¼ δμ y in both planes for all the particles. Note that the perturbation is autonomous in this model. Therefore, the IOTA toy-model can be a good test-bed for comparing REM against FMA on autonomous systems. To visualize the comparison of the two indicators, we will use dynamic aperture plots.
Throughout this and later sections, we use the following procedure to produce dynamic aperture plots. We start by preparing particles on a rectangular grid of points bounded by the circle x 2 n þ y 2 n < 1 in the transverse normal coordinate plane. The transverse normal momentum is set to zero. Then, we numerically track particles for T forward in time. This is followed by tracking an additional time forward T (for FMA) or backward −T (for REM). Then, each particle's initial condition in the transverse coordinate plane is mapped to a color based on the indicators Δ FMA or Δ REM , so that the brighter color indicates chaotic orbits while the darker color indicates regular orbits.
First, we perform reliability tests of REM by estimating the dynamic aperture of the map associated with the IOTA toy-model. Figure 14 shows one such test with various values of the tune advance error δμ. Note that dynamic aperture plots obtained using FMA have finer and richer structure than the REM results generally. However, such microstructures tend to vanish as we sample the data for a longer time, as shown in Fig. 16. Apart from the differences in such fine-scale structures, the agreement in the overall structures of the regular regions are apparent for various values of δμ in Fig. 16. The reliability can be further checked by tracking for a much longer time as shown in Fig. 15. Comparing the figures, we see that the majority of particles indicated to be on a chaotic orbit using the shorttime tracking data as shown in Fig. 14 are eventually lost after long time tracking as in Fig. 15.
Second, we perform tests of numerical convergence with the length of the tracking interval. Figure 16 shows the dynamic aperture prediction for the IOTA toy-model at ν ¼ 0.3, τ ¼ −0.4, and δμ ¼ 0.01 obtained by FMA and REM using data for various numbers of turns T per interval (Recall that two successive forward intervals are used for FMA, while one forward and one backward interval are used for REM). FMA is a powerful tool for revealing resonance structures in the phase-space [34], as it measures the tune diffusion over time. This is illustrated in Fig. 17.
Here, contour lines of ν x =ν y corresponding to low-order resonances are drawn over the background dynamic aperture plot, where ν x and ν y are the horizontal and vertical fractional tunes computed using the NAFF algorithm. The background is the same FMA dynamic aperture plot shown in Fig. 16, except that it is for T ¼ 128 instead of T ¼ 256, and the color scale was modified to improve visibility by increasing contrast with the resonant lines. The lines indicate resonances of the ideal (unperturbed) system. That is, these are resonances of the map N of Eq. (12), defined by n x ν x þ n y ν y ¼ 0 with n x ; n y ∈ Z. The resonant lines largely coincide with the chaotic regions indicated by FMA. However, in Fig. 16, much of the finer and richer structure (including the resonant lines) appearing in the FMA dynamic aperture plot at T ¼ 256 tends to vanish at T ¼ 2048. This suggests that such fine microstructures FIG. 14. Dynamic aperture plots for the IOTA toy-model at ν ¼ 0.3, τ ¼ −0.4 and various δμ with a fixed number of turns T ¼ 1024 for each interval. Two successive forward intervals are used for FMA, while one forward and one backward interval are used for REM. The color maps in the left and right columns are Δ FMA , and Δ REM respectively. From top to bottom δμ is increasing from 0.0 to 0.04 in increments of 0.01. Red x-marks denote the singular points of the nonlinear insert potential. The red curve is where I ¼ 2H. We draw the red curve only on the x n > 0 side for dynamic aperture visibility, but the I ¼ 2H curve should also be present on the side x n < 0, due to the reflection symmetry of H and I.  [34] for further details. The Henon-Heiles potential example presented in Appendix A further supports this argument. We believe that this occurs because the accuracy of FMA is very sensitive to the precision of the frequency measurement that can be obtained using a time series of finite length. In the limit of continuous-time data, it can be shown that the numerical error of the frequency vector obtained using the NAFF algorithm scales as an inverse power of the time length of the data, i.e., ∝ T −n for n > 0 [13]. Practically, for discrete-time data, the convergence can be worse. The sensitivity of FMA to the time length is especially high near resonances, as reported in [35]. On the other hand, the REM results in Fig. 16 show that the structure of the dynamic aperture plots is already converging with a time length as short as T ¼ 256, verifying the relatively rapid convergence of REM. Finally, we analyze in more detail the dynamic aperture of the IOTA toy-model at the nominal parameters used to produce Fig. 14. We emphasize that these figures illustrate the breakdown of integrability in the presence of a perturbation. The perturbation is given by an error in tune advance δμ over the arc section, which is motivated by the space-charge induced tune depression in the IOTA ring. When δμ ¼ 0, the system is perfectly integrable everywhere except at the singular points ðx n ; y n Þ ¼ ðAE1; 0Þ inherent to the potential in Eq. (8). This may be seen clearly in the uppermost figure in Fig. 14. However, numerical errors due to the finite integration step size and finite precision can break some of the invariant tori. When the invariant tori are broken by such small errors, it is likely that these tori will also be broken when physical perturbations including the tune error over the arc δμ are considered. Among the broken tori are those containing initial conditions that lie along the red curve shown in Fig. 14, which corresponds to the locus of points I ¼ 2H in the x n − y n plane (with p x;n ¼ p y;n ¼ 0, y n ≠ 0), where H are I are given in Eq. (8) and Eq. (10). It can be shown that the invariant level sets corresponding to these initial conditions contain unstable critical periodic orbits, and lie on a separatrix-like structure separating distinct dynamical regions in the phase space [31]. Note that the I ¼ 2H curves also form the innermost boundary of the dynamic aperture at the various values of δμ in Fig. 14.

V. SELF-CONSISTENT SPACE-CHARGE SIMULATION OF REALISTIC IOTA LATTICE
In this section, we apply FMA and REM to simulated particle orbits obtained in the presence of space charge, using a realistic IOTA lattice with a 2D symplectic spacecharge solver [36] implemented in IMPACT-Z [37]. We choose a 2D solver: (1) in order to model the nominal case of a long beam that nearly fills the ring, and (2) in order to isolate transverse space-charge effects from off-momentum chromatic effects, which are enhanced in the presence of a longitudinal space-charge force. The nonlinear integrability is preserved only for on-momentum particles unless a special chromatic correction is applied [38]. In this study, we choose a small beam current I ¼ 0.41 mA targeting a test proton operation of IOTA before ramping up to high current. The nonlinear effects over the arc section including geometric nonlinearity of dipoles, nonlinear kinetics, and nonlinear fringe fields (but linear fringe field effects are considered) are ignored so that weak space-charge contribution can be solely explored. Lattice parameters were adjusted accordingly such that δμ x ¼ δμ y ¼ 0 for particles near the beam center [39]. In order words, the phasesadvance of particles of zero amplitude limit is integer multiple of π. The nominal settings for the nonlinear insert are used: For the initial beam, we use a waterbag distribution, defined by the distribution function: where ΘðHÞ ¼ 1 for H ≤ 0 and ΘðHÞ ¼ 0 for H > 0. One million macroparticles are sampled from Eq. (15) with H 0 ¼ 0.06 (which corresponds to the nominal beam size) to represent the beam. In addition to the particles of the beam, we also prepared test particles of zero charge over the disk x 2 n þ y 2 n < 1 in the x n − y n plane, with zero normal momentum, as we did in Sec. IV, to obtain a plot of dynamic aperture for test particles moving in the combined space charge and external fields. Figure 18 shows the resulting dynamic aperture plots obtained using FMA and REM. The initial boundary of the waterbag beam is indicated by the white contour. None of the particles in the beam are lost and the evolution of the beam density is shown in Fig 19. Note that the FMA dynamic aperture plot looks noisy. This can likely be attributed to the time dependence of space-charge effects, including fluctuations of the beam density and numerical errors due to finite resolution of the space-charge fields. On the other hand, the REM result looks less sensitive to such fluctuations. Nevertheless, both indicators show that the I ¼ 2H curve is the major boundary of the dynamic aperture, as was the case for the IOTA toy-model. Note also that these figures are similar to the δμ ¼ 0.02 and δμ ¼ 0.03 cases shown in Fig. 14. This can be understood in the following way. Each particle is affected by spacecharge tune depression, and the particles farther from the beam center experience less tune depression [We observed the beam density decrease monotonically for the waterbag distribution Eq. (15)]. Since the tune error near the beam center is nearly zero, this means that the outer particles have positive tune error, making the dynamic aperture resemble the toy-model dynamic aperture at δμ ≤ 0.027. (The linear tune depression is 0.027 and 0.016 for the horizontal and vertical tunes, respectively. See Appendix C.) The agreement with the toy-model can be further improved by assigning the toy-model to have the different tune error parameters δμ x and δμ y for different particles. To illustrate, the tune advance over the arc section in Fig. 13 is supposed to be an integer. However, in the presence of space-charge, different test particles have different values of tune error over the arc section. We calculate the tune errors δμ x , δμ y in the presence of space-charge by tracking the beam and test particles over the arc section during the first turn. Although the beam has not yet reached an equilibrium state after the first turn, the computed values of δμ x , δμ y obtained using the first turn differ only slightly from those obtained over later turns, as the beam is exactly matched to the lattice at zero current, and the current is weak. In other words, the initial beam is already close to (but not in perfect) equilibrium as shown in Fig. 19. Figure 20 shows the dynamic aperture obtained using this modified toy-model, in which each test particle's tune errors δμ x , δμ y are obtained as described above. The result by REM agrees well with Fig. 18 at T ¼ 256.
In order to further study the stability of the particles in the beam, we applied REM to the orbits of the 1 million macroparticles sampled from Eq. (15). In order to visualize Δ REM of the 1 million particles, we divided the x n − y n plane into cells of a rectangular grid. Figure 21 shows the maximum value of Δ REM among the particles whose x n − y n initial conditions are in the corresponding cell. Comparing to the Δ REM values appearing in Fig 18, we can say that the particle orbits in the beam is quite regular.

VI. DYNAMIC APERTURE AT VARIOUS NONLINEAR MAGNET PHASE ADVANCE
In this section, we study the dependence of the dynamic aperture of the IOTA toy-model to the parameter ν appearing in Eq. (12) (the tune advance across the nonlinear magnet). In Section IV, we observed that in the case ν ¼ 0.3, the geometric structure of the dynamic aperture is determined primarily by the curve I ¼ 2H (that served as the innermost boundary of the dynamic aperture), which does not change much as the tune error δμ is varied. This makes sense, as the orbits with I ¼ 2H lie on a separatrixlike structure in the phase space, containing unstable periodic orbits that result in homoclinic points (and chaos) when a small perturbation is applied.
The geometric structure of the I ¼ 2H curve depends on the nonlinear magnet strength τ but not on the nonlinear magnet tune advance ν, as can be seen from Eq. (8) and Eq. (10). However, in the presence of a noninteger tune error δμ, the dynamics of the map N R can no longer be described by the Hamiltonian in Eq. (8) alone, so it is possible that the geometric structure of the dynamic aperture may change as we vary ν. In this section, we vary ν at a fixed perturbation δμ ¼ 0.01, to see if the innermost boundary of dynamic aperture changes. For this purpose, we divide the unit disk in the x n − y n plane into two geometric regions, as shown in Fig. 22. We refer to the blue region as region I, and we refer to the orange region as region II. So far, we have observed (at ν ¼ 0.3, in the presence of various perturbations, including space-charge) that most of the particles starting from region I lie on regular orbits, while many of the particles starting from region II lie on chaotic orbits. It can be shown that the particles starting in region I encircle the stable fixed point at the origin, while the particles starting in region II do not cross the x n axis, and thus are confined to the upper or lower half of the plane [31]. Therefore, the y n time series for each particle starting in region II has a large direct current (DC) component in the frequency domain. This fact can be used to detect a change in the structure of the separatrix as ν is varied. Figure 23 shows the dynamic aperture of the IOTA toymodel at various ν and fixed δμ ¼ 0.01. Notice that the boundary of the dynamic aperture is significantly changed in the two cases ν ¼ 0. 36  working point for the nonlinear magnet tune advance roughly in ν ∈ ½0.3; 0.34.
In the case ν ¼ 0.36, the dynamic aperture is drastically reduced when compared to the red curve. The chaotic band along the red curve bifurcates, and the bifurcated curve forms a bulge toward the origin, so that chaotic orbits appear in region I. To see if this is associated with a change in the geometric structure of the separatrix, we collected all of those particles whose orbital time series data y n has a dc frequency component larger than the dominant nonzero frequency component. Figure 24 shows the initial conditions of all such particles. Note that these particles paint the band of the borderline of the bulge in Fig. 23. This suggests that the geometric structure of the dynamics, as defined by the separatrix, is changed by a small perturbation δμ ¼ 0.01 when ν ¼ 0.36.
In the case ν ¼ 0.28, a new chaotic arc appears inside region I that is not connected to the I ¼ 2H curve. Figure 25(a) shows those particles whose orbital time series data y n has a dc frequency component larger than the dominant nonzero frequency component. Note that there are no such particles visible in region I. Figure 25(b) shows a tune footprint for all particles with initial conditions in region I. Note that many of the chaotic particles in region I appear to lie on resonance lines. This suggests that the new chaotic arc appearing in region I when ν ¼ 0.28 is not due to deformation of the separatrix, but is instead due to the presence of low-order resonances.

VII. CONCLUSION
The chaos indicator known as the reversibility error method (REM) was compared with the use of frequency map analysis (FMA) for studying the IOTA ring dynamic aperture in the presence of weak space-charge. The reliability of FMA is questionable when the Hamiltonian possesses nonperiodic time dependence (such as spacecharge). This does not present a difficulty for REM, which is better-suited for such systems. In addition, the REM indicator converged more quickly with respect to the length of the numerical tracking interval than FMA, making it a suitable application for numerically heavy simulations. The two indicators agreed well when fully converged for an autonomous system. Application to studies of the IOTA dynamic aperture revealed that a separatrixlike structure played a major role in limiting the dynamic aperture for various perturbation strengths. In addition, a simplified toymodel was used to represent the effects of space-charge, and we observed reasonable agreement with the realistic space-charge simulation. The toy-model is further used to study the effect of the nonlinear insert section phaseadvance on the dynamic aperture in the presence of a small perturbation.

ACKNOWLEDGMENTS
We thank the anonymous referees for their insightful comments that helped improved the paper. We also thank Ji Qiang for the support with computer resources. This work was supported by the U.S. Department of Energy under Contract No. DE-AC02-05CH11231 and used computer resources at the National Energy Research Scientific Computing Center.

APPENDIX A: CONVERGENCE OF CHAOS INDICATORS ON THE HENON-HEILES PROBLEM
In this section, we present a numerical convergence test of REM for comparison against FMA on the Henon-Heiles problem in Eq. (3). Figures 26 and 27 show dynamic aperture plots corresponding to orbits with E ¼ 1=12 and E ¼ 1=9, given for two values of the tracking time interval T. The orbits with E ¼ 1=12 are regular everywhere except on a set whose phase-space volume is zero [22]. In other words, there should be no visible chaotic band on the dynamic aperture plot. However, this fact is not well reflected by FMA in comparison against REM when the Here ν x and ν y are the horizontal and vertical fractional tunes.
T ¼ 512. In the case of E ¼ 1=9, chaotic bands appear on the dynamic aperture plot. In this case also, FMA tends to largely overestimate the chaotic region at T ¼ 512 while REM tends to slightly underestimate the chaotic region at T ¼ 512. Note that some chaotic resonant islands reported by FMA for E ¼ 1=9 and T ¼ 512 are vanishing for T ¼ 2048. In both energies, both indicators converged and agreed with each other as T becomes long. However, the convergence looks relatively faster for REM.

APPENDIX B: COMPARISON OF CHAOS INDICATORS ON A PERIODIC FOCUSING SYSTEM NEAR RESONANCE
Here, we illustrate the sensitivity of FMA to errors when used in a system near resonance. Consider a system described by the Hamiltonian Hðx; p; tÞ ¼ with period of length We interpret (B1)-(B2) as representing a storage ring near the quarter-integer resonance. We consider two cases by inserting thin and thick octupoles: (i) integrable case: A thick octupole field, of length L oct ¼ 2πð1=4 − ϵÞ, is superposed on top of the constant focusing field such that the map outside of the thick octupole field is an identity map. The potential is where Θ is the Heaviside step function. (ii) nonintegrable case: A thin octupole is added at t ¼ 0, given by the potential: where δ is the Dirac delta function. Figure 28 shows the dynamic aperture plot for the thick octupole case (B3). Since the system is integrable, the entire phase-space area is filled with regular orbits. However, FMA indicates orbits near the quarter-integer resonance (forming a bright ring on the x-p plane) to be chaotic compared to the surrounding phase-space area. This is because the accuracy of the NAFF algorithm is degraded near the resonance. (This phenomenon is illustrated in Ref. [34].) The accuracy can be improved by sampling data for a longer time interval. It is not shown in the plot, but we observed that the width of the bright ring appearing in the FMA dynamic aperture plot becomes narrower as the time interval T becomes longer. Although it is often informative to indicate the resonant lines, this can confuse interpretation, blurring the distinction between chaotic and regular resonant orbits in the dynamic aperture plot. Figure 29 shows the dynamic aperture plot for the thin octupole case (B4) obtained using FMA and REM. The resonance islands and the chaotic band near the separatrix are well captured by both indicators.

APPENDIX C: ESTIMATES OF SPACE-CHARGE TUNE DEPRESSION IN THE IOTA RING
It can be shown that the linear space-charge tune shifts in a proton ring are given by [12,40]: where s is the longitudinal coordinate along the reference orbit, β x;y and σ x;y denote the betatron amplitude(s) and rms beam size(s), and K SC ¼ 2nr 0 β 2 0 γ 3 0 . Here r 0 ¼ 1 4πϵ 0 q 2 mc 2 is the classical proton radius and n is the line density of number of particles. Using the waterbag distribution Eq. (15), we find that the linear space-charge tune shifts are ξ x ¼ −0.027 and ξ y ¼ −0.016. We also find that the ratios of tune depression over the arc section and nonlinear insert section are R insert dsβ x =σ x ðσ x þ σ y Þ R arc dsβ x =σ x ðσ x þ σ y Þ ¼ 0.036 ðC2Þ R insert dsβ y =σ y ðσ x þ σ y Þ R arc dsβ y =σ y ðσ x þ σ y Þ which verifies our claim that the tune depression over the arc dominates the tune depression over the nonlinear insert section.