Weak-strong beam-beam simulations for the Large Hadron Collider

A weak-strong simulation code is used to study the single-particle stability in the presence of triplet field errors, head-on collisions, and long-range beam-beam interactions at the Large Hadron Collider. We present the dependence of the simulated transverse diffusion rate on various parameters, such as starting amplitude, working point in tune diagram, crossing angle, beta function at the interaction points (IPs), beam current, triplet nonlinearities, tune modulation, and a transverse offset at one of two IPs. For several examples, we perform a frequency map analysis à la Laskar, to obtain tune footprints and the tune variation in time. A cursory look at the effect of a Möbius lattice is also reported.


I. INTRODUCTION
The beam-beam interaction has limited the performance of most storage ring colliders (see [1][2][3][4] for past and recent experience). It also imposes serious constraints on the design and operation of the Large Hadron Collider (LHC), a proton-proton collider under construction at CERN, which will bring into collision, at a cm energy of 14 TeV, two counterrotating beams, each one consisting of 2835 particle bunches. The LHC will enter a new regime of beam-beam interaction, in that there are almost 20 long-range beam-beam encounters on either side of each interaction point, before the beams are separated into two independent rings.
These long-range encounters can introduce a strong perturbation at amplitudes where particles come close to the center of the opposing beam. Using this criterion, a diffusive or chaotic aperture could be expected at where u c denotes the full crossing angle between the two beams, u x the rms divergence, and s x the rms beam size at the interaction point (IP). For the nominal LHC parameters, listed in Table I, x da amounts to about 9.5s x . Simulation studies conducted by Irwin for the Superconducting Super Collider (SSC) showed that, in the presence of tune modulation, the diffusive aperture of the SSC was actually about 2.5s x smaller, a reduction which was found to be independent of the crossing angle [6]. Hence, extrapolating from these results, we would predict the LHC diffusive aperture at about 7s x . However, the LHC and SSC parameters differ substantially, e.g., the SSC bunch population was only 7.3 3 10 9 , 15 times less than that of the LHC. Thus, a simple extrapolation may not be valid. In addition, since the beam-beam interaction is the most important limitation of the LHC performance at top energy, it is important to study the effect on the beam stability of various related parameters, such as the beam current, the beams' crossing angle, or the IPs beta function, the values of which have a direct impact on the luminosity.
Here, we report the results of a dedicated weak-strong simulation study for the LHC. Different from previous beam-beam dynamic aperture studies [7], in this work we evaluate several indicators of diffusion, such as the long term action variance or the short term tune variation. In the latter case, we use Laskar's method of frequency map analysis [8][9][10], for the first time in this context.
The article is organized as follows. In Sec. II, we describe the employed weak-strong model. Section III presents the simulation results, i.e., tune footprints, diffusion maps, and the dependence of the action diffusion on different beam-beam parameters. We summarize and conclude our work in Sec. IV.

II. MODEL
The simulation study was performed following the recipe given by Irwin [6], and using the LHC parameters of Table I. The simulation is four-dimensional: the horizontal and vertical motion of single particles is calculated under the influence of the field of the opposing beam. We treat two IPs, symmetrically spaced around the ring, one with a horizontal crossing angle, the other with a vertical crossing angle, so that the linear tune shifts induced by the long-range interactions cancel between the IPs [11]. At each IP, we apply a series of five kicks: the first represents the lumped nonlinear effect due to the superconducting quadrupole triplet, the next corresponds to the long-range beam-beam interactions on the incoming side, then a kick for the head-on collision effect is applied, another for the long-range interaction on the outgoing side, and finally a kick due to the triplet nonlinearities on the outgoing side.
An optional tune modulation can also be selected. These tracking elements are described in more detail below.

A. Arcs
Between the two IPs, denoted by subindices i and j, we perform a linear rotation of the form 0 with the rotation matrix The coordinates are those at the successive IPs and primes denote the particle slopes. Unless noted otherwise, we assume that the two rotation matrices, from IP 1 to IP 2 and from IP 2 to IP 1, are identical (i.e., R 12 R 21 ), and also that the beta functions at the interaction points, b x and b y , are the same in both planes and at both IPs and equal to the LHC design value b x b y 0.5 m. The bare half-ring phase advances were chosen as f x 2 3 p 3 31.655 and f y 2 3 p 3 29.66, corresponding to the nominal working point of LHC optics version 5.

