ARTICLE IN PRESS Nuclear Instruments and Methods in Physics Research A 606 (2009) 708–712
Contents lists available at ScienceDirect
Nuclear Instruments and Methods in Physics Research A journal homepage: www.elsevier.com/locate/nima
BDSIM: A particle tracking code for accelerator beam-line simulations including particle–matter interactions I. Agapov a,, G.A. Blair b, S. Malton b, L. Deacon b a b
CERN, CH-1211 Geneva, Switzerland John Adams Institute at Royal Holloway University of London, Egham, Surrey TW20 0EX, UK
a r t i c l e in f o
a b s t r a c t
Article history: Received 6 February 2009 Received in revised form 24 April 2009 Accepted 24 April 2009 Available online 8 May 2009
BDSIM, a code for particle–matter interaction studies in accelerators, is presented. BDSIM is a toolkit for simulation of particle transport in accelerator beam-lines. It provides a collection of classes representing typical accelerator components, a collection of physics processes for fast tracking, procedures of runtime geometry construction, analysis routines and tools to interface it to other particle tracking codes. & 2009 Elsevier B.V. All rights reserved.
Keywords: Particle tracking code Accelerator
1. Introduction During the design and operation phases of a particle accelerator extensive simulations are required. In the first approximation one calculates linear beam envelope functions (or the socalled Twiss parameters) and other optics parameters to match the beam requirements. This will give a rough machine layout and focusing properties. Single-particle tracking can be used to further optimise the beam-line parameters. General-purpose beam optics codes to perform such calculations include MADX [1], BMAD [2], and DIMAD [3]. Since linear (and largely nonlinear) beam transport theory has been well established the focus of these codes is largely on effective beam parameter matching to facilitate beam-line design. A particle accelerator can be several kilometres long and it is typically required to calculate a particle trajectory through it to a very high precision. Calculating particle trajectories by straightforward Runge–Kutta methods for such extended systems is not computationally effective. Instead, particle beam dynamics are usually treated in a perturbative fashion so that in the first approximation the particles perform linear oscillations about a design orbit. Also, for each magnet a linear or nonlinear map can be computed which corresponds to transporting a particle from its entrance to its exit. Particle tracking then reduces to computing superpositions of these maps, so that fast calculations can be performed to any desired order. Other codes are used to simulate dynamic effects such as ground motion and beam jitter; one such code is PLACET [4]. There are a number of codes
Corresponding author.
E-mail address:
[email protected] (I. Agapov). 0168-9002/$ - see front matter & 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.nima.2009.04.040
that treat self-consistent electromagnetic interactions inside particle beams (collective effects) [5]. A number of codes have been developed to simulate interaction of particle beams with materials, for example, Geant4 [6], FLUKA [7], MARS [8] and MNCP [9]. They are typically used for the design of particle detectors or radiation shielding applications. In an advanced particle transport code such as Geant4 one can in principle describe and simulate single-particle beam dynamics in electromagnetic fields. However, inefficient numerical methods and inflexible input (requiring the geometry description to be provided as a piece of C þ þ code) can make it impractical for application to a large-scale system. For computations requiring both tracking particles through long beam-lines and simulating particle interactions with the accelerator components, one often uses several codes. A fast tracking code is used to record positions where particles hit a collimator or the beam-pipe, and these positions are then input into a detailed Monte-Carlo radiation transport code. This approach, however, has several drawbacks: the geometry description capabilities of such codes are rarely equivalent (beam optics codes normally have very limited geometry descriptions), there is a need to transfer the data between two codes and so on. BDSIM [10,11] was developed to overcome these bottlenecks. It is based on Geant4, a Monte-Carlo framework, thus giving access to many electromagnetic and hadronic interaction models as well as a powerful geometry description framework. On top of this, fast particle tracking routines and some additional physics processes are introduced, and a high level geometry description language GMAD [12] is added. Since GMAD is an extension of MAD, a standard for beam optics description, this allows complex accelerator descriptions to be loaded from existing repositories with just a few modifications.
ARTICLE IN PRESS I. Agapov et al. / Nuclear Instruments and Methods in Physics Research A 606 (2009) 708–712
BDSIM has been applied to the simulations of a future linear electron–positron collider [13–17] and to the forward detector systems at the Large Hadron Collider [18]. It has also been used in the detailed simulation of a laser-wire beam diagnostics system [19,20]. BDSIM has been benchmarked against DIMAD and PLACET for particle tracking in vacuum [15,21] giving a good agreement for optical functions and tracking over several kilometres. Because the tracking is not yet symplectic, the code cannot be used for long-term tracking in a large circular accelerator. For electromagnetic and hadronic physics a set of cross-checks was performed, although here we rely on the Geant4 benchmarking effort. In the following, some details of the code design and implementation are given. Possible extensions and further applications are also discussed.
709
Compton scattering and photo-electric effect; synchrotron radiation; muon production and transport; neutron and proton elastic and inelastic scattering; neutron capture; fission; radioactive decay.
Certain methods are provided to increase the speed of the calculations:
biasing schemes are introduced for processes with small crosssections such as muon pair production;
controlling production cuts for various reactions in different regions of the geometry is possible;
‘‘fast’’ physics processes are sometimes provided by BDSIM to optionally replace the default Geant4 routines.
2. Architecture The architecture of BDSIM is sketched in Fig. 1. Geant4 provides general run management framework, a set of physics processes and a geometry construction scheme. BDSIM provides the GMAD parser and additional geometry drivers (currently a Mokka [22] database driver with SQL syntax, with LCDD [23] and GDML, a new Geant4 geometry description standard, drivers in development); these are used to build internal Geant4 geometry representations. This is done with a set of C þ þ classes that implement accelerator components: quadrupole and dipole magnets, collimators, beam pipes, and so on. Each element has a ‘‘stepper’’ associated with it which implements the particle transportation inside this element. C þ þ classes that implement various physics processes additional to Geant4 are also provided.
3. Physics BDSIM can exploit all physics processes that come with Geant4. Depending on the specific problem being addressed, different sets of physics processes can be turned on. These processes are grouped into the so-called ‘‘physics lists’’. The following processes are generally available:
multiple scattering, ionisation, bremsstrahlung, and positron annihilation;
gamma conversion;
Fig. 1. BDSIM design chart.
The accelerator geometry is composed of ‘‘components’’. These could be standard entities such as quadrupoles and bending magnets, or user-defined, such as detector components and collimators with non-standard apertures. Every element is made of several ‘‘logical volumes’’. For instance, a quadrupole can have inner vacuum, beam pipe, and iron yoke volumes. Each logical volume is associated with a certain material, a set of physics processes, particle production cuts, computation methods and several other parameters. The method by which particle transportation is calculated depends on the logical volume it is in and on the physics processes that are invoked. When no additional processes are present, a particle is transported in steps that end at boundaries between logical volumes. One step can be produced through a series of smaller sub-steps to ensure the required accuracy. The sub-step size can be controlled by setting a global Geant4 parameter (deltaChord, as listed later in the example in Section 4) that specifies the upper bound on the estimated error due to track curvature. In specific cases, however, it can be more efficient to avoid global changes in accuracy and, instead, to increase local accuracy by dividing up a single longer element into a set of shorter ones. A further Geant4 tracking parameter (deltaIntersection in the example in Section 4) defines the accuracy of intersection between volumes. If an element is a quadrupole, drift, or a bending magnet, a step for a primary particle with energy close to nominal inside the beam pipe is a mapping given by an appropriate transfer matrix (as described, e.g., in Ref. [5]). For sextupoles and octupoles a step is accomplished in a series of helical steps to obtain the required precision. For multipoles of higher order as well as for fields specified by value maps and for highly off-momentum or off-axis particles a Runge–Kutta method of fourth order is used. This procedure makes tracking of the primary beam faster, which is essential for accelerator applications. Detailed tracking comparisons have been performed [24] where it is shown that, for the CLIC beam delivery system, a final vertical spot-size of 0.716 nm was achieved with BDSIM, agreeing to within 0.002 nm with that of conventional machine tracking codes. Excellent agreement was also obtained when synchrotron radiation and dispersive effects were included. Additional processes such as synchrotron radiation may be present. Within Geant4 the processes are categorised as either ‘‘continuous’’ or ‘‘discrete’’. For a continuous process (which typically models a continuous energy loss such as multiple scattering) the step length is determined in the same way as before, but at the end of the step the particle momentum is altered. For a discrete process (e.g., electron–positron pair
ARTICLE IN PRESS 710
I. Agapov et al. / Nuclear Instruments and Methods in Physics Research A 606 (2009) 708–712
production) the step length is determined additionally by a random Monte-Carlo trial which computes the mean free path. At the end of the step the momentum change is computed and secondary particles produced. The processes present in beam–matter interactions (bremsstrahlung, etc.) are invoked only if the logical volume has an associated material other than vacuum. A comparison with other full-simulation codes [17] has shown good agreement; the fact that BDSIM runs on top of Geant4 means that the code benefits from continued improvements and benchmarking of the Geant4 physics lists. 4. Geometry description, interfacing to other codes There are two common ways to introduce geometry description to particle tracking codes. One is to provide a set of libraries in a high level programming language such as C þ þ, python, or similar, and let the user describe the geometry in this language. Another approach is to provide a separate input language. Programs such as MADX [1] adopt the latter approach. The MAD description language dates back to the 1980 and was written in FORTRAN; it is now more or less an industrial standard for accelerator description, at least in Europe. There have been recent developments to create a standardised XML-based description language called AML [25]. For large accelerators, accurate geometry description becomes a complicated problem. For an accelerator such as the LHC currently being commissioned at CERN, thousands of elements are installed in the tunnel and each must have its parameters—position, misalignment, aperture, magnetic field, etc—accurately documented and maintained. This information is typically stored across several databases; if MADX is used for optics modelling then a corresponding file repository is also provided. We developed an extension of the MAD language which allows us to introduce both geometry and material information, as well as parameters relevant to a Monte-Carlo radiation transport simulation such as energy cut-offs for secondary particle production. An example of an input file for simulating a beam-line containing a laser-wire detector is given below: !some beam-line parameters option; beampipeRadius ¼ 3:8 cm, boxSize = 2.0*m, tunnelRadius = 2.0*m, beampipeThickness = 0.1*cm, ! Transportation parameters chordStepMinimum = 0.0000000001*m, deltaIntersection = 0.00000001*m, deltaChord = 0.001*m, lengthSafety = 0.00001*m, thresholdCutCharged = 1*MeV, thresholdCutPhotons = 1*MeV; ! beam-line definiton laser1 : laser, l=10*cm, x=1, y=0, z=0, ! set direction of laser to wrt beam axis waveLength=532e-9*m, ! wavelength of laser 532 nm spotsize=1e-6*m, intensity=10e6, position=f 0,0,0g ; mk1: marker; cryoModule: element, l=12.4784*m, geometry=‘‘mokka:cryo.sql’’;
Fig. 2. Visualisation of a sample beam-line in BDSIM.
Predefined materials can be used: vcol: rcol, l=0.3, xsize=1.0, ysize=vcolgap, material=‘‘Tungsten’’;
or material definitions can be created on the fly: CarbonDioxide : matdef, density=1e-14, components=[‘‘C’’,‘‘O’’], componentsWeightsf 1,2g
BDSIM is a single-particle tracking code, where each particle is tracked from start to end before the next begins tracking. This makes calculating collective effects such as wake-fields impossible, as this requires a description of the entire bunch at the same time. Equally, it is necessary to track the bunch core simultaneously with the halo, as the forces acting upon the halo are dependent on the halo particles’ positions relative to the bunch core. To overcome this, BDSIM has been interfaced to PLACET [4]. The beam halo is tracked up to the first collimator—or any other position in the beam-line where the wake-fields first become significant—using BDSIM. Simultaneously, the bunch core is tracked to the same location in PLACET. BDSIM then passes the halo description through a FIFO to PLACET. PLACET tracks both the halo and the bunch core through the collimator, calculating the wake-field-induced effect of the bunch core on the halo. It is assumed that this effect is well approximated by a discrete particle momentum change or ‘‘kick’’. The halo description after the kicks have been applied is then passed back through the FIFO to BDSIM, which applies the kicks to the bunch distribution and then continues tracking either to the next collimator, where this process is repeated, or to the end of the beam-line if no more collimators are found. This procedure allowed the beam halo propagation to be simulated through a linear collider beam delivery system taking wake-fields and energy deposition into account [21].
5. Visualisation, data analysis, deployment on the grid BDSIM uses the standard visualisation frameworks provided with Geant4, (an example using OpenGL is shown in Fig. 2). The output of a simulation consists of energy deposition histograms and of particle distributions sampled at specified locations in the beam-line. This information can be provided either as an ASCII file or in the ROOT [26] format. A set of ROOT scripts for result analysis is provided. An example is given in the following section.
rb1:sbend, l=0.5*m, angle = 2e-6, outR=1*m; modb1: line = (cryoModule,rb1); test2 : line = (laser1,modb1,mk1); use,period=test2; sample, range=mk1; option, ngenerate=1000, physicsList=‘‘lw’’,randomSeed=-1; beam, particle=‘‘e-’’, energy=250, sigmaE=1e-4, nparticles=1e+3, distrType=‘‘gauss’’;
6. Example BDSIM beam-line construction and particle data analysis Here we demonstrate a simple beam-line constructed in BDSIM. We construct a beam-line made from five cells—an example of one cell is shown in Fig. 2. Each cell consists of one focusing and one defocusing quadrupole and two dipole magnets,
ARTICLE IN PRESS I. Agapov et al. / Nuclear Instruments and Methods in Physics Research A 606 (2009) 708–712
7. Planned future extensions of BDSIM Several directions in which the code can be extended are now mentioned. When LHC experiments at CERN become operational and achieve full luminosity more attention will be drawn to problems of beam-induced background. Such studies will involve simulating particles hitting the interaction region, scattering of the particles by the collimator jaws and simulating the processes which lead to beam halo formation. Many of these effects can be readily studied with BDSIM. The study of halo formation should
0
20
40
60
80
100
0
20
40
60 z (m)
80
100
0.15 0.1 y (m)
0.05 0 -0.05 -0.1 -0.15
Fig. 3. Vertical coordinates for an ensemble of particles with respect to nominal trajectory for a sample beam-line.
0
20
40
60
80
0
20
40
60 z (m)
80
100
Energy deposited (GeV)
3 2.5 2 1.5 1 0.5 0 100
Fig. 4. Position of energy deposits due to the interaction of the beam particles with the beam-line elements.
0.003 0.002 0.001 y’ (rad)
separated by short lengths of beam-pipe. This structure is then repeated, with sampler volumes between each cell to record particle data. The initial particle distribution is given by a Gaussian beam with sx ¼ sy ¼ 10 mm and sx0 ¼ sy0 ¼ 1 mrad; the geometry includes a beam-pipe with radius 14.5 cm and thickness 1 mm. The values listed are the widths of the distributions for transverse position (x, y), and for the angle of the particle horizontal and vertical trajectories with respect to the beam axis (x0 , y0 ). The initial longitudinal position (z) and particle energy (E) are fixed values of 0 mm and 1 GeV, respectively. Fig. 3 shows the vertical displacement of each particle relative to the nominal beam axis. The position is recorded at the end of each transportation step; in this case, since we have no physics processes in the vacuum, steps are between the boundaries of the logical volumes. Individual trajectories can be discerned from the data, and the alternating focusing and defocusing can be clearly seen. As described above, where tracked particles interact with the beam-line secondary particles are produced and energy deposits are recorded. An example of the energy deposition histogram is shown in Fig. 4. Further discrimination is also possible, with the histograms of energy deposited by each type of particle. A display utility is included in the BDSIM package which visualises the beam-line on each plot. This can be seen in Figs. 3 and 4. In this display, blocks to either side of the central line represent focusing or defocusing quadrupoles, while blocks which cross the line are bending magnets. Other types of element will also be displayed depending on the composition of the beam-line. Particle distributions can be sampled at marker locations along the beam-line. Fig. 5 shows a typical depiction of the particles’ vertical positions (y) and angles (y0 ). Primary electrons—input particles to the simulation—are shown in red, while secondary photons created by the interaction of the primary particles with the beam-line are drawn in green.
711
0 -0.001 -0.002 -0.003 -0.004 -400
×103 -300
-200
-100 0 y (µm)
100
200
300
Fig. 5. Distribution y and y0 with respect to the nominal beam trajectory at the end of the beam-line. The photons (green) have been created by the interaction of particles with the beam-line elements. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
include effects like beam diffusion due to RF noise or nonlinear imperfections, effects which are beyond the capabilities of BDSIM; however, combining BDSIM simulations with such studies could lead to a comprehensive background simulation. Unlike in a linear accelerator, simulating the energy deposition due to losses in a proton storage ring requires the tracking of particles over a large number of turns. To achieve the necessary accuracy would require symplectic integration schemes [27] which are not yet present. For electron accelerators in the energy range of several GeV such as those used as synchrotron light sources, beam dynamics is dominated by collective particle interactions. In Ref. [21], an example of exploiting the particle–matter interaction capabilities of BDSIM in conjunction with the collective effects calculated with another code (PLACET) is given; additional such interfaces are possible. Intensity-induced effects play a very important role in many accelerators. Enabling self-consistent bunch dynamics would be a major step forward and further extend the range of applications of BDSIM. 8. Conclusion BDSIM has been successfully applied to a number of problems, such as the collimation system design for a future linear collider and laser diagnostics system studies. Applications and corre-
ARTICLE IN PRESS 712
I. Agapov et al. / Nuclear Instruments and Methods in Physics Research A 606 (2009) 708–712
sponding code modifications relevant for the LHC at CERN and its injector chain are under way. Several improvements which would extend the application range for both high energy accelerators such as the LHC and lower energy ones such as synchrotron light sources are planned.
Acknowledgements We acknowledge the input and applications of R. Appleby, J. Carter, and O. Dadoun. References [1] H. Grote, F. Schmidt, MAD-X: an upgrade from MAD8, in: Proceedings of Particle Accelerator Conference (PAC 03), Portland, Oregon, 12–16 May, 2003. [2] D. Sagan, Nucl. Instr. and Meth. A 558 (2006). [3] R.V. Servranckx, TRIUMF design note TRI-DN-93-K233, July 1993. [4] A. Latina, et al., Recent improvements in the tracking code PLACET, CERN-AB2008-017. [5] A. Chao, M. Tigner, Handbook of Accelerator Physics and Engineering, Word Scientific, Singapore, 1999. [6] J. Allison, et al., IEEE Trans. Nucl. Sci. NS-53 (2006). [7] G. Battistoni, et al., AIP Conf. Proc. 896 (2007) 31. [8] N. Mokhov, Status of MARS code, in: FERMILAB-CONF-03-053, Presented at the 6th Meeting on Shielding Aspects of Accelerators, Targets and Irradiation Facilities (SATIF-6), Menlo Park, CA, 10–12 April, 2002. [9] J. Briesmeister (Ed.), MCNP: A General Monte Carlo N-Particle Transport Code, LA-12625-M, 1993. [10] I. Agapov, et al., The BDSIM toolkit, EUROTEV-REPORT-2006-014. [11] I. Agapov, et al., BDSIM: beamline simulation toolkit based on GEANT4, EUROTEV-REPORT-2006-035. [12] I. Agapov, GMAD accelerator description language, EUROTeV-Memo-2006-002-1.
[13] G.A. Blair, Simulation of the CLIC beam delivery system using BDSIM, CERNOPEN-2002-057. [14] G.A. Blair, Simulation of the beam delivery system, interaction region and backgrounds using BDSIM, in: Proceedings of 2005 International Linear Collider Workshop (LCWS 2005), Stanford, CA, 18–22 March, 2005. [15] R. Appleby, et al., Particle tracking in the ILC extraction lines with DIMAD and BDSIM, LAL-RT-06-05. [16] R. Appleby, et al., The 2-mrad crossing angle interaction region and extraction line, in: Proceedings of European Particle Accelerator Conference (EPAC 06), Edinburgh, Scotland, 26–30 June, 2006. [17] J. Carter, et al., Simulation of the ILC collimation system using BDSIM, MARS15 and STRUCT, SLAC-PUB-11951. [18] F. Roncarolo, et al., Machine induced backgrounds for FP420, in: EPAC08TUPP062, Proceedings of 11th European Particle Accelerator Conference (EPAC 08), Magazzini del Cotone, Genoa, Italy, 23–27 June, 2008. [19] G.A. Blair, Simulation of laser-wires at CLIC using BDSIM, in: 26th Advanced ICFA Beam Dynamics Workshop on Nanometer Size Colliding Beams (Nanobeam 2002), Lausanne, Switzerland, 2–6 September, 2002. [20] L. Deacon, G.A. Blair, Laser-wire simulation in the ILC beam delivery system, EUROTeV-Report-2008-018. [21] G. Blair, et al., Simulation of beam halo in CLIC collimation system, in: Proceedings of 11th European Particle Accelerator Conference (EPAC 08), Magazzini del Cotone, Genoa, Italy, 23–27 June, 2008. [22] P. Mora de Freitas, H. Videau, Detector simulation with MOKKA/GEANT4: present and future, in: International Workshop on Linear Colliders (LCWS 2002), Jeju Island, Korea, 26–30 August, 2002. [23] J. McCormick, Full detector simulation using SLIC and LCDD, in: 2005 International Linear Collider Workshop (LCWS 2005), Stanford, CA, 18–22 March, 2005. [24] S. Redaelli, et al., Comparison of different tracking codes for beam delivery systems of linear colliders, CERN-SL-2002-041-AP. [25] D. Sagan, The accelerator markup language and the universal accelerator parser, in: SLAC-PUB-11918, European Particle Accelerator Conference (EPAC 06), Edinburgh, Scotland, 26–30 June, 2006. [26] R. Brun, F. Rademakers, Nucl. Instr. and Meth. A 389 (1997) 81. [27] E. Forest, Beam Dynamics: A New Attitude and Framework, Harwood, Sidney, 1998 (p. 463).