Space charge compensation studies of hydrogen ion beams in a drift section

Demands of high intensity beams are increasing for material science and engineering applications. To control the high intensity beam dynamics, the nonlinear space charge forces, produced by the beam and the residual plasma at low energy, have to be carefully studied. Interactions between the beam and the residual gas tend to neutralize the beam space charge. The temporal and spatial evolutions of these interactions complicate the beam dynamics, even into a drift section. A 1.5D (xy for motion, r for space charge) PIC code, CARTAGO, has been written in order to simulate the beam dynamics in a low energy beam transport line including the space charge compensation effects. This paper relates the structure and the numerical methods used in the code and shows studies of the space charge compensation regime for H and H beams in a drift section.


I. INTRODUCTION
In the low energy part of an accelerator, high intensity beams undergo important space charge (SC) forces. Such beams are usually transported in a space charge compensation (SCC) regime induced by the secondary charged particles produced by ionization of the residual gas. This phenomenon occurs naturally and favors the reduction of the defocusing SC forces. Depending on the beam charge sign, the produced ions and electrons are repelled or trapped until a steady state is reached. If no clearing mechanism is present, the neutralizing secondary particles accumulate up until the SC becomes negligible. The beam SC is then assumed fully neutralized [1].
Analytical and experimental studies [2 -4] show that the beam charge in a low energy beam transport line is not always fully neutralized. Several studies have also recorded that the SCC depends on time [5,6]. These conditions complicate considerably the beam dynamics in a space charge compensation regime. In Ref. [7], a sophisticated numerical model is shown to investigate neutralization in case of the heavy ion beam transport to a reactor chamber. This development is dedicated to the study of transport line configurations which are specific for inertial fusion applications and do not apply directly to the low energy light ion beam transport line for accelerators dedicated to nuclear or particle physics (ADNPP). For instance, the neutralization in [7] is produced by an external source of gas (FLiBe) compared to ADNPP where the residual gas comes from the ion source. This induced different time scales which produced different effects on the beam. In this paper, we propose to use transport line parameters for ADNPP, and derive useful basics rules for ADNPP designers.
We propose to study such a SCC regime in the case of a simple drift. The study assumes that no longitudinal component of the electrical field and no external magnetic field (quadrupole or solenoid lenses) are present in order to simplify the modeling and to only illustrate the transverse behavior of the particles in such regime. We also neglect the secondary electrons produced when particles hit the wall, the validity of this hypothesis is discussed below.
First, Sec. II of this paper describes tools to define, quantify, and study the space charge compensation mechanism: a definition of the space charge compensation degree (SCCD) is given, a characteristic time constant for SCC is also defined. It will help to normalize the speed of this process during its transient.
This section also detailed the dedicated numerical code CARTAGO which has been written to study the 1.5D (xy for motion and r for the Poisson solver) behavior of the charged particles in a SCC regime.
In Sec. III, the SCCD and the transient evolution of the compensating plasma are studied for various cases of continuous (DC) rigid beams. The influence of the residual gas pressure on the time building of the SCC process is also studied for H and H ÿ beams.
In Sec. IV, the influence of the SCC on the beam dynamics in a drift section is shown. This study is performed for both positive and negative hydrogen ion beams. Comparisons with experimental observations are discussed.

A. The space charge compensation degree
The evolution of the beam rms sizer is given by the envelope equation [8]: where hAi is the average value of function A over the beam, " r is radial rms emittance, and F r is the total applied force acting on the beam. When the external focusing force is linear and the only other force is the space charge force F sc;r , this equation becomes where k r0 is the wave number induced by the external focusing. Whatever the particle distribution in the beam, the behavior of the rms size can be modeled with (2) with a very good approximation using the linear space charge force of a uniform beam with the same rms size. When the beam is neutralized, the space charge force is reduced and even canceled by the charges trapped in the beam potential well. Many definitions of space charge neutralization have been used to quantify this compensation degree. In Refs. [1,9], the ratio between the beam potential well depths with and without compensation is used. In Refs. [6,7], the SCCD is quantified with the ratio of the number of trapped particles over the number of beam particles. A drawback of these definitions is that they cannot be used directly in (2) to predict the beam transport in a SCC regime.
To reach such a goal, we propose to use the following parameter: hr F comp;r i hr F sc0;r i ; where hr F comp;r i is the linearized compensated space charge force (induced by all ionized particles) and hr F sc0;r i the linearized beam space charge force at the beginning of the process. With this definition, one has 0 ! no compensation 1 ! full compensation: The usual way of using the SCCD is simply to decrease the beam current (or charge) by the factor I ssc 1 ÿ I 0 : The envelope equation with space charge compensation is then d 2r ds 2 k 2 r0 r ÿ " 2 r r 3 ÿ 1 ÿ hr F sc0;r ĩ r 0: The SCCD of a cylindrical beam at a given time t is then defined as t R 1 0 r F comp;r r; t r; t dr R 1 0 r F sc0;r r; t r; t dr : The beam space charge force has one electric and one magnetic component. In the paraxial approximation, the radial space charge magnetic force F sc0;r;B partially compensates the radial space charge electric force F sc0;r;E [6]: If one considers that the field induced by the trapped particles is essentially electric, one can easily prove that where E comp;r is the electric field induced by ionized particles and E sc0;r the one induced by the beam. At low energy, the compensation is only due to the compensation of the electric force. As at high energy, the magnetic force reduces partially the electric force, a lower compensation electric field is needed to reach a complete neutralization.

