Kinetic description of neutralized drift compression and transverse focusing of intense ion charge bunches

1098-4402= A kinetic model based on the Vlasov equation is used to describe the axial drift compression and transverse focusing of an intense ion charge bunch propagating along the axis of a solenoidal focusing field Bsol x . The space charge and current of the ion charge bunch are assumed to be completely neutralized by the electrons provided by a dense background plasma. In the absence of self-field forces, the Vlasov equation is solved exactly for general initial distribution function fb x; v; 0 , using the method of characteristics. It is shown that the Vlasov equation possesses a class of exact, dynamically evolving solutions fb W?; Wz , where W? and Wz are transverse and longitudinal constants of the motion. Detailed dynamical properties of the charge bunch are calculated during axial compression and transverse focusing for several choices of distribution function fb W?;Wz .


I. INTRODUCTION
Intense charged particle beams have a wide range of applications ranging from basic research in high energy and nuclear physics, to spallation neutron sources, to heavy-ion-driven high energy density physics and fusion, to nuclear waste transmutation, to mention a few examples [1][2][3][4][5][6][7][8][9].Of considerable practical importance for heavy ion beam applications to high energy density physics and fusion is the axial compression and transverse focusing of the (initially long) charge bunch to a small spot size at the target location.Indeed, considerable technical progress has been made in analytical and numerical studies of the axial drift compression and transverse focusing of long, unneutralized charged bunches [10 -22], including the recent development and application of a self-consistent warm-fluid model, optimization of the axial pulse shape, and the lattice design for transverse focusing to a small spot size [18,19,21,22].At the very high space-charge intensities characteristic of intense heavy ion beams, one of the major challenges is compressing the charge bunch axially (against space-charge forces) to the very short pulse length ideal for maximizing the beam intensity (number density) focused on the target.To facilitate the compression, one approach now under active theoretical and experimental investigation [23][24][25] is called neutralized drift compression.In neutralized drift compression, the space charge and current of the intense ion charge bunch are neutralized by the electrons from a dense background plasma.When the charge bunch enters the plasma with an initial axial velocity tilt, the absence of axial spacecharge force greatly facilitates the axial compression of the charge bunch to a short axial length, and also facilitates the transverse focusing of the charge bunch to a small transverse spot size.In this paper, we make use of a kinetic description based on the Vlasov equation [1] to develop a self-consistent dynamical model that describes the axial compression and transverse focusing of the charge bunch for the case where the background plasma completely neutralizes the space charge and current of the charge bunch.
In this paper, we consider the axial drift compression and transverse focusing of an intense charge bunch made up to positively charged ions for the case where the space charge and current of the charge bunch are completely neutralized by a background plasma.The charge bunch propagates in the z direction along the axis of an applied solenoidal magnetic field B sol x, and we make use of the Vlasov equation [Eq.(1)] to describe the evolution of the beam distribution function f b x; v x ; y; v y ; z; v z ; t in the six-dimensional phase space x; v x ; y; v y ; z; v z , neglecting self-field effects E s 0 B s because of charge and current neutralization.The analysis is carried out in the beam frame in the paraxial (thin-beam) approximation, and the axial force on the particle motion is treated as negligibly small F z ' 0. The theoretical model and assumptions in the present analysis are described in Sec.II.In Sec.III, we first consider the case where the full Vlasov equation ( 1) is integrated over the transverse phase space x; v x ; y; v y , and the evolution of the longitudinal distribution function F b z; v z; t R dxdydv x dv y f b x; v x ; y; v y ; z; v z ; t is determined.A wide variety of distribution functions is considered that demonstrate the flexibility of the model, and the dynamic axial compression of the charge bunch is investigated in detail for the class of distribution functions F b W z that have a negative velocity tilt at t 0. Here, W z is the constant of the motion defined in Eq. (39).In Sec.IV, for general initial distribution function f b x; v x ; y; v y ; z; v z ; t 0 the full Vlasov equation (1) is solved exactly using the method of characteristics.We then specialize to the class of distribution functions f b W ? ; W z , where W z and W ? are the constants of the motion defined in Eqs.(39) and (100), and the time-PHYSICAL REVIEW SPECIAL TOPICS -ACCELERATORS AND BEAMS 8, 064201 (2005) 1098-4402=05=8(6)=064201 (15) 064201-1  2005 The American Physical Society dependent axial and radial scale factors, z b t and r b t, solve Eqs. ( 40) and (101), respectively.As demonstrated in Sec.IV, the dynamical properties of the axial compression and transverse focusing of the charge bunch can be calculated for a wide variety of distribution functions f b W ? ; W z .What's most remarkable about the present analysis is the enormous flexibility there is in the choices of distribution function f b W ? ; W z that undergo axial drift compression and transverse focusing, and correspondingly have different profiles for number density, line density, longitudinal and transverse pressures, etc.This important feature, of course, has been facilitated by the assumption of complete neutralization of the charge and current of the charge bunch by the background plasma.The detailed optimization of the solenoidal field configuration, velocity tilts, etc., that provide simultaneous axial drift compression and transverse focusing to a small, localized spot size is being investigated numerically, and will be subject of a future publication.
An important motivation for the present paper is its application to the high-intensity ion beams characteristic of ion-beam-driven high energy density physics and fusion.Such beams tend to be space-charge dominated (prior to charge and current neutralization), and are characterized by high volume number density and low transverse emittance.It is important to recognize, however, that the present analysis is quite general, and applies equally well to low-intensity and moderate-intensity charge bunches, characterized by lower number density and higher transverse emittance.The main assumption here is that the background plasma completely charge neutralizes and current neutralizes the charge bunch, and that the selfelectric and self-magnetic fields are treated as negligibly small.

