Beamstrahlung-enhanced disruption in beam-beam interaction

The radiation reaction (beamstrahlung) effect on particle dynamics during interaction of oppositely charged beams is studied. It is shown that the beam focusing can be strongly enhanced due to beamstrahlung. An approximate analytical solution of the motion equation including the radiation reaction force is derived. The disruption parameter is calculated for classical and quantum regime of beamstrahlung. The analytical model is verified by QED-PIC simulations. The model for head-on collision of long beams undergoing a number of betatron oscillation during interaction is also developed. It is demonstrated that the beamstrahlung-enhanced disruption effect can play a significant role in future lepton colliders with high-current particle beams.


Introduction
Interaction of particle flows is a fundamental problem of plasma sciences and high-energy physics. On one hand, it plays a key role in many astrophysical processes, e.g. the relativistic plasma jets are associated with gamma-ray bursts, tidal disruption events, active galactic nuclei and blazars. In the collapsar model of gamma-ray bursts [1,2] the jet interacts with star envelope formed after star collapse. Another phenomenon is development of quantum-electrodynamical (QED) cascades near the polar gaps of neutron stars [3]. The electron and positron flows produced as result of QED cascading can interact with each other in the neutron star magnetosphere [4]. On the other hand, particle colliders, which are the main research tool of high-energy physics, are based on head-on collisions of high-energy charged beams. There are several projects aiming at constructing high-energy lepton colliders with record parameters such as ILC [5] and CLIC [6]. Recently plasma-based acceleration methods have been considered as a promising direction for development of compact high-gradient linear colliders [7]. Intense electromagnetic fields can be generated at the interaction point in such colliders thereby making possible manifestation of some strong-field phenomena such as disruption [8][9][10], beamstrahlung [11][12][13], electron-positron pair production [14,15] or even effects of non-perturbative strong field quantum electrodynamics [16][17][18]. Nowadays, a primary way for exploring strong-field physics is concerned with using multi-PW laser facilities such as ELI [19], SULF [20], Apollon [21], and in future 100 PW lasers such XCELS [22], SEL [23], etc. Although, reaching large intensities poses increasingly more harsh requirements for the pulse contrast, stability, quality, etc [24]. In this connection high-current high-energy colliders which feature high beam quality and stability can provide an attractive 'laser-less' alternative for strong-field experiments. For example, the FACET-II project devoted to study of plasma-based acceleration has been discussed in this context [16,25].
In ultrarelativistic regime the dynamics of a beam particle is governed mainly by the field of the counter-propagating beam, while the field of the own beam can be neglected [26,27]. In this approximation the Lorentz force acting on the beam particle can be written as where v is the particle velocity, c is the speed of light, E and B are the electric and the magnetic fields of the counter-propagating beam, respectively, q = ±e is the particle charge, r is the distance from the particle to the beam axis, ω 2 b = 4πe 2 n b /m is the electron (positron) plasma frequency squared, n b is the density of the counter-propagating beam, m and e > 0 are the positron mass and charge, respectively. The positive sign in equation (1) refers to the case of e − e − or e + e + collisions, where the net force causes a defocusing of the two colliding beams. On the other hand, in the case of e − e + collisions the beam particles undergo transverse betatron oscillations with frequency ω b / √ γ [10,28], where γ is the Lorentz factor of the beam particle. The disruption (or pinching) time can be introduced as the time it takes the particle to reach the beam axis, and it can be estimated (up to numerical factor) as follows If the beam length σ z fulfills the condition σ z /c > T D , then the beam radii are significantly changed during the interaction. The beam distortion in the interaction region can be quantified by the so-called disruption parameter which reads for a uniform charge distribution of the beam with length σ z and radius r b [8]. Note that it is π −1/2 2 3/2 ≈ 1.6 times greater than the disruption parameter for a beam with a Gaussian charge distribution, that has the same total charge, rms length equal to σ z and rms radius equal to r b [8,28]. The expression for D can be generalized to other beam charge distributions and can be also used to characterize e − e − beam interactions. While large disruption can be desirable for increasing luminosity [29], at the same time it may lead to the increase of the background noise and interfere precision measurements. For that reason, the regime of D 1 is desirable for an experimental probing of non-perturbative strong-field QED [16]. The bending of the particle trajectory at the interaction point is accompanied by synchrotron radiation, which is known under the term beamstrahlung in the collider physics community [12,28]. The total power of the photon emission losses depends on the electron (positron) dynamical QED parameter χ [30][31][32] where α = e 2 /( c) is the fine structure constant, is the Planck constant, E S = m 2 c 3 /(e ) is the critical Sauter-Schwinger field [31], K ν (z) and Γ(z) are the modified Bessel functions of the second kind function and Gamma function, respectively [33]. In both the classical (χ 1) and the strong QED (χ 1) limits equation (4) can be reduced to simple power-law expressions If the beam length is so small, that only few photons are emitted by a single beam particle during the interaction, then the quantum nature of the synchrotron radiation should be taken into account even in the limit χ 1. In addition to beamstrahlung, other quantum effects are possible at the interaction point such as electron-positron pair production from beamstrahlung photons in strong electromagnetic fields, trident process, etc [14,34]. In general, the interplay between radiation of hard-photons and pair production can lead to a very fast growth of the total number of particles-an effect known as QED cascade, which has recently drawn a lot of attention [35][36][37][38][39][40][41][42][43][44][45][46][47][48][49][50][51]. Such QED cascades can also develop in colliding beam scenarios. It is thereby likely that the cascade development is similar to that observed in laser-solid interactions at scales much less than the laser wavelength [50,51] due to similar field configuration in the cascading regions. All in all, beamstrahlung and QED cascades may cause beam distortion due to energy depletion and in general play a negative role on clean collider operation. In the context of particle physics, colliders are therefore usually designed to mitigate beamstrahlung as much as possible. Understanding the collective effects at the interaction point is nevertheless crucial not only for optimal collider operation but also for high energy density physics. Here, the regime of beam-beam interactions with strong beamstrahlung can be exploited, for example, to produce bright gamma-ray sources or to explore strong-field QED experimentally [17,25,52].
Up to now, analytical models of the beam-beam interaction considered both disruption and beamstrahlung independently. In this paper we advance these works by focusing our study on synergic disruption-beamstrahlung effects, and we find modified expressions for the disruption parameter accounting radiation reaction in both the weak-(χ 1) and the strong-field (χ 1) limit. The coupling of these two processes comes from the fact that beamstrahlung causes loss of the particle energy and, since the disruption time is proportional to √ γ (see equation (2)), leads to a reduction of the disruption time and thus an increase of the disruption parameter D. It is important to estimate strength of this effect not only qualitatively but also quantitatively, whether small or large disruption is required by the experiment of interest.
The paper is organized as follows. In section 2 we formulate the basic equations describing the physics of the beam-beam interaction. In section 3 we adopt several assumptions to approximately solve the equations of motion and find an analytical estimate for the beamstrahlung-enhanced disruption parameter. The interaction of long beams, which may be relevant for the interaction of beams with significantly different energy densities in the center-of-mass reference frame, is discussed in section 4. The results of QED-PIC simulations are presented in section 5. They are compared with the model predictions. Applicability of our model to future colliders is discussed in section 6. Section 7 contains conclusions.
Throughout the manuscript, the equations will be given in normalized units. The normalization is mainly determined by the initial plasma frequency of the beam ω b . Then, time is measured in 1/ω b , lengths in c/ω b , momenta in mc, electromagnetic fields in mcω b /e, and power in mc 2 ω b .

