Emittance preservation in plasma-based accelerators with ion motion

In a plasma-accelerator-based linear collider, the density of matched, low-emittance, high-energy particle bunches required for collider applications can be orders of magnitude above the background ion density, leading to ion motion, perturbation of the focusing fields, and, hence, to beam emittance growth. By analyzing the response of the background ions to an ultrahigh density beam, analytical expressions, valid for nonrelativistic ion motion, are derived for the transverse wakefield and for the final (i.e., after saturation) bunch emittance. Analytical results are validated against numerical modeling. Initial beam distributions are derived that are equilibrium solutions, which require head-to-tail bunch shaping, enabling emittance preservation with ion motion.


I. INTRODUCTION
Plasma accelerators (PAs) have received substantial theoretical and experimental interest because of their ability to produce large accelerating gradients, enabling compact accelerating structures [1]. In a PA, a short and intense laser pulse or charged-particle beam propagating in a plasma drives (via the ponderomotive force in the case of the laser, or via the transverse space-charge field in the case of the beam) an electron plasma wave (or wakefield). The plasma wave has a relativistic phase velocity (of the order of the driver velocity), and can support large longitudinal and transverse fields, suitable for accelerating and focusing a particle beam properly delayed with respect to the driver. PAs have demonstrated the production of high-quality (mono-energetic) electron bunches, and the generation of accelerating gradients in the range of 10s to 100 GV=m, several orders of magnitude larger than that obtained in conventional accelerators, which are presently limited to ∼100 MV=m by material breakdown [2][3][4][5].
The rapid development and the properties of PAs make them interesting candidates for applications to future highenergy linear colliders (LCs) [6][7][8][9][10]. Even though the properties of future LCs will be determined by high-energy physics experiments that are currently underway, it has been anticipated that a center-of-mass energy ≳1 TeV and a luminosity ≳10 34 cm −2 s −1 will be required [11,12]. Typically, this implies using electrons and positron bunches with N b ∼ 10 10 particles and normalized horizontal and vertical emittances such that ðε n;x ε n;y Þ 1=2 < 100 nm in order to guarantee a sufficiently small bunch size at the interaction point [7,10].
Preservation of the (small) bunch emittance during the acceleration process, which usually requires cascading several PA stages, is of fundamental importance to ensure the viability of a PA-based LC. Emittance preservation is ensured by matching the bunch in the plasma wave. This is possible because the transverse wakefield in a beam-driven PA operating in the blowout regime [13], or in a laserdriven PA operating in the quasilinear regime using a nearhollow plasma channel [14], varies linearly with the transverse position. Denoting by W ⊥ the amplitude of the transverse wake (the force experienced by a relativistic electron is F ⊥ ≃ −eW ⊥ , where W x ¼E x −B y , W y ¼ E y þ B x , and where E x;y and B x;y are, respectively, the transverse components of the electric and magnetic fields in the wake), we have that (an axisymmetric driver is assumed) here r ¼ ðx; yÞ represents the transverse coordinates, k p ¼ ω p =c, and E 0 ¼ mcω p =e, where ω 2 p ¼ 4πn 0 e 2 =m is the plasma frequency, n 0 is the background electron plasma density (density in the channel for the near-hollow channel case), c is the speed of light, and e (m) is the electron charge (mass). Using Eq. (1), the matched rms bunch sizes are σ x½y ¼ k −1=2 p ð2ε 2 n;x½y =γ b Þ 1=4 , where γ b ≫ 1 is the beam energy normalized to mc 2 .
The linear dependence of the transverse wakefield on the transverse coordinates relies on the assumptions that the background ion distribution is uniform and stationary. However, as the bunch accelerates and the matched beam sizes adiabatically decrease, the bunch peak density, n b;0 , increases, and so does the amplitude of the bunch spacecharge fields. As was initially pointed out in Ref. [15], when the fields of the accelerating beam become large enough so that the background ions move significantly during the bunch transit, the transverse wake can be strongly perturbed (e.g., the wakefield strength acquires a nonlinear dependence from the transverse coordinate and changes slice-by-slice along the bunch), resulting in a potentially severe degradation of the bunch emittance. This is anticipated to occur for high bunch densities, longer bunches, and lighter background ions, and can be quantified by the parameter where Z i is the ion charge state, M i the ion mass, and L b the bunch length. Ion motion and the related bunch emittance growth is potentially a serious issue for future PA-based LCs. For instance, for the LC design presented in Ref. [10], the bunch parameters are N b ¼ 10 10 , L b ≃ 20 μm, ε n;x ¼ 10 μm, ε n;y ¼ 35 nm, and the bunch energy in the first PA stage (n 0 ≃ 10 17 cm −3 , nonlinear blowout regime) is 25 GeV. The matched rms bunch sizes are then σ x ≃ 1 μm, σ y ≃ 60 nm, yielding n b;0 =n 0 ≃ 12 000, and so Γ ≃ 10 for a Hydrogen ion background. We expect ion motion to be important in this case. Ion motion is reduced for heavier background ions. However, in this case, the multiple ionization states available could lead to uncontrolled plasma formation inside the bunch, resulting in the degradation of the wakefield [15,16].
In Ref. [17], the emittance growth problem was addressed by considering an adiabatic matching section, located at the entrance of the PA, consisting of a short plasma section with a decreasing ion mass to allow for the beam to remain matched to the focusing wakefield.
In this article we analyze the interaction of a highly relativistic bunch propagating in a stationary ion background that provides focusing for the bunch. We compute the ion response to the bunch transverse space-charge field by solving the fluid equations describing the ion distribution together with Maxwell's equations for the associated electromagnetic fields. We provide, for any arbitrary bunch shape, analytic expressions for the transverse wake in the bunch region, and we compute the final (i.e., after phasespace mixing) bunch emittance. Analytical results, valid in the nonrelativistic ion motion regime where the induced ion density perturbation, δn i ¼ n i − n i;0 (n i;0 ¼ n 0 =Z i being the unperturbed ion density and n i the ion density including ion motion effects), satisfies δn i =n i;0 ≲ 1, are compared to fully nonlinear and self-consistent particle-in-cell (PIC) simulations performed with the code INF&RNO [18,19]. Finally, we derived a class of initial beam distributions that are exact equilibrium solutions, enabling ion motion without emittance growth.
This article is organized as follows. In Sec. II we derive an analytic expression for the transverse wakefield in presence of ion motion. The ion-motion-induced bunch emittance growth is discussed in Sec. III. Analytic expressions for the final bunch emittance are presented in Sec. IV. In Sec. V we discuss the equilibrium bunch distributions that enable ion motion without emittance growth. Conclusions are presented in Sec. VI.