II. THEORETICAL MODEL AND ASSUMPTIONS
We consider an intense ion charge bunch with characteristic diameter 2r b and length 2z b that propagates with average axial velocity V b e z relative to the laboratory frame through a static solenoidal magnetic field B sol x l B z z l e z ÿ 1=2@B z z l =@z l x l e x y l e y in the thinbeam approximation.Here, x l ; y l ; z l are configurationspace variables in the laboratory frame, and x l e x y l e y is the transverse displacement from the beam axis.For present purposes, the space charge and current of the ion charge bunch are assumed to be fully charge neutralized and current neutralized by background plasma E s 0 B s , and the particle motion in the beam frame is assumed to be nonrelativistic.Denoting the phase-space and time variables in the beam frame by x; v x ; y; v y ; z; v z ; t, where v p=m b is the (nonrelativistic) velocity, the Vlasov equation for the one-particle distribution function f b x; v x ; y; v y ; z; v z ; t in the beam frame can be approximated by [1,26,27] In Eq. ( 1), b 1 ÿ V 2 b =c 2 ÿ1=2 is the relativistic mass factor, e b and m b are the ion charge and rest mass, respectively, c is the speed of light in vacuo, and B z and B 0 z are evaluated at z l b z V b t [27], i.e., Furthermore, in obtaining Eq. ( 1) we have neglected selffield forces E s 0 B s because of the assumptions of charge and current neutralization, and we have neglected the z component of magnetic force F sol z in comparison with F sol x and F sol y , which is a valid approximation provided v 2 ? =V 2 b 1, i.e., provided the directed axial velocity is large in comparison with the transverse particle speed.
Equation ( 1) constitutes the final form of the Vlasov equation to be used in the subsequent analysis.Because of the assumptions of charge and current neutralization, note that Eq. ( 1) is a linear equation for the beam distribution function f b , and is not nonlinearly coupled to Maxwell's equations.Indeed, because Eq. ( 1) is linear, it can be readily solved exactly for general initial distribution function f b x; v x ; y; v y ; z; v z ; t 0 [1,26].Furthermore, because Eq. (1) includes the effects of a z-dependent solenoidal field, it can be used to investigate the transverse focusing of the charge bunch as well as the axial drift compression of the charge bunch produced by an initial axial velocity tilt.Here, an axial ''velocity tilt'' is introduced prior to charge neutralization by applying a voltage pulse over the length of the charge bunch, which causes the front portion of the bunch to move more slowly than the rear portion, thereby causing the bunch to compress axially. For The integration domains in Eqs. ( 3)-( 6) cover the full available phase space, and it is assumed that f b 0 for v x ! 1, v y ! 1, v z ! 1, or z ! 1, and that f b 0 for jxj, jyj > r b .Here it is assumed that there are no beam particles beyond some transverse distance r b from the axis of the charge bunch.The Vlasov equation (1) of course conserves the total number of particles N b in the charge bunch.Operating on Eq. (1) with R dxdydzdv x dv y dv z and integrating by parts with respect to x, v x , etc., readily gives as expected.

III. ONE-DIMENSIONAL MODEL OF AXIAL DRIFT COMPRESSION
Before analyzing the three-dimensional case, we consider the implications of the Vlasov equation when Eq. ( 1) is integrated over the transverse phase-space variables x; v x ; y; v y .Introducing the one-dimensional longitudinal distribution function F b z; v z; t defined by we note that Eq. ( 7) represents the projection of the full distribution function f b onto z; v z ; t space.Operating on Eq. (1) with R dxdydv x dv y , integrating by parts with respect to x, y, v x , and v y , and making use of R dxdydv x dv y v x @f b =@x 0, R dxdydv x dv y v y @f b =@y 0, etc., it readily follows that @ @t v z @ @z F b z; v z ; t 0: From Eqs. (4), (6), and (7), note that the normalization on where b z; t is the line density and N b is the number of particles in the charge bunch.Equation ( 8) of course has the exact solution for general initial distribution function at t 0.
The striking feature of Eqs. ( 1) and ( 8) is that no matter how complicated the x-y dynamics and transverse focusing of the charge bunch in the solenoidal field, the longitudinal distribution function F b z; v z ; t solves the simple onedimensional kinetic equation (8).Of course, this is a consequence of the assumptions of complete charge and current neutralization E s 0 B s which neglect the effects of collective interactions with the background plasma, and the assumption that v 2 ?=V 2 b 1, which allowed us to neglect the small axial force F sol z ' 0 associated with the solenoidal focusing field in Eq. (1).The simple form of Eq. ( 8) [and Eq. ( 10)] permits a detailed examination of the longitudinal dynamics of the charge bunch during neutralized drift compression.For future reference, it is useful to introduce the axial flow velocity V zb z; t and the effective particle pressure P b z; t in the beam frame defined by and Note from Eq. ( 11) that an initial velocity tilt with V zb z; 0 Þ 0 in the beam frame corresponds to an initial distribution function f b z; v z ; 0 that is asymmetrically skewed in v z space, e.g., by application of a voltage pulse over the length of the charge bunch prior to entering the neutralizing background plasma.Note also from Eq. ( 12) that the pressure P b z; t is a measure of the axial velocity spread of the particles composing the charge bunch, and that can be identified with the effective longitudinal temperature.
As noted earlier, the exact solution to Eq. ( 8) is given by Eq. (10) for general initial distribution function F b z; v z ; 0. To illustrate the rich physics content of Eq. ( 8), we now examine the rate equations for the evolution of global averages, such as root-mean-square bunch dimension, axial emittance, etc. (Sec.III A), investigate a class of exact dynamical solutions to Eq. ( 8) undergoing axial drift compression (Sec.III B), and develop a density inversion theorem to determine the distribution function F b for specified line density profile b (Sec.III C).