B. The SCC time constant
To study the time evolution of the SCCD, it is convenient to define characteristic rise time, T SCC . This parameter can be used to normalize the time evolution in order to compare different cases even when, for instance, the pressure changes. This time constant, depending on gas nature and pressure, nature, and velocity of the beam, gives the average time needed for one particle to produce one compensating particle [6]: where i (m 2 ), n g (m ÿ3 ), v f (m s ÿ1 ) are, respectively, the ionization cross section, the gas density, and the beam velocity.

C. The CARTAGO code
To calculate the nonlinear transport in such a regime, PIC simulations are one possible candidate. However, if the modeling of the SCC is performed with the only help of a weighting of the fields computed with a Poisson solver which embeds only the beam, it is supposed that the compensation force is distributed in space like the selfspace-charge force of the beam. It means that each beam particle is at the same location of one neutralizing particle.
To study the behavior of the different particle distributions and to check this improbable statement, a dedicated PIC code, CARTAGO, has been developed. CARTAGO is a beam dynamics simulation code which includes the effects of the nonlinear space charge compensation. The time is used for the integration of the motion. This ensures a temporal coupling between the computed space charge and the macroparticle dynamics. The algorithm to simulate the beam and the plasma dynamics is composed of four basic parts (Fig. 1). Part 1.-Every time step, electrons and ions are created by the background gas ionization with the initial position given by colliding beam particle and the velocity given by the differential ionization cross section [10]. The plasma population increases with the ionization and decreases with particle losses. The ionization process is the only source for secondary particles, no wall emission is assumed in the simulation as the beam losses are assumed negligible, and the residual pressure is relatively high.
This hypothesis has been validated by comparing measured and simulated neutralization rise time. This macroparticle system is stored in the computer memory with the method of chained lists. This technique allows a convenient dynamical evolution of the number of particles. Having a rough idea of the plasma rise time (deduced from T SCC ), and choosing the number of macro-ions and electrons created during T SCC , the charge and the number of macroparticles generated every time step can be estimated. This method simplifies significantly the read access memory management during the computation.
Part 2.-The charge distributions, deduced from the positions of electrons, ions, and beam particles, are weighted in a 1D mesh (r). The Gauss theorem is used to solve the Poisson equation.
Part 3.-The SC forces are deduced from the potential and applied to the particles via the step by step ''leap frog'' scheme [11].

III. TRANSIENT SCC IN A DRIFT SECTION
The purpose of this section is to simulate the transient of the SCC process for DC parallel and rigid H and H ÿ beams which propagate in a medium with different residual gas pressures. The only difference between these two beams is the charge of the particle. A similar ionization cross section is considered for both beams. As the study is restricted to a transverse slide, the longitudinal velocity component of compensating particles is neglected. The beams (95 keV, 100 mA, 3 cm radius) are continuous, axisymmetric, and uniform. The only source of secondary charges taken into account is the ionization of the H 2 residual gas molecules at 300 K.

