Tune compensation in non-linear Fixed Field Alternating Gradient accelerators

In this paper, we investigate the stability of the particle trajectories in Fixed Field Alternating Gradient (FFAG) accelerators in the presence of field errors: The emphasis is on the scaling FFAG type: a collaboration work is on-going in view of better understanding the properties of the 150 MeV scaling FFAG at KURRI in Japan, and progress towards high intensity operation. Analysis of certain types of field imperfections revealed some interesting features about this machine that explain some of the experimental results and generalize the concept of a scaling FFAG to a non-scaling one for which the tune variations obey a well defined law. A compensation scheme of tune variations in imperfect scaling FFAGs is presented. This is the cornerstone of a novel concept of a non-linear non-scaling fixed tune FFAG that we present and discuss in details in the last part of this paper.

In this paper, we investigate the stability of the particle trajectories in Fixed Field Alternating Gradient (FFAG) accelerators in the presence of field errors: The emphasis is on the scaling FFAG type: a collaboration work is on-going in view of better understanding the properties of the 150 MeV scaling FFAG at KURRI in Japan, and progress towards high intensity operation. Analysis of certain types of field imperfections revealed some interesting features about this machine that explain some of the experimental results and generalize the concept of a scaling FFAG to a nonscaling one for which the tune variations obey a well defined law. A compensation scheme of tune variations in imperfect scaling FFAGs is presented. This is the cornerstone of a novel concept of a non-linear non-scaling fixed tune FFAG that we present and discuss in details in the last part of this paper.

I. INTRODUCTION
Scaling Fixed Field Alternating Gradient (FFAG) accelerator is a concept that was invented in the 1950s almost independently in the US, Japan and Russia [1] [2]. Several electron machines were built in the US. However, it was not until the 1990s that the interest for FFAGs was revived in Japan. Several machines were built among which a 150 MeV proton machine at Kyoto University Research Reactor Institute (KURRI). One main feature of this machine is its potential for high power applications, hence its use as a proton driver for an Accelerator Driven Sub-critical Reactor (ADSR) [3] at KURRI. A large dynamic acceptance can be achieved since the crossing of the betatron resonances is avoided in this concept. This is achieved by introducing a large field increase with the radius resulting in the beam experiencing the same focusing throughout the acceleration, therefore keeping the tunes constant. The magnetic field profile allowing this is written in the form B = B 0 (R/R 0 ) k where B is the vertical component of the magnetic field in the median plane, R is the radial coordinate with respect to the center of the ring, B 0 the reference field at R = R 0 and k is the average field index of the magnets, sometimes also referred to as the scaling factor, defined by k = R/B × dB/dR and designed to be a constant value everywhere in the ring. One shall insist here that this form of the magnetic field is a sufficient but non-necessary condition in order to obtain a fixed tune FFAG. In the scaling FFAG concept, the orbit excursion is uniquely determined by the field index. Given that all built machines in the past opted for a phase advance per cell below 180 • , the orbit shift is therefore large, of the order of 1 meter which makes the magnet larger than typical synchrotron magnets. Besides, due to the complexity of the field profile and flutter of the magnets, field imperfections can be problematic * hajtaham@gmail.com and difficult to cure since the orbits move outwards from the centre of the machine with increasing energies. This can lead to the crossing of several betatron resonances at low speed and to beam deterioration in consequence. For instance, this is a strong contribution to the overall low beam transmission in the KURRI FFAG [4] [5].
In the original work of Symon [1], it was established that the scaling FFAG satisfies the two conditions: where p is the particle momentum and n is the field index of the magnet defined by n = −ρ/B × dB/dx Here B is the vertical component of the magnetic field in the median plane and ρ is the bending radius of the particle trajectory. These conditions are known as the cardinal conditions: the first one expresses the constancy of the field index with respect to the momentum at any given azimuth, while the second one expresses the similarity between the different closed orbits. In other words, any closed orbit can be obtained from a photographic enlargement of any other closed orbit. These conditions lead to the following approximate values of the ring tunes ν x,y : where F is the magnetic flutter.
In what follows, we will discuss the validity of the Symon model by taking into account the field imperfections. One objective of this paper is to derive general expressions of the ring tunes that take into account the non-scaling of the orbits due to field imperfections and compare with the measured as well as the simulated values. This will be a crucial result in order to establish a correction scheme to minimize the tune variations in imperfect scaling FFAG.

