Vlasov description of the beam response to noise in the presence of wakefields in high-energy synchrotrons: beam transfer function, diffusion, and loss of Landau damping

Noise can have severe impacts on particle beams in high-energy synchrotrons. In particular, it has recently been discovered that noise combined with wakefields can cause a diffusion that leads to a loss of Landau damping after a latency. Such instabilities have been observed in the Large Hadron Collider. This paper, therefore, studies the beam response to noise in the presence of wakefields, within the framework of the Vlasov equation. First, a wakefield beam eigenmode transfer function (MTF) is derived, quantifying the amplitude of a wakefield eigenmode when excited by noise. Then, the MTFs of all the wakefield eigenmodes are combined to derive the beam transfer function (BTF) including the impact of wakefields. It is found to agree excellently with multi-particle tracking simulations. Finally, the MTFs are also used to derive the single-particle diffusion driven by the wakefield eigenmodes. This new Vlasov-based theory for the diffusion driven by noise-excited wakefields is found to be superior to an existing theory by comparing to multi-particle tracking simulations. Through sophisticated simulations that self-consistently model the evolution of the distribution and the stability diagram, the diffusion is found to lead to a loss of Landau damping after a latency. The most important technique to extend the latency and thereby mitigate these instabilities is to operate the synchrotron with a stability margin in detuning strength relative to the amount of detuning required to barely stabilize the beam with its initial distribution.


Vlasov equation
This section introduces the Vlasov formalism and notation used throughout this paper. It follows closely the instructive explanation of direct Vlasov solvers in Ref. [5], which itself is based on Ref. [9].
The beam will in general be described by a distribution function Ψ (x, x , y, y , z, δ; t), and the single-particle dynamics are governed by an effective Hamiltonian H(x, x , y, y , z, δ; t), both dependent on the phase space coordinates (x, x ) in the horizontal plane, (y, y ) in the vertical plane, and (z, δ) in the longitudinal plane. The vertical and longitudinal phase space coordinates can be expressed in terms of the corresponding action angle coordinates as y = 2J y β y cos(θ y ), y = 2J y β y sin(θ y ), where β y is an effective vertical beta function, corresponding to R/Q y0 in Ref. [5], β z is an effective longitudinal beta function, corresponding to η s R/Q s = η s v/ω s in Ref. [5], p 0 = γ m 0 v is the longitudinal momentum of the synchronous particle, R is the average radius of the circular machine, Q y0 is the unperturbed vertical tune, Q s is the synchrotron tune, ω s = ω 0 Q s is the synchrotron frequency, ω 0 is the (angular) rotation frequency of the beam around the machine, η s is the slippage factor, v is the speed of the synchronous particle, γ is the relativistic factor, and m 0 is the mass of the synchronous particle. By convention, the transverse angles are written as θ x and θ y , while the longitudinal angle is written φ. The horizontal phase space coordinates (x, x ) can be expressed as in Eq. (1) with all instances of y changed to x. The evolution of the distribution Ψ due to the dynamics described by H is governed by the Vlasov equation where the time evolution for this Hamiltonian system is governed by Hamilton's equationṡ The Vlasov equation is typically solved using perturbation theory, assuming that the Hamiltonian H = H 0 + ΔH can be written as the sum of an unperturbed part H 0 , which here will express the focusing of the particles around the design orbit, and a first-order perturbation ΔH, which expresses the weak forces due to, e.g., wakefields or noise H 0 = ω 0 (Q x0 + Q x δ)J x + ω 0 (Q y0 + Q y δ)J y − ω s J z (6) = ω 0 Q x J x + ω 0 Q y J y − ω s J z , where Q y0 and Q y are the unperturbed betatron tune and linear chromaticity, respectively, in the vertical plane, and Q x0 and Q x are the equivalents in the horizontal plane. Here, it has been expressly stated that the vertical coherent force F coh y only depends on the longitudinal position and time. Hence, this formalism cannot study effects due to, e.g., beam-beam interactions, but it is sufficient for our study of dipolar wakefields and noise. Note that our derivation considers a vertical coherent force, while all results will focus on noise and wakefields in the horizontal plane.
Within such a perturbation formalism, also the distribution Ψ = Ψ 0 + ΔΨ can be treated as the sum of an initial equilibrium part Ψ 0 and a first-order perturbation ΔΨ . The constant equilibrium distribution can be given as dependent only on the invariants of the unperturbed motion, being the actions. The dependence of H 0 on φ through Q δ is negligible [5]. Due to negligible coupling between the transverse planes and the longitudinal plane, the equilibrium distribution is separated in the transverse and longitudinal distribution functions f 0 and g 0 , respectively. The equilibrium distribution is here normalized such that 2π 0 dφ ∞ 0 d J z g 0 (J z ) = 1.
The goal of this derivation is to arrive at an expression for the distribution perturbation ΔΨ . Since ΔH depends on y and z, ΔΨ depends on θ y and φ, as well as all the actions, but not θ x . The linearized Vlasov equation can now be found by only keeping the terms in Eq. (3) that are linear in the perturbations of the distribution and Hamiltonian as where it has been assumed that the derivatives of ΔH with respect to the longitudinal action and angle are negligible. Because the coherent force is only in the vertical plane, the horizontal coordinates do not enter at this point. To solve for the distribution perturbation, we first assume that it consists of a single oscillation mode with a complex angular frequency Ω close to ω 0 Q y0 . By making a Fourier expansion in φ, it can be shown that the distribution perturbation is given by (Eq. (68) in Ref. [5]) where Q y z/Q s β z is the headtail phase factor due to chromaticity. The transverse part of the distribution perturbation is fully specified. The longitudinal distribution perturbation consist of the Fourier sum over the so far unspecified functions R l (J z ), which can be found by solving the linearized Vlasov equation that now takes form (Eq. (69) in Ref. [5]) In the following, we will interchangeably use the radial coordinate r z = √ 2J z β z to simplify the notation and analytical calculations. This change of coordinate is not used when expressing Hamilton's equations; hence, it does not need to be a canonical transformation [10].