Influence of the gas pressure
The SCC mechanism of the H beam has been simulated for different gas pressure during 5 T SCC . According to Eq. (10), the neutralizing time T SCC goes like the inverse of the residual gas pressure. For instance, it is around 10 s for 4 10 ÿ5 hPa. Figure 2 illustrates the time behavior of the SCC for different residual gas pressures (4 10 ÿ6 hPa, 4 10 ÿ5 hPa, and 4 10 ÿ4 hPa). The pipe has been set to 3.3 times the beam radius.
At low pressure, the SCCD is about 90% after one T SCC and tends to 100% after 5 T SCC for all cases. The numerical results of the Fig. 2 are in a good agreement with our expectation: if the secondary interactions between the neutralizing electrons and the other particles are neglected, 100% space charge compensation must be observed in the stationary state [1]. Measurements of SCCD close to 100% have been already recorded for high intensity H beams at low energy [12]. However, for a good prediction of the steady state SCCD, the secondary interactions have to be taken into account.
A particular behavior of the SCC time evolution can be noticed at high pressure (4 10 ÿ4 hPa). This case shows 3 consecutive stages detailed in Fig. 3.
The first stage is visible only at high pressure. At the beginning of the process, the secondary ions are not immediately repelled (heavy particles created with a very low energy). They need around 100 ns to exit the beam. This mechanism slows down the initial SCCD evolution.
We propose to introduce the parameter to quantify this stage in the SCC process. It is defined as the ratio between the average time which is required for an ion to exit the beam exit and the compensation time with hdi the average distance between the ion creation position and the beam edge and hvi the average ion velocity to go over hdi.
To give an order of magnitude, the average distance can be calculated with where R b is the beam radius. For the average velocity hvi, we can use where V b is the potential difference between the beam center and its edge, m i and q i are the ion mass and charge, respectively. Equations (11)-(13) lead to In our simulation conditions, one finds that exit is equal to 100 ns. A detailed calculation presented in Appendix A shows that about 75% of the ions have exited a uniform beam after exit .
The parameter for each curve has been added in Fig. 2. It shows that stage 1 can be ignored if 10 ÿ3 (at low pressure).
During the second stage, the ions are lost at a constant rate when the electrons continue to be trapped. The total space charge force decreases linearly with time.
In the last stage, the gradual reduction of the SC potential well slows down the repulsion of the ions and favors loss of some electrons which have been created with an initial kinetic and/or potential energy. The SCCD evolution is therefore more and more slowed down.
The presence of the first stage has also consequences on the steady state particle distribution in the beam. The particle distributions and corresponding electric fields are plotted on Fig. 4 when SCCD reaches 90% for two studied gas pressures (Ne, Ni, and Nf are, respectively, electron, ion, and beam distribution).
When is much smaller than 10 ÿ3 , ions can exit the beam before a significant compensation occurs. The initial time evolution of the SCCD is linear (stage 1 is very short) and a very small fraction of the ions are accumulated in the beam when reaching 90% compensation. With these conditions, the effect of the ions on the compensation process can be neglected. At the opposite, for larger , the SC is rapidly compensated and the repulsive forces for the secondary ions (H 2 ) are significantly reduced. An important accumulation of H 2 in the beam occurs.

Nonlinearities
The nonlinear aspect of the SC in the transient regime can also be seen on Fig. 4. The beam core is not neutralized as the edge. This is induced by the nonuniformity of the compensating particle distribution. When the electrons are created with an initial kinetic energy, they oscillate around the central part of the beam on an elliptical trajectory with a minimum and a maximum radius (kinetic orbital moment conservation).
These elliptical trajectories reduce the electron density at the beam center (see Fig. 4). The impact of this nonlinear field on the beam is discussed in Sec. IV.

Influence of the beam and the pipe radii
In this section, b1 indicates the ratio between the pipe radius R p and the beam radius R b : For a uniform beam, the expression of the potential difference, V total , between the pipe edge and the beam center can written has a function of b1: where I b and v b are, respectively, the intensity and the velocity of the beam. Equation (16) shows that the higher the parameter b1 the higher the required energy for electrons to be lost on the pipe. At the opposite, when the b1 parameter is reduced, the time evolution of the SCC process may slow down with the electron losses because their initial energy can be significantly greater than the beam potential well. Figure 5 shows the SCC time evolution of the H beam for different b1 to illustrate this mechanism.

Influence of the beam current
The dynamics of the SCC is linked to the electron population evolution and especially to their initial velocity spectrum. The behavior of this electron population should depend on their average initial energy over the potential well ratio: More, the reduction of the beam intensity increases also the time exit needed for ions to exit the beam (at 4 10 ÿ5 Pa and with a beam current of 5 mA, exit is equal to 1:8 s). The time evolution of the SCC process is then slowed down when the parameter b2 is increased.
This phenomenon is also present when b1 is varying but it is less sensitive because the potential is proportional to the current and the log of b1.

