Beam loss studies in high-intensity heavy-ion linacs

090101-1 The proposed Rare Isotope Accelerator (RIA) Facility, an innovative exotic-beam facility for the production of high-quality beams of short-lived isotopes, consists of a fully superconducting 1.4 GV driver linac and a 140 MV postaccelerator. To produce sufficient intensities of secondary beams the driver linac will provide 400 kW primary beams of any ion from hydrogen to uranium. Because of the high intensity of the primary beams the beam losses must be minimized to avoid radioactivation of the accelerator equipment. To keep the power deposited by the particles lost on the accelerator structures below 1 W=m, the relative beam losses per unit length should be less than 10 5, especially along the high-energy section of the linac. A new beam dynamics simulation code TRACK has been developed and used for beam loss studies in the RIA driver linac. In the TRACK code, ions are tracked through the three-dimensional electromagnetic fields of every element of the linac starting from the electron cyclotron resonance (ECR) ion source to the production target. The simulation starts with a multicomponent dc ion beam extracted from the ECR. The space charge forces are included in the simulations. They are especially important in the front end of the driver linac. Beam losses are studied by tracking a large number of particles (up to 106) through the whole linac considering all sources of error such us element misalignments, rf field errors, and stripper thickness fluctuations. For each configuration of the linac, multiple sets of error values have been randomly generated and used in the calculations. The results are then combined to calculate important beam parameters, estimate beam losses, and characterize the corresponding linac configuration. To track a large number of particles for a comprehensive number of error sets (up to 500), the code TRACK was parallelized and run on the Jazz computer cluster at ANL.


I. INTRODUCTION
The Rare Isotope Accelerator (RIA) driver linac will deliver 400 kW beams to the production targets [1].To avoid extra radioactivation of the accelerator equipment due to beam losses, the sources of beam halo formation must be carefully studied and measures to avoid uncontrolled beam losses should be taken.The driver linac is based on about 400 individually phased superconducting (SC) resonators [2].As was shown in previous studies [3], the most challenging dynamics in the driver linac are those of a uranium beam due to the simultaneous acceleration of multiple-charge states (multi-q) and the double stripping.In accelerators one can separate beam losses into controlled and uncontrolled losses.The controlled losses are the particles collected by intercepting any beam halo with specially designed collimators located at designated areas along the accelerator.All other losses anywhere along the linac are considered as uncontrolled losses.In this paper we do not study any losses related to the faults or failures of the accelerator systems.The main source of beam losses is the formation of beam halo either in the longitudinal or transverse phase planes.In the simulation, a particle is considered lost if it is outside the physical transverse aperture of the linac.In order to study beam losses at the level of 10 ÿ5 per meter we have used the specially developed beam dynamics code TRACK.The following reasons motivated the development of the new simulation code: (i) There was no available code to track multicomponent heavy-ion beams in threedimensional external and space charge fields.(ii) No existing simulation codes included the passage of heavy ions through strippers.(iii) There was no code capable of performing end-to-end simulations of the driver linac.
In addition, for reliable beam loss studies the code must run 10 6 particles and many randomly seeded accelerators with misalignments and errors.Therefore the code TRACK has been parallelized to run on a multiprocessor environment.Extensive use of the code allowed us to study many problems related to the design of highpower accelerators.
In the next section the newly developed beam dynamics code TRACK is described with an emphasis on the new features of multiple-charge-state acceleration and beam stripping.In Sec.III, the results of end-to-end beam dynamics simulations of two options for the RIA driver linac are presented.The sensitivity of each design option to different sources of error is also discussed.In Sec.IV, beam losses are studied for both design options and different values of the most critical errors.Conclusions and ideas for future studies are given in the last section.

II. SIMULATION CODE A. Prerequisites for the development of a new code
Detailed design and beam dynamics simulations in the RIA accelerators required the development of a new computer code to track multicomponent heavy-ion beams directly from the electron cyclotron resonance (ECR) ion source to the production targets.The code had to include a wide range of electromagnetic elements such as high-voltage accelerating tubes, rf resonators, electrostatic and magnetostatic multipoles with full threedimensional fields, etc.In addition, the interaction of heavy-ion beams with material media such as strippers has to be included.For careful localization of beam losses the code had to calculate particle trajectories in the close vicinity of the apertures of electromagnetic devices.A medium-energy SC linac contains a large number of accelerating, focusing, and bending elements.Hence the beam dynamics become quite sensitive to complex sets of errors and misalignments of these elements.The new code should be used for the statistical analysis of beam parameters obtained by tracking large numbers of particles through hundreds of accelerators each with randomly seeded errors.We have started developing the code TRACK [4] due to the lack of any code capable of simulating multiple-charge-state heavy-ion beam dynamics through such a wide range of 3D elements.
There are a number of available codes such as TRANSPORT [5], TRACE3D [6], COSY [7], and GIOS [8] based on matrix formalism for the design and study of beam-optics systems.Some of them include detailed beam space charge calculations [9].However, none of these codes provide a full description of beam dynamics in all the elements of the RIA driver linac.A specially developed raytracing code such as TRACK is necessary for detailed studies of the accelerator and beam transport systems and for finalizing the design.There are several obvious unique properties of tracking codes such as (i) The external field can be represented accurately within the physical aperture of the device while in a matrix code the external field is represented as a Taylor's expansion.Frequently, the convergence of the Taylor's expansion is not obvious and the field approximation is not valid in the vicinity of the aperture.
(ii) Particle coordinates are known at any point in the space.
(iii) Nonelectromagnetic elements, such as beam degraders or strippers, can be accurately simulated.
(iv) Beam space charge fields, especially for multicomponent ion beams, can be calculated with high accuracy.
(v) Beam losses can be calculated in the presence of complex sets of field errors and device misalignments.