General approximations
In general, the motion equations including the radiation reaction force for the ultrarelativistic electron are dp dτ where P refers to the power in normalized units. These equations describe the classical motion of the electron in an electromagnetic field with the radiation reaction effect taken into account in semiclassical approximation (with QED corrections reducing total radiation power for large values of χ) [53][54][55][56][57].
Although stochastic nature of beamstrahlung leads to individual particles being focused differently, the beam disruption relates to the diameter of the whole beam, and thus should be calculated by averaging the distance from the axis over all particles. Such kind of averaging leads exactly to the equations of motion in semiclassical approach. Below we show that results of QED-PIC simulations, where the beam collision is modeled in a self-consistent way with account of stochastic nature of quantum processes, coincide fairly good with results of our analytical model, which also justifies application of this approach. In order to analytically explore the beam disruption effect during head-on collisions of electron and positron beams, we make additional assumptions. First, as mentioned in section 1, the self-force generated by an ultrarelativistic beam can be neglected in equation (8) since it scales with γ −2 [26,27]. Second, it is sufficient to concentrate the study on the transverse dynamics of the particles located at the beam front, since they start to feel the collective force from the counter-propagating beam earlier than other particles. And third, we further restrict our analysis to the particles at the periphery of the beam, i.e. the particles which experience the largest force and thus are more likely to emit photons. As beamstrahlung leads to a decrease of the energy and thus inertia of the particles, it is exactly those particles at the periphery and at the front of the beams that are expected to experience the strongest disruption. The analysis of the motion of such particles is greatly simplified due to the fact that their dynamics is only affected by the unperturbed part of the opposite beam. Finally, we assume that the electron and positron beams have identical initial parameters, in which case the beams evolve symmetrically. In addition, the beams are considered to have cylindrical symmetry. In that case, one can write the beam density distribution as n(ξ ± , r) = n 0 η z (ξ ± )η r (r), where n 0 is the maximum beam density, ξ ± = z ± τ describes the longitudinal coordinate for beams that travel at the speed of light, and the functions 0 η r,z 1 determine the shape of the density distribution 3 . The electric field generated by such a beam is mostly transverse and can be expressed using Gauss's law as where ρ = r/r b is the transverse coordinate measured relative to the distance, r b , from the beam axis at which the electric field reaches its maximum. For electrons with v z = const = c interacting with the counter-propagating beam, ξ + is decisive and one finds ξ + = 2τ . With all that in mind and defining η(τ ) ≡ η z (2τ ), the governing equations of motion reduce to Here, a S = eE S /(mcω b ) = mc 2 /( ω b ) represents the Sauter-Schwinger field in normalized units. In the derivation of these equations we assumed that electric and magnetic components of the Lorentz force acting on the particle are almost equal to each other (hence, the factor 1/2 in equation (10) is canceled), which is valid if v z c v r and γ 1. This also allows us to assume that the radiation friction force acts mostly along the z-axis. Thus, it is not explicitly present in the equation for the transverse coordinate ρ.
As mentioned above we will be interested in particles experiencing the largest fields, i.e. particles for which the initial displacement r 0 from the beam axis equals r b and thus ρ 0 ≡ ρ(τ = 0) = 1.