B. SCC of H ÿ beam
Negative beams are neutralized by heavy positive ions compared to positive beams which are neutralized by light electrons. Moreover, the initial ion energy is negligible compared to the electron case. This induces significant changes for the SCC mechanism. This section details the SCC process of H ÿ beams for the same conditions of the H simulations. The SCCD evolution is plotted for different gas pressures on Fig. 7. The simulations of the SCCD show a full neutralization at the steady state. SCCD exceeding 100% are observed in all cases during the tran-sient. This SC overcompensation of H ÿ beams has been recorded in several experimental studies [9,13]. It is justified by the slow reaction of the neutralizing secondary ions when the full compensation is exceeded.
Stage 1 observed with H does not appear with H ÿ as the mobility and initial velocity of electrons are very high. The electrons can then exit the beam much faster than the compensation build up. The required pressure to record phase 1 for H ÿ beams is unused in an accelerator.
We will limit our SCC mechanism study at low gas pressure (4 10 ÿ5 hPa) for the following simulations. As ions are created with a negligible kinetic energy compared to their potential energy, they oscillate in the beam. At the opposite, electrons are quickly repelled to the vacuum chamber. Assuming a uniform beam distribution, the potential well is parabolic and the ion radial density nir can be easily calculated with Eq. (18) (see Appendix B): where R is the beam radius and ni 0 a normalization con-  stant. The ion radial distribution which neutralizes the H ÿ beams (plotted on Fig. 8) is very sharp on the axis because each ion trajectory, created with zero initial velocity, crosses the beam center. This behavior is then very different compared to the electron distribution for positively charged beams. The analysis of the SCC at 4 10 ÿ5 hPa (T SCC 10 s) shows that the initial ion radial distribution after 0:8 s (Fig. 9) is very close to the expected one (Fig. 8). A large concentration of secondary ions is observed on the beam axis. This leads to a very fast compensation of the beam center. A positive field appears when the secondary ion density exceeds the beam density. The local force repels ions and reduces their concentration near the beam axis. In addition, this force enables some electron trapping (see Fig. 9). The steady state starts when the positive and the total negative density are similar (Fig. 9). We can notice that the nonlinear SC field is recorded during the whole transient phase, especially close to the axis. At high gas pressure (4 10 ÿ4 hPa), Fig. 7 shows some fluctuations with time of the SCC evolution in the case of a uniform beam. This behavior is also observed when the time step is refined. During the SCC transient regime, the H 2 ions are trapped by the remaining negative field and repulsed by the induced positive field near the beam axis.
When the beam is uniform, the initial linear SC field brings back all the produced H 2 on the beam axis at the same time. Thereafter, they are simultaneously repulsed. This reaction is less significant with a low rate of H 2 production and especially with long time constants (case of low residual gas pressure). The described SCC fluctuations are reduced when nonuniform beam distributions are used [14].

Influence of the beam and the pipe radii
We showed in the case of the H beams that the variation of the b1 parameter [the ratio between the pipe radius and the beam radius, of the Eq. (15)] can increase the mechanism of losses of the energetic neutralizing particles (electrons). For the H ÿ beams, because the H 2 ions have a negligible initial energy compared to their potential energy, they remain in the space of the beam and are insensitive to b1 and especially when the beam radius changes. The simulated SCC time evolution, for several values of b1, is plotted in Fig. 10. These results confirm the insensitivity of the SCC evolution to the b1 parameter.

Influence of the beam current
Because the ion initial kinetic energy is negligible and the quickly repelled electrons have no influence on the SCC, it is straightforward to consider that the beam current  has almost no influence on the starting phase of the SCC. For H ÿ beams, overcompensation is recorded during the last phase of the transient when the H 2 ions begin to be repelled and start to escape the beam. We showed that the escape time, exit , of these particles is proportional to the inverse of the square root of the potential well depth [see Eq. (14)]. Thus, for beams with a low intensity, the overcompensation degree is enhanced. This is in very good agreement with the simulation results of Fig. 11.

IV. SCC IMPACT ON THE BEAM DYNAMICS IN A DRIFT
In this section, the dynamics of H and H ÿ beams in a 25 cm drift section is studied. The pressure of the H 2 residual gas is 4 10 ÿ4 hPa (T SCC 1 s). At the drift entrance, the beam (95 keV, 100 mA) is continuous, axisymmetric, and uniform. Stationary optical parameters are used. They have been selected to be close to those of the H beams at the exit of the electron cyclotron resonance source Source d'Ions Légers de Haute Intensité [15]. The horizontal rms normalized emittance and the Twiss parameters are respectively: " x rms 0:2 mm mrad, ÿ2 and 1:74 mm= mrad. The beam radius is 10 mm.
To study the impact of the SCC on the beam dynamics during the transient phase and the steady state, the integration of the motion of the beam itself is added. The time evolutions of the beam emittance and the Twiss parameters at the drift exit are given. Because of the computation method of SC solver (1D), the longitudinal effects are not taken into account.