B. Head-on collision
For the beam-beam interaction, we assume round Gaussian beam profiles. The effect of head-on collisions is then given by with s ϵ s x s y the rms beam size at the IP, r p x 2 1 y 2 the radial distance to the origin, r p the classical proton radius, g the Lorentz factor, and N b the bunch population.

C. Long-range interactions
For the long-range interactions, we lump together the effect of all n par parasitic collisions on each side of the IP. Since they occur at a betatron phase advance close to p͞2 from the IP [6,11], the kick can be approximately expressed as a change in the IP coordinate, while the trajectory slope at the IP stays unchanged. In case of a horizontal crossing we have where u t ϵ ͓͑x 0 1 u c ͒ 2 1 y 0 2 ͔ 1͞2 and u x,y is the rms beam divergence at the IP. The effective number of parasitic crossings per side of one IP, n par , is about 18 [12]. In fact, assuming that we have a parasitic crossing every 3.75 m for the nominal bunch spacing of 25 ns, there are around 16 long-range collisions up to the separation dipole. The two remaining crossings occur inside the dipole before the beams get completely separated. In the first part of Eq. (5), we subtract the average horizontal dipole kick on the bunch, since its effect would be canceled by steering correctors. Note that the kick is the same on both sides of the IP, because the betatron phase advance of about 180 ± compensates for the opposite direction of the beam-beam separation. The vertical crossing is treated in complete analogy.

D. Triplet nonlinearities
The integrated effect of the higher-order multipoles in the low-b quadrupoles can be written in complex form. Under some simplifying assumptions (equal b functions in the two transverse planes), the nonlinear kick on the incoming side of the IP with horizontal crossing is given by where the complex coefficient G n represents the effective strength of the nth order multipole kick, through the sum ͑a n,k 1 ib n,k ͒b n͞2 k ! over the four low-beta quadrupoles on one side of the IP. The latter expression is evaluated prior to the tracking, with b k taken to be the geometric mean of the actual horizontal and vertical beta function at the center of the kth quadrupole, and a n,k and b n,k the skew and normal multipole components relative to the main quadrupole field at the same reference radius r 0 (see Tables II and III).
where l quad denotes the quadrupole length (l quad ഠ 5 m) and K 1 the nonintegrated quadrupole gradient (K 1 ഠ 0.01 m 22 ). As before, x 0 and y 0 are the trajectory slopes at the IP. Note that the dipole kick, as well as the static quadrupole and skew quadrupole components induced by the crossing angle, are subtracted, because we assume that in the real machine the changes in the closed orbit, tunes, and coupling, due to the field errors, will all be corrected by standard tuning methods. This could be a difference with respect to previous applications of similar kick-map codes for the SSC and the LHC [6,7].
Due to the simplification of equal b functions in both planes and for both sides of the IP, the geometric mean for the outgoing side is 2G n , where the opposite sign reflects the asymmetry in the polarity of the triplets before and after the IP. Taking into account the 180 ± phase advance difference between the two sides of the IP, and the opposite direction of the beam-beam separation, for the outgoing side we use the same formulas, but without the "2" sign in ͑ 1 r 0 ͒ n21 . This means that the net effect of the systematic field errors of even order n would cancel if there were no head-on collisions at the IP. Finally, the case of the vertical crossing is treated identically except that u c ͞2 is added to y 0 instead of x 0 .
In the present study, we consider a single random seed for the errors calculated according to Tables II and III. Let us finally point out that we include the errors of only one type of quadrupole (Fermilab and KEK design) for each study, in order to compare their effect, even if in the actual LHC design it is foreseen to mix the two magnet families [13].

E. Tune modulation
As a further ingredient in our simulation, a tune modulation can be added. It is described by a linear transport matrix of the form TABLE II. Harmonic multipole content in low-b quadrupoles (Fermilab design), after application of tuning shim correction [13,14]. The harmonic values are quoted in units of 10 24 of the main quadrupole field, for a reference radius r 0 17 mm. The uncertainty in the mean, as well as the estimated standard deviation, is also listed.
Order Mean Uncertainty rms Mean Uncertainty rms ͑n͒ ͑b n,M ͒ ͑ b n,U ͒ ͑ b n,rms ͒ ͑ a n,M ͒ ͑ a n,U ͒ ͑ a n,rms ͒  where Df y 2pDQ y sin͑2pf y t͒ .
Here, DQ x,y and f x,y denote the modulation amplitude and frequency, respectively, and t is the time. For instance, synchrotron oscillations and residual chromaticity result in a modulation of the betatron tune for off-energy particles. Assuming Q 0 ഠ 1, a particle at 1s d experiences an effective tune modulation of amplitude 10 24 at the 22 Hz synchrotron frequency. Ground motion, magnet vibrations, and power supply ripple may also induce tune modulation at similar amplitudes and frequencies.