B. Brief description of the code TRACK
In the code TRACK, the transport of a charged particle is described by the equation of motion, d p=dt q Ẽ ṽ B; (1) where p is the particle momentum and q is its charge, Ẽ Ẽext Ẽint and B Bext Bint are the sums of external and internal electric and magnetic fields, and ṽ is the particle velocity.The structure of the TRACK code and its ideology is close to the RAYTRACE code [10] except that the TRACK code has many additional features and capabilities.Unlike RAYTRACE, TRACK integrates the equations of motion of all tracked particles for a short distance and calculates space charge fields.In the TRACK code particle motion is generally described in three coordinate systems: two rectangular Cartesian coordinate systems related to both the entrance and the exit of each ion-optics device and a third Cartesian coordinate system used for the definition of the electromagnetic field distribution of the device.The particle trajectory is calculated in the input coordinate system for all types of devices.
Only the transformation from the field-distribution coordinate system to the input coordinate system of the device is needed for the description of the external fields in the equation of motion.The transformation from the input to the output coordinate system is later performed to prepare for tracking through the next device.
The fourth-order Runge-Kutta method is used for the integration of particle trajectory through an ion-optics device.TRACK uses the independent variable z for the tracking of phase-space coordinates of the particles (x; x 0 dx=dz, y; y 0 dy=dz, v=c; ), where v j ṽj, is the particle phase with respect to the rf field at the given section of the accelerator.For static ion-optics devices represents the time difference between the particle being tracked and the reference particle (RP) which provides the beam pulse transformation through the device.In the different accelerator sections represents the particle phase with respect to the current frequency of the rf resonators.For direct current (dc) beams extracted from the ion source the coordinate is uniformly distributed within at the frequency of the downstream rf resonator f 0 .
We acknowledge that the time-dependent integration of the equation of motion including space charge fields is the most natural and physically correct method.The development of the TRACK code was mainly motivated by the dynamics of multicomponent ion beams in the linac with negligible space charge effects but noticeable nonlinear components of the external fields in the presence of strippers, field errors, and device misalignments.Therefore the use of the Cartesian z coordinate as an independent variable was justified.Most features of the code were already developed before the implementation of the space charge routine.For the applications being considered in this paper, as will be seen from the examples given in the following sections, the TRACK code produces results which are consistent with a time-dependent integration code.Currently the TRACK code is being modified to support time-dependent integration for the simulation of space charge dominated beams.
The following procedure is applied for particle tracking through an ion-optics device.The same integration routine is applied either to the RP or any other particle.First, the RP is tracked from the beginning to the end of the device.The phase and velocity of the RP and the length of the central trajectory are calculated.The central trajectory is equidistantly divided into M 1 points P m (m 0; M).In each point P m the rectangular coordinate system is formed with the unit vectors ñx m ṽm =v, ñy m ñy 0, ñz m ñx m ñy m, where ṽm is the RP velocity at the point P m .The coordinate system at the point P 0 , f ñx 0; ñy 0; ñz 0g, coincides with the input coordinate system while the coordinate system at the point P M , f ñx M; ñy M; ñz Mg, corresponds to the output coordinate system of the ion-optics device.For ionoptics devices with straight central trajectory the same integration step is applied for all particles to move from one plane to the next.For beam-optics devices with curvilinear central trajectory different particles must move with different step size to reach the next plane which is achieved by the help of a simple iterative procedure.
The calculations start with the beam spatially placed on the plane f ñx 0; ñy 0g of the input coordinate system.The trajectories of all beam particles are integrated stepby-step from the plane f ñx m; ñy mg to the plane f ñx m 1; ñy m 1g using the z coordinate of the input coordinate system as an independent variable.The tracking continues until all particles reach the plane f ñx M; ñy Mg of the output coordinate system.
As in the code RAYTRACE [10], when sufficiently small step sizes are used, the accuracy is limited only by the uncertainties in our knowledge of the magnetic and electric fields.The calculation of the space charge fields is performed in the curvilinear coordinate system of the RP.Appropriate transformation of space charge forces from the curvilinear coordinate system to the input coordinate system is required.The set of equations used for the stepby-step integration routine is where 1= 1 ÿ 2 p , h 1= 1 x 02 y 02 p , jej=m a c, A is the mass number, m a is the atomic mass unit, Q q=jej, and B x , B y , B z , E x , E x , E x are the components of the magnetic and electric fields.The trajectory equations (2) are directly derived from the Lagrangian of a charged particle in any time-dependent electromagnetic field (see, for example, Ref. [11]) and do not include any simplifications.
Depending on the geometry and the device type, external fields in the code can be defined in one of the following formats.
(i) Three-dimensional tables of Ẽ and B in the input rectangular coordinate system obtained by the help of external codes.For the calculation of the field value at the particle location a quadratic interpolation routine is used.
(ii) Two-dimensional tables in the plane fr; zg for the elements with axial symmetry such as solenoids or Einzel lenses.
(iii) Two-dimensional tables of the B y component of the magnetic field in the median plane fx; zg for rectangular dipole magnets.Off-median component B y and components B x and B z are evaluated using the method described in [10].
(iv) The fringe field falloff for dipole and multipole elements is described by a six-parameter Enge function [12], Fz 1 1 expa 0 a 1 z=D 1 a 5 z=D 5  ; where z is the distance along the line which is perpendicular to the effective field boundary and D is the full air gap of the element.Space charge fields of multicomponent ion beams are defined as the result of solving the corresponding Poisson equation.The charge distribution in the mesh points of the rectangular mesh is defined using the ''clouds-incell'' method.The code calculates both two-dimensional (for dc beams) and three-dimensional (for bunched beams) space charge fields.The three-dimensional Poisson equation is solved with rectangular boundary conditions in the transverse direction and periodic conditions along the direction of beam propagation.The detailed description of the space charge routine can be found in Ref. [13].A special routine has been written for the integration of multicomponent dc ion beams through bending magnets.
PRST-AB 7 P. N. OSTROUMOV, V. N. ASEEV, AND B. MUSTAPHA 090101 (2004) 090101-3 090101-3 Three-dimensional interpolation of the electromagnetic fields is performed within the full aperture of ionoptics devices.For most common elements the fields can be described by several methods which makes the code TRACK compatible with many existing beam transport codes.This feature is valuable for the extraction of transfer matrices and the comparison of the simulation results with existing codes.
The code TRACK supports the following electromagnetic elements for acceleration, transport and focusing of multicomponent ion beams: (i) any type of rf accelerating resonator with realistic 3D fields; (ii) static ion-optics devices with both electric and magnetic realistic threedimensional fields; (iii) radio frequency quadrupoles (RFQ); (iv) solenoids with fringing fields; (v) bending magnets with fringing fields; (vi) electrostatic and magnetic multipoles (quadrupoles, sextupoles, etc.) with fringing fields; (vii) multiharmonic bunchers; (viii)axial-symmetric electrostatic lenses; (ix) entrance and exit of high-voltage decks; (x) accelerating tubes with dc distributed voltage; (xii) transverse beam steering elements; (xiii) stripping foils or films; (xiv) horizontal and vertical jaw slits for beam collimation.