II. TRANSVERSE WAKEFIELD IN PRESENCE OF ION MOTION
We consider a relativistic bunch with a density profile parametrized as where ζ ¼ z − ct is the comoving longitudinal coordinate (z is the longitudinal coordinate, t is the time), g ∥ ðζÞ and g ⊥ ðr; ζÞ describe, respectively, the longitudinal and the ζ-dependent transverse profile of the bunch. We assume that the bunch head is located at ζ ¼ 0 and that the bunch extends for ζ < 0. We require that, for any longitudinal slice, so that the bunch current density profile only depends on the choice of g ∥ ðζÞ, and this can be arbitrary. The ion dynamics may be described using a cold fluid model [1]. By making use of the quasistatic approximation [20], the equations for the ion density, n i , and for the normalized fluid momentum, u i ¼ γ i β i (cβ i being the ion fluid velocity, and γ i the associated relativistic factor), are, respectively, and E and B being the total fields due to the ions and the bunch. The fluid equations for the plasma are coupled to Maxwell's equations, describing the evolution of the electromagnetic fields in the wake. In particular, the amplitude of the transverse wakefield in presence of ion motion is given by The transverse current density in the regime of interest is due solely to ion motion (the contribution due to the relativistic bunch is negligible), J ⊥ ¼ Z i en i cβ i;⊥ . This current provides the coupling between ion fluid equations and Maxwell's equations.
In the nonrelativistic limit for the ion motion, i.e., jβ i j≪1 and u i ≃ β i , the equation for the transverse ion fluid momentum becomes ∂ ζ β i;⊥ ≃ −ðZ i e=M i c 2 ÞE ⊥ , the current is J ⊥ ≃ Z i en i;0 cβ i;⊥ ¼ en 0 cβ i;⊥ , and the transverse wakefield equation is where at the head of the bunch (ζ ¼ 0) the wakefield is from the uniform ions [Eq. (1)].
In general, both the ion column density and the bunch density contribute, via Gauss's law, to the transverse electric field E ⊥ in Eq. (8). However, as described in Ref. [15], ion motion is important when n b;0 =ðZ i n i;0 Þ ≫ 1, and so the bunch contribution to E ⊥ is the dominant one. Furthermore, since the bunch is (highly) relativistic, its space-charge field is mainly transverse (i.e., jE z j ≪ jE ⊥ j) and so, we have, Analytical solutions to this equation to be inserted in Eq. (8) can be obtained for a specific beam distribution (e.g., round or flat bunch distribution).
In the case of a round bunch (σ x ¼ σ y ), the transverse quantities depend only on the radial coordinate r ¼ jrj ¼ ðx 2 þ y 2 Þ 1=2 , W ⊥ ¼ W rr , and E ⊥ ¼ E rr ,r being the transverse radial versor. The bunch field is E r ¼ −½4πen b;0 g ∥ ðζÞ=r R r 0 g ⊥ ðr 0 ; ζÞr 0 dr 0 , and Eq. (8) becomes For a bunch with a longitudinally flat-top distribution of length L b (i.e., g ∥ ðζÞ ¼ 1 for −L b ≤ ζ ≤ 0, and zero elsewhere), and a transverse Gaussian profile, g ⊥ ¼ expð−r 2 =2σ 2 x Þ, the expression for the wakefield Eq. (10) becomes where HðqÞ ¼ ð1 − e −q Þ=q. The wake strength increases quadratically in ζ going from the head towards the tail of the bunch. Furthermore, as a result of the nonuniform transverse bunch profile, the wake also acquires a nonlinear dependence from the transverse coordinate. In the case of an uniform transverse profile, g ⊥ ¼ 1 for r < R ¼ 2σ x , the wake is linear in r within the bunch cross section, i.e., HðqÞ ¼ 1 in Eq. (11). Figure 1 shows the transverse lineouts of the wakefield W r obtained for the case of a bunch with a longitudinal flattop distribution with k p L b ¼ 1, and a Gaussian transverse distribution with k p σ x ¼ 0.015 (round beam) propagating in an ion channel (Hydrogen ions) for different values of the bunch peak density. Lineouts are taken at ζ ¼ −L b . Colored lines refer to (fully nonlinear) numerical modeling performed with the code INF&RNO (quasistatic modality), where the bunch and the ion distribution are modeled using a particle-in-cell (PIC) approach. Dashed black lines are the theoretical values given by Eq. (11). The analytical solution is in good agreement with modeling for n b;0 =n 0 ≲ 2000 (corresponding to Γ ≲ 1). For higher bunch densities the space-charge field of the bunch is strong enough that the ion motion becomes relativistic during the bunch transit time, violating one of the assumptions of Eq. (11). In this regime, the ion distribution collapses toward the center of the bunch, generating a high-density filament with a characteristic size much smaller than the bunch size, resulting in a major disruption of the linear confining wakefield toward the bunch core [21]. This is shown in Fig. 1 where we plot (orange dashed line) a transverse lineout of the normalized ion density (n i =n 0 ) at ζ ¼ −L b for n b;0 =n 0 ¼ 8000.
In the case of a flat bunch (σ x ≫ σ y ) and within the bulk of the bunch cross section, we have ∂ x E x ≪ ∂ y E y , and so the components of the fields of the bunch are given by E x ≃ 0, and E y ≃ −4πe R y 0 n b ðx; y 0 Þdy 0 . In this case, the wakefield in the y direction is affected by ion motion, while the one in the x direction is essentially unperturbed, namely W x =E 0 ¼ k p x=2. From Eq. (8), assuming, as before, a longitudinal flat-top distribution and a transverse Gaussian distribution, where . We see that in the nonsymmetric case the degree of perturbation of the wakefield is approximately twice as large compared to the symmetric case Eq. (11).