F. Möbius insertion
In addition, a Möbius transformation [16] may be inserted in order to study the stabilizing, or destabilizing, properties of such a scheme. The Möbius twist is of the simple form 0 where the indices label the coordinates before and after the application of the twist. Because of the additional symmetry of a Möbius lattice, there is only one independent tune value. The two tune lines, Q 6 , are placed symmetrically above and below 0.25, where Q x and Q y denote the nominal tunes without the Möbius twist.

A. Tune footprints and diffusion maps
Frequency map analyses have long been used in celestial mechanics [8] and recently in accelerator models [9,10]. In this section we present results from an ap-plication of this technique to the beam-beam interaction, through the kick-map model discussed previously. Figure 1 presents tune footprints à la Laskar [9,10], obtained by tracking single particles over 1000 turns and, subsequently, frequency analyzing the tracking data with the SUSSIX program [17]. Through this analysis, we compute, with a very high precision [18], the fundamental frequencies of motion, for a large number of starting conditions, with initial horizontal and vertical amplitudes varying from 0 to 10s x,y , and initial slopes set to zero. By plotting the fundamental frequencies in tune space, phase-space distortions representing resonances or chaotic regions become visible.
The dramatic effect of the long-range collisions is revealed through the comparison of Figs. 1b and 1a, which show footprints obtained with and without the longrange kicks. Up to initial particle amplitudes of around 6s x,y , the effect of the head-on collisions dominates. Then, the long-range effect takes over and the frequency footprint flips, as the tune shift with amplitude changes direction. This nonmonotonic dependence of the tune with respect to the amplitude is potentially dangerous for the stability of particles beyond this limit.
The additional detrimental influence of the triplet nonlinearities can be observed in the plots Fig. 1c and 1d, where we include the effect of the errors in the quadrupoles designed in FNAL and KEK, respectively. The tune shift with amplitude gets larger and the difference resonances of order 5 and above are getting more pronounced. From this representation, however, it is difficult to say which magnet family is more harmful with respect to the beam stability.
The two remaining images [Figs. 1e and 1f] show the pure effect of the KEK triplet errors and the combined influence of triplet errors and parasitic collisions, respectively. It is clear that the long-range effect is dominant. Indeed, the tune spread induced for amplitudes up to 10s x,y is of similar size as with the head-on collisions included (compare Figs. 1d and 1f).
In all the cases, except the one with head-on collisions only, some particles diffuse out to the ͑1, 21͒ resonance, as our working point is quite close to the frequency space diagonal. An additional outcome of the frequency map analysis is displayed in the plots of Fig. 2, where we depict the variation of the betatron tunes jDQj that occurs between the first and second sets of 500 turns, as a function of the starting amplitude [9,10]. The different colors correspond to different amounts of tune variation on a logarithmic scale, extending from jDQj # 10 27 to jDQj . 10 22 . The gray regions correspond to particles with a tune variation less than or equal to the precision of the tune calculation for this number of turns. Thus, their tune variation is consistent with no variation at all, and they may be considered as completely stable. The two types of green areas are weakly unstable. We speculate that the blue, magenta, and brown regions are strongly chaotic, and that particles in the black areas also might be lost, after a larger number of turns. In these plots we can observe directly the traces of resonances which limit the region of stability. The conclusions of the previous paragraph regarding the dominant destabilizing role of the long-range collisions are also confirmed here. The additional effect of head-on collisions and triplet nonlinear fields is negligible.
A further frequency analysis was performed for a model with a reduced crossing angle of 200 mrad (Fig. 3). The two images can be directly compared to Figs. 2d and 1d, which correspond to the same model, but at the nominal Vertical Position crossing angle (300 mrad). For the smaller crossing angle, the diffusive aperture should be diminished, according to Eq. (1). Indeed, the detuning is further increased (Fig. 3a) and the particle motion is heavily perturbed at amplitudes beyond about 4s x,y (Fig. 3b). The influence of the application of a Möbius twist (8) to the phase space of the system is presented in the diffusion map of Fig. 4. In that case, tune footprints cannot be provided, as there is only one independent tune. The diffusion map now presents features which are symmetric with respect to the diagonal of the initial amplitude space. It seems that by this twist, instead of stabilizing the system, we mirror instabilities in other parts of the phase space. This was also reported in a recent experimental study [19]. An illustrative comparison of the phase space stability for all the previous cases is given in Fig. 5. We plot the average tune difference in logarithmic scale over all the tracked particles. The two curves correspond to a simple average of the tune change and another average of the tune variation normalized by the particles' initial amplitudes. The head-on case presents the smallest tune diffusion coefficient. Adding the long-range force increases the coefficient by 2 orders of magnitude. The   The blue curve corresponds to a simple average of the tune change over all tracked particles and the red one to a normalized average by dividing the tune variation of every particle with its initial amplitude.
can even distinguish a small difference between the case with FNAL and KEK design triplets. Interestingly enough, the case without head-on but with long-range and triplet errors seems more unstable: the linear tune shift due to the head-on effect puts the tune to a position further away from the dangerous ͑0, 3͒ skew resonance. As we have seen before, the situation deteriorates for a reduced crossing angle of 200 mrad, and also when a Möbius twist is included [19].