C. Stripper simulation
When accelerating heavy ions, stripping the beam increases the charge state and makes the acceleration more effective to reach higher beam energies.Scattering off the atoms of a stripper foil, a given ion not only loses electrons but also loses energy and deviates from its original direction.Therefore, the stripping process changes the properties of the beam depending on its energy and the stripper thickness.Strippers are usually thin and chosen to keep acceptable beam properties and minimal losses by inelastic interactions.Knowing the distributions of the beam energy loss and angle scattering after a stripper is important for beam dynamics simulations.For calculations where ions are tracked individually through the entire accelerator it is also important to be able to calculate or generate the energy loss and scattering angle ion by ion.
Full Monte Carlo calculations of ion's energy loss and angular scattering by codes such as SRIM [14] are useful but usually slow and not easy to couple to a beam dynamics code like TRACK.The calculation time is proportional to the stripper thickness because an incident ion is tracked through the different atomic layers of the stripper foil.In the case of the second stripper of the proposed RIA driver linac, the SRIM calculation time is comparable to tracking the same number of ions through the entire linac which would double the beam dynamics calculation time.
As an alternative, tables or formulas of energy loss, energy straggling, and angular straggling could be used for fast generation of the energy loss and scattering-angle distributions for most beam-stripper combinations.Such formulas reproduce reasonably well the peak of the distribution but they usually ignore the tails which could be responsible for beam halo formation and eventual beam losses after the stripper.Another problem of these formulas is that the energy loss and angle distributions were parametrized independently and lack any energy-angle correlation which could be critical in the acceleration process especially when accelerating multiple charge states simultaneously.
For the reasons mentioned above we opt for the parametrization of the correlated energy loss and angular distributions calculated with the code SRIM for a given beam-stripper combination.Once the parameters are found they are used to regenerate the calculated distributions using the Monte Carlo method.In this way the energy-angle generator based on the parametrization is invoked every time an ion passes through the stripper.Note that this parametrization is local and has to be performed for every beam-stripper combination to be used.Two different parametrizations were performed and used for the two strippers of the RIA driver.It is also important to note that after a stripper foil the beam is a mixture of several charge states and that this parametrization is charge-state independent.
To describe the parametrization procedure we consider the case of the second stripper in the RIA driver linac for which we assumed a 15 mg=cm 2 carbon foil.The beam is an 85 MeV=u 238 U beam.SRIM calculations were performed for 10 6 incident ions.The corresponding results are shown in Fig. 1.Observing Fig. 1(a), we notice that the energy distribution has a peak at about 81:7 MeV=u (energy loss of 3:3 MeV=u) with a tail extending toward lower energies (higher energy loss).The scattering-angle distribution of Fig. 1(b) peaks at about 1 mrad with a tail extending to more than 20 mrad exceeding the acceptance of the high-energy section which can be limited to 7 mrad [15].Therefore some of the scattered ions will be stopped at the collimators placed behind the stripper.Figure 1(c) shows a strong correlation between the energy loss and the scattering angle.This correlation becomes linear in the energy-angle squared plane; see Fig. 1(d).This correlation is consistent with the kinematics of single elastic scattering.
Studying the energy and angle distributions and their correlations, we noticed that for a selected small interval in angle squared 2 the corresponding energy distribution E has a Gaussian-like shape; see Fig. 2(c).The 2 distribution itself is a portion of a decaying exponential; see Fig. 2(b).We noticed also that for different intervals in 2 , the slope of the exponential as well as the center and width of the corresponding energy Gaussian are different.The linear correlation between E and 2 indicates that the center of the energy Gaussian should be a linear function of 2 .Figure 2 illustrates the method to use for the parametrization of ion energy and angle distributions after a stripper foil.First, the center and width of energy Gaussians are parametrized as functions of the angle squared 2 .Second, the 2 distribution is subdivided into an appropriate number of regions represented by different slopes of the exponential and different event fractions or probabilities.For the normalization of the exponential we used the event fraction of the corresponding region.The normalized integral over all regions should give 1.In this way the statistical importance of each region is conserved.More details on the parameterization procedure could be found in [16].
The code SRIM assumes that the thickness of the stripper is well defined and it does not also consider statistical variations of the ion's charge state as it traverses the foil [17].The combination of these two effects leads to an increased energy spread of the ion beam downstream of the stripper foil.To mock up these effects, simulations were performed for various foil thicknesses.SRIM calculations performed for 5% and 10% deviations from the nominal thickness of the second stripper showed that a 5% thickness fluctuation would increase the beam energy dispersion by a factor of 3 and a 10% fluctuation would increase it by a factor of 5 relative to the single thickness case (0% fluctuation).
In order to consider the effect of thickness fluctuation and/or charge-state variations in our parametrization we have repeated the procedure described above for the cases of 5% and 10% thickness deviations.The five points corresponding to the 0%, 5%, and 10% deviations from the original thickness have been used to study and characterize the dependence of the parameters on the stripper thickness.
The energy-angle squared correlation seems to be identical for the five thicknesses.More detailed analysis [16] showed that the slopes of the energy loss and straggling as function of 2 are almost the same for the 0%, 5%, and 10% thickness deviations.The average values are used for these slopes.The values at 0 of both the energy loss and straggling show a linear dependence on the thickness as is expected from formulas of energy loss.Concerning the parameters of the angle distribution we have studied the exponential slope as well as the event fraction for each of the regions in 2 as functions of the thickness.They all seem to fit reasonably well with linear functions with different slopes.For the exponential slope, the first regions showed a linear dependence on the thickness which is washed out for regions of larger angles due to the lower statistics.For these regions the slopes are considered to be constant with thickness and the average value is used.At this level the parametrization of the energy loss and scattering-angle distributions is complete including straggling due to fluctuations in stripper thickness and/or variations of the ion's charge state within the foil.In order to use this parametrization for the generation of the energy loss and the scattering angle ion by ion we will need to know the thickness distribution of the stripper foil.Experimental measurements [18] suggest that the thickness distribution of a stripper foil could be represented by a Gaussian distribution.The center of which is the average thickness and its full width at half maximum (FWHM) measures the thickness fluctuation.Since the minimal and maximal possible thicknesses are finite the Gaussian should be truncated at a certain limit.We choose to truncate the Gaussian at FWHM so that the full width is 2 FWHM.
For an ion incident on the stripper foil, first a thickness is generated from the appropriate Gaussian thickness distribution.Using the generated thickness the values at 0 of the energy loss and energy straggling are calculated from their parametrization as function of thickness.The exponential slopes and the event fractions for all regions in 2 are also determined for the current thickness.At this stage most parameters are determined and ready to use for the generation of the energy loss and the scattering angle for the incident ion.The generator was first tested with 0% thickness fluctuation in order to check if it reproduces SRIM results for the original thickness.The comparison is shown in Fig. 3.We notice that the generator reproduces reasonably well the results of SRIM.In order to check how well the thickness fluctuation is represented in the parametrization we have tested it against experimental data.Recent measurements at the AGS-BNL where a carbon foil is used to strip a 100 MeV=u Au beam from charge state 31 to 77 showed that even within the beam size these fluctuations could reach 10% FWHM [18].The measurement was performed using a 4 MeV proton beam on a 0.005 in.carbon foil.Figure 4 shows the comparison of the energy spectrum of the protons between the data, SRIM calculation (0% thickness fluctuation), and the parametrization of SRIM results considering a 10% FWHM thickness fluctuation.We clearly notice that the experimental data are very well reproduced by the parametrization with a 10% FWHM thickness fluctuation.For high-intensity heavy-ion beam the stripping foil experiences high thermal load, especially the first stripper of the RIA driver linac.A liquid lithium film is under development [1] for the first stripper to accommodate the thermal load and to produce higher charge states.For the second stripper and due to the higher beam energy inelastic interactions become more probable which could result in the contamination of the beam by radioactive products of nuclear reactions.For uranium at 85 MeV=u, 0.25% of the beam interacts inelastically in the 15 mg=cm 2 carbon foil.A preliminary study using intranuclear cascade and fission-evaporation codes showed that about 5%-10% of the products could be accepted in the high-energy section of the linac producing a possible beam contamination of about 10 ÿ4 .A more detailed study of the eventual contamination by radioactive products and their dynamics in the accelerator is a special topic and will be presented in a future publication.