Timescales
Before solving equations (12) and (13), it is useful to estimate characteristic normalized timescales present in the problem, i.e. the timescale of the electron trajectory τ D 0 and the timescale of the energy loss due to beamstrahlung τ BS where χ 0 = r b γ 0 E (ρ 0 ) /a S and γ 0 = γ(τ = 0) are the initial values of the χ parameter and Lorentz-factor of the particles, respectively. Let us also introduce a parameter κ in the following way, This parameter determines the regime of the beam-beam interaction. In case κ 1, which is explored in section 3, significant energy losses due to beamstrahlung occur on a timescale much shorter than the time it takes the particle to reach the beam axis. In the opposite limit κ 1, which is considered in section 4, it takes many betatron periods for beamstrahlung to significantly decrease the beam energy.
Utilizing the relation between γ 0 and χ 0 , and noting that P(χ) ≡ αa S Φ(χ), the parameter κ can be also expressed as follows It means that the beamstrahlung effect is determined by two initial parameters of the interaction: the beam radius r b and the parameter χ 0 . It will be shown below that these two parameters are enough to calculate the relative change of the disruption parameter caused by beamstrahlung. In the classical and QED regime equation (18) can be rewritten as follows 3. Beamstrahlung dominated regime

Constant force approximation
Obtaining the general solution of equations (12) and (13) seems infeasible, thus to make some analytical estimates first we resort to constant force approximation which corresponds to substituting the electron coordinate ρ in the rhs of equation (12) with its initial value ρ 0 = 1. In that case equations (12) and (13) take the form According to equations (6) and (7) in both classical (χ 1) and QED (χ 1) limits, the function P can be approximated as a power function of χ In that case we can obtain the solution in quadratures where ν = 2 for the classical regime and ν = 2/3 for the QED regime, P 0 = P (χ 0 ),ρ 0 =ρ(τ = 0). Let us analyze the obtained solution for a uniform beam η z = η r = η = 1 for which E(ρ) = ρ. In this case all the integrals can be calculated explicitly. In particular, we get the following solutions for γ and ρ whereρ 0 = 0 is assumed. It is interesting to note that the dependence of the electron energy on time is identical in terms of τ D 0 for both classical and QED regimes at the beginning of the interaction (0 < τ If we introduce the time τ γ after which the electron energy halves because of beamstrahlung, then this time is about 1.6 times smaller in the classical regime than in the QED regime, This is the expected result as the beamstrahlung losses according to the classical expression are greater than that according to the quantum one. This fact shows that beam collision in quantum regime is desirable if we want to suppress energy losses due to radiation (see also reference [58] for detail). The beamstrahlung-affected disruption time can be found from the condition ρ(τ = τ D ) = 0. When we use the relation D ∝ τ −2 D between disruption parameter and disruption time, and when we further consider the case that beamstrahlung sets the timescale of the interaction, κ 1, then the expression for the disruption parameter including beamstrahlung takes the form In virtue of equation (18) we can rewrite equation (3) in terms of r b and χ 0 as follows (32) Figure 1 shows that while both solutions in quantum and classical regimes describe energy loss quite well, the particle trajectory according to this solution diverges from the real trajectory quite strongly and overestimates beam disruption, thus this simple model can serve only for rough estimates of the disruption parameter, which can be sufficient in cases when only its order of magnitude is of interest.

Corrections to the model
The accuracy of the analytical model can be greatly improved by two modifications. First, we use the mean transverse coordinate in the rhs of equation (12) instead of its initial value ρ 0 , And second, we stitch the solutions in the QED and classical regimes at some time instance τ 1 at which the particle χ parameter reaches some threshold value χ 1 ∼ 1, if its initial value was large enough, i.e. χ 0 > χ 1 .
So for τ < τ 1 the equations of motion have the following solution Please note that the variableκ 0 (and alsoκ 1 , see the next equations) includes an additional factor 1/3 as compared with the definition of κ in equation (17). This re-definition simplifies later calculations. The time instance τ 1 is found from the condition From this, the ratio γ 1 /γ 0 can be found as follows where we introduced ζ = χ 1 /χ 0 . For τ > τ 1 , the classical formulas are used so that the solution of the equation of motion reads Now the disruption time can be calculated from the condition ρ C (τ D ) = 0 which can be explicitly solved, but the expression is too cumbersome to include it here. Instead, it can be estimated with a bit less bulky expression (see appendix A for derivation) For the case χ 0 < χ 1 , τ 1 ≡ 0 and χ 1 has to be substituted with χ 0 in τ 2 . Analogously to equation (18) the parameter governing significance of beamstrahlung can be expressed as followsκ = α 2 9 r b a S × ⎧ ⎨ ⎩ (μχ 0 ) 2/3 , χ 0 < χ 1 , and it can be shown thatκ 1 can be expressed in terms ofκ. The beamstrahlung-affected disruption parameter in the corrected model is (51) Figure 2 shows that calculating disruption parameter using this corrected model is significantly superior to using simple expressions (31). According to equations (47) and (48), D can be expressed as being explicitly dependent only on two initial parameters: the beam radius r b and the value of χ 0 . This allows us to scan over only a two-dimensional map of parameters to calculate the value of D from results of full 3D QED-PIC simulations.
Note that although μ should be calculated in a self-consistent way from the solution obtained above, numeric analysis shows that the value of μ is close to 0.5. So to actually find an analytical solution, we treat μ as a free parameter which we set to 0.5. This is also justified by the fact that varying μ in the range 0.3-0.7 does not significantly alter the final value of the disruption parameter.

Interaction of long beams
In this section we discuss the interaction of long uniform beams of oppositely charged particles when the number of the betatron oscillations is large σ z /(cτ D 0 ) 1 and radiative losses due to beamstrahlung are insignificant during a single period which corresponds to the limit κ 1. Such a configuration may correspond to an interaction of the beams with significantly different energy densities in a center-of-mass reference frame, resulting from either larger mass, Lorentz-factor or density of the particles of one beam compared to the other one in a laboratory reference frame. This may be the case e.g. in an electron-proton collision. In this scenario the characteristic time-scale of the evolution of a more energetic beam is much larger than that of the opposite one, thus all the particles of the opposite beam feel almost unperturbed field (in contrast to only the particles at the front of the beam, which was considered in the previous section). For simplicity, we consider the interaction of uniform beams η z = η r = η = 1, for which It is convenient to introduce the following variables where a and φ are the amplitude and the phase of the betatron oscillations (ρ = a cos φ). In the new variables, equation (21) takes the form To calculate slowly the varying component of the betatron amplitude, A = a , we average equation (54) over φ and neglect the contribution of the fast varying component of a and γ which is valid under made assumptions whereγ = γ . Introducingχ = χ = χ 0 Aγ/γ 0 one obtains the system describing the electron dynamics averaged over the betatron oscillations The system has the constant of motion In the classical limit (χ 1), one has f 2 (v) = 4f 1 (v) = P C (v)/2 and the constant of motion takes the formγ It follows from equations (60) and (63) that In the QED regime (χ 1), one has f 2 (v) = (8/3)f 1 (v) = Γ(5/6)Γ −1 (4/3)π −1/2 P Q (v) and the constant of motion takes the formγ Equations (60) and (67) then yield Figure 3 demonstrates the numerical solution of equations (12)- (14) and the analytical result given by equations (64) and (68) for γ(τ ) which are in a good agreement.

QED-PIC simulations
To confirm the prediction of the model developed in section 3. We performed 3D QED-PIC simulations using the QUILL code [59], which enables modeling of the QED effects stochastically via the Monte-Carlo  (64) and (65) and (68)   method. Choosing z as the axis of beam propagation, the simulation parameters were Δt = 0.6Δz, Δx = Δy = 2.5Δz = r b /20. For all performed simulations the resulting time-step Δt was much smaller than the average delay between consecutive QED processes, i.e. WΔt 1 where W is the total probability of some QED process (emission of the gamma-quant or birth of the electron-positron pair). A hybrid FDTD scheme [60] was used for the numerical solution of Maxwell's equations and the Vay pusher [61] was used to push the particles. Simulations were also performed using the VLPL code [62][63][64] in combination with the dispersionless RIP solver [65]. Differences between the results of the simulations using two different codes were insignificant. Figure 4 shows an example of such a simulation (see supplemental material (https://stacks.iop.org/NJP/23/103040/mmedia) [66] for a corresponding video). It shows that at χ 0 = 10 abundant creation of secondary electrons and positrons occurs, which is an evidence of QED cascading. As this process does not affect motion of the beam particles at the front, formation and development of such cascade is not discussed in detail. Also note that development of transverse kink instability is triggered in simulations with QED processes taken into account. This is probably due to the fact that QED processes are stochastic and thus lead to perturbation of the initially symmetrical particles distribution acting as a seed for kink instability.
We performed a set of simulations with varying initial radius r b and χ 0 of the beam particles. The length of the beam was chosen in such a way that for each simulation the uncorrected disruption parameter, i.e. D 0 , was equal to 10. This was done in order to have a clear way of determining pinching time used in  (12) with constant value of γ, which corresponds to neglecting beamstrahlung, and blue line represents approximate analytical solution (36) and (37) and (41) and (42). Different rows correspond to different initial parameters r b and χ 0 . calculation of the disruption parameter. In this regard our QED-PIC simulations do not represent any particular strong-field experiment as the latter requires a very short interaction time in order to mitigate large radiative losses [16]. Instead, QED-PIC simulations were used as a mean to solve the particles motion equations in a self-consistently calculated electromagnetic field and with account of stochastic nature of QED processes. For each simulation we were tracking several hundreds of the particles located at the front and periphery of the electron beam to calculate the mean time of crossing the beam axis. Examples of individual tracks, numerical solution of the equations (12) and (13) and approximate analytical solution are shown in figure 5. For each pair of values r b and χ 0 we performed a simulation with QED processes (Breit-Wheeler and nonlinear Compton scattering) and a reference simulation in which these processes were turned off artificially. By comparing the mean disruption times in these two simulations we are able to calculate the disruption parameter for a wide range of parameters. A map of the value of D/D 0 obtained from QED-PIC simulations is given in figure 6, together with the estimate from equations (46)-(48) (in which we used μ = 0.5, χ 1 = 1) and with the result from the numerical solution of equations (12) and (13).
It is important to note that for large values of r b and χ 0 we introduced a different numeric criterion for the calculation of disruption. The region of parameters where this criterion was used is marked with the red border in figure 6. This was done due to the fact that in such simulations energy loss due to beamstrahlung is so strong that after some time, beam particles are no longer relativistic, and their longitudinal velocity becomes comparable to their transverse velocity so eventually the particles stop their directional motion and start spinning without crossing the beam axis (see supplemental material [66]). In such cases instead of the time of reaching the beam axis, we used the mean time it takes the particles to decelerate down to 0.5c. As our analytical model assumes that longitudinal velocity is always larger than the transverse one it cannot be applied in these cases. Exploring beam collisions in this region of parameters requires a separate dedicated study. Our analytical model predicts that the ratio D/D 0 does not explicitly depend on the energy of the particles. We performed several QED-PIC simulations with different particle energies but the same values of χ 0 and r b . The simulation results show, that as long as the particle energy is large enough to remain ultrarelativistic before reaching the beam axis, the resulting ratio D/D 0 does not significantly depend on the particle energy.
We did not perform QED-PIC simulations of the beam-beam interaction in the regime when beamstrahlung takes many betatron oscillations to significantly decrease particles energy (κ 1) due to several reasons. First, such simulations would take significantly more time. And second, as this regime is mostly related to the interaction of beams with significantly different energy densities (in center-of-mass reference frame) during which the more energetic one does not deform much, such interaction can be sufficiently simulated by a single particle dynamics in an undisturbed field of the slowly evolving (energetic) bunch, which was done in section 4.

Discussion
The scheme, discussed in reference [16], requires a controlled beam collision, corresponding to disruption parameter being very small. This requirement poses tight restriction on the beam length and diameter. However, as we have shown, a common way of calculating disruption does not account for radiation reaction. At the same time it is useful to also have an analytical estimate to verify, whether a certain set of beam parameters meets the requirement of small D. Our analytical model shows that for beam parameters necessary to reach value of χ ∼ 1600 at the beam energy 125 GeV, total charge 1.4 nC and radius 0.01 μm, the disruption enhancement due to beamstrahlung reaches 55%. For upcoming colliders CLIC and ILC, on the contrary, radiation reaction may slightly soften requirements for the beam parameters to reach desired luminosity at the interaction point, as it is partially achieved via beam flatting.
While in the derivation of our analytical estimate for the ratio D/D 0 we considered cylindrical beams, the results can be applied to flat beam configurations, proposed for using at CLIC and ILC facilities. For particles aligned with the axes of the flat beam, their motion remains planar and thus equations (12) and (13) remain relevant. Thus, by calculating values of χ 0 and r b with respect to the beam charge distribution (see e.g. reference [66] for distribution of the electromagnetic fields for elliptical Gaussian charge distribution), we can calculate disruption parameter for a specific axis using our analytical model. Carrying out this process indicates that using a flat beam reduces effect of radiation reaction on the beam disruption, which is the main reason of using such beams. In particular for the expected beam parameters at CLIC we obtain that the disruption parameter increases around 35% for the longer axis and only 5% for the shorter one. For the round beam with the same total charge and equal area of the cross-section the beamstrahlung enhancement is around 35% for both axes. For ILC parameters disruption enhancement does not exceed 5% for both axes due to quite a low value of χ. These calculations serve as a sanity check of our analytical model. Figure 7 demonstrates how radiation reaction affects disruption for different beam parameters. In particular, it shows that collision of beams with quite a large total charge (>10 nC) with small radii can be significantly altered by beamstrahlung. Another interesting trend that can be observed is that although increasing particle energy and/or decreasing beam length (while preserving the same total charge) decreases disruption parameter, at the same time it increases significance of beamstrahlung. Overall the ratio D/D 0 can be used to determine whether radiation effects are significant or not.

Conclusion
The beam dynamics during the interaction of beams with opposite charges is studied in the beamstrahlung dominated regime. It turns out that beam radius r b and χ 0 are the key parameters determining enhancement of the beams focusing or disruption. For a uniform beam, the parameter χ 0 can be calculated from the beam density n e or total beam charge Q as follows where ε b is the beam particle energy. According to the constant force approximation model the ratio of the beamstrahlung-affected disruption parameter to the beamstrahlung-free disruption parameter can be also expressed in terms of the beam parameters. In the classical regime (χ 0 1), it reads and in the QED regime (χ 0 1) In above expressions we used the following expression for the beamstrahlung-free disruption parameter D 0 D 0 ≈ 1.8 × 10 −4 Q (nC)σ z (μm) ε b (100 GeV) r b (μm) 2 , ≈ 0.9 × 10 −4 n e (10 21 cm −3 ) σ z (μm) 2 ε b (100 GeV) . (76) A more accurate value of D can be calculated from the corrected model presented in section 3.2. The analytical model describing the interaction of long oppositely charged beams in a weak beamstrahlung regime is also presented. It is assumed that the beam particles perform a number of betatron oscillations along the beam axis. The dependence of the particle energy and its amplitude of betatron oscillations are calculated in the classical regime as well as in the QED regime of beamstrahlung. Predictions of models for both weak and strong regimes of beamstrahlung are in good agreement with results of numerical simulations.
In order to study disruption effect the dynamics of the electrons located only at the beam front were examined, although the beam-beam interaction is determined by all the particles of the beam. This is especially the case for intense regime of interaction accompanied by QED cascading. Abundant production of the secondary particles in such a cascade observed in numerical simulations leads to disturbance of the initially symmetrical distribution, which serves as a seed for development of the kink instability. This effect requires further study as it may be as limiting for collider operation as the disruption effect itself. Although to accurately describe dynamics of the beams as a whole one would need to self-consistently calculate electromagnetic field distribution, which makes the problem much more complicated for analysis.
Finally, an approximate solution of equation (A1) can be found as x = min{x 1 , x 2 }. To find τ 1 , we perform the following substitution x → τ 1 τ D 0 ,

A2. Estimation of τ 2
To estimate the disruption time from the condition ρ C (τ D ) = 0 (ρ C is defined in equation (42)), let us consider the following equation on x For large values k 3 , a rough estimate for solving this equation is For smaller values k 3 , we can first find a solution by setting k 3 = 0, i.e.
The above equation has the solution By assuming that the solution of equation (A7) is only slightly different from x , i.e. x = x + x , we can expand this equation in x and keep only linear terms: From this we obtain that x can be approximated in the following way And finally we choose the smallest of x 1,2 , i.e. x = min{x 1 , x 2 }. To find τ 2 , we perform the following substitution x → τ 2 τ D 1 , k 1 → ρ 1 , k 2 →ρ 1 , k 3 →κ 1 . (A13)