B. Action diffusion
The diffusion rate of particles in the chaotic part of the phase space quantifies the detrimental effect of the considered perturbations. Following standard ideas popularized by Chirikov and co-workers [20,21], a diffusion coefficient can be estimated by calculating the variance of the unperturbed actions for a large number of turns. This approach has already been followed by Irwin, for the study of the beam-beam interaction in the case of the SSC [6].
First, we study the evolution of the rms action spread for a group of particles, launched with random phases at the same values of transverse action of the unperturbed mo-tion. To suppress short-time fluctuations, e.g., caused by deformations of the invariant tori in phase space due to resonances, we compute the rms spread of the action average over 1000 turns. The evolution of the rms spread, for particle trajectories of 10 5 turns, with a starting amplitude of 5s x,y in both transverse planes is illustrated in Fig. 6, comparing the situation of head-on collisions only, with the cases of long-range collisions and triplet errors. In all pictures, the diffusion is quite limited. Thus, an amplitude of 5s x,y seems to lie inside the diffusive aperture. Figure 7 shows a similar picture for a starting amplitude of 6s x,y . While the head-on case looks comparable, the action spread shows notably larger variation when the long-range collisions are present. If the triplet errors are also switched on, some of the tracked particles experience a rapid growth in amplitude, leading to a steep growth in the calculated spread of action values. When the first particle is lost, we stop the simulation. This accounts for the much shorter time scale on the two bottom pictures. Finally, compared with the FNAL case, the action variance seems bigger when we include the KEK triplets and particle losses occur a little earlier, for this random seed.  More systematically, we can compute the diffusion at many different amplitudes and, for each case, compute the average increase in the rms action spread per turn. The result is illustrated in Fig. 8, which compares different combinations of head-on collisions, long-range interactions, triplet field errors, tune modulation, and Möbius twist. The tune modulation has little, if any, effect on the diffusion. The Möbius twist appears to increase the diffusion at low amplitudes, while the diffusion at larger amplitudes remains unaffected. In all cases, including long-range collisions, there is a well defined diffusive aperture, between 5.5 and 6s x,y , beyond which the motion is unstable. In the case of the SSC, where both the long-range bunch separation (7.5s x,y ) and the bunch population were smaller, the equivalent limit without tune modulation was found to be between 6 and 7s x,y [6], a similar number. Figure 9 shows the dominating effect of the long-range collisions: it compares the diffusion generated by head-on plus long-range beam-beam interactions with that due to the uncorrected triplet field errors alone. The latter give rise to strong diffusion at an amplitude of about 8s x,y , which is 2s larger than for the beam-beam effect. The differences between FNAL and KEK field errors appear to be marginal.
We can use the tune variation of the tracked particles in order to confirm the established diffusive aperture thresholds. In Fig. 10, we plot the tune difference versus the amplitude, averaged over all initial x 2 y amplitude ratios (Fig. 10a), and for a fixed ratio f ϵ arctany͞x 45 ± corresponding to particles with equal initial positions (Fig. 10b), for the same cases as in Fig. 1. We also mark two thresholds corresponding to the precision boundary and to an empirical crude loss boundary for tune changes bigger than 10 24 . For all the cases where long-range collisions and triplet field errors are included, the loss boundary is located at the same point, around 5.5s x,y . For the case where the triplet field errors are not added to the beam-beam effect, the threshold is reached a little further, around 6s x,y . The case with only the KEK triplets is clearly more stable, but indeed there is still a visible effect for larger initial amplitudes. No effect whatsoever can be observed for the case with only the head-on effect included, where the tune variation is very close to the precision limit of the method.
The complementary picture (Fig. 10b) for a fixed initial amplitude ratio gives the same qualitative information regarding the dynamical influence of the various perturbations included in the model. The fluctuation of the tune variation with the initial amplitudes is due to the presence of some high order resonances which are identified and indicated in the plot.
In the following, we study the dependence of the diffusion on various parameters. Unless stated differently, we use a model which includes head-on and longrange effects in addition to the KEK triplet errors. For simplicity, we choose a fixed launch amplitude of 5s x,y in both transverse directions, which is close to the limit of stability, for this model. and 4th integer resonances. Figure 11b shows the result of another tune scan performed parallel to the tune space diagonal, i.e., at a constant distance to the coupling resonance. This scan indicates that the nominal working point is close to optimal.