II. BEAM STABILITY ANALYSIS
Technically, it is impossible to make a field which corresponds exactly to the designed one. Therefore, it is important to understand the effect of small imperfections of the field on the beam dynamics. The general equations of motion including non-linear terms and imperfections are defined by: where x and y represent the deviation around the closed orbit, s the curvilinear coordinate and P and Q are real analytic functions of x, y and s. In the following, we investigate the stability of the particle trajectories that can arise with different field errors. We use three different approaches to investigate the beam stability due to field errors: the first approach is based on the hard edge model: we assume that the magnetic field drops abruptly to zero at the edges of the magnet which greatly simplifies the analysis and allows the development of a geometrical model of the lattice. The second model is an analytical model based on the Bogoliubov method of averages which we use to seek an approximate solution of the particle equations of motion. The last and third model is the ZGOUBI [6] tracking model of the magnet, which is the most accurate one. The ZGOUBI model solves the non-linear equation of motion using field maps or user-implemented analytical models. To conclude, we establish a comparison between the different results and comment on the outcome of this study.

A. Field imperfections in scaling FFAGs
In order to carry out parametric studies of the field defects in scaling FFAGs, one decided to extend the original definition of the average field index of the magnet k to take into account the azimuthal and radial field variations. Thus, one introduces the average field index of the focusing (F ), defocusing (D) magnet and the drift (drif t) in the following way: where B i is the vertical component of the magnetic field in the median plane of the FFAG magnet that is averaged over the width of the element. Since the drift space between the magnets is likely to contain the fringe fields, it is important to assign a field index to it to determine its effect on the beam dynamics. However, in the ideal case, k drif t = 0, which we will assume in the future for simplicity, yet without loss of generality. Now, assuming that the k-values have no radial dependence (a complete derivation of the expression of the field when k is Rdependent can be found in the Appendix A), Eq. (5) can be integrated and the magnetic field expressed in cylindrical coordinates: where F F and F D are the fringe field factors (or flutter functions) that describe the azimuthal variation of the field in the F and D magnets respectively, and a is a scale factor that allows to vary the FD ratio of the magnet (a ≥ 0). It is important to note that the field is a separable function in radial and azimuthal coordinates since the fringe fields merge to zero in our model between the magnets as can be seen in Fig. 1, thus F F (θ)F D (θ) = 0. Also, note that if k F = k D , the field writes in the standard form of a scaling FFAG. The lattice considered for this study is a radial sector KURRI-like DFD triplet [7]. B. Analytical expression of the betatron wave numbers

Hard edge model
In this approach, the analysis is restricted to a linear motion around the closed orbits. Therefore Eqs. (4) simplify to: which is the set of linearized equations with periodic coefficients valid in a close neighbourhood of the equilibrium orbit. This equation is a Hill's equation. For each closed orbit, there is a different set of linearized equations, i.e. different pairs (K x , K y ) so that a solution of the above equations cannot be generalized to the whole system. However, in the ideal case of a scaling FFAG, there exists only one pair (K x , K y ), which guarantees that the betatron wave number remains constant for all energies. One shall recall that the linear system is only a first order approach to the problem so that the following stability analysis is mainly valid in the vicinity of the closed orbits. However, for small perturbations, the stability of a non-linear system may often be inferred from the stability of its associated linear system. In order to solve the above equations (7), one used a piecewise method that consists in following the particle trajectories through a series of elements defined by their transfer matrices. The main assumption in this model is that the field within each magnet is constant as a function of the azimuthal angle, yet drops abruptly to zero at the edges. Thus, the particle trajectory within each magnet is described by an arc of a circle with a constant bending radius as shown in Fig. 2. The radial dependence of the field is maintained. This allows the development of a geometrical model of the ring. To begin with, one re-writes the expression of the field index of the magnet by neglecting the scalloping of the non-equilibrium orbits: In reality though, the second term containing the azimuthal derivative of the field yields an additional focusing effect, known as Thomas focusing, or edge focusing that will be taken into account later on in this paper. The relation among the various geometrical quantities can be found in the appendix B. Now, based on the Symon formula (Eq. 3), one makes the following conjecture: where (x 1 , x 2 , x 3 , x 4 , x 5 ) and (y 1 , y 2 , y 3 ) are constants that are strictly related to the geometry of the cell. If the conjecture is true, then solving a system of 5 nonlinear equations in the horizontal plane, and a system of 3 linear equations in the vertical plane, will provide for a given cell the horizontal and vertical tunes, respectively. For instance, we solve the following system in the vertical plane, The results are finally summarized in Figs. 3 and 4 which show good agreement between the hard edge model and the conjectured formula.
In reality though, the similarity of the orbits is no longer valid in presence of a difference of the scaling factor of the F and D magnet. To show this, one can proceed using a reductio ad absurdum: let's suppose that ρ i /R i = const. Then, from Eq (B2) to (B5), it results that ρ F /ρ D = const and R F /R D = const. However, if we write, this yields, which is R-dependent if κ = 0. This is in contradiction with ρ F /ρ D = const, quod erat demonstrandum.