III. EMITTANCE GROWTH FROM ION MOTION
For a bunch initially matched in the linear (unperturbed) wakefield, the modification of the transverse wake due to ion motion results in the growth of the total projected rms bunch emittance. Figure 2(a) (black curve) shows the evolution, as a function of the (normalized) propagation distance k β z [k β ¼ k p =ð2γ b Þ 1=2 being the betatron wave number], of the normalized emittance for a bunch in an ion column (no acceleration) with parameters taken from the TeV-class LC presented in Ref. [10], namely N b ¼ 10 10 , L b ≃ 20 μm (we use a flat longitudinal current profile), and γ b ¼ 49 000 (first PA stage). For simplicity, we consider a symmetrized bunch with initial emittances ε n;x ¼ ε n;y ¼ ε n;0 ¼ 0.6 μm (this value is the geometric average of the emittances in Ref. [10]), yielding a matched rms bunch size σ x ¼ σ y ¼ 0.245 μm. The transverse bunch profile is Gaussian. The background Hydrogen ion density is n i;0 ¼ n 0 ¼ 10 17 cm −3 , yielding n b;0 =n 0 ≃ 12 000, and so Γ ≃ 10. The modeling shows that in a propagation distance k β z ≃ 4 (corresponding to ≃2 cm), the projected bunch emittance increases by ∼20%. Furthermore, the slice-dependent nature of the wake perturbation causes the final (i.e., after saturation) bunch emittance to be slicedependent (i.e., we expect the final emittance of a bunch slice located towards the tail of the bunch to be higher than that of a slice towards the head of the bunch). This is shown in Fig. 2(b) [black line], where we plot the bunch slice emittance for k β z ¼ 12. A slice-dependent emittance can pose major practical issues for bunch transport and manipulation (e.g., transport between PA stages, final focus). For instance, staging of PAs to reach high energy [22] requires coupling the bunch extracted from a depleted PA stage to the entrance of the subsequent PA stage, using some focusing optic [23]. Conventional focusing optics, which, for any fixed distance, are designed to focus bunches with a given value of the emittance and energy, a slice-dependent emittance can result in an imperfect coupling (mismatch) of the bunch to the entrance of the PA stage.