Beam transfer function
A transfer function of a system represents the relationship between the output signal and the input signal. The transverse beam transfer function (BTF) typically represents the relationship between the transverse oscillation amplitude of the center of mass of the beam (output) and the excitation signal provided by magnets (input) [6][7][8]. The excitation signal can either be broadband or consists of a single frequency.

Beam transfer function with chromaticity
First, we want to derive an expression for the BTF including chromaticity. We consider the coherent force F coh y to be an external harmonic driving force of amplitude A y and frequency Ω Putting it in Eq. (13) and using the Jacobi-Anger expansion (Eq. (8.511.4) in Ref. [11]) give where J l (·) are the Bessel functions of order l. By equating terms with e −ilφ , one gets an expression for the longitudinal modes Inserting this into the expression for the distribution perturbation in Eq. (12) gives where we have introduced a weak transverse detuning that is a function of the transverse actions as in Ref. [6], such as the detuning driven by octupole magnets [12] Q y (J x , J y ) = Q 0y + a y J y + b y J x .
We are now in a position to calculate the BTF, defined as the beam response to the driving force. By realizing that only terms with k = l do not vanish when integrating ΔΨ over φ and that the equilibrium distribution is normalized as in Eqs. (9) and (10), one finds where y dip is the dipolar moment of the bunch. The dipolar moment is scaled by the beta function β y so that it is given in terms of y , as seen in Eq. (1). The force is scaled by the momentum p 0 , as in Eq. (34) in Ref. [5], to be the change of y . The BTF can be put in the form with the familiar dispersion integral [6] and a weight that depends on chromaticity In the absence of chromaticity, w 0l is 1 for l = 0 and 0 otherwise, returning T 0 as the BTF. For a Gaussian distribution with variation σ 2 z = J z β z and use of Eq. (6.633) in Ref. [11], we recover the result in Ref. [8] where I l (·) are the modified Bessel functions of the first kind.