Bogoliubov method of averages
In the previous section, it was established that the effect of scaling imperfections is to introduce a non-scaling of the orbits which is a source of tune variations. In order to determine the effect of the field defects on the tunes, one used the Bogoliubov-Krilov-Mitropolsky (BKM)'s method of averages [8] to compute approximately the frequencies of the betatron oscillations and their dependence on the average field index of the F and D magnets. This method is characterized by the expansion of the dependent variables in power series of a small parameter as detailed below. We start from Eq. (7) which is re-written into a cylin-drical coordinate system in the form: where µ(R, θ) = R(θ)/ρ(θ) is a geometrical parameter that is independent of the radius in the case of a scaling FFAG, x and y are the normalized transverse particle coordinates in the Serret-Frenet coordinates system and R is a parameter characterizing the scalloping of the orbits (see Appendix C). Now Eqs. (13) can be re-written in the standard form where x stands for x or y, R is the local radius of the closed orbit at a given azimuth θ and g x (R, θ) is periodic in θ with period 2π/N , N being the number of sectors in the ring. Then, applying the averaging method, we find that the betatron wave number per turn is given by: where the first order term g 1 is the mean value of the forcing term in Eq. (14). Now, the preceding formalism applies by replacing g(R, θ) by Next, we introduce the azimuthal variation of the scaling factor in the following way Starting from the vertical plane, we calculate the first order term: where E is the kinetic energy of the particle and α i (E) is what we will refer to in the future as the 1 st order index of similarity of the closed orbits: The azimuthal variation of the field yields an additional positive edge focusing effect (see [9]), so that, to the first order, the vertical tune writes: where ξ is the spiral angle of the magnet and F is the magnetic flutter (not the flutter function F ) defined by: The same steps are repeated in the horizontal plane by calculating the first order term where β i (E) is what we will refer to in the future as the 2 nd order index of similarity of the closed orbits: Thus, to the first order in the method of averages, the betatron wave number in the horizontal plane can be approximated by: where the azimuthal field variations yield an additional term included here ( [9]) although often neglected for a large number of sectors.
The assumption on the similarity of the closed orbits is withdrawn here given that α i is allowed to change with the radius. Also, note that for the F magnet, R/ρ¿0, so that α F < 0. However, for the D magnet, α D > 0. Therefore, by means of Eq. (18), one can predict that, to the first order in the method of averages, the D magnet increases the vertical focusing, while the F magnet plays the opposite role. Also, the drift is included in this formula. Since the drift space between the magnets is susceptible to contain the fringe fields, it may be important in some cases to assign a scaling factor to it to determine its effect on the beam dynamics. However, in the ideal case, ρ drif t → ∞ so that α drif t ≈ 0. Now, let's compute analytically the expression of the index of similarity of the closed orbits. A closed orbit curve can be defined in a cylindrical coordinates system in the following way: where h(N θ) has period 2π in N θ and average value 0. The radius of curvature can thus be computed [10] and the index of similarity for a given energy E writes in the following way: Inserting Eq.(25) into (26) and integrating, it results that, for any closed orbit: where α is a parameter measuring the orbit scalloping and that is energy-dependent (0 α < 1). Similarly, one can calculate the index of similarity β. Note that, in the case where the scalloping is neglected, α = 0, β = 1. Inserting Eq. (27) into (18), one obtains: which represents the changes of the square number of betatron oscillations due to the average magnetic gradient encountered by the particle. This is a generalization of the hard edge model conjecture where one assumed that all orbits remain similar, i.e. α F (E) = const. Now, one evaluates α and β as a function of the scalloping amplitude: one assumes that h(N θ) = cos(N θ). Inserting Eq.(25) into (26) and integrating it yields the results shown in Fig. 5. The index of similarity β grows rapidly with the scalloping amplitude f . In addition, the changes of β with the non-scaling of the orbits, i.e. the variation of the scalloping amplitude, are dominant in comparison to α. This is expected due to the fact that the orbits are closed so that when calculating α, the effect of the scalloping in the F and D magnet cancel each other. Note however that, in the Symon model [1], R/ρ is defined as a periodic function of θ only, so that any effect of the non-scaling of the orbits cannot be taken into account. Now, assuming that k F = k D and neglecting the scalloping of the orbits, one recovers Eq.(3).