D. General features of the code
Depending on the task, the code TRACK can be compiled for the simulation of up to  contains about 400 rf resonators operating at frequencies from 57.5 to 805 MHz, 90 solenoids, 16 bending magnets, and about 100 quadrupoles.The simulation of 10 4 particles on a 3 GHz PC requires about 27 min.To study beam halo formation and eventual beam losses the simulation of large number of particles through multiple accelerators with randomly seeded errors is required.Considering the number of error sets to use for each design (up to 500) and that a new simulation is needed for every variation in the accelerator design, a fast multiprocessor machine is absolutely needed.For these reasons, the code was parallelized [19] using the message passing interface (MPI) and run successfully on the new Jazz cluster at the Laboratory Computing Resource Center of Argonne National Laboratory.Jazz is a new teraflop machine (10 12 floating point operations per second) with 350 ÿ2:4 GHz processors.
Prior to the simulation with the code TRACK three main steps have to be performed: (1) The elements of the accelerator system must be defined using first and, if necessary, higher order optimization codes.
(2) Preparation of field tables.This step requires the use of external codes such as Microwave Studio (MWS) and Electromagnetic Studio (EMS) [20].A special preprocessor code has been written for the transformation of the calculated fields to the rectangular mesh with uniform step which is required for the code TRACK.
(3) At the first passage through the accelerator the code TRACK sets phases in the SC resonators according to the user defined values.The optimization of phases using TRACK or other codes is required prior to phase setting.The phase setting should maximize the longitudinal acceptance of the linac sections and minimize the effective emittance of multi-q beams at the location of the strippers.
One of the useful features of the code TRACK is the ability to calculate the transverse and the longitudinal acceptance along the accelerator system by simulating large phase-space volumes and tracking survived particles.The phase-space image of the survived particles at the entrance of a given section of the accelerator defines the corresponding acceptance.An extensive number of beam parameters is calculated after each element and written into an output file.The list of parameters includes beam centroids, rms values and maximum beam sizes in x; y and the phase, and three types of normalized emittances in the three phase-space planes.These are the rms emittances and the emittances containing p% and 100% of particles, where p is a percentage to specify as input, for example p 99:5.The emittance is calculated as the area of the corresponding ellipse divided by .The ellipse orientation is defined from the second moments of the particle distribution (see, for example, [21]).The Twiss parameters and the number of lost and survived particles are also determined.The output data is usually treated using a scientific software for statistical analysis.