Wakefield beam eigenmode transfer function
In Sect. 3.1, we considered the coherent driving force in Eq. (14) with beam-independent amplitude. Now, we will study the driving by the coherent wake force F wake y , which depends on the transverse beam motion. In the following, we adopt the weak wakefield approximation, assuming that the beam modes of different azimuthal number l are linearly independent. To distinguish between different modes, we introduce the subscript m, such that mode m has azimuthal number l m and so forth. Note that different modes m = n can have the same azimuthal mode number l m = l n , but a particular mode has only a single l m . For ease of notation, the longitudinal polar radius r z = √ 2J z β z is used instead of the longitudinal action. The wake force due to the passings of a beam eigenmode at all turns k can be given as (Eq. (85) in Ref. [5]) where Z y is the impedance function, i.e. the Fourier transform of the wake function, and N is the number of particles in the bunch. With this expression for the coherent force, Sacherer's integral equation [13] takes the form (Eq. (88) in Ref. [5]) where it has been used that the modes have a single azimuthal mode number l m . This set of equations can be solved numerically with direct Vlasov solvers, such as DELPHI [14], to find the transverse eigenmodes ΔΨ m and their frequencies Ω m . These mode details can also be found using codes that implement the circulant matrix model, such as BimBim [15]. Both these codes have been used in this paper. The wake force due to a dipolar eigenmode is often expressed as proportional to the dipolar transverse offset, F wake y / p 0 = −2ΔΩ y /β y [6], where ΔΩ is the tune shift driven by the wakefields. The factor β y / p 0 occurs here as in Eq. (19) to properly scale the force to the change of position, as in Ref. [5]. Here, we aim at generalizing this proportionality of the wake force for an eigenmode with transverse offset dependent on the longitudinal phase space. To that end, we calculate the transverse offset of the mode as a function of the longitudinal coordinates where it has been used that the average offset due to the equilibrium distribution Ψ 0 in the numerator is 0 and the integral over the wakefield eigenmode ΔΨ m in the denominator is 0. In general, there does not exist a scalar proportionality constant between F wake y m (z; t) and y(r z , φ; t) m at any given instant of time, since the former depends only on z, while the latter depends on both r z and φ. However, by combining Eqs. (24)-(27), one can find that where the time average of the ratio is taken over one synchrotron period τ s = 1/ω s . This time averaging corresponds to an averaging over the longitudinal phase φ, which also enters in z as given by Eq. (2). Therefore, we suggest instead an effective force, as in Ref. [3], that is proportional to the transverse offset, including its dependence on the longitudinal coordinates A noise force acting on the beam is now introduced in addition to the wake force Since the noise does not depend on the beam, it is not expected to change the shape of the wakefield beam eigenmodes ΔΨ m . However, it will excite more than one mode at once. Hence, we can no longer consider a single wakefield beam eigenmode. Instead, we get where we have introduced as in Ref. [3] the time-dependent mode amplitude χ m (t) and normalized constant eigenmode shape functions m m (r z , φ), which by comparison to Eq. (27) can be found to be where A m is a normalization constant such that where the superscripted * implies a complex conjugation and it has been used that g 0 is normalized, as given by Eq. (10). Note that the average over the longitudinal coordinates, denoted by · z , includes the unperturbed longitudinal distribution g 0 in the integrand. These excited modes are not necessarily orthogonal, but are still assumed to be linearly independent and can therefore be treated independently.
The noise can be modeled as a single kick per particle per turn. The noise signal ξ(z; t) can be decomposed in orthonormal noise shape functions Ξ i (z), which are polynomials of order i, being Ξ 0 = 1, Ξ 1 = z/σ z , and higher-order polynomials for a longitudinally Gaussian distribution. These noise shape functions can themselves be decomposed in the eigenmode shape functions m m where the noise signals ξ i (t) are assumed to only have power in positive frequencies ω ∈ [0, ω 0 /2) and the noise-mode moments η i m are the projection coefficients of the noise shape function Ξ i on the mode shape function m m . The noise shape functions are normalized over the initial equilibrium distribution function as where δ i, j is the Kronecker delta andη j m ≡ m * m Ξ j z is equal to η j m if the modes are orthogonal. If not, η j m can be found by solving a linear system of equations. Most noise sources produce noise that is constant across a bunch, being proportional to the constant Ξ 0 (z) = 1. This type of noise will be referred to as rigid or dipolar noise. The shape of crab cavity amplitude noise, on the other hand, is proportional to the longitudinal position, Ξ 1 (z) ∝ z. To the authors' knowledge, there are no significant noise sources for which higher-order terms in z are necessary.
The full coherent force can now be decomposed as: where the coherent force F coh y m proportional to mode m has been defined. To get the transfer function of mode m, one must now calculate the expected excitation due to noise at a given driving frequency Ω. The mode is then forced to oscillate at the driving frequency, even if it is different from the natural eigenfrequency of a mode Ω m . In general, the closer the driving frequency is to the eigenfrequency, the larger the amplitude of the modes motion is expected to be. By combining Eqs. (12) and (13) with the assumption of a single azimuthal mode number l m per linearly independent mode, meaning that the sums over l are exchanged for a single term, one gets the expression where we have again introduced a weak transverse detuning. The average offset of mode m at frequency Ω can be calculated with the first line of Eqs. (27) and (38) where the dispersion integral T l m (Ω) in Eq. (21) has reappeared. By rearranging the terms, we get an expression for the wakefield beam eigenmode transfer function (MTF) due to the excitation by noise at frequency Ω where the notationf is an important quantity in beam stability, often referred to as the stability diagram. Examples of the stability diagram are illustrated in Fig. 6. A mode is considered to be stabilized by Landau damping when its complex tune is inside (below) the stability diagram [6]. Note also the appearance of η i m in the summand in the denominator on  37) can with this MTF be expressed as:

Beam transfer function including wakefields
The MTF for a single linearly independent beam mode with a single azimuthal mode number l m is given in Eq. (40). To express the measurable dipolar BTF, the MTF of all modes must be combined. Since the BTF typically is measured by controlled excitation by a dipolar noise source, which excites all particles in the bunch equally, we will in this section assume that ξ 0 (t) is the only noise signal of nonzero amplitude. The dipolar mode moment can be calculated as where it was used that Ξ 0 (z) = 1. By combining with Eq. (40), and assuming only dipolar noise, we get the desired expression for the BTF including wakefields where F noise y is defined in Eq. (30) and T wake l m (Ω) is defined in Eq. (40). In the absence of chromaticity and wakefields, Eq. (43) returns the standard BTF, T 0 (Ω). In the limit of negligible wakefields, but still including chromaticity, such that T wake l m = T l m , we can by comparison to Eq. (20) calculate the weight of sideband l as The weight w MTF 0l is, in the limit of negligible wakefields, equivalent to the weight w 0l in Eq. (22), which was derived in the presence of chromaticity but completely without wakefields.
So far, the derivation has not made any assumptions as to the orthogonality of the modes, leading to the distinction betweenη and η in Eq. (44). If the modes are orthogonal, one hasη * 0 m η 0 m = |η 0 m | 2 . This has been found to be true within the numerical accuracy of BimBim for the LHC.
The weights w 0l in Eq. (22) and the weights w MTF 0l in Eq. (44), assuming orthogonal modes and using negligible wakefields in DELPHI, are compared in Fig. 1a. This figure shows how the dipolar moment of modes with nonzero azimuthal mode number l m = 0 increases with the chromaticity, a behavior that is important for the mechanisms in this paper. How the weights are calculated in DELPHI is explained in Appendix A.
One can alternatively imagine a headtail BTF (BTF ht ), excited by headtail noise ∝ Ξ 1 whereupon the headtail motion is measured. This is further explained and derived in Appendix B. The theoretical weights w 1l in Eq. (B.7) are in Fig. 1b   The expressions derived for the BTF in Eqs. (20) and (43) will here be benchmarked against multi-particle tracking simulations run with COMBI [16,17]. Calculations made with BimBim, using the impedance model for the LHC at top energy in 2018, have been used to get the mode details required by Eq. (43), i.e. the noise-mode-moments η 0m and tune shifts ΔQ coh m . The bunch intensity, synchrotron tune, and root-mean-squared (rms) particle momentum spread were the same as given by Table 1, see below. The tracking simulations were run with 10 7 macro-particles for 10 5 turns, from which the BTF was found by calculating the power spectral density of the beam signal and noise using the Welch algorithm. The horizontal noise used was approximately white with an rms kick per turn of 5 × 10 −4 σ x , in units of the beam size.
The new theories and the tracking simulations, both with and without wakefields, show excellent agreement in Fig. 2. It is clearly shown that it is the wakefields in the LHC that cause the shift of the main peak (and side peaks as well) of the BTF, as has been hypothesized based on measurements in the LHC (Fig. 7 in Ref. [7]). Although it is not shown here visually, the theory in Eq. (20) is also able to explain the appearance of loops at the sidebands in the stability diagram reconstructed based on the BTF [18].

Loss of Landau damping driven by diffusion
Instabilities have been observed in the LHC that occur after a long latency, in the order of tens of minutes, after reaching top energy. The mechanism behind these instabilities has just recently been understood [3]: the interplay between noise and wakefields drives a diffusion that eventually renders the beam unstable. In Ref. [3], the expression for the diffusion was based on the MTF for a mode modeled as a single underdamped stochastic harmonic oscillator (USHO). In this section, we will derive the diffusion by use of the MTF in Eq. (40), which is based on the linear Vlasov theory. The new derivation was inspired by the second part of Ref. [19], yet relaxing strong assumptions that prevented quantitative estimates for realistic machines. Most important was to relax the rigid bunch assumption and the simplified wakefield model, which corresponded to only considering mode shape functions m m that were independent of the longitudinal coordinates r z and φ. In addition, we have extended from using a one-dimensional tune spread, setting b y in Eq. (18) to 0, to using a two-dimensional tune-spread. The new Vlasov theory will be compared to the existing USHO theory.