Property
In this section, we establish the following property: In presence of scaling imperfections, the number of betatron oscillations per turn increases (resp. decreases) with the energy if κ > 0 (resp. κ < 0). Besides, the variation of the tune squared is, to the first order, From what preceded, it can be predicted that the tune variations are imposed by the β term in the horizontal plane and the magnetic flutter in the vertical plane. Let's start from the horizontal plane: This means that, in order to explain the tune variations, one has to focus on explaining the changes of R/ρ F with the energy of the particle: Let's consider the case of a scaling FFAG, for which k F = k D = k and the associated non-scaling FFAG for which k D is perturbed. Writing the expression of the average magnetic field of the non-scaling FFAG, one obtains: where ns (resp. s) denotes the non-scaling (resp. scaling) case and A 3 is the FD ratio of the magnet, defined by: Defining the average radius R E of a closed orbit of energy E such that B R E = p/q where p is the particle momentum, this yields: which is an increasing (resp. decreasing) function of the energy if κ > 0 (resp. κ < 0). Knowing that ρ F satisfies Given that β s F = const, this proves the first part of the property. The result can be interpreted in the following way: when the radial field increase of the D-magnet is faster than that of the F-magnet, the scalloping of the orbits increases as well, therefore increasing the β term as shown in Fig. 5 above, which is a measure of the horizontal restoring force. Now, if we write R = R 0 + ∆R where R is the maximum radius of the closed orbits and R 0 is the injection radius, then a Taylor expansion of Eq.(34) yields: which proves the property above in the horizontal plane. Now, we establish the result in the vertical plane by cal-culating the magnetic flutter. From Eq. (6), it results that: Also, Thus, given that F F (θ)F D (θ) = 0, and setting F F and F D such that F F (θ) = − F D (θ) = 1, one obtains which is an increasing (resp. decreasing) function of the energy if κ > 0 (resp. κ < 0). Thus the vertical tune is an increasing (resp. decreasing) function of the energy if κ > 0 (resp κ < 0). Now, the vertical tune excursion writes which proves the property above in the vertical plane. It is important to highlight that the validity of the above result fails when the magnet contains an overlapping region of the field of the F and D magnet. However, the above calculation can be extended to take this effect into account. Equation (41) shows that reducing the FD ratio helps reduce the tune variations. This is expected since, in our model, the D magnet is the source of the field defect. Furthermore, increasing the alternation of the gradient increases the sensitivity of the tunes to the field imperfections via the second order moments of the flutter functions, i.e. F 2 D (θ) and F 2 F (θ) . This is further discussed in the next section. In addition, it is shown that the effect of the scaling imperfections on the tune variations grows linearly with the radial excursion of the orbits in both horizontal and vertical planes. This shows the advantages of having an FDF triplet configuration (rather than a DFD one) with much larger field index, since then the orbit excursion can be reduced by a factor of 5 or more [11] leading to much lower tune variations due to field errors of the type described above.

C. Tracking simulations
Now, we demonstrate that our findings with the previous analysis are reinforced by numerical simulations and provide an application example which is the KURRI 150 MeV scaling FFAG. From the hard edge model and the Bogoliubov method of averages, it was found that the tune is sensitive to the average field index k F and k D of the F and D magnet respectively. In other words, breaking the scaling law, although a major source of imperfection in scaling FFAG, can also be utilized in order to control the tune path in FFAG. In order to quantify the source of imperfection, we introduce two new quantities in the calculation: the average value of the tunes ν m x = ν x and ν m y = ν y over the range of energies to quantify the average focusing strength of the applied forces on the beam, and the rms value of the tunes ν rms x = σ νx and ν rms y = σ νy to quantify the scaling imperfections in terms of tune variations. One could instead use the |max − min| value of the tunes to account for the oscillations. However, the rms quantities have the merit to be average quantities, thus more appealing to use in order to obtain smooth variations of the described quantities.

Benchmarking work
Following the FFAG14 workshop held at BNL[12], a simulation campaign was established to benchmark several simulation codes. The main objective is to provide reliable modelling tools for FFAG type of accelerators and to better explain the results of the experiments at the KURRI 150 MeV scaling FFAG [5]. The first benchmarking test was carried out for the calculation of the betatron tune as a function of the momentum and shows excellent agreement between the different codes as shown in Fig. 6 as the details of the simulation can be found in [13].