A. Dynamical evolution of statistical averages
As is customary, we define the statistical average hi of a phase-space function z; v z by [1,28,29] where F b evolves according to Eq. ( 8), and Several interesting and exact dynamical properties of the charge bunch follow from Eqs. ( 16)- (18).For example, if the charge bunch is initially centered at z 0 with hzit 0 0, and the z-averaged axial flow velocity is initially zero with hv z it 0 0, then it follows from Eqs. ( 16) and ( 17) that hzit 0 and hv z it 0 at all subsequent times t.Furthermore, it is convenient to introduce the root-mean-square bunch length z b t and (dimensional) axial emittance " z t defined by and where use have been made of Eqs. ( 16)- (18).Equation ( 21) implies that the axial emittance is constant with as expected during axial drift compression described by the present simple model.We denote initial values of quantities at t 0 by a subscript 0, i.e., z b0 z b t 0, etc.It then follows that b0 can be expressed as where is the characteristic initial thermal speed squared.Substituting Eqs. ( 22) and ( 23) into Eq.( 20) and solving for z 2 b t, we obtain Equation ( 24) is a closed expression for the mean-square bunch length z 2 b t, valid for arbitrary initial distribution function F b z; v z ; 0. As such, it has a wide range of applicability no matter how complicated the evolution of F b z; v z ; t in phase space according to Eq. ( 8).Note that _ z b0 is related to an average of the initial axial flow velocity of the charge bunch.Making use of z b _ z b hzv z i [Eq.( 16)], it follows from Eq. ( 14) that where V zb z; t 0 is the initial axial flow velocity defined by Eq. ( 11).If, for example, b z; t 0 is an even function of z centered at z 0, and V zb z; t 0 is an odd function of z with V zb z > 0; t 0 < 0, then Eq. ( 25) gives z b0 _ z b0 < 0, corresponding to an initial axial compression of the charge bunch.
The solution for z b t in Eq. ( 24) is illustrated in Fig. 1 for the case where _ z b0 < 0: Note that z b t contracts from z b0 at t 0 to a minimum value z bf at time t t f where and The values of t f and z bf z b t t f are calculated by setting _ z b t tt f 0 in Eq. ( 24) and solving for t f .For t > t f , note from Eq. ( 24) and Fig. 1 As would be expected, the smaller the initial thermal speed v T0 in comparison with _ z b0 , the smaller the compressed bunch length z bf in comparison with z b0 .Indeed, for v 2 T0 _ z 2 b0 , Eqs. ( 26) and ( 27

B. Class of exact dynamical solutions for F b z; v z ; t
Equation ( 10) can be used to determine the exact solution for F b z; v z ; t for a wide range of choices of initial distribution function F b z; v z ; 0. For present purposes, in this section we investigate a class of exact dynamical solutions to Eq. ( 8) undergoing axial drift compression with a linear velocity tilt, i.e., V zb z; t is linearly proportional to z over the axial extent of the charge bunch.
As an application of Eq. ( 10), we consider an example where the initial distribution function is given by where N b , z b0 , _ z b0 , and v 2 T0 2T 0 =m b are constants.From Eqs. ( 9) and ( 28), note that the initial line density has the Gaussian profile where b0 N b =z 2 b0 1=2 , and is the number of particles in the charge bunch.Furthermore, from Eqs. ( 11), (28), and ( 29), the initial average flow velocity profile is which is a linear function of z for the shifted Maxwellian in Eq. ( 28).
To determine the distribution function F b z; v z ; t at subsequent times, we now make use of F b z; v z ; t F b z ÿ v z t; v z ; 0. Substituting Eq. ( 28) into Eq.( 10), some straightforward algebraic manipulation shows that where z 2 b t is defined by z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 [see also Eq. ( 24)].Equation ( 32) reduces to Eq. ( 28) at t 0, but otherwise has a similar functional form to Eq. ( 28) in appropriately scaled variables.Substituting Eq. ( 32) into Eqs.( 9) and ( 11)-( 13), it is readily shown that the dynamical expressions for various macroscopic properties of the charge bunch are given by b z; and For negative initial velocity tilt with _ z b0 < 0; it follows from z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 and the discussion in Sec.III A that z b t continues to compress from its initial value z b0 at t 0 to the minimum value z bf at t t f [see Eqs. ( 26) and ( 27 Moreover, because dz b =dt tt f 0, the instantaneous value of flow velocity at peak compression is V zb z; t t f 0, which follows from Eq. (34).While Eq. (32) provides an interesting example of an exact solution to Eq. ( 8) for the specified initial condition in Eq. ( 28), it is also suggestive of the existence of a FIG. 1.For negative initial velocity tilt _ z b0 < 0 the rms bunch length z b t [Eq.( 24)] compresses to a minimum size z bf at time t t f [Eqs.( 26) and ( 27)].
general class of exact dynamical solutions to Eq. ( 8) undergoing drift compression.In this regard, we consider the class of distribution functions F b z; v z ; t of the form where W z z; v z ; t is defined by and F b W z is a smooth differentiable function.In Eq. (39), it is assumed that z b t solves which has the exact solution z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 [see also Eq. ( 24)].Here, z b0 z b t 0 and _ z b0 dz b =dt t0 are initial conditions, and v 2 T0 is a constant with dimensions of thermal speed-squared.We now demonstrate that the class of distribution function in Eq. ( 38) exactly solves the Vlasov equation (8).Operating on Making use of the definition of W z in Eq. ( 39), some algebraic manipulation and rearrangement of terms shows that where use has been made of Eq. (40).It follows from Eq. ( 41) and (42) that which demonstrates that F b W z represents a class of exact dynamical solution to Eq. ( 8).That W z defined in Eq. (39) plays the role of a constant of the motion in the Vlasov equation ( 43) can also be demonstrated by a standard Courant-Snyder analysis [1,30].Because F sol z ' 0 is assumed in the present analysis, the axial motion is described by v z dz=dt and dv z =dt 0, or equivalently , where A z and z0 are constant, it is readily shown that d 2 z b =dt 2 z 2 b0 v 2 T0 =z 3 b , which is identical to Eq. (40).It can also be shown that the constant A z is exactly related to W z defined in Eq. (39) by W z A 2 z =z b0 v T0 const.Introducing the dimensionless variables Z z=z b t and 9) can be expressed in the simple form Furthermore, it readily follows from Eqs. ( 11), (38), and (39) that the axial flow velocity is given exactly by for the entire class of distribution functions in Eq. (38).In obtaining the result in Eq. ( 45), use has been made of v z z=z b dz b =dt z b0 v T0 =z b V z in integrating over V z .Note from Eq. ( 45) that V zb z; t has a velocity tilt that is linear in z and proportional to z ÿ1 b dz b =dt over the entire length of the charge bunch.Other interesting macroscopic properties of the charge bunch include the line density b z; t and pressure P b z; t defined in Eqs. ( 9) and (12).For the class of distribution functions described by Eqs. ( 38) and (39), some straightforward algebraic manipulation shows that for general choice of F b W z .Note from Eqs. ( 46) and (47) that the line density and pressure depend on z and t through the scaled variable z=z b t, where z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 solves Eq. ( 40).There are of course many possible choices of the functional form for the distribution function F b W z .For example, the distribution function in Eq. ( 32) is given simply by which has the equilibrium properties described earlier in this section.For example, Eq. ( 48) has the Gaussian line density profile for b z; t given by Eq. (33).As a class of distribution functions that has a sharp spatial edge to the line density profile, we consider the case where Here, n 0; 1; 2; . . . is an integer, W z z; v z ; t is defined in Eq. (39), and the normalization of , as required [see Eq. ( 44)].As illustrated in Fig. 2, in the phase space z; v z the distribution function in Eq. (49) has a sharp elliptical boundary described by W z 1, i.e., where Note that the phase-space density F b is nonzero inside the ellipse 0 W z < 1 and zero outside the ellipse W z > 1 in Eq. ( 50) and Fig. 2. Indeed, it follows from Eq. (49) that for the special case where n 0, F b W z const inside the ellipse W z 1, whereas for n 1; 2; . . ., F b W z decreases monotonically to zero as W z ! 1.As would be expected, the corresponding line density profiles and pressure profiles will be different for different choices of n 0; 1; 2; . . .(see discussion later in this section).
In Eq. (50) and Fig. 2, it is assumed that _ z b0 dz b =dt t0 < 0. At t 0, the ellipse has full length 2z b0 projected along the v z axis.For 0 t < t f ÿz b0 _ z b0 = _ z 2 b0 v 2 T0 , the ellipse rotates counterclockwise, and then reverses direction at t t f .Here, t f is calculated from the condition dz b =dt tt f 0. At time t t f , the full length projection of the ellipse along the z axis is 2z bf (peak compression), and the full width projection along the v z axis is 2v zf , where For t > t f , the ellipse in Fig. 2 and Eq.(50) rotates clockwise, and time asymptotically t t f becomes a long, thin ellipse with full length 2z b t t f 2 _ z 2 b0 v 2 T0 1=2 t projected along the z axis.During the dynamic rotation and evolution of the ellipse in Eq. (50) and Fig. 2, the area of the ellipse enclosed by the curve W z 1 in z; v z phase space of course remains constant and equal to z b0 v T0 .
We now return to the problem of evaluating closed expressions for the line density b z; t and pressure P b z; t for the choice of distribution function in Eq. ( 49).Substituting Eq. (49) into Eq.( 46), the line density can be expressed as for 0 z 2 < z 2 b t, and b z; t 0 for z 2 > z 2 b t.Changing variables to V z 1 ÿ z 2 =z 2 b 1=2 sin in Eq. ( 51) and making use of R =2 0 dcos 2n1 2n!!=2n 1!!, the integration can be carried out to give where the constant b0 b z 0; t 0 is defined by For n 0, which corresponds to uniform phase-space density with F b W z const for 0 W z 1 in Eq. ( 4), the line density profile in Eq. ( 52) is proportional to 1 ÿ z 2 =z 2 b 1=2 .For n 1; 2; . . ., however, the density profile in Eq. ( 52) becomes increasingly peaked at smaller values of z 2 < z 2 b t.Substituting Eq. ( 49) into Eq.( 47), it can similarly be shown that the pressure P b z; t can be expressed as where P b0 P b z 0; t 0 is defined by (55) and b0 is defined in Eq. (53).From Eq. (54), the pressure profile also becomes increasingly peaked for increasing values of n 0; 1; 2; . . . .Moreover, substituting Eqs. ( 52) and (54) into Eq.( 13), the effective longitudinal temperature T b z; t P b z; t= b z; t can be expressed over the interval 0 z 2 < z 2 b t as for the entire class of distribution functions F b W z defined in Eq. (49) for n 0; 1; 2; . . . .Finally, for W z defined by Eq. (39), where z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 , it should be noted that z 2 b t generally differs from the mean-square bunch length hz 2 i by a constant scale factor for the class of distribution function F b W z in Eq. (38).In particular, it follows from Eqs. ( 14), (38), and (39) that FIG. 2. For _ z b0 < 0, the ellipse W z 1 rotates counterclockwise in the phase space z; v z from its initial position at t 0. The ellipse becomes vertically aligned and slows to a stop at t t f (maximum axial compression), and then rotates clockwise for t > t f .For t t f , the ellipse expands indefinitely along the z axis.Throughout its dynamics the ellipse W z 1 remains centered at z; v z 0; 0 and has constant area z b0 v T0 .
where the scale factor g is defined by For example, for the choices of F b W z in Eqs. ( 48) and (49), it follows that g 1=2 and g 1=2n 4, respectively.

C. Density inversion theorem
For the class of distribution functions in Eq. ( 38) the detailed properties of the charge bunch can be calculated for a wide variety of choices of F b W z .For example, the line density b z; t is given by Eq. ( 46), where W z V 2 z Z 2 , and Z z=z b t and V 2 z =V z , and making use of W z V 2 z Z 2 , Eq. ( 46) can be rewritten as Because the dependence of F b W z on V 2 z and Z 2 is so intertwined through the variable W z V 2 z Z 2 , a knowledge of the line density profile b Z 2 is actually sufficient information to reconstruct the distribution function F b W z .Indeed, inverting Eq. ( 59), it can be shown that the solution for F b W z is given by (60) Equation ( 60) can be verified by direct substitution.Substituting Eq. ( 60) on the right-hand side of Eq. ( 59) gives The integration in square brackets in Eq. ( 61) can be carried out exactly to give (for Here, we have changed variables to wz W z ÿ Z 2 , and then to wz Z 02 ÿ Z 2 sin 2 , to obtain the value of the integral in Eq. (62).Substituting Eq. (62) into Eq.(61) then gives where b Z 2 ! 1 0 has been assumed.The result in Eq. ( 63) is identical to Eq. ( 59), which completes the proof of the density inversion theorem.Of course the density inversion theorem in Eq. ( 60) can be used to determine the distribution function F b W z self-consistently for a wide variety of choices of line density profile b Z 2 .

A. Exact solution by the method of characteristics
To determine the effects of transverse focusing as well as longitudinal drift compression, we make use of the Vlasov equation ( 1) for the distribution function f b x; v x ; y; v y ; z; v z ; t in the full six-dimensional phase space.The simplest approach for solving Eq. ( 1) is based on the method of characteristics [26], which makes use of the fact that the distribution function f b is a constant following the exact particle motion in the field configuration acting on the particles.In particular, we denote the particle trajectories by the ''primed'' orbits x 0 t 0 and v 0 t 0 , where t 0 is the running time variable.We also choose ''initial'' conditions such that x 0 t 0 ; v 0 t 0 passes through the phase-space point x; v at time t 0 t, i.e., x 0 t 0 t x; v 0 t 0 t v: It then follows from d=dt 0 f b x 0 t 0 ; v 0 t 0 ; t 0 0 that f b x 0 t 0 ; v 0 t 0 ; t 0 const (independent of t 0 ), which gives f b x 0 t 0 ; v 0 t 0 ; t 0 j t 0 t f b x 0 t 0 ; v 0 t 0 ; t 0 j t 0 0 : (65) Substituting Eq. (64) into Eq.( 65) then gives f b x; v; t f b x 0 t 0 0; v 0 t 0 0; t 0 0: (66) Therefore, from Eq. (66), a determination of the particle orbits x 0 t 0 and v 0 t 0 in the field configuration acting on the particles, and enforcement of the ''initial'' conditions in Eq. ( 64), in principle permits a determination of the distribution function f b x; v; t for general choice of initial distribution f b x; v; 0. If constants of the motion exist for the primed orbits x 0 t 0 and v 0 t 0 , then the distribution function f b can also be constructed from such constants.The equations of motion for x 0 t 0 and v 0 t 0 take on a particularly simple form because F z ' 0 and self-field effects have been neglected in Eq. ( 1).Introducing the Larmor frequency L defined by it follows from Eq. ( 1) that the orbit equations for x 0 t 0 and v 0 t 0 solve In Eq. ( 68), L b z 0 V b t 0 is defined in Eq. ( 67).The equations for z 0 t 0 and v 0 z t 0 in Eq. ( 68) are readily integrated to give where z; v z are phase-space coordinates, and use has been made of z 0 t 0 t z and v 0 z t 0 t v z .As expected, the axial motion in Eq. ( 69) is free streaming, which is a consequence of the assumption F z ' 0 in Eq. ( 1).
Referring to the transverse orbit equations in Eq. ( 68), we note from Eq. (69) that the Larmor frequency occurs with Therefore, from Eq. ( 70), it follows that which can be used to simplify the transverse equations of motion in Eq. ( 68).It is convenient to introduce the timedependent phase L t 0 ÿ t; t for the Larmor motion defined by z @ L =@ b z 0 , where use has been of Eqs. ( 71) and (72).Substituting into Eq.( 68), the transverse orbit equations for x 0 t 0 and y 0 t 0 become The form of the transverse orbit equations in Eqs. ( 73) and ( 74) is particularly amenable to direct calculation.Combining Eqs. ( 73) and (74), it is readily shown that x 0 t 0 iy 0 t 0 solves It is convenient to transform the orbit Eq. ( 75) for x 0 iy 0 in the laboratory frame to an orbit equation for X 0 iY 0 rotating at the Larmor frequency relative to the nonrotating beam frame.That is, we introduce X 0 t 0 iY 0 t 0 defined in terms of x 0 t 0 iy 0 t 0 by [31] x 0 iy 0 X 0 iY 0 expi L ; (76) or equivalently where L t 0 ÿ t; t is defined in Eq. ( 72).Denoting d=dt 0 by a superdot, Eq. (76) gives Substituting Eqs. ( 76), (78), and (79) into Eq.( 75), some straightforward algebra shows that or equivalently, In Eqs. ( 81) and ( 82), X 0 t 0 and Y 0 t 0 are the transverse particle orbits in the rotating frame and the Larmor fre- The advantage of transforming the orbit equations to the Larmor frame is evident.First, the X 0 and Y 0 motions in Eqs. ( 81) and (82) are uncoupled.Second, the X 0 and Y 0 forces are restoring, and linearly proportional to X 0 and Y 0 , respectively.Once the solutions for X 0 t 0 and Y 0 t 0 are determined from Eqs. ( 81) and (82) (Sec.IV B), the back transformation to the nonrotating beam frame is straightforward.In particular, it follows directly from Eqs. ( 76 and (78) that where _ L d L =dt 0 L , x 0 x 0 t 0 , _ x 0 dx 0 t 0 =dt 0 v 0 x t 0 , X 0 X 0 t 0 , etc. Making use of Eqs. ( 69) and ( 83)-( 86), the formal solution for the distribution function f b x; v x ; y; v y ; z; v z ; t is then given by where use has been made of Eq. (66).Although the primed orbit x 0 t 0 ; y 0 t 0 ; . . ., pass through the phase-space point x; y; . . ., at time t 0 t [Eq.( 64)], it should be kept in mind that x 0 t 0 0, y 0 t 0 0; . . ., do depend explicitly on the time variable t [e.g., z 0 t 0 0 z ÿ v z t in Eq. ( 69)].

B. Dynamically focused solutions depending on constants of the motion
The procedure described in Sec.IVA, based on the method of characteristics, can be used to determine the distribution function f b x; v; t for a wide variety of choices of initial distribution function f b x; v; 0 by solving the orbit equations ( 81) and (82) and making use of Eqs.(83)-(87).For the present purpose, we examine the class of distribution functions f b constructed from the constant of motion W z [Secs.III B and III C], and the two constants of the motion A 2 x and A 2 y [31], which are exact consequences of the transverse orbit equations (81) and (82).For example, expressing X 0 t 0 A x wt 0 cos R t 0 t dt 00 =w 2 t 00 x , where A x and x are constants, it is readily shown from Eq. (81) that X 0 t 0 and dX 0 t 0 =dt 0 exactly satisfy [31] X 02 w 2 t 0 w 2 t 0 dX 0 dt 0 ÿ Here, the constant A 2 x is independent of t 0 , and the envelope function wt 0 solves [31] where the Larmor frequency L f b z V b t 0 ÿ t v z t 0 ÿ t V b tg is defined by Eqs. ( 67) and (70).A similar conservation constraint pertains to the Y 0 motion described by Eq. ( 82), with where wt 0 solves Eq. (89).Although A 2 x and A 2 y are independently conserved, for present purposes in the remainder of Sec.IV, we consider the class of distribution function f b that depend explicitly on A 2 x and A 2 y through the linear combination A 2 x A 2 y const, which have the desirable feature that they are axisymmetric @=@ 0 about the z axis of the charge bunch.For dimensional convenience, we further define w 2 t 0 r 2 b t 0 =r b0 v ?0 , where r b t 0 has dimension of length, and r b0 and v ?0 are constant scale factors with dimensions of length and perpendicular speed, respectively.It then follows that multiplying Eqs. ( 88) and (90) by r b0 v ?0 ÿ1 , and adding the two resulting equations gives the exact constant of the motion W ? defined by In Eq. ( 91), the (dimensionless) constant of the motion W ? is independent of t 0 , and from Eq. ( 89) r b t 0 r b0 v ?0 1=2 wt 0 solves where L f b z V b t 0 ÿ t v z t 0 ÿ t V b tg is defined by Eqs.(67) and (70).In Eq. (92), r b t 0 has the dimensions of length, and the constant r 2 b0 v 2 ?0 plays the role of an effective transverse emittance.As will become apparent in the subsequent analysis, r b t 0 differs from the rms bunch radius by a (constant) scale factor that depends on the specific choice of distribution function f b and its dependence on the constant of the motion W ? .
We now transform the constant of the motion W ? defined by Eq. (91) in the Larmor frame back to the nonrotating beam frame where the orbits are x 0 t 0 , _ x 0 t 0 dx 0 t 0 =dt 0 ; . . . .From Eqs. ( 76) and (78), it is readily shown that the inverse transformation to Eqs. ( 83)-( 86) is given exactly by  93)-(96), some straightforward algebra shows that X 0 t 0 , _ X 0 t 0 dX 0 t 0 =dt 0 ; . . .are related to x 0 t 0 , _ x 0 t 0 dx 0 t 0 =dt 0 ; . . .by (97) Substituting Eqs. ( 97) and (98) into Eq.(91) gives where Equation ( 99), where r b t 0 solves Eq. ( 92), provides an exact representation of the constant of the motion W ? in terms of the nonrotating beam-frame orbits x 0 t 0 and y 0 t 0 .Equation ( 99) is valid for general t 0 .Evaluating Eq. (99) at time t 0 t, where x 0 ; dx 0 =dt 0 ; . . .passes through the phase-space point x; v x ; y; v y , then gives Equation (100) provides an exact expression for the constant of the motion W ? in terms of the nonrotating beamframe phase-space variables x; v x ; y; v y .In Eq. (100), evaluating r b t 0 and L at t 0 t, it follows from Eq. (92) that r b t solves where L L t 0 t is defined by (102) and use has been made of Eqs.(67) and (70).Note that Eq. (101) plays the role of an effective envelope equation for the evolution of r b t including the effects of the solenoidal focusing through the term proportional to 2 L r b .Moreover, because L depends on both z and t through the combination b z V b t in Eq. ( 102), the quantity r b t similarly depends on z and t.
Because W z [Eq.( 39)] and W ? [Eq.( 100)] are both constants of the motion, the Vlasov equation ( 1) in beam-frame variables supports the class of exact solutions Depending on the specific choice of f b W ?; W z , Eq. ( 103) can be used to investigate the detailed dynamic focusing of the charge bunch, both axially and transversely, for a wide range of choices of distribution function and solenoidal focusing lattice.Detailed properties of the dynamic evolution of the charge bunch can be calculated for a variety of choices of distribution function f b W ?; W z .Because of the functional forms of W z and W ? defined in Eqs. ( 39) and (100), respectively, several macroscopic properties of the charge bunch are readily calculated.For example, the average flow velocity b zdz b =dt.In cylindrical polar coordinates r; ; z, where x r cos and y r sin, Eq. ( 104) can also be expressed as where 0 dV ?V ? . . ., Eq. ( 109) can then be expressed as Note from Eq. (113) that n b x; t generally depends on the scaled variables R 2 r 2 =r 2 b t and Z 2 z 2 =z 2 b t, and that the integral in Eq. ( 112) is multiplied by the timedependent factor r 2 b tz b t ÿ1 .Similarly, without presenting algebraic details, some straightforward algebra show that the pressure tensor defined by [1] Here, the effective transverse and longitudinal temperatures are defined by The specific functional form of n b x; t, T ?b x; t, T zb x; t, and other properties of the charge bunch of course depend on the choice of distribution function f b W ?; W z .As an example of a charge bunch with diffuse spatial boundaries, consider the case where is the number of particles in the charge bunch, and W ? and W z are defined by Eqs. ( 110)-(112).Substituting Eq. (118) into Eq.( 113) readily gives for the number density In Eq. (119), z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 , and r b t solves Eq. (101).Note from Eq. (119) that during the dynamic evolution of the charge bunch, the number density is a maximum at r; z 0; 0 and that n b is constant on spheroidal surfaces with r 2 =r 2 b t z 2 =z 2 b t const.Furthermore, substituting Eq. (118) into Eqs.( 116) and (117), we find that the temperature profiles are isothermal (independent of x) but depend on time t according to Other choices of distribution function are of course possible.As a second example, where the charge bunch has a sharp radial boundary, we consider the case where f b W ?; W z is specified by where the dependence on W ? in Eq. ( 122) corresponds to the Kapchinskij-Vladimirskij distribution [32].Substituting Eq. (122) into Eq.( 113), we obtain for the number density profile For the choice of distribution function in Eq. ( 122), note from Eq. ( 123) that the number density profile has a sharp radial boundary at r r b t.Indeed, for each z, the density profile in Eq. ( 12) is uniform (independent of r) out to radius r r b t, and equal to zero for r > r b t.By contrast, for the choice of distribution function in Eq. ( 118), the number density profile has a Gaussian radial profile proportional to expÿr 2 =r 2 b t.However, because the dependence of f b W ?; W z on W z is the same in both Eqs. ( 118) and ( 122), the corresponding line density profiles b z; t calculated from Eqs. ( 119) and ( 123) are identical, with A detailed calculation of the dynamic compression of the charge bunch for the class of distributions f b W ?; W z of course requires that r b t be calculated self-consistently from Eq. ( 101), which requires a specification of the solenoidal field lattice occurring in the definition of L .This will be the subject of a future investigation in which Eq. ( 101) is integrated numerically and the lattice configuration is optimized to assure maximum transverse focusing of the charge bunch in combination with the longitudinal drift compression.For our purposes here, we consider the case where the charge bunch has passed through the final focus magnet at t 0 with _ r b0 dr b =dt t0 < 0, corresponding to an inward transverse focusing of the charge bunch.Setting L 0 (corresponding to B z 0 in Eq. (101) then gives r b r 2 b0 v 2 ?0 =r 3 b , which can be integrated to give which is similar in form to the solution to Eq. (40) given by z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 .For _ r b0 < 0, the solution to Eq. ( 125) decreases to a minimum value r 2 bf at time t t f? , where _ r b tt f? 0, and t f? and r 2 bf are defined by For distribution functions f b W ?; W z of practical interest, the number density n b x; t in Eq. (113) assumes its maximum value at the center of the charge bunch corre-sponding to r; z 0; 0. Denoting n b r 0; z 0; t nb t, and defining the density compression ratio at time t by c d t nb t= nb t 0, it follows from Eq. ( 113) that As r b t and z b t decrease, it is clear from Eq. ( 128) that the density compression ratio c d t continues to increase.Moreover, from z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 , the axial compression decreases to z 2 bf z 2 b0 v 2 T0 = _ z 2 b0 v 2 T0 at time t t f ÿz b0 _ z b0 = _ z 2 b0 v 2 T0 .To achieve the maximum number density compression c d max , which is required for the most effective interaction of the charge bunch with a solid target, _ r b0 , r b0 , _ z b0 , and z b0 should be chosen such that t f? t f , i.e., such that the maximum transverse and axial compressions occur simultaneously.In this case, the charge bunch number density compresses to a maximum compression ratio c d max r 2 b0 z b0 r 2 bf z bf (129) at time t t f? t f .

V. CONCLUSIONS
In this paper, we considered the axial drift compression and transverse focusing of an intense charge bunch made up of positively charged ions for the case where the space charge and current of the charge bunch are completely neutralized by a background plasma.The charge bunch propagates in the z direction along the axis of an applied solenoidal magnetic field B sol x, and we made use of the Vlasov equation (1) to describe the evolution of the beam distribution function f b x; v x ; y; v y ; z; v z ; t in the sixdimensional phase space x; v x ; y; v y ; z; v z , neglecting self-field effects E s 0 B s because of charge and current neutralization.The analysis was carried out in the beam frame in the paraxial (thin-beam) approximation, and the axial force on the particle motion was treated as negligibly small F z ' 0. The theoretical model and assumptions in the present analysis were described in Sec.II.In Sec.III, we first considered the case where the full Vlasov equation bility of the model, and the dynamic axial compression of the charge bunch was investigated in detail for the class of distribution functions F b W z that have a negative velocity tilt at t 0. Here, W z is the constant of the motion defined in Eq. (39).In Sec.IV, for general initial distribution function f b x; v x ; y; v y ; z; v z ; t 0 the full Vlasov equation (1) was solved exactly using the method of characteristics.We then specialized to the class of distribution functions f b W ? ; W z , where W z and W ? are the constants of the motion defined in Eqs. ( 39) and (100), and the timedependent axial and radial scale factors, z b t and r b t, solve Eqs. ( 40) and (101), respectively.As demonstrated in Sec.IV, the dynamic properties of the axial compression and transverse focusing of the charge bunch can be calculated for a wide variety of distribution functions f b W ? ; W z .What's most remarkable about the present analysis is the enormous flexibility there is in the choices of distribution function f b W ? ; W z that undergo axial drift compression and transverse focusing, and correspondingly have different profiles for number density, line density, longitudinal and transverse pressures, etc.This important feature, of course, has been facilitated by the assumption of complete neutralization of the charge and current of the charge bunch by the background plasma.The detailed optimization of the solenoidal field configuration, velocity tilts, etc., that provide simultaneous axial drift compression and transverse focusing to a small, localized spot size is being investigated numerically, and will be the subject of a future publication.
Finally, it is important to reiterate the limitations of the present model, which has assumed complete charge and current neutralization, with zero self-generated electric and magnetic fields, E s 0 B s .As the beam-plasma system evolves dynamically, however, one would expect perturbations in charge density, P jb;e;i e j R d 3 pf j , and current density, J P jb;e;i e j R d 3 pvf j , to develop, where P jb;e;i denotes species summation over the beam ions j b and background plasma electrons and ions j e; i.These perturbations would of course produce perturbed self fields, E s and B s , which in turn would cause changes in the distribution functions f j x; v; t for j b; e; i.Therefore, the analysis in the present paper should be viewed as a first step in a more complete Vlasov-Maxwell description of the evolution of the beam-plasma system, in which we have determined the class of dynamically evolving beam equilibria in the absence of perturbations in charge density and current density.

2 W
? const: (100) vf b in the nonrotating beam frame is readily shown to beV b x; t ÿ e z ; (104)where L is defined in Eq. (102),z 2 b t z b0 _ z b0 t 2 v 2 T0 t 2 ,and r b t solves Eq. (101).Equation (104) readily follows from Eq. (103) because, W ? and W z are even functions of v x L y ÿ x=r b dr b =dt, v y ÿ L x ÿ y=r b dr b =dt, and v z ÿ z ÿ1 ) is diagonal, with P P b P ?b x; te x e x e y e y P zb x; te z e z n b x; tT ?b x; te x e x e y e y T zb x; te z e z : RONALD C. DAVIDSON AND HONG QIN Phys.Rev. ST Accel.Beams 8, 064201 (2005) 064201-12

( 1 )
was integrated over the transverse phase space x; v x ; y; v y , and the evolution of the longitudinal distribution function F b z; v z; t R dxdydv x dv y f b x; v x ; y; v y ; z; v z ; t was determined.A wide variety of distribution functions was considered that demonstrate the flexi-KINETIC DESCRIPTION OF NEUTRALIZED DRIFT . . .Phys.Rev. ST Accel.Beams 8, 064201 (2005) 064201-13 future reference, it is useful to introduce the volume number density n b x; y; z; t and axial line density b z; t (5); v x ; y; v y ; z; v z ; t:(4)Denoting the total number of particles in a charge bunch by N b Z dxdydzdv x dv y dv z f b x; v x ; y; v y ; z; v z ; t;(5) b e V zb e z ;Here, use has been made of e r cose x sine y , and e cose y ÿ sine x , where e r and e are unit vector in the r and directions.Equation (107) corresponds to a rotation of the charge bunch at the Larmor frequency L defined in Eq. (102), as expected, whereas Eq. (106) shows that the radial flow velocity is proportional to r=r b dr b =dt, which corresponds to transverse focusing of the charge bunch whenever dr b =dt < 0.