Diffusion driven by wakefield eigenmodes
If the vertical noise and wakefields are sufficiently weak to be modeled as a perturbation, they drive a second-order diffusion of the equilibrium distribution. This requires a limited change of the distribution over the correlation period of the coherent force. The diffusion in the vertical plane can be modeled analytically by [20,21] where the equilibrium distribution has been renamed Ψ 0 → Ψ eq since it will evolve with time. The autocorrelation of the change of action is averaged over the longitudinal coordinates, as how the transverse diffusion depends on the longitudinal action is irrelevant. An equivalent horizontal diffusion term can be added to Eq. (45) to model the diffusion driven by horizontal noise and wakefields. There is assumed no coupling between the two transverse planes in terms of diffusion.
Since the wake force in Eq. (24) is complex, it will give a complex expression forJ y by use of Eq. (4). Therefore, we calculate a real expression for the change of action aṡ where θ y is the angle of the single particle. It can be shown as in Ref.
[3] that the vertical single particle angle is Since F wake y m ∝ exp(iΩt), with Ω ∼ Ω m , and F noise is assumed to only have power in positive frequencies, only the close-to-resonant term exp(iθ y ) ∝ exp(−iω 0 Q y t) will give a nonnegligible integrated change of action. The headtail phase factor Q y z/Q s β z in the single particle angle θ y in Eq. (47) will then cancel the headtail phase factor of the coherent force in Eq. (37), contained within the functions m m in (33). Note that as we are no longer focused on the behavior of the modes, but rather on their consequences on the individual particles, we are now referring to ω y ∼ ω 0 Q y0 , the frequencies of the individual particles, instead of Ω m ∼ ω 0 Q y0 + l m ω s , the frequency of the mode.
The integrand in Eq. (46), neglecting the averaging over z, becomeṡ where ω y = ω y (J x , J y ) = ω 0 Q y (J x , J y ) is the transverse frequencies of the individual particles, which depends on the transverse actions of the particles, the subscript z = 0 in F coh y m,z=0 means that the headtail factor is canceled out, both c.c. and the superscripted * imply complex conjugation, O(0) refers to terms that later in the calculation would become 0, and the time parameter s has been exchanged for u = t − s. In the first line, O(0) refers to terms that later will become 0 due to the lack of noise in negative frequencies. In the second line, O(0) additionally refers to non-resonant terms proportional to exp[iω y (s + t)] that will average to zero. For this reason, the factors exp(iθ 0 y ) cancel. In the second line, the dependence on the particles longitudinal phase φ + ω s t, including how it evolves with time, is also taken out of the coherent force.
To calculate the diffusion coefficient driven by noise and noise-excited wakefields, as modeled by the Vlasov theory derived in this paper, the last line in Eq. (49) must be treated as given by Eq. (46), with s exchanged for u = t − s. Remember from Eq. (34) that · z means an averaging over both φ and J z . The averaging over φ gives that l m = l n , i.e., that modes close to different sidebands do not interact to drive a common diffusion since they are orthogonal. The calculation then becomes the double sum of the Fourier transform of the cross-correlation function of F coh y m and F coh y n , which is the product of their Fourier transforms where a factor 2 has come from the inclusion of the complex conjugate in Eq. (49). Note that this expression is real, as it contains the sum over pairs of complex conjugated expressions. Looking at the expression forF coh y m in Eq. (41), we find that orthogonal modes with m m m * n = 0 will not contribute to the diffusion coefficient. In the unperturbed case, all modes are orthogonal. Since we have adopted the hypothesis of weak wakefields, we assume that the contribution of the terms with m = n remains negligible compared to the diagonal terms. This hypothesis was confirmed numerically in all the cases studied. Therefore, we only keep the dominant terms with m = n. The diffusion coefficient can finally be expressed as a function of the transverse actions as The term l m ω s in the argument to ΔΩ SD l m will cancel the same term in Eq. (21). Hence, a mode with eigenfrequency Ω m close to a sideband ω 0 Q y0 + l m ω s can still drive diffusion of particles with betatron frequency ω y ∼ ω 0 Q y0 , which is far away from the sidebands. That is often the case in the LHC. Note that the noise signals ξ i are given in units of the beam size, σ y = 2J y β y , so the factors β y will cancel.
We will now compare the new Vlasov theory to the USHO theory derived in Ref. [3]. The USHO theory assumes that the motion of the unstable mode with complex frequency shift ΔΩ m (with positive imaginary part) can be modeled as a single underdamped harmonic oscillator with complex frequency shift changed to ΔΩ LD m (with negative imaginary part) when stabilized by Landau damping. If Im ΔΩ LD m > 0, the Landau damping is not sufficient to stabilize the mode. The interpretation of ΔΩ LD m is further explained in Ref. [3]. Assuming small tune shifts, the MTF in the USHO theory iŝ where Δω y = ω y − ω 0 Q y0 is the frequency shift of the individual particles. Equation (53) can be compared to Eq. (40). All consequential differences originate in the different expressions for the MTFs. Assuming orthonormal mode shapes m m , the diffusion coefficient due to both noise and wakefields according to the USHO theory (only diffusion driven by wakefields was treated in Ref. [3]) can in the notation and framework used in this paper be written as Note that the noise was assumed dipolar and white in Ref. [3], giving |ξ 0 | 2 = σ 2 ξ 0 f rev for ω ∈ [0, ω 0 /2) and zero for other frequencies and noise signals. Here, σ ξ 0 is the rms amplitude of the noise kick per turn in units of the beam size. In addition, normalized units are used in Ref. [3], corresponding to β y = 1. A factor β y has therefore been added in the denominator of Eq. (54) to cancel that of the noise signal squared.
The derivation modeling the mode as an USHO made stronger assumptions on the dynamics and is therefore assumed to be less accurate than the new Vlasov-based theory. Even so, the only difference between the expressions is that the MTF in the Vlasov theory depends on the distance from the stability diagram to the undamped tune shift, ΔΩ m − ΔΩ SD l m , while the MTF in the USHO theory depends on the distance from the single-particle betatron frequency to the damped mode frequency, ΔΩ LD m − Δω y . Separately, a great challenge with using the USHO theory, was the need to calculate ΔΩ LD m , which in Ref. [3] introduced additional assumptions and inaccuracies through the use of a Taylor expansion. Hence, the MTF and diffusion coefficient derived in this paper are both assumed closer to reality and do not require further assumptions to be calculated.