ZGOUBI tracking model
Based on the successful benchmarking test, we carry out parametric studies based on the ZGOUBI tracking code. The ZGOUBI code solves the non-linear equation of motion using truncated Taylor expansions of the field and its derivatives up to the 5 th order. Thus, it is more accurate than the linear approach. Given that the energy gain per turn is small in scaling FFAGs (typically 2 keV per turn in the KURRI machine), one can reasonably assume that the accelerated orbit trajectory for any given energy is quasi the same as the closed orbit trajectory. Thus, the procedure employed for the calculation of the betatron wave number is based on the closed orbits formalism described below: First, one generates a median plane field map for a given (x, k F , k D ) as shown in Fig. 1. The field fall-off at the end of the magnets is obtained by using an Enge-type fringe field model [14]. Extrapolation off the median plane is then achieved by means of Taylor series: for that, the median plane symmetry is assumed and the Maxwell equations are accommodated. This yields results in excellent agreement with the 3D field map calculation.
Second, search for NCO closed orbits between injection and extraction using the built-in fitting routines in ZGOUBI. NCO was chosen to be 30 in order to have good statistics and ensure the convergence of the calculated quantities. A typical example of 4 closed orbits search is illustrated in Fig. 7.

FIG. 7. Example of several closed orbits for a scaling FFAG .
Lastly, for each closed orbit, the betatron wave number is calculated in both planes. Fig. 8 shows the stability diagram obtained by varying the average field index of the magnets (κ = 0) as well as their FD ratio, therefore the scale factor x in Eq. (6). One can observe that, on the top and bottom right, the stability limits are set by the horizontal and vertical cell tunes, respectively. On the left side, the physical size of the magnets (here we choose a radial excursion limited to 10m) determines the boundary limits. Now, we choose to focus our analysis on the average field index of the magnets. For that, we fix the FD ratio. We choose a = 1 which corresponds to the green line in Fig. 8. A scan on k F and k D provides the stability diagram of the DFD triplet in the transverse plane (see Fig. 9). Qualitatively, it shows that, in the case where k F = k D = k, the average cell tune exhibits the expected behavior predicted by the Symon formula i.e. Eq. (3): increasing k increases the horizontal tune and decreases the vertical one. One can also observe that for large k values, the stability diagram shrinks, thus any design imperfection will make the orbits quickly unstable. This is explained in the following way: on the right side of the stability diagram, i.e. when κ < 0, the stability limit is set by the condition that ν 2 y is to remain positive (Floquet resonance), given that the tunes decrease with the energy. Using Eq. (41), the lower limit is set by: where ν y is the tune of the unperturbed lattice and A is given by: Now, on the left side of the stability diagram, i.e. when κ > 0, the stability limit is set by the radial π-mode stopband resonance, given that the tunes increase with the energy. Using Eq. (36), the upper limit is set by: where B is given by: This shows that the maximum allowed scaling imperfection |κ| decreases linearly as a function of k for a lattice with fixed FD ratio. Note that a second stability island exists for larger k values [11]. However, we restrict our analysis to a KURRI type FFAG for which the design value of k ≈ 7.6. Now, calculating the RMS tune variations shows that the latter exhibit the expected behaviour in the vicinity of the line k F = k D where they become negligible. This is shown in Fig. 10. When field imperfections such that κ = 0 are introduced, one can observe that the tune variations increase with |κ| as demonstrated earlier.  Based on all the above, we compare the tracking results with those obtained from the analytical formula established in the previous section. This is shown in Figs. (11) and (12) for the horizontal and vertical plane, respectively. The blue points are the simulation results while the red points represent the analytical formula: for κ > 0 (upward-pointing triangle), the tunes increase with the energy in both planes, while for κ < 0 (down-pointing triangle), the tunes exhibit the opposite behaviour. This confirms the findings of the previous section. Kinetic Energy (MeV) (k F ,k D )=(7.6,7.0) (k F ,k D )=(7.6,7.6) (k F ,k D )=(7.6,8.2) (k F ,k D )=(7.6,7.0) (k F ,k D )=(7.6,7.6) (k F ,k D )= (7.6,8.2) FIG. 11. Example of tune calculation using the tracking code ZGOUBI (blue) and comparison with the 1st order analytical formula 24 (red). To the first order, good agreement (< 4%) is obtained in the cases (kF , kD) = (7.6, 7.0) and (kF , kD) = (7.6, 7.6). However, for the last case, i.e (kF , kD) = (7.6, 8.2) the large discrepancy observed ( 16%) is due to the fact that the beam is located nearby the boundary of stability. FIG. 12. Example of tune calculation using the tracking code ZGOUBI (blue) and comparison with the 1st order analytical formula 20 (red). To the first order, good agreement (< 4%) is obtained in the cases (kF , kD) = (7.6, 7.0) and (kF , kD) = (7.6, 7.6). However, for the last case, i.e (kF , kD) = (7.6, 8.2) the large discrepancy observed ( 16%) is due to the fact that the beam is located nearby the boundary of stability.
The main finding of the analytical model is that scaling imperfections produce an orbit distortion that manifests through a radial dependence of the index of similarity of the orbits as well as the magnetic flutter. The tracking simulations confirm our findings: as shown in Fig.  (13), the magnetic flutter and the index of similarity β explain the monotonic behavior of the tune as a function of the energy for various κ values. Figs. 14 and 15 show a comparison of the analytical formula (35) and (36) with the tracking results: in our example, the FD ratio gates the impact of the field profile on the tune variations of the FFAG in presence of scaling imperfections. The idea is to maintain a fixed FD ratio of the magnet (A 3 fixed) in order to have the same closed orbits for different lattices. By adjusting the field profile of the D magnet, one aims at investigating the sensitivity of the lattice to scaling imperfections. As shown in Eq.(41), reducing the mean square value of the flutter function of the D magnet helps reduce the tune variations due to magnetic field gradient errors. This is confirmed through numerical simulations using the analytical model "FFAG" [15] in ZGOUBI where one varied the flutter functions (Fig. 16). The results, in agreement with the formula are illustrated in Fig. 17 below: as expected, reducing the width of the D-magnet and increasing the amplitude of its field (in order to maintain the same FD ratio) helps increase the focusing strength in both planes. However, the sensitivity to the gradient errors increases accordingly. This is explained by the fact that the azimuthal field variations increase when reducing the width of the magnet while maintaining the same average bending angle of the particle trajectory, thus contributing to the increase of the orbit scalloping. The variation of the orbit scalloping with the radius increases when gradient field errors are present so that the variation of the Thomas focusing increases accordingly.