E. Benchmarking of the TRACK code
For all the elements supported by the code TRACK the first-and second-order transformation matrices have been calculated numerically and compared with the ones used by the codes TRANSPORT, COSY, TRACE3D, and GIOS.Complete agreement of the matrix elements has been found.The front-end section consisting of the low-energy beam transport (LEBT), RFQ, and medium-energy beam transport (MEBT) has been compared with the DYNAMION [22] and PARMTEQ codes [23].The different linac sections have been compared with LANA simulations [24].An effort is currently underway to modify the codes PARMTEQ and IMPACT [9] for the simulation of heavy-ion beams in RIA-type accelerators.Initial comparisons of energy gain and beam second moments in the high-energy section of the RIA driver linac show very good agreement [25].

III. BEAM DYNAMICS IN THE RIA DRIVER LINAC
A detailed configuration of the 1.4 GV RIA driver linac was described in Ref. [3].The linac consist of a front end and three sections of SC linac: low-, medium-, and high-sections.The front end includes an ECR ion source, a LEBT system, a multiharmonic buncher (MHB), a RFQ, and a MEBT system.The three sections of the linac are separated by two stripper areas each with a stripper foil or film and a poststripper magnetic transport system (MTS).Beam dynamics in the driver linac are the most challenging for the multiple-charge-state uranium beam [3].The baseline design of the driver linac has been optimized for simultaneous acceleration of two charge states (28 and 29 ) in the front end and the first section of the linac up to the first stripper, five charge states between the two strippers (average charge state is 74 ) and five charge states in the high-section (average charge state is 88 ).The acceleration of multi-q beams not only increases the total intensity but also reduces significantly the power to dump at the strippers.For example, the five charge states after the second stripper represent 98% of the total intensity.End-to-end simulation of the linac is an essential part of the beam dynamics design and precedes the beam loss studies.The simulation starts with the multicomponent ion beam at the exit of the ECR ion source.Beam parameters from the ion source were adopted from the recent experimental data and extrapolation to higher intensities based on estimates by the VENUS Ion Source Group [26].Because of the cw regime of the driver linac, the space charge effects are negligible in all accelerator sections except the ECR source and following extraction optics and beamline.The required beam intensity in the LEBT is 500 A for protons and 250 A for uranium to produce 400 kW accelerated beams in cw mode.After separation and selection of ion species the beam optics becomes emittance dominated and the space charge effects produce small perturbations with respect to the ''zero-current'' beam optics.