Comparison to multi-particle tracking simulations
The expressions derived for the diffusion coefficients will here be benchmarked against multi-particle tracking simulation run with COMBI. Since diffusion is a second-order effect, more simulations are needed to get an accurate estimate of the diffusion coefficient that depends on the spread of action ΔJ 2 . 24 tracking simulations have been run with 10 7 macro-particles for each configuration. Only horizontal rigid white noise ∝ Ξ 0 is modeled. The diffusion is presented in units of D 0 = σ 2 ξ 0 /2β x , the diffusion expected due to noise in the absence of wakefields and a feedback system.
First, the diffusion due to a single mode with η 0m = 1 is studied. The wakefields are in COMBI modeled by an anti-damper, i.e., a kick that is constant over the bunch length and that is proportional to the average position offset. This simplified model allows for an easy control of the corresponding tune shift in the numerical simulations. The numerical diffusion coefficient is in Fig. 3 compared to the Vlasov theory in Eq. (52) and the USHO theory in Eq. (54). In Fig. 3a, the tested mode has a tune shift qualitatively similar to that of the most problematic modes in the LHC, with a significant negative real part relative to the imaginary part. The two theories are close to identical and agree well with the simulations. In Fig. 3b, the tested mode has a tune shift such that the Taylor . The stability margin in octupole current is 25% in both configurations. The red shaded area corresponds to 1 standard deviation from the solid red line expansion used to get ΔQ LD m in the USHO theory is inaccurate. The Vlasov theory is found to be superior to the USHO theory for this mode. For this setup with one rigid mode, the Vlasov theory in Eq. (52) is identical to its inspiration in Ref. [19] when using a two-dimensional tune spread.
Second, the diffusion due to the full LHC wakefields will be studied. Calculations made with BimBim, using the impedance model for the LHC at top energy in 2018, have again been used to get the mode details. The bunch intensity, synchrotron tune, and rms particle momentum spread were the same as given in Table 1. The numerical diffusion coefficient is in Fig. 4 compared to the Vlasov theory in Eq. (51). There is good agreement both with and without a feedback. Note that the diffusion coefficient with a feedback is the sum of a multitude of D Vlasov m , driven by individual modes, which are shown individually in Fig. 5a. The sums for various stability margins are shown in Fig. 5b.

Distribution and stability evolution
The diffusion driven by noise and wakefields causes an evolution of the distribution, which leads to a change of the stability diagram since the latter is inversely proportional to the dispersion integral in Eq. (21). We consider here an illustrative example for the horizontal plane in the LHC at top energy in 2018 that is also considered in Ref. [22], with most parameters given in Table 1, a ratio between the detuning coefficients of b x /a x = −0.7, and a detuning coefficient a x that is 25% larger than the initial stability threshold for a bunch with a Gaussian transverse distribution. This is the same case as is tested in Fig. 4b. The diffusion coefficient is a sum of a multitude of D Vlasov m that are shown individually in Fig. 5a. The modes that are inherently stable also without Landau damping, in part due to the feedback, are bowl-shaped, i.e., zero at a center tune and increasing on both sides. This is qualitatively the same shape that is found in Ref. [23] for a bunch driven by noise and damped by a feedback system, in the absence of wakefields. It is found in Ref. [23] that diffusion coefficients of this shape are not critical for the loss of Landau damping that is observed in the (a) (b) Fig. 5 Diffusion coefficient with the LHC wakefields, using the parameters in Table 1, according to the Vlasov theory in Eqs. (51) and (52). (a) displays the largest diffusion coefficients driven by individual wakefield-driven modes with a stability margin of 25%. The stability margin is a global measure, being the margin for the least stable mode, which corresponds to the most peaked line in (a). (b) displays the total diffusion coefficient for various stability margins. The green curve with 25% stability margin is equal to the blue curve in Fig. 4b and the sum of the curves in (a)   [3], based on the USHO theory. b Equal assumptions as when deriving Eq. (50) in Ref. [3] LHC. On the other hand, the modes that would be unstable if not for Landau damping are peaked and narrow as the ones in Fig. 3. It is found in Ref. [3] that diffusion coefficients of this shape lead to a local flattening of the distribution in action space and a drilling of a hole in the stability diagram until the mode is unstable. This is the cause for loss of Landau damping after a latency that has been observed in the LHC. The distribution and stability evolution according to the Vlasov theory, both driven by the diffusion of only the least stable mode and by the diffusion of all the modes, are presented in Fig. 6. These evolutions have been calculated using the code PyRADISE  Table 2, which also include the latencies found with the USHO theory numerically and analytically. The relative change of the distribution driven by only the least stable mode in Fig. 6a appears quite different from that driven by all the modes in Fig. 6c. However, what matters most is the local flattening of the distribution that is visible in both figures as a line across which particles have moved from the dark blue area on the left to the dark red area on the right. The stability diagram therefore evolves similarly due to only the least stable mode in Fig. 6b as due to all the modes in Fig. 6d. Since the additional modes drive a measurable diffusion that in this case partly counteracts the distribution evolution driven by the least stable mode, their inclusion increases the latency from 24 s to 37 s. This increase by about 50% is considered small on the vast scale of possible latencies [3], as will be exemplified in Sect. 4.4.