Application to the KURRI 150 MeV scaling FFAG
We will benchmark the analytical formula against the simulated values obtained from particle tracking using 3D field maps and compare with the measurement. According to the previous results, the number of betatron oscillations per turn in the KURRI DFD triplet are given by: The tracking results provide the closed orbits for various particle energies which are then exploited to calculate the horizontal and vertical tunes by applying the 1st order Bogoliubov method of averages (see Eqs. (46) and (47)). The hard edge model formula that was established earlier is also shown (see Eqs. (9) and (10)). The main difference between the BKM method and the hard edge model formula is that the only parameters that may vary in the latter are the average field indices k F and k D . The results are also compared to the 1st order Symon formula (ν 2 x = k +1 and ν 2 y = −k +F 2 ) in which the scaling factor k as well as the magnetic flutter F 2 are allowed to change with the radius and were obtained from the closed orbits calculation by Zgoubi. Comparison between the different results are summarized in Figs. 18 and 19. To begin with, the Symon formula is less in agreement with the benchmarking results than the others. The main reason is the missing index of similarity β (Eq. (29)) of the orbits in the horizontal plane which results from the mean restoring force of the scalloped orbit. The latter, along with the variation of the mean field index, explain the monotonic behavior of the tunes as a consequence of κ > 0.
Although the scalloping of the orbits in the KURRI machine is only of the order of 3%, this is enough to explain the discrepancy between the Symon model and the tracking results. In the vertical plane, the discrepancy is due to the cross-talk between the F and D magnets which creates an additional focusing effect in the drif t space as well as produces a variation of the magnetic flutter. In summary, the Bogoliubov method of averages yielded good agreement with the tracking results as well as the measured values.

III. CORRECTION SCHEME AND ADVANCED FFAG CONCEPT
Practically, it is difficult to correct the orbit and optics distortion in fixed field accelerators for the entire momentum range since the beam moves outward radially during acceleration. Therefore, a dedicated correction system should be implemented along the radius of the magnet to produce the desired field profile. From the point of view of cyclotrons, this consists in the implementation of trim coils to correct the isochronism, i.e. the revolution time of the orbits. From the point of view of a scaling FFAG, the main target is to fix the betatron wave number in both planes which will allow to avoid the crossing of transverse resonances and maximize the overall beam transmission from injection to extraction. From what preceded, one was able to explain the monotonic behaviour of the tunes as a function of the energy as well as the amplitude of the tune variations. This is a crucial result if one aims to reduce the tune excursion. One major outcome of this study is that gradient errors in the FFAG magnet yield a non-scaling of the orbits and lead to a change of the average focusing forces applied on the beam. This means that fixing the field defect of the FFAG magnet by aiming to produce a constant average field index k (by considering the average field over the entire circumference) is not sufficient since the azimuthal variations of this quantity yield a non-scaling of the orbits. In what follows, we present a novel scheme to correct the field errors in FFAG that relaxes the constraint of having the scaling law valid everywhere in the ring.