A. Beam parameters at the exit of the front end
The front end includes different types of ion-optics devices including the MHB and the RFQ.Detailed design, optimization, and simulation of the front end is extremely important to produce realistic beam distribution in the six-dimensional phase space at the entrance of the SC linac.Several publications have been devoted to this problem [27][28][29].The simulation by the code TRACK includes transport in the LEBT of multicomponent ion beams from the ECR to the MHB, bunching and acceleration in the RFQ, and transport in the MEBT to the entrance of the SC linac.After the final optimization of the front-end system, the simulations including space charge of multicomponent ion beams have been carried out with different numbers of particles from 2 10 3 to 10 6 .Figure 5 shows the rms envelopes of a dual chargestate uranium beam along the MHB-RFQ-MEBT section obtained from the simulation of 10 6 particles.Each charge state carries 4.4 particle A (p A) current.The RFQ transports 100% of injected particles including 17.7% nonaccelerated particles.The latter can be intercepted by several collimators between the quadrupole lenses.In Fig. 5, the beam size downstream of the RFQ is shown only for accelerated particles.The longitudinal phase-space plots shown in Fig. 6 correspond to the entrance of the SC linac.Figure 7 shows the fraction of particles 1 ÿ N=N 0 outside of a given longitudinal emittance.The simulation of larger numbers of particles reveals an increased beam halo.The MHB forms an extremely low longitudinal emittance of 1:6 keV=u ns containing 99% of particles.As was shown in Ref. [28], it is not possible to avoid the formation of halo in the longitudinal phase space due to the transverselongitudinal coupling in the RFQ.As seen from Fig. 7, the total emittance for the accelerated particles (8:23 10 5 ) could reach 12 keV=u ns.
The beam dynamics including space charge effects in the LEBT and the RFQ has been simulated using both the t-dependent code DYNAMION [22] and the z-dependent code TRACK.The dual charge-state uranium beam distributions in the longitudinal phase space obtained for 2 10 5 particles by both DYNAMION and TRACK are given in Fig. 7(b).The total uranium beam current is 250 A. One can notice the good agreement between the DYNAMION and TRACK codes.Both DYNAMION and TRACK codes show a 100% transmission through the RFQ including a 17.7% fraction of nonaccelerated particles.

B. Two options of the linac
Massive parallel-computer simulations have been performed for two options of the SC driver linac in order to investigate possible beam losses and determine the exact location of these eventual losses.The first option is the baseline design of the driver linac which was described in Refs.[2,3].Since these publications the baseline design  (ii) The accelerating lattice and the phase setting of the resonators were optimized to minimize the effective emittance of the multi-q beam at the location of the strippers.
(iii) The accelerating and focusing lattice of the 805 MHz section were modified to increase the longitudinal acceptance.For this purpose the length of the cryostat was minimized [30] and four G 0:49 elliptical cavities are placed per cryostat.
(iv) The focusing field of SC solenoids has been restricted to 9 T to reduce the cost of the linac.
The second option of the driver linac is based on triplespoke resonators in the high-section of the linac [31].In what follows we refer to the first option as the ellipticalcell linac (ECL) and the second option as the triple-spoke linac (TSL).A focusing system with room-temperature quadrupole doublets located outside of the cryostats is used in the ECL design and could be used for the TSL design.We have chosen, however, to use 50 cm long, 9 T SC solenoids in the TSL design.The solenoids provide strong transverse focusing for all ion beams and can be incorporated into the cryostats, economizing space along the beamline.As was mentioned in Ref. [31] the obvious advantage of the TSL option is a significantly larger longitudinal acceptance compared to the ECL option.

C. Beam collimation
A strong correlation between the ion's energy and scattering angle is clearly observed after a stripper foil; see Fig. 1(c).Ions with lower energies are scattered to larger angles.This correlation suggests a simple way to remove the low-energy halo by a system of collimators appropriately located along the poststripper MTS system as de-scribed in Ref. [15].The main collimator is located in a highly dispersive area to dump all unwanted charge states.Five other collimators are designed to clean the beam halo in the transverse phase planes.Both the ECL and TSL are designed to accept five charge states of uranium beam in the high-section, therefore only 2% of the initial intensity, which is about 2 kW, has to be dumped after the second stripper.The transverse acceptance of the MTS calculated at the location of the second stripper with collimator openings of 7:5 mm is shown in Fig. 8.The acceptance of the MTS with collimators is about 10 times smaller than the acceptance of the subsequent high-section.In the linac sections there are no uncontrollable mechanisms for beam halo formation.The main source of halo is the strippers.Therefore appropriate beam collimation in the poststripper MTS is the key solution to avoid or minimize beam losses associated with the beam dynamics.

D. End-to-end simulation of the linac without errors
The particle distribution at the front-end exit has been used as an initial distribution for the simulation of the whole SC linac.All presimulation procedures mentioned in Sec.II C have been applied to both the ECL and TSL.Particularly, the setting of the accelerator which includes (i) beam matching between the different focusing periods; (ii) providing minimal beam size at the stripper location; (iii) transformation of the multi-q beams in the six-dimensional phase-space by the MTS after the stripper; (iv) setting of the reference phases of the resonators to provide minimal effective emittance of the multi-q beam at the location of the strippers; (iv) adjustments of the collimator openings in the MTS.
Figures 9 and 10 show the evolution of beam envelopes and emittances along the baseline linac obtained from the simulation of 82.3% of the 10 6 particles injected into the front-end system.In this calculation, the stripper thickness fluctuation was set at 5% FWHM.After the first and second strippers 0.3% and 0.2% of particles are inter-