Dependence of latency on the stability margin
The latency depends on multiple key parameters, such as the noise and the noise-mode moments, but on none as critically as the stability margin, which is how much more detuning is used than necessary to stabilize a Gaussian beam. Here, this dependence will be illustrated by calculating the latency for the same test case as in Sect. 4.3, but with different stability margins.
The latencies calculated with the same five approaches as in Table 2 are shown in Fig. 7. There are many key results. First of all, the latency increases by more than 4 orders of magnitude by increasing the stability margin by 1 order of magnitude. Secondly, the analytical latency agrees well with the simplified simulation based on the USHO theory for all stability margins. This is expected, as the two approaches make the same assumptions. In general, the analytical latency gives an acceptable estimate also compared to the less simplified simulations, but it is slightly longer. The Vlasov and USHO theories with a single mode agree well, especially at small margins. Note that the Vlasov theory with a single mode has a numerical challenge at large stability margins that will be discussed in Sect. 4.5.
The model most true to the actual physics are the simulations based on the Vlasov theory including all the modes. This is also the most complex mechanism. The initial diffusion coefficient for different stability margins is displayed in Fig. 5b. At low margins, the peak due to the least stable mode is clearly visible, and it is this mode that eventually goes unstable. However, the diffusion due to the other modes slows down the destabilizing process. The latency with the Vlasov theory with all the modes is longer than the latency with only the least stable mode, and the difference increases with the stability margin, for small margins. Then, starting at a margin of 50% the peak due to the least stable mode is no longer visible in Fig. 5b and the total diffusion drills a wider hole. Eventually, it is actually a different mode that goes unstable, which is centered more in the middle of the peak of the total diffusion coefficient. This transition can also be seen in Fig. 7.  Fig. 7 Latency as a function of the stability margin for the LHC case otherwise described by the parameters in Table 1. The different curves correspond to different approaches (the same as in Table 2) to calculate the latency

Discussion
Several aspects of this mechanism deserve further discussion. Here, we will discuss the physics, the differences between the USHO and the Vlasov theories, and various numerical challenges with both.
It is only the Vlasov theory that can accurately model the full wake-driven diffusion. The USHO theory can model the impact of additional almost unstable modes, but not the stable ones that cause the bowl-shaped diffusion coefficients exemplified in Fig. 5a. As also discussed in Ref. [3], it is challenging to predict the impact of the secondary modes in addition to the least stable one. In some cases, as the one tested here, they will extend the latency, but they may in principle also shorten it. If this level of understanding is required, it must be evaluated on a case-by-case basis.
The Vlasov theory for the MTF has been derived in this paper. It is illustrated in Fig. 3 that the Vlasov theory is more accurate than the more simplified USHO theory derived in Ref. [3]. However, the fact that the USHO theory allows for the derivation of an analytical equation for the latency (Eq. (50) in Ref. [3]) is quite useful. As shown in Fig. 7, the analytical estimate is in general within short range of the simulated latencies, in comparison with the many orders of magnitude over which the latency varies with the stability margin. Therefore, we suggest to use the analytical latency equation to estimate how much stability margin is required to achieve a desired latency. The required margin depends in particular on the noise and wakefield eigenmodes.
An aspect, which has not been focused on here, is the fact that the diffusion modeled by both theories is the expected motion as the number of turns goes to infinity. As seen by the variability between the separate simulations in Figs. 3 and 4, the expected diffusion coefficient is clearly not achieved each turn. This was discussed in more detail in Sect. IV.C in Ref. [3]. In fact, it was found that the wakefields driven by a mode particularly close to the stability limit drove a single-particle motion more comparable to resonant motion than diffusion.
A numerical issue has been discovered with the Vlasov theory. Loops start to develop at the edge of the drilled hole in the stability diagram when there is a large initial stability margin, especially when it models only a single mode. It is visible already in Fig. 6. These loops continue to grow so much that it interrupts the drilling of the initial hole, in addition to making the physical interpretation of the stability diagram challenging. We have found that this behavior is caused by initially insignificant details in the diffusion coefficient, which begin to enhance themselves. Based on the discussion in the previous paragraph, about the diffusion being the expected distribution evolution and not the actual evolution per turn, we believe this enhancement of minor details to be a non-physical consequence of the theory. Furthermore, this behavior is not reproduced with the USHO theory. For this reason, the latency is found by extrapolating the drilling from when the initial stability margin is reduced by 75%.
A final numerical challenge that we would like to point out is the convergence with numerical parameters. By reducing the discretization step in time and action space, the latency does converge towards a single value, as expected. Similarly, the stability diagram is calculated as the one-sided limit when the small parameter → 0 + [6]. However, for → 0 − and = 0, one gets different results. In other words, the dispersion integral in Eq. (21) is discontinuous in the imaginary part of Ω, which causes numerical challenges for small positive values of . In the PyRADISE calculations presented here, we have used = |a x | /200, where a x is the in-plane detuning coefficient due to octupole magnets. With this value of , the latency and stability diagram calculation appear to have converged, but it is quite close to a threshold below which the results are inaccurate due to the discontinuity of the dispersion integral.

