Generalized threshold of longitudinal multi-bunch instability in synchrotrons

Beam stability is an essential requirement for particle accelerators. Longitudinal coupled-bunch instabilities (CBI) are driven by beam interaction with long-range wakefields induced in the resonant structures with narrow-band impedance. Single-bunch loss of Landau damping (LLD) is mainly determined by short-range wakefields excited at any geometry change of the beam pipe (broadband impedance) and leads to undamped bunch oscillations. Up to now, to define the threshold beam intensity or impedance, these two effects were evaluated separately. We developed an approach to numerically solve the stability problem in a more general case and derived a new analytical threshold. We have shown that LLD can modify the mechanism of multi-bunch instability and reduce the CBI threshold below the LLD threshold. This effect explains the existing observations in the CERN SPS and should be considered for future accelerators, such as HL-LHC, EIC, FCC, and others.

The interaction of charged particles with the accelerator environment (impedance) can result in the development of undamped or exponentially growing bunch oscillations.In the absence of synchrotron radiation damping, an accelerator design relies either on beam stabilization only by the natural frequency spread of individual particles, called Landau damping [1], or on active damping systems in addition.
It is common in beam stability analysis to separate single-and multi-bunch collective effects since they are driven by short-and long-range wakefields, respectively [2][3][4][5][6][7][8].For the longitudinal beam motion considered below, many impedance sources, related to the wake functions via the Fourier transform, can be modeled by the resonances with a shunt impedance R, quality factor Q, and resonant frequency ω r = 2πf r .The cases of Q ∼ 1 and Q ≫ 1 correspond to broadband (BB) and narrowband (NB) impedances, or short-and long-range wakefields, respectively.The synchrotron oscillations inside the bunch can be described as van Kampen modes [9][10][11][12].Their frequencies are modified by the presence of a BB impedance.At a certain intensity, the coherent mode moves outside the band of incoherent oscillation frequencies of individual particles leading to a loss of Landau damping (LLD), see Fig. 1, left.Only the reactive (imaginary) part of the accelerator impedance is responsible for this effect, but the resistive impedance is needed to drive bunched-beam instabilities.A single-bunch instability resulting from the coupling of different coherent modes usually appears at bunch intensities significantly higher than the LLD threshold and is not discussed in this Letter.Coupledbunch instability (CBI), leading to the coherent coupled motion of individual bunches, is possible if the decay time of the wakefield 2Q/ω r is longer than the bunch spacing.This instability takes place even if its coherent frequency lies inside the incoherent frequency band (Fig. 1, center).
The threshold of CBI driven by a NB impedance was accurately calculated [4] using the matrix equation derived by Lebedev in 1968 [2].The LLD threshold in the longitudinal plane was first found in 1973 [3] and it was revised recently [13], using the Lebedev equation together with the approach [12] developed originally for the analysis of single-bunch instability [14].So far, the thresholds for these two effects were calculated separately, with few examples when the CBI growth rates were found in the presence of two impedance sources [15,16].
In this work, we evaluate the CBI threshold for two impedance types using the Lebedev equation.First, we consider them separately, and then we show that, if the LLD threshold is comparable to the CBI threshold, the instability in the presence of two impedances has a significantly reduced threshold and can even emerge below the LLD threshold.The mechanism of this multibunch instability is modified as well (Fig. 1, right).We also derive an analytical expression for the general instability threshold and compare it with the results of self-consistent numerical calculations using code MELODY [17] and macroparticle tracking simulations with code BLonD [18,19].
Let us consider a beam of M identical equidistant bunches, each containing N p particles with a charge q.We will use a coordinate system relative to the synchronous particle of the first bunch with the design energy E 0 and the rf phase ϕ s0 , so that ∆E and ϕ are the particle energy and rf phase deviations.They are connected via equation, e.g., [6], where ω 0 = 2πf 0 , f 0 is the revolution frequency, β the normalized particle velocity, h the harmonic number, η = 1/γ 2 tr − 1/γ 2 the slip factor, γ the Lorentz factor with value γ tr at transition energy.In this work, we analyze a single-rf case V rf (ϕ) = V 0 sin(ϕ s0 + ϕ), with V 0 being the rf voltage amplitude, but it can be extended to other rf waveforms.Each particle performs synchrotron  6).Left: LLD with BB impedance, no instability.Center: CBI with NB impedance.Right: CBI with BB and NB impedances.Thresholds found from exact Eq. ( 5) are shown with dashed lines while corresponding approximate solutions (15), (14), and ( 17) with dotted lines.Examples for nine bunches in the scaled LHC ring with parameters from Table I and k oscillations with energy E and phase ψ , where ω s0 is the angular oscillation frequency of smallamplitude synchrotron oscillations in a bare single-rf potential U rf and ω s (E) -oscillation frequency in total potential U t , modified by intensity effects.Evaluation of beam stability is usually done in two steps.First, a stationary solution should be found, e.g., by the iterative procedure [12] for some distribution function F = F(E), bunch intensity, and impedance model.Below we consider the binomial family of the stationary distribution functions µ , where E max and ϕ max are the maximum energy and phase of the synchrotron oscillations, and µ defines the bunch shape.As a result, one obtains the line density The total voltage V t (ϕ) = V rf (ϕ) + V ind (ϕ), contains contribution from the stationary beam-induced voltage with k = f /f 0 and λ k being the Fourier harmonics of the line density normalized by M N p which are non-zero only for k = pM, p = 0, ±1, ... The synchronous phase shift due to intensity effects ∆ϕ s satisfies the relation To proceed with the second step, we assume a small perturbation F(E, ψ, t) = F(E, ψ, Ω)e iΩt to the stationary function F(E).For M equidistant bunches, the perturbation Fl is characterized by a coupled-bunch mode number l (l = 0, 1, ..., M − 1), defining the phase advance e −i2πl/M between consecutive bunches, and it should satisfy the linearized Vlasov equation (see, e.g., [6]) Here the perturbed induced potential Ũ l ind (ϕ, t) is defined similarly to Eq. ( 2).The solution of Eq. ( 4) must be periodic function of the phase ψ, and can be expanded in the harmonics e −imψ (m ̸ = 0).Then the corresponding harmonics of the line-density perturbation λl k (Ω), defined similarly to Eq. ( 3), are connected by the Lebedev equation [2]: with Z k (Ω) = Z(kω 0 + Ω) and the intensity parameter Note that λl k are non-zero only for k = pM + l.The detailed derivations in variables (E, ψ) can be found, e.g., in [13].Below we consider a stationary case with ϕ s0 = π.The matrix elements G k ′ k (Ω) in Eq. ( 5) are and they contain functions The coherent frequency Ω is the solution of Eq. ( 5) when its determinant is zero.The beam is unstable if ImΩ < 0. The exact solution in a general case has to be obtained numerically, e.g., with code MELODY.Approximate analytical solutions can be also found under certain conditions.Using a general matrix property where tr(X) is the trace of a square matrix X, and taking into account that X(ε) = X(0) + ε(dX/dε)(0) + ..., the expansion up to the first order of ε yields: Then equating the determinant of Eq. ( 5) to zero, we obtain, similarly to [13], the generalized instability threshold For the solution Ω = Ω g , valid for any impedance model, the imaginary part of the sum is zero since ζ is a real quantity.Here we will consider the sum of NB and BB resonator impedances, Z = Z nb + Z bb , described by Eq. ( 1).The small parameter ε ∝ ζ th (Ω g ), while its dependence on other parameters can be defined from comparison with the exact numerical solution [13].
Let us start with a more simple case of the narrowband impedance.If Z bb = 0 and where ∆ω nb = ω r,nb /2Q nb is the resonator bandwidth, all elements in Eq. ( 5), except with k nb = ⌊ω r,nb /ω 0 ⌋, can be neglected (⌊x⌋ denotes the rounding of x to the nearest integer).Then the instability threshold is given by Eq. ( 10) with only one element in the sum.Under these assumptions, the solution is exact and can be found by analyzing the equation for ImΩ → −0, ReΩ = Ω nb = mω s ( Ẽm ), 0 < Ẽm < E max .The approximate analytic expression can be also obtained for a certain particle distribution from the stability diagrams [4].In a short-bunch approximation, when ϕ max ≪ π, we have ), and the functions I mk (E) can be replaced by Bessel functions J m (x) of the first kind and the order m The instability threshold for µ > 1 is the lowest for the dipole mode m = 1 [20] and it is mainly defined by the value of R nb /k nb , since ImZ nb = 0 at ω r .
For numerical calculations, we used the parameters of the Large Hadron Collider (LHC) at injection energy from Table I and binomial distribution function with µ = 2 and zero-intensity ϕ max = 1.3.Due to the very large number of LHC bunches (3564), a direct comparison of macroparticle simulations with the model predictions is computationally too expensive.So we simulated nine equidistant bunches by reducing the harmonic number to nine and scaling other parameters to keep ζ and ω s0 /ω 0 unchanged.
For nine bunches and a single NB resonator with k nb = 11, a coupled-bunch mode l = 2 with m = 1 should become unstable above a certain threshold.Figure 1 (center) shows the dipole bunch oscillation mode (m = 1) as a function of intensity obtained using the Oide-Yokoya method [14].A mode inside the synchrotron frequency spread becomes unstable for ζ > ζ nb th , and the instability growth rate increases with intensity.The instability threshold found as the exact solution of Eq. ( 12) coincides with the emergence of the unstable mode.The solution ( 14) obtained in the short-bunch approximation gives about a 10% higher threshold, which is not surprising since the bunch with ϕ max = 1.3 is not so short.
Assuming only BB impedance, the LLD threshold can be also derived for the same distribution functions in the short-bunch approximation.Evaluating the matrix elements G kk for Ω = Ω bb = ω s (0) and performing summation in Eq. ( 10), one obtains [13] The function   (16) with the effective cutoff-frequency number k eff which maximizes the cumulative sum in the nominator.For a broadband resonator k eff = k bb = ω r,bb /ω 0 .The results of MELODY calculations for a pure BB resonator impedance with k bb /h = 5 are shown in Fig. 1 (left).The LLD mode emerges at ζ = ζ bb th above the maximum incoherent frequency ω s (0) leading to undamped but stable oscillations.Similarly to the previous case, the LLD threshold computed from the Lebedev equation agrees with the Oide-Yokoya method.The analytical approximation overestimates it by ∼ 30%.
Finally, if we combine the BB and NB impedance contributions, the threshold (10) can be written in the form The coherent mode Ω g differs from Ω nb and Ω bb found for each impedance separately.Nevertheless, for a first estimate, one can use ζ nb th (Ω nb ) and ζ bb th (Ω bb ) in Eq. ( 17).
, Ω g ≈ Ω bb and the actual CBI threshold is reduced, as shown in Fig. 1 (right), where it is even below the LLD threshold.In the opposite case, the BB impedance has a negligible impact and Ω g ≈ Ω nb .The relative role of each contribution in Eq. ( 17) can also be seen in Fig. 2 from instability thresholds numerically obtained for different (ImZ/k) eff and k bb using Eq. ( 5).As expected, a larger BB impedance leads to a lower CBI threshold, except for the case when ζ bb th (Ω bb ) ≫ ζ nb th (Ω nb ).BLonD simulations follow closely the growth rates obtained using the Oide-Yokoya method and the direct solution of the Lebedev equation confirms the thresholds.
The impact of k eff is shown in Fig 2 (bottom).The CBI threshold is reduced for larger k bb since the LLD threshold is also reduced, and its effect on the overall instability threshold is increased.We see, however, that the growth rates above a certain intensity become smaller for larger k bb .For this new instability mechanism, the mode is localized in the bunch center (Fig 3,right), similar to a single-bunch LLD (Fig 3,left).For larger k bb , the mode spectrum shifts towards higher frequencies as for LLD mode [13] and interacts weaker with NB impedance since k nb < k bb .For Z bb = 0, the perturbation looks very different and involves mainly high-amplitude particles for the considered k nb (Fig 3,center).
Considering LHC parameters without scaling (Table I), in the full ring, the lowest CBI threshold due to a higher-order mode (HOM) of the future crab cavities [21] with R nb = 280 kOhm, f r,nb = 582 MHz, and Q nb = 1360 is three times smaller with BB impedance [23].The SPS, being the LHC injector, also provides beams for the fixed-target physics.For these beams filling the whole ring, CBI is driven by a higher-order mode (HOM) in the main 200 MHz rf system with f r ≈ 914 MHz [24].The CBI thresholds found for the full SPS impedance model [25,26] and for only HOM impedance practically coincide (see Fig. 4) since the LLD threshold is higher by an order of magnitude.If the beam fills every 5th rf bucket (LHC-type), the HOM-driven CBI threshold is higher and the resulting threshold is more affected by the BB part of the SPS impedance.This explains why the LHC-type bunch train in the SPS has a very low CBI threshold, which also weakly depends on the number of bunches in a train [23].To deliver high-intensity LHC beams, an additional 800-MHz rf system leading to the raised CBI threshold is routinely deployed [22].
To summarize, we proposed the approach to analyze beam stability in the presence of both broadband and narrowband impedance sources.The broadband impedance can significantly reduce the threshold Growth rates of the most unstable mode (solid lines) versus bunch intensity found with MELODY for SPS with M = h (blue), M = h/5 (orange) and LHC train of 72 bunches spaced by 5/f rf (green) with corresponding CBI thresholds driven by only HOM at 914 MHz (dashed lines).The LLD threshold is shown by a black dotted line.The SPS parameters from Table I, zero-intensity ϕmax = 0.9, and µ = 1.5.
of coupled-bunch instability driven by the narrowband impedance, and there is a new instability mechanism associated with it.The derived generalized analytical expression shows the key role of LLD and demonstrates how two impedance contributions add up.For the LLDdominated case, the values of the effective BB impedance and its cutoff frequency are important.The main conclusions are verified by macroparticle simulations and they are consistent with beam observations in the SPS.This understanding can help in finding mitigation measures aimed at increasing Landau damping for existing highcurrent synchrotrons.The discovered effect should be also taken into account in the design of future rings.