E. Simulations considering different sources of error
We can classify the possible sources of error into three groups: (i) misalignment errors: affect all the elements of the accelerator system: accelerating cavities, quadrupoles, solenoids, etc.; (ii) rotation errors: affect mainly quadrupoles, multipoles and bending magnets; (iii) rf field er- rors: affect the field level as well as the phase of an accelerating cavity.
For heavy ions that need to be stripped, another important source of error is to be considered, it is the fluctuation in the thickness of the stripper foil or film.The errors are two types: static and dynamic.Misalignments of accelerator elements are considered as static errors.A jitter of rf and focusing fields is an example of dynamic errors.The phase and amplitude setting of the accelerating cavities when first tuning the accelerator or restoring a tune from the computer memory is also a source of static errors.In a real machine, the effect of static errors can be partially corrected using beam measurements.We are currently developing a modification to the parallelized version of the TRACK code that can provide automatic steering of the beam position along the linac using beam position monitors [32].As was found from preliminary studies, the transverse steering can be effectively provided.In the current version of the code a simple correction procedure is used.The beam centroid position after every other cryostat is shifted back to the center of the phase planes by ''kicks'' with residual errors of 0:1 mm for the position and 0:2 mrad for the angle.Similar procedures can be developed for the correction of static errors of rf fields and phase settings.In the current simulations we do not distinguish between static and dynamic errors of rf fields.The purpose of this paper is to study the effect of different types of errors on beam losses.
Table I lists the errors to be considered as well as their typical amplitudes.The current errors are to be generated randomly according to the corresponding distribution.The Gaussian distributions are truncated at 3 where is the rms value except for the stripper thickness fluctuation which is truncated at FWHM.The uniform distributions are generated between the extreme values max.The displacement errors are applied to the x and y positions of element ends.The rotation errors are applied around the z axis (beam axis).Figure 11 shows the effect of these errors on the transverse and longitudinal emittances of the beam at the accelerator exit for both the ECL and TSL designs.The emittance growth "=" " ÿ " 0 =" 0 relative to the case with no errors (" 0 ) is plotted as function of error numbers from Table I.For each case 50 sets of errors are generated and applied to the corresponding elements of the accelerator system.For each set 2 10 5 particles are tracked making a total of  particles.The emittance " is 4 times the rms emittance averaged over all error sets using the following expression: where N is the number of error sets, here N 50.We clearly notice that the transverse emittances are mainly affected by solenoid displacements (error 2).The longitudinal emittance is mainly affected by the rf field and phase errors (errors 5 and 6) and by the fluctuation in the stripper thickness (error 7).
Comparing the sensitivity of both the ECL and TSL designs of the driver linac to the different errors, no significant difference is seen for the transverse emittances.For the longitudinal emittance, the ECL design shows a much larger sensitivity to rf errors (error 5: field; error 6: phase) compared to the TSL design.This applies also to the stripper thickness fluctuation (error 7).This makes the baseline design less tolerable to rf errors and energy spread induced by the strippers.Consequently the amplitudes of these errors should be kept as small as possible to avoid or limit beam losses.
Figure 12 presents particle coordinates in the three planes of the phase space for both designs and selected errors from Table I.The acceleration of multiple charge states causes oscillations in the effective emittance as the beam propagates along the linac which result in an effective emittance growth in the longitudinal phase space.This emittance growth can be minimized by selecting larger phase angles for the reference (synchronous) particle.Because of the large acceptance of the TSL section there is no need to minimize the longitudinal emittance and the synchronous phase can be set at ÿ25 .Therefore the longitudinal phase-space plots in Fig. 12 show more nonelliptical structure for the TSL than for the ECL which has ÿ30 for synchronous phase and higher frequency of small longitudinal oscillations.

IV. BEAM LOSS ANALYSIS
In the previous section we studied the effect of each source of error separately and we showed that the most critical sources of error are the rf errors (field and phase) and the fluctuation in the stripper thickness.In this section we apply all the errors simultaneously and we study the beam dynamics and the eventual beam losses.Different combinations of amplitudes for rf errors and thickness fluctuation are used as shown in Table II Both the baseline and triple-spoke designs of the accelerator were simulated using 200 sets of errors for each combination.Figure 13 shows particle coordinates at the exit of the accelerator for the different error combinations.From these plots we notice that, while the transverse size of the beam stays unchanged for the triplespoke design, it is increasing for the baseline design signaling a growth in the transverse emittance.It is seen more clearly in Fig. 14 where the relative transverse and longitudinal emittance growth "=" are plotted for  than increasing the thickness fluctuation from 5% FWHM to 10% FWHM (combinations 3 and 4) while keeping the rf errors at (0.5%, 0.5 ).This signals a more important beam halo formation for rf errors compared to stripper thickness fluctuation.Therefore rf errors should be kept as small as possible.The error on the thickness fluctuation should also be limited to about 5% FWHM or less in the baseline design.
To study the stability of theses losses, we have calculated the case of error combination 4 for the ECL design with up to 500 sets of errors (500 seeds of random numbers).Figure 15 shows the average lost fractions in MTS-1, MTS-2, and the high-section as functions of the number of seeds.For MTS-1 and MTS-2 we notice that the losses are quite stable.In the high-sections, we notice some periodic fluctuations which averages out at about 200 seeds which justifies our choice to use 200 error sets (200 seeds) for this study.We also have performed simulations with 5 times more particles (10 6 ) and 200 seeds and obtain the same beam losses.Reducing the number of accelerated charge states in the high-energy section from five to three produced the same beam losses in the ECL design.II and both the baseline and triple-spoke designs of the driver linac.The losses at the collimators following both strippers (MTS-1-2) are controlled losses.They are in the 2-4 10 ÿ3 range and vary slightly with the error combination.

Combination
Baseline Triple-spoke  For hands-on maintenance the limit on power lost on accelerator structures is about 1 W=m. Figure 16 shows the power lost per meter for both the baseline and triplespoke designs and the combination of errors in Table II.The horizontal line shows the 1 W=m limit.For the ECL design, uncontrolled losses in the high-energy section are either absent or small for combinations 1-2, approaching the 1 W=m limit for combinations 3-4, and exceeding 10 W=m for combinations 5-6, which have the largest errors.No uncontrolled losses are observed for the TSL design.The 1 W=m limit does not apply to losses at the strippers because they are controlled.
From this study it is obvious that the current baseline design has more limitations concerning beam losses, whereas the triple-spoke design shows more tolerance to errors.