A. Transport of H beam with SCC
When the SCC is just starting, the radial defocusing force of the uniform beam is still linear. The consequence is that the beam size grows in the drift whereas the emittance remains constant. Figure 12(a) shows the first outgoing beam at the drift exit (after 58 ns) in phase space.
The resulting x-x 0 distribution confirms the linear behavior due to the initial SC, whereas the beam radius reaches 18 mm.
During the transient regime and in agreement with the study performed in Sec. III, the radial force produced by the total SC is nonlinear. The beam emittance at the drift exit grows up rapidly and reaches, for instance,  The beam SC is almost compensated inside the drift after 3 T SCC (corresponding to 3 s). The emittance is conserved after this phase. The use of a nonuniform beam distribution, at the drift input, emphasizes this emittance growth. In the case of a uniform beam, the nonlinearity of the SC induces an emittance growth and changes the beam distribution inside the drift.
The beam Twiss parameters at the drift exit vary also considerably during the transient regime (Fig. 13). This complicates the beam matching of the front end for pulsed beams. This mismatch may lead to a beam emittance growth and a halo formation in the accelerator. The effect of the transient nonlinearities is increased when the beam propagates in a longer drift. A 0:35 mm mrad maximum rms normalized emittance has been reached at a 40 cm drift exit. The transient dynamics, described above, must not be neglected when the SCC time scale is not negligible compared to the beam pulsed duration.

B. Transport of H ÿ beam with SCC
The same simulations have been carried out for H ÿ beams. Except the beam charge sign, all parameters are similar compared to the preceding section. The H ÿ beam size also grows when it propagates through the drift whereas the emittance remains constant when the SCC is just starting. The nonlinear force produced by the SCC may be at the same time focusing (close to the axis) and defocusing (close to the edge) compared to the proton beam case. Figure 9 illustrates this phenomenon. To show the effects of this particular nonlinearity, Fig. 14(a) shows the beam in the phase space after 0:7 s. At the drift exit, the normalized beam rms emittance reached 0:28 mm mrad.
We saw in Fig. 7 that the H ÿ SC is almost completely compensated after 1 T SCC (here 2 s). The initial emittance is then preserved in the drift after this time. Because of Twiss parameter variations at the drift exit, a beam mismatch is also expected in the case of H ÿ beams. Figure 15 shows the detailed behaviors of the beam emittance growth and Twiss parameters normalized to their steady state values.
In agreement with the study of the Sec. III, the same fluctuations of the SCCD (Fig. 7) are observed on the emittance evolution. However, the influence of the SCC on the H ÿ beam dynamics is very close to the H beam case.
The main difference comes from the steady state which can be very different for positive or negative beams as the mobility and the initial kinetic energy of the compensating particles are significantly different.

Different sizes
Other studies of the H ÿ beam transport showed that the evolution of the emittance and the SCC process is very dependent on the initial beam parameters. This behavior is rather different when, for instance, the size of the H ÿ beam is changed. The SCCD of a rigid H ÿ beam reached 130% at the transient regime when the initial radius is 50 mm. Intense focusing and defocusing forces are then produced during this overcompensation and the resulting normalized emittance reaches 0:42 mm mrad at the exit of a 25 cm drift.

V. CONCLUSIONS
This paper detailed analytical and numerical studies of the space charge compensation in a drift section. They show that at low energy and when the ionization is the only effective process, the beam space charge is fully neutralized for continuous H or H ÿ beams at the steady state. The corresponding rise times for different gas pressures and for various parameters of hydrogen ion beams are given. The observed effects during the transient regime in the numerical simulations are in a good agreement with the analytical expectations.
The transport of H or H ÿ continuous beams into a drift section allowed us to estimate the emittance growths during the transient regime. These simulations confirm that a beam mismatching occurs in the front end during this phase. This may lead to halo formation in the high energy part of the accelerator. To refine these results at the steady state, a special effort is required to understand and to introduce other collision processes like heating.

ACKNOWLEDGMENTS
Thanks are due to all the members of the CEA SILHI group at Saclay for the experimental assistance and the helpful discussions. This work is supported by the European Community-Research Infrastructure Activity under the FP6 ''Structuring the European Research Area'' program (CARE, Contract No. RII3-CT-2003-506395).

APPENDIX A: TIME TO EXIT A UNIFORM BEAM FOR A GIVEN FRACTION OF IONS
In this Appendix, we look for the time for a given fraction of ionized ions to exit the uniform beam in which they have been created.