IV. EXPRESSION FOR THE ION-MOTION INDUCED EMITTANCE GROWTH AT SATURATION
Analytical expressions for the final (saturated) bunch emittance can be obtained by analyzing the transverse dynamics of the particles in the wakefield given by Eq. (8). For simplicity, we will neglect acceleration (i.e., γ b is (b) Slice emittance measured for k β z ¼ 12 (saturation) for the three types of bunches in (a). The bunch parameters are N b ¼ 10 10 , L b ≃ 20 μm (flat longitudinal current profile), γ b ¼ 49 000, ε n;x ¼ ε n;y ¼ ε n;0 ¼ 0.6 μm. The background Hydrogen ion density is n i;0 ¼ n 0 ¼ 10 17 cm −3 . constant), and we assume that the transverse motion satisfies juj ≪ γ b , u ¼ ðu x ; u y Þ being the transverse electron momentum normalized to mc. We will consider a bunch with a uniform longitudinal current profile of length L b , and with an initial Gaussian transverse phase-space distribution uniform in ζ, namely where σ u x [σ u y ] describes the rms horizontal [vertical] momentum distribution. Generalization to other longitudinal current profiles is straightforward. The bunch is initially matched in the unperturbed (linear) wakefield if σ u x ½u y ¼ ðγ b =2Þ 1=2 k p σ x½y , and the initial bunch emittances are then ε n;x½y ¼ σ x½y σ u x ½u y ¼ ðγ b =2Þ 1=2 k p σ 2 x½y . We also make the assumption that, during the ion motion-induced phase mixing, the structure of the transverse wake (and so the bunch shape) does not change significantly. Under these assumptions, the transverse equations of motion for an electron in the bunch are with the propagation distance z ¼ ct the independent variable. These equations can be derived from the Hamiltonian where the potential U ⊥ ðrÞ satisfies W ⊥ ¼ ∇ ⊥ U ⊥ . We recall that all the quantities are, in general, ζ-dependent. We also note that since the electrons in the bunch are highly relativistic, there is no longitudinal slippage, and so the motion in each longitudinal slice of the bunch is decoupled from the motion in other slices. Computing the final bunch emittance requires evaluating the second order transverse phase-space moments of the particle distribution at saturation, namely, for the x-plane, hx 2 i and hu 2 x i, where the average h·i is done over all the particles in the bunch. We recall that, at saturation, hxu x i ¼ 0, and so the final equilibrium emittance is given by ε Ã n;x ¼ ðhx 2 ihu 2 x iÞ 1=2 . Similar definitions apply to the moments in the y-plane.
In general, if the transverse wakefield has a nonlinear dependence from the transverse coordinates, the motion in the x and y planes is coupled and difficult to analyze analytically. However, since we are interested in the average properties of the transverse dynamics, and we made the assumption that the structure of the transverse wake does not change significantly during the phase-space mixing, we can simplify the problem by performing a dimensionality reduction. We can derive a 1D dynamical system describing the (approximate) dynamics of the bunch in the x-plane by averaging the 2D Hamiltonian H ⊥ given by Eq. (15)  does not depend explicitly on time (i.e., propagation distance), it can be used to analyze the particle orbits in the ðx; u x Þ phase plane. Given a particle with the initial condition ðx 0 ; u x0 Þ, and defining h 0 ¼ H ð1DÞ x ðx 0 ; u x0 Þ, then, in any point of the orbit identified by the chosen initial condition, the relationship between particle momentum and position is u x ðxÞ ¼ AEf2γ b ½h 0 − eU ð1DÞ ðxÞ=mc 2 g 1=2 . Using the equation of motion dx=dt ¼ cu x ðxÞ=γ b , we find that the time spent by a particle moving between the positions x and x þ dx is dt ¼ ½γ b =u x ðxÞdx=c. From this, we can evaluate the period of the closed orbit (i.e., the betatron period), The time average over the closed orbit (i.e., the equilibrium value after phasespace mixing) of any dynamical variable Qðx; u x Þ is then where we assume that Qðx; u x Þ is an even function of x and u x . In particular, for Q ¼ x 2 and Q ¼ u 2 x , we obtain, respectively,x 2 ðx 0 ; u x0 Þ ¼ ∂ h P 2 =∂ h P 0 j h¼h 0 ¼H ð1DÞ x l ½h − eU ð1DÞ ðxÞ=mc 2 1=2 dx. The second order transverse phase space moments at saturation in a given longitudinal slice along the bunch are obtained by averagingx 2 ðx 0 ; u x0 Þ andū 2 x ðx 0 ; u x0 Þ over the initial phase space distribution in the x-plane, namelȳ A similar expression holds forū 2 x . Finally, the bunch moments are obtained by averaging in ζ the slice moments, A similar expression holds for hu 2 x i. For a symmetric bunch (σ x ¼ σ y , ϵ n;x ¼ ϵ n;y ¼ ϵ n;0 ), where the transverse wake is given by Eq. (11), and U ⊥ ðrÞ ¼ R r 0 W r ðr 0 Þdr 0 , the second order moments at saturation in a given longitudinal slice arex 2 ðζÞ=σ 2 The final slice emittance is then ϵ n ðζÞ=ε n;0 ≃ 1 þ 0.005Γðζ=L b Þ 2 þ 0.0017Γ 2 ðζ=L b Þ 4 . By averaging longitudinally the slicedependent moments we obtain the saturated rms bunch moments that read hx 2 i=σ 2 x ≃ 1 − 0.039 Γ þ 0.003 Γ 2 , and hu 2 From this, we find that the final saturated bunch emittance is As expected, the ion-induced emittance growth effects increase with the parameter Γ (i.e., the final saturated emittance is higher for higher bunch densities, longer bunches, and lighter background ions). In Fig. 3 we show (red dots) the final bunch emittance (ε Ã n;x =ε n;0 − 1) obtained via fully non-linear modeling performed with INF&RNO, plotted as a function of the bunch peak density for a bunch propagating in an Hydrogen ion column. The bunch longitudinal and transverse distribution, energy, emittance, as well as the other plasma parameters are the same as in Fig. 2. The dashed black line is the theoretical result Eq. (19). The theoretical result is in good agreement with modeling for n b;0 =n 0 ≲ 2000 (corresponding to Γ ≲ 1). For higher bunch densities, where ion motion is relativistic, the expression for the wakefield Eq. (11) used in deriving Eq. (19) is no longer accurate.
In the strongly asymmetric case (σ x ≫ σ y , ϵ n;x ≫ ϵ n;y ), ion motion affects mainly the wakefield in the vertical y-plane, while the wakefield in the horizontal x-plane is left essentially unperturbed. In this case, the horizontal emittance is preserved (ε Ã n;x ≃ ε n;x ), while the saturated vertical emittance is ε Ã n;y =ε n;y ≃ 1 þ 0.0027 Γ þ 0.0053 Γ 2 . As before, we expect this result to be accurate as long as Γ ≲ 1. We note that the saturated emittance in the vertical plane is approximately twice as large compared to a round bunch.