V. CONCLUSION
A new tracking code TRACK has been developed to include up-to-date multicomponent heavy-ion beam physics in all types of focusing and accelerating elements being used in linear accelerators.For reliable beam loss studies in the presence of field errors and element misalignments the code has been parallelized to run up to 10 6 particles in hundreds of randomly seeded accelerators.
Detailed end-to-end beam dynamics simulations for multiple-charge-state uranium beam using the TRACK code have been performed in the intermediateenergy heavy-ion linac.These simulations have been iterated repeatedly with the design of the overall linac architecture.After the finalization of the linac design, massive multiprocessor simulations including misalignments of focusing and accelerating elements, random errors of the rf fields, and thickness fluctuation of the stripping films have been carried out.The concept of beam halo collimation in designated areas of the poststripper transport system has been developed and applied to minimize beam losses in the high-energy part of the accelerator.
Beam losses have been studied for two options of the RIA driver linac.The first option is based on the ECL and the second option is the TSL.The main sources of beam halo formation are rf field errors and the beam energy spread caused by the stripping films.Some particles become first unstable in the longitudinal phase space and eventually hit the aperture due to the loss of transverse stability.The studies show that the ECL design has more limitations concerning beam losses.The error in rf field amplitude and phase should be kept below 0.5% and 0.5 , respectively.In addition, the thickness fluctuation of the stripper should be limited to less than 5% FWHM.The losses in the ECL are extremely sensitive to the longitudinal tuning of the linac in terms of phase setting to produce the lowest possible effective emittance at the location of the strippers.The triple-spoke design is more tolerant of errors without producing any uncontrolled beam losses in wide range of rf errors, thickness fluctuations of the strippers, and overall longitudinal tuning of the linac.
Further beam loss studies will include tracking of radioactive ions after the second stripper and consider cavity-dependent field levels which is inherent to a linac comprising a large number of SC resonators.

FIG. 1 .
FIG. 1. Results of SRIM calculations for 10 6 238 U ions at 85 MeV=u incident on 15 mg=cm 2 carbon stripper foil.(a) Energy distribution, (b) distribution of scattering angles, (c) strong energy-angle correlation which becomes linear on the energy-angle squared plane (d).

FIG. 2 .
FIG.2.Basis of the parameterization.(a) Linear correlation between the energy and angle squared, (b) projection on angle squared of the selected small interval on (a) showing an exponential distribution, and (c) projection on energy of the same angle interval suggesting a Gaussian distribution for every small angle interval.

FIG. 4 .
FIG. 4. (Color)Effect of thickness fluctuation.Comparison of the experimental data with SRIM calculation (0% thickness fluctuation) and the parametrization using 10% FWHM thickness fluctuation.The data were measured using a 4 MeV proton beam on a 0.005 in.carbon foil[18].

FIG. 3 .
FIG. 3. (Color) Comparison of the energy and angle distributions and their correlation generated by the generator described in the text to the original ones calculated using SRIM.

FIG. 7 .
FIG. 7. (Color) Fraction of particles outside a longitudinal emittance given by the abscissa, (a) at the entrance of the SC linac obtained for different numbers of simulated particles and (b) at the exit of the RFQ obtained by the DYNAMION and TRACK codes.Space-charge effects are included in the simulations.

FIG. 11
FIG. 11. (Color) Comparison of the sensitivity of the baseline and triple-spoke designs to the different source of errors of TableI.The plots show the increase in transverse and longitudinal beam emittances relative to the case without errors.The blue curves are for the baseline design and the red ones are for the triple-spoke design.

FIG. 12
FIG.12.(Color) Phase-space plots for the two accelerator designs and selected errors from TableI.The top plots are for the case with no errors.The other three are for the most critical errors-5, 6, and 7, respectively.Errors 1-4 (not shown) show no significant difference from the case without errors.

FIG. 14
FIG. 14. (Color) Transverse and longitudinal emittance growth for the different error combinations of Table II relative to the case without errors.The blue curves are for the baseline design and the red ones are for the triple-spoke design.
FIG. 15. (Color) Evolution of beam losses as a function of the number of error sets or seeds of random numbers.Calculated for the ECL design and the error combination 4 of TableII.

FIG. 16
FIG. 16. (Color) Beam power lost along the accelerator in W=m for both the baseline and triple-spoke design and the error combinations of TableII.The plots are for combination 1 (top) to combination 6 (bottom), respectively.The horizontal line shows the 1 W=m limit to not exceed for hands-on maintenance.The first two peaks on each figure correspond to the losses at the two strippers, which are controlled losses.
).The peak surface electric field in all drift-tube SC resonators is assumed to be 20 MV=m except for the first seven 4-gap quarter wave resonators where 16 MV=m is used.

TABLE I .
Different sources of error and their typical values.For solenoid displacement, the error depends on the length of the solenoid.Smaller errors are applied to shorter solenoids.
. Other errors have the same values of Table I and are kept unchanged.

TABLE II .
Combinations of rf errors and stripper thickness fluctuations used to study beam dynamics and beam losses for both the baseline and triple-spoke designs of the driver linac.Other errors are kept unchanged at the values given in TableI.
FIG. 13. (Color) Phase-space plots for both accelerator designs and the error combinations of Table II.The plots are for combination 1 (top) to combination 6 (bottom), respectively.The logarithmic density isolines are represented by different colors.

TABLE III .
Beam losses for the error combinations of Table