FIG. 1 .
FIG.1.Evolution of coherent modes (black and colored circles) and incoherent frequency bands (gray) obtained with semianalytical code MELODY versus normalized intensity parameter ζ defined in Eq. (6).Left: LLD with BB impedance, no instability.Center: CBI with NB impedance.Right: CBI with BB and NB impedances.Thresholds found from exact Eq. (5) are shown with dashed lines while corresponding approximate solutions (15),(14), and (17) with dotted lines.Examples for nine bunches in the scaled LHC ring with parameters from Table I and k nb /h = 11/9, Q nb = 100, R nb = 37 kOhm, k bb /h = 5, Q bb = 1, R bb = 30 kOhm.

3 FIG. 2 .
FIG.2.The growth rate of the most unstable mode versus the normalized intensity parameter ζ for various strengths (top) and cutoff frequencies (bottom) of the BB impedance.Instability and LLD thresholds found from exact Eq. (5) are shown with vertical dashed and dotted lines, respectively.BLonD simulations are marked with crosses.The parameters from Table I scaled to h = 9; k nb /h = 11/9, Q nb = 100, R nb = 37 kOhm, Q bb = 1, and R bb = 3A bb × (k bb /h) kOhm.Top: k bb /h = 5, bottom: A bb = 2.

FIG. 3 .
FIG. 3. Perturbed particle distribution function F(ϕ, φ) in phase space for LLD mode (BB impedance) at ζ/ζ nb th = 0.5 (left), unstable mode driven by NB impedance at ζ/ζ nb th = 1.1 (center), and for unstable mode with NB and BB impedances at ζ/ζ nb th = 0.5 (right).The dotted line is the outermost particle trajectory in the bunch.Other parameters are as in Fig. 1.