V. EQUILIBRIUM BUNCH DISTRIBUTIONS ENABLING ION MOTION WITHOUT EMITTANCE GROWTH
Preserving the emittance of a bunch in presence of ion motion is of fundamental importance for the design of future PA-based LCs. A strategy to mitigate emittance growth based on the use of an adiabatic matching section at the entrance of the PA stage was described in Ref. [17].
Here we propose to longitudinally tailor the transverse bunch profile (and phase space) in such a way that, even though ion motion is enabled and the transverse wake is perturbed, the bunch transverse distribution at each longitudinal location is an equilibrium solution (i.e., the bunch properties do not evolve), and so the bunch remains matched at all times. Achieving this requires that, at any longitudinal location, the transverse phase space distribution of the bunch is a stationary (equilibrium) solution of the Vlasov equation. This solution allows for constant slice emittance, greatly simplifying the transport between PA stages. For simplicity, we will illustrate the construction of the matched solution in the case of a symmetric bunch with equal emittances in the x and y planes. However, this technique is also applicable in the strongly asymmetric case. Also, we consider, as before, a uniform longitudinal current profile of length L b . Generalization to arbitrary longitudinal profiles is straightforward.
The Vlasov equation governing the evolution of the transverse phase-space bunch distribution at any ζ location along the bunch, f 0;⊥ ðr; u; ζ; zÞ, is where f·; ·g are the Poisson brackets [24]. Stationary solutions to this equation are, by construction, distributions of the form f 0;⊥ ðr; u; ζÞ ∝ F½H ⊥ ðr; u; ζÞ=H 0 ðζÞ; where F is any positively defined function such that N ðζÞ ¼ R R f 0;⊥ ðr; u; ζÞd 2 rd 2 u < þ∞ (i.e., the distribution is normalizable), and H 0 ðζÞ is a ζ-dependent scale FIG. 3. Emittance growth at saturation for a Gaussian bunch initially matched in the unperturbed wakefield as function of the bunch peak density. Red dots are the simulation results, the black dashed line is the theoretical prediction Eq. (19). All the other bunch and plasma parameters are the same as in Fig. 2. parameter used to control locally the properties of the distribution. The second order phase-space moments at any slice along the bunch arē and the rms slice emittance is ϵ 2 n ðζÞ ¼x 2 ðζÞū 2 x ðζÞ. We require that the slice emittance is constant along the bunch, namely, for any ζ,x 2 ðζÞū 2 x ðζÞ ¼x 2 ðζ ¼ 0Þū 2 x ðζ ¼ 0Þ. This can be enforced by properly choosing the scale parameter H 0 ðζÞ. We note that since the knowledge of the phasespace distribution relies on the knowledge of the Hamiltonian, which depends on the wake potential that, in turns, is determined by the bunch distribution and ion dynamics via Maxwell's equations, obtaining an explicit expression for the bunch and wake quantities requires solving self-consistently the coupled set of Maxwell-Vlasov equations.
As an example, we will study the matched solution obtained taking FðxÞ ¼ expð−xÞ. We also define H 0 ðζÞ ¼ ðk p σ x Þ 2 =2ηðζÞ, where the slice-dependent parameter ηðζÞ satisfies ηðζ ¼ 0Þ ¼ 1. The transverse bunch profile including the proper normalization [see Eq. (4)] reads with To set the normalization we considered that for ζ ¼0 the wakefield is unperturbed [from Eq. (1) we have U ⊥ ðr;ζ¼0Þ¼ E 0 k p r 2 =4], and so the argument of the exponential function in Eq. (23) is ½2=ðk p σ x Þ 2 ðeU ⊥ =mc 2 Þ ¼ r 2 =2σ 2 x . We see that the transverse bunch profile at the head of the bunch (ζ ¼ 0) is Gaussian with rms bunch size σ x . For ζ < 0 the transverse bunch distribution Eq. (23) is, in general, non-Gaussian. From Eq. (23) we find that the second order spatial moment along the bunch is From Eq. (21) we have that the transverse momentum distribution at any slice along the bunch is Gaussian, and the second order moment of the distribution is The slice emittance for ζ ¼ 0 is ϵ n ¼ ðγ b =2Þ 1=2 k p σ 2 x . Imposing that, slice-by-slice, the emittance is constant (and equal to ϵ n ) we obtainx and this, using Eq. (25), can be rewritten as Equation (28) determines the values of ηðζÞ defining the desired matched solution with constant slice emittance. In general, solutions to Eq. (28) can be found numerically and require implementing an iterative procedure. Starting from an ansatz for the bunch shape, we compute the associated wake potential including ion motion, then we solve Eq. (28) for ηðζÞ, and subsequently we use ηðζÞ together with the wake potential to compute a new transverse density profile using Eq. (23). The cycle is repeated until convergence is reached. The numerical solutions is particularly fast in the nonrelativistic regime (Γ ≲ 1), where the expression for the wakefield is given by Eq. (10). We note that for the matched bunch solution the projected rms emittance differs from the slice emittance owing to the fact that the second order phase space moments are slice-dependent. In fact, the rms moments of the matched bunch are and hu 2 The projected rms bunch emittance is then ε n;x ¼ ðγ b =2Þ 1=2 k p σ 2 x ½hηih1=ηi 1=2 . In Fig. 4(a) we show a plot of the normalized rms bunch size along the bunch, ½x 2 ðζÞ 1=2 =σ x , for the matched equilibrium solution obtained for different values of the bunch peak density. The bunch length, energy, slice emittance, and the background ions parameters are the same as in the example discussed in Fig. 2. We see that the matched bunch size shrinks towards the tail of the bunch in response to the increase of the strength of the transverse wakefield. Figure 4(b) shows a 2D ðζ; xÞ snapshot of the bunch density for the equilibrium solution obtained for n b;0 =n 0 ≃ 12 000. The projected rms emittance for the matched solution is ε n;x ¼ 0.63 μm. The red curve in Fig. 2(a) shows the emittance evolution for the bunch depicted in Fig. 4(b). In Fig. 2(b) we show (red curve) the corresponding slice emittance measured for k β z ¼ 12. No emittance growth is observed in this case, and the slice emittance is constant along the bunch.
The experimental realization of the matched equilibrium solution requires an exact preparation of the initial transverse phase space of the bunch, and this can be difficult to achieve in practice. For instance, the transverse bunch shape given by Eq. (23) changes continuously along the bunch and is, in general, non-Gaussian. We investigated the behavior of an approximate equilibrium solution with a bunch structure much simpler compared to the exact equilibrium case. The approximate solution is such that, at any longitudinal location, the transverse density profile is Gaussian with the same rms properties of the exact equilibrium. The blue line in Fig. 2(a) shows the emittance evolution for the approximate equilibrium bunch. In Fig. 2 (b) we show (blue curve) the corresponding slice emittance measured for k β z ¼ 12. As expected, we observe some emittance growth in this case. However, the final emittance increase is very moderate (∼3%), and the slice emittance shows a maximum increase ≲10%, much smaller compared to the nontailored Gaussian case (black curve). The possibility of generating a tapered, low-emittance electron bunch via ionization injection [25] using a laser pulse with a properly tailored evolution is currently under investigation.

VI. CONCLUSIONS
In this article we analyzed ion motion and the associated emittance growth of a relativistic bunch propagating in an ion column for parameters relevant to the design of next generation PA-based LCs. The model derived here well describes the transverse dynamics of a bunch in the nonlinear wakefield generated by a particle or a laser driver in a PA. We provided analytical expressions, valid in the regime of nonrelativistic ion motion, for the structure of the transverse wake, and for the corresponding ionmotion-induced bunch emittance growth. The analytical results are in good agreement with numerical modeling. We proposed and analyzed a solution that completely eliminates ion-motion-induced emittance growth. This solution requires a head-to-tail shaping of the bunch distribution. Controlling emittance growth is critical to high-energy physics applications of plasma accelerators.