Conclusion
In this paper, we have studied the excitation of beams in high-energy synchrotrons by noise. Within the framework of the linearized Vlasov equation, the MTF has been derived, modeling how a wakefield beam eigenmode is excited by noise. The initial excitation by noise depends on the shape of the eigenmode and the correlation of the noise along the bunch, modeled through the noise-mode moments η i m . The future evolution of the mode, after being excited by noise, depends on the strength of the wakefields and the details of the detuning. It has here been assumed that the wakefields are weak, so that modes of different azimuthal mode numbers are linearly independent.
The BTF including chromaticity and wakefields has been derived based on the MTF. By comparing to multi-particle tracking simulations, it has been found that the new semi-analytical expression accurately models the impact of wakefields on the BTF. In particular, the shift of the BTF amplitude peaks in tune space has been thoroughly explained. The new theory can, in the limit of negligible wakefields, also be used to get an expression for the BTF including chromaticity. This expression also agrees excellently with multi-particle tracking simulations.
The main motivation for this work was to better understand the loss of Landau damping driven by noise and wakefields, which has been observed in the LHC. The MTF based on the linear Vlasov equation has been used to derive an expression for the diffusion driven by the interplay between noise and wakefields. The Vlasov-based theory can both model the diffusion driven by almost unstable modes, for which the diffusion coefficients are peaked in tune space, and the diffusion driven by inherently stable modes, for which the diffusion coefficients are bowl shaped in tune space. The peaked diffusion coefficients driven by the almost unstable modes can cause a local flattening of the distribution, leading to the drilling of a hole in the stability diagram. This mechanism can cause a loss of Landau damping. By comparing to multi-particle tracking simulations, the Vlasov theory is found superior to the previously derived USHO theory, which modeled the eigenmodes as underdamped stochastic harmonic oscillators. Nevertheless, by use of the PDE-solver PyRADISE, both theories produce similar latencies within about a factor 2, which is fairly close compared to the many orders of magnitude that the latency varies over. The main technique suggested to increase the latency and thereby mitigate the loss of Landau damping is to increase the stability margin, i.e., operating with more detuning than predicted to barely stabilize a beam with a Gaussian transverse distribution. How much margin that is needed to achieve a certain latency depends in particular on the noise and wakefield eigenmodes and can be estimated with the analytical latency equation derived in Ref. [3].
Going forward, it will be of particular interest to study operational conditions of real accelerators with the new diffusion theory for the loss of Landau damping. For instance, the impact of crab cavity amplitude noise must be studied in detail to estimate the accepted noise levels in the HL-LHC. Such a study may also find alternative operational configurations, by e.g., optimizing the linear chromaticity and gain of the transverse feedback system so that the latency is maximized. We note that the relation (8.971.5) is not valid for n = 0. Nevertheless, it can be shown that the expression obtained remains valid for n = 0 starting from Eq. (A.5), using the integral (6.631.1) in Ref. [11], identifying the relations with the confluent hyper-geometric function (9.215.1) and (9.212.1) in Ref. [11], exploiting the link between the confluent hyper-geometric function and Laguerre polynomials (8.972.1) in Ref. [11], and finally expressing the Laguerre polynomials with Eq. (8.973.2) in Ref. [11].

Appendix B Headtail beam transfer function with chromaticity
In this appendix, we will derive an alternative BTF to the one in Sect. 3.1. Here, both the external driving force and the measured response will have a linear headtail dependence, i.e., being proportional to z = r z cos(φ). The coherent force F coh y is now F coh y (z, t) = r z cos(φ)A y e iΩt , (B.1) which is qualitatively equivalent to crab cavity amplitude noise. The rms kick strength is rms(F coh y ) = σ z A y e iΩt , where σ z is the rms longitudinal spread of the particles in the bunch. Putting it in Eq. (13) and using the Jacobi-Anger expansion (Eq. (8.511.4) in Ref. [11]) give