A. Alternating scaling imperfections
In the context of the present study and for the sake of simplicity, the field defect is due to one of the magnets, either the focusing (F) or the defocusing (D) one such that κ = 0. Without any loss of generality, we assume that the D magnet is the source of imperfection. In order to minimize the tune variations, one way is to reduce the FD ratio of the DFD triplet. This can be achieved by reducing its excitation current or by sliding it to outer radii so that the average field encountered by the particle at any radius is lower. However, this approach leads to the loss of focusing in both planes as shown earlier (see Fig. 8). Based on the property established earlier, the tune variation with the energy exhibits antagonistic behavior based on the sign of κ. Therefore, the idea of the following correction scheme is to introduce a perturbation of the field every P sectors in order to counteract the already existing imperfections. For instance, if we choose P = 2, then a 12-fold symmetry machine is replaced by a 6-fold symmetry in the following way: let's note D i (resp F i ) the Defocusing (resp Focusing) magnet with scaling factor k Di (resp k F i ). The original design 12 × (D 0 F 0 D 0 ) is replaced by 6 × (D 0 F 0 D 0 D 1 F 0 D 1 ) in the following way: if k D0 > k F 0 then k D1 < k F 0 and viceversa. Thus, instead of aiming to fix the design imperfections by correcting the field profile to match with the ideal one for every magnet and make the orbits scale at every azimuthal position, one can fix the scaling of the orbits on an average sense by creating alternation of the difference of the average field index of the magnets. This has a major advantage of reducing the cost of the cor-rection system since then only 12 D-magnets (out of 36 magnets) will require trimming coils to be implemented. The number can be further reduced if increasing P . This is the cornerstone of the fixed tune non-linear non-scaling FFAG concept that is discussed in detail in what follows.

B. Fixed tune non-scaling FFAG
In what follows we shall examine the characteristics of this concept which incorporates the alternation of the difference of the gradients scheme to FFAG. The orbital trajectories in the median plane are shown in Fig. 20 where one compares the orbits of a 12-fold symmetry scaling FFAG 12 × (k F 0 , k F 0 ) with those of a 6-fold symmetry novel concept 6 × (k F 0 , k D0 ), (k F 0 , k D1 ): the orbit scalloping in the non-scaling case differs mainly in the dominant F-magnet. If κ < 0, the radial field increase in the F-magnet is faster than in the D-magnet so that the equilibrium orbit (pink) in the F-magnet is at lower radii compared to the scaling case. The opposite effect occurs when κ > 0 (light blue orbit). As a consequence FIG. 20. Closed orbits in a scaling FFAG (black) and in a fixed tune non-scaling FFAG with alternating κ (pink and light blue). For the sake of clarity, the distance between the closed orbits of the scaling and the non-scaling FFAG is amplified.
of the alternation of κ, the monotonic behavior of the phase advance per cell is alternating (increasing if κ > 0 and decreasing if κ < 0). This is illustrated in Fig. 21 where one can see that the combination of the two cells yields a fixed average tune per cell. This results from the alternation of the monotonic behavior of the horizontal restoring force as well as the Thomas focusing effect which is due to the azimuthal change of the orbit scalloping. The magnetic flutter (Eq. (21)) calculation (k F ,k D )=(7.6,7.8) (k F ,k D )=(7.6,7.8), (7.6,7.4) (k F ,k D )=(7.6,7.4) FIG. 21. Tune variations as a function of the energy before and after correction: the corrected scheme is shown in green where (kF 0, kD0) = (7.6, 7.8) and (kF 1, kD1) = (7.6, 7.4). These results are obtained from multi-particle tracking assuming a Gaussian distribution of the particles in the transverse plane. The errorbars represent the overall tune variation from the distribution.
within each cell as well as for their combination is shown in Fig. 22 and is in agreement with the property established above. Similar result is obtained for the β term.