Phase advance
We have calculated the rms action variance as a function of a split in the horizontal phase advance between the two half arcs, with the nominal horizontal tune, while the vertical phase advance is kept constant. Figure 12 demonstrates that a difference in horizontal phase advance has little effect on the diffusion rate. Thus, the weak-strong beam-beam interaction does not constrain this parameter, which might be adjusted for optimum chromatic correction or for minimizing strongstrong beam-beam effects.

Offset
We have also investigated the effect of a transverse offset between the two colliding beams. Figure 13 presents the change of the action spread per turn as a function of the horizontal amplitude for different values of horizontal separation at one of the two IPs. There is no noticeable effect for an offset varying from 0.2s x,y up to 2s x,y , which is an indication that the long-range effect remains dominant. initial amplitude for four different bunch populations. It is obvious that bigger currents will reduce the diffusive aperture. The dependence of the difference between the diffusive aperture and the beam separation at the parasitic collision points on current follows the square root law found for the SSC, which is also expected from a simple scaling argument for a long-range force of the form 1͞r [6]. Figure 15 shows a scan of the diffusion rate versus the crossing angle. For crossing angles smaller than 300 mrad, the diffusion increases by many orders of magnitude. The design crossing angle thus appears to be a good compromise.

Crossing angle
The tune variation averaged over initial amplitude ratios versus the initial amplitude for the two crossing angles of 300 and 200 mrad is displayed in Fig. 16. There is more than 1s difference between the locations of the empirical loss boundary for the two cases.

IV. CONCLUSIONS
We have performed a series of weak-strong beambeam simulations for the LHC. The simulation model is similar to the approach followed by Irwin for the SSC [6], studying the diffusion in the action variable of a group of particles launched at the same transverse amplitude with random betatron phase. We added some new features, such as a Möbius twist element, and the tracking data was further processed by a frequency map analysis.
Preliminary simulation results indicate that the stability of particle motion is completely determined by the longrange beam-beam interaction, which causes substantial diffusion at amplitudes beyond about 6s x,y . If triplet nonlinearities are also taken into account, unstable particles at these amplitudes can be lost within a few 10 000 turns, while without triplet errors no particle loss is observed within the first 10 5 turns, due to the nature of the beambeam force, which decreases at large amplitude. In the presence of long-range collisions, the simulation results for the FNAL and the KEK triplet errors are almost indistinguishable in the action and small in the average tune changes. The uncorrected triplet field errors alone cause a strong diffusion at a threshold amplitude of about 8s x,y , i.e., 2s larger than for the long-range collisions. Compared with both long-range interactions and triplet errors, the effect of the head-on collisions is negligible. A tune modulation of amplitude 10 24 at 22 Hz has only marginal effect on the diffusion. Equally small effect has a transverse offset between the two beams in one of the head-on collisions. With long-range collisions present, the diffusive aperture of about 6s x,y is 3.5s smaller than the beam-beam separation at the parasitic collision points. The distance between the diffusive aperture and the beambeam separation increases as the square root of the bunch population, in accordance with previous studies [6].
The present nominal working point in the tune diagram corresponds to a broad local minimum in the diffusion rate calculated for a starting amplitude of 5s x,y . The diffusion rate is unaffected by a difference in the horizontal phase advance between the two IPs. For crossing angles below 300 mrad, the diffusion at 5s x,y increases by many orders of magnitude. For large crossing angles, the diffusion is roughly constant. Therefore, the design crossing angle of 300 mrad appears to be optimal. Finally, the IP beta function could be squeezed down to about 0.35 m without appreciable increase in the diffusion rate at 5s x,y . For even smaller values of b x,y , the 5s x,y diffusion rate increases dramatically. At b x,y 0.25 m there is strong diffusion below 4s, if the triplet errors remain uncorrected.