C. Dynamic Acceptance
Although strong non-linearities of the field are inherent to the scaling FFAG, large Dynamic Acceptance (DA) is obtained with this concept. Therefore, one main question to answer is how the DA of the fixed tune non-scaling FFAG compares to that of the scaling FFAG. In our analysis, the DA is defined as the maximum transverse invariant value that the beam can have without loss due to single particle dynamics effects. Particle tracking at fixed energy is employed for our analysis. A particle with original displacement from the closed orbit is defined as stable if it survives 1000 turns. Given that the vertical aperture in fixed field accelerators is the limiting factor due to the small gap size, we focus our analysis on the horizontal plane. However, in our simulations, it is noted that a vertical beam size up to 1 cm at injection can be transported without any losses and that the horizontal DA is insensitive to the vertical beam size. The main idea is to generate two lattices that have the same tunes in both planes. This is obtained by first generating a non-scaling fixed tune lattice with κ ≈ 0.3, then matching its tunes by finely tuning the average field index as well as the FD ratio of the scaling lattice. The results are shown in Fig. 23 where one can see that both results match very well. Comparison of the calculated DA in both cases shows (Fig. 24) that, for the same tunes, the horizontal beam acceptance is the same even though the orbits do not scale. This is only valid if resonance crossing is not occurring: in the non-scaling case discussed here, and for P = 2, the resonance population is doubled in comparison with the scaling case. Thus, one can expect that the resonance crossing problem is more severe in the non-scaling case. This requires further investigation. Comparison of the phase space trajectories between the two cases is finally shown in Fig. 25. The trajectories in both cells are symmetric with respect to X = 0.

IV. CONCLUSION
In this paper, one analyzed the stability of the particle trajectories due to field errors. Several approaches to the problem were developed and analyzed. Comparison of the results showed that the linear approach is only valid in the vicinity of the closed orbits: all parameters of the hard edge model evolve due to the non-scaling of the orbits. Therefore, one relied on the non-linear approach based on tracking simulations. A crucial result was to establish a relationship between the betatron wave number and the field defects. A key parameter to measure the amplitude of the defects is the κ-value defined as the difference of the average field index of the focusing and defocusing magnets. Furthermore, analysis of the stability diagram ( Fig. 9) showed that the tolerance to scaling imperfections becomes lower when increasing the average field index of the magnets. The limit is set by equations (42) and (44). Based on these results, a new scheme to remediate the variation of the betatron oscillations with the energy was proposed. The main idea consists in alternating the κ-values of the magnets, every two (or more) sectors. This lead to the new concept of the fixed tune non-scaling FFAG that one developed in section 4: in addition to the fact that this demonstrates that the conditions of scaling are non necessary to obtain a fixed tune FFAG, the newly developed concept is easier to implement by means of trim coils that can be adjusted to find the condition of minimum tune excursion and avoid the crossing of harmful resonances. Analysis of the DA showed that the lattice with alternating κ-values has the same DA as the equivalent scaling FFAG case. Given that the alternating-κ FFAG reduces the number of super-periods in the accelerator, therefore doubling the resonance population, one can expect that the impact of the resonance crossing is more severe than the scaling FFAG case. This needs further investigation.

ACKNOWLEDGMENTS
The authors would like to express their special gratitude to the members of the international Kyoto University Research Reactor Institute collaboration. The research leading to these results has received funding from Brookhaven Science Associates, LLC under Contract No. de-sc0012704 with the U.S. Department of Energy.
Appendix A: Analytical expression of the magnetic field to account for radial defects In order to obtain the radial dependence of the field when the mean field index k is R-dependent, let's assume that k can be fitted with an n-order polynomial. Thus, k writes in the following way: Then by equating Eq. (5) and Eq. (A1), one obtains: which gives after integration: so that the general form of the magnetic field becomes, B(R, θ) = B(R)F (θ), i.e., Appendix B: Geometry of the equilibrium orbits in the hard edge model The geometrical model of the lattice is shown in Fig. 2. All the geometrical relations are expressed as a function of the bending radius ρ F of the focusing magnet.
Applying the sinus law in the triangle OO F B yields: Applying the sinus law in the triangles OBI, O F BI and OBC yields: Applying the sinus law in the triangle OBC yields: Applying the sinus law in the triangles OCD and O D CD yields: Applying the sinus law in the triangles ODE and OCD yields: The cardinal condition given by Eq.
(2) implies that the ratio ρ i /R i is constant and shall be the same for different closed orbits. This implies that the bending angles θ F and θ D are constant and all